3
$\begingroup$

In a previous post, I was asking how to force a fine spatial meshing in NDSolve using the FiniteElement method. The solution was

NDSolveValue[{op == 0, ics}, u, {t, 0, 2}, {x} \[Element] reg, 
Method -> {"FiniteElement", "MeshOptions" -> {"MaxCellMeasure" -> 0.005}}]

giving the correct result below

Correct with the Method option specifying the MaxCellMeasure sub-sub-option

instead of the incorrect one obtained (without any error or warning message) when using simply

NDSolveValue[{op == 0, ics}, u, {t, 0, 2}, {x} \[Element] reg]

Incorrect with default NDSolve parameters

But this solution is not satisfactory because in the general case, (1) it is not so obvious to detect that a solution is False and (2) it is long and boring to run two subsequent NDSolve integrations with two different meshings in order to check whether it gives twice the same result, and (3) the user can be tempted to use a too fine meshing and consequently too much time in order to be sure to get a correct result. One could expect NDSolve to have an optional algorithm that automatically refine the spatial meshing if needed.

1) Does such an option exist?

2) Is it compatible with the FiniteElement method?

Thanks for your help.

Denis

Improve this question
$\endgroup$
6
  • $\begingroup$ “...it is long and boring to run...” is one of the reasons society has invested in the development of computers and programming. What I’m getting at is that one could write a program to carry out adaptive refinement. I think ToElementMesh[] probably has the tools/features needed. It would probably not be slow on 1D problems (ODEs) but higher dimensional PDEs might be more challenging. I’ve not done it nor do I know if there are error estimation schemes other than what you describe above. $\endgroup$
    – Michael E2
    Commented Feb 14, 2019 at 13:28
  • $\begingroup$ It is necessary to publish all the code, not two pictures. $\endgroup$
    – Alex Trounev
    Commented Feb 14, 2019 at 15:12
  • $\begingroup$ Here's an example for the method of lines, which is easier than what I had in mind for FEM since a message is emitted when further reduction of step size is needed. It uses recursion instead of a loop, but I assumed it would converge after just a few iterations. mathematica.stackexchange.com/a/73848/4999 $\endgroup$
    – Michael E2
    Commented Feb 15, 2019 at 3:58
  • $\begingroup$ Sorry Alex, the code is in the previous post linked at the top. I did not want to repeat it here because the question is different. $\endgroup$
    – user3650925
    Commented Feb 15, 2019 at 13:08
  • $\begingroup$ Thank you Michael. I'll see what I can do along this line, but note that in my example, no error or warning is generated by NDSolve when the integration happens to be completely wrong $\endgroup$
    – user3650925
    Commented Feb 15, 2019 at 13:24

2 Answers 2

Reset to default
6
$\begingroup$

No, adaptive mesh refinement is not implemented in version 11.3 for PDE solving with the FEM. NIntegrate with the FEM uses an adaptive mesh refinement.

Improve this answer
$\endgroup$
0
$\begingroup$

Here is a brute-force solution consisting in integrating the problem with stepsizes divided by 2, 4, 8, 16... until the absolute difference between two consecutive solutions goes below a tolerance threshold:

reg = Line[{{0}, {1}}];
shape = Cos[16 π (x - 0.5)] D[0.125 Erf[(x - 0.5)/0.125], x];
op = D[u[t, x], {x, 2}] - D[u[t, x], {t, 2}];
ics = {u[0, x] == shape, Derivative[1, 0][u][0, x] == 0};
sol[maxCellMeasure_] := NDSolveValue[{op == 0, ics}, u, {t, 0, 2}, {x} \[Element] reg, 
      Method -> {"FiniteElement",
                 "MeshOptions" -> {"MaxCellMeasure" -> maxCellMeasure_}}];

decreaseStep[{maxCellMeasure_, so_}, print_: False] := 
  Block[{}, 
   If[print, Print["maxCellMeasure=" <> ToString[maxCellMeasure/2]]];
   {maxCellMeasure/2, sol[maxCellMeasure/2]}]

compare[sol1_, sol2_] := 
  Block[{points = sol1[Grid][[-1]]},
    ( (sol2[Sequence @@ #] & /@ points)
    - (sol1[Sequence @@ #] & /@ points)) // Abs // Max]

solu[initialMaxCellMeasure_, tol_, maxSteps_: 10, print_: False] :=
  NestWhileList[decreaseStep[#, print] &,
                {initialMaxCellMeasure, sol[initialMaxCellMeasure]}, 
                (compare[#1[[2]], #2[[2]]] > tol)&,
                2, maxSteps][[-1]]
solu0 = solu[0.05, 0.001, 10, True]

enter image description here

Any more elegant solution would be appreciated.

Denis

Improve this answer
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.