6
$\begingroup$

So I wish to stop integrating my PDE when the spatial gradient is larger than some value, let's say 10. That is to say, my WhenEvent condition wants to be stop when the the maximum gradient over all x at a given t exceeds a threshold.

However, when I try implement the code I get "The function value is not True or False errors".

Here's my attempt:

NDSolve[{Derivative[0, 1][h][x, t] + 2000 h[x, t]*Derivative[1, 0][h][x, t] == 
  1/4 D[(1 - 4 h[x, t]^2)^3 Derivative[1, 0][h][x, t], x], 
  h[-1/2, t] == 0, h[1/2, t] == 0, h[x, 0] == 0.1 Sin[8 π (x - 1/2)], 
  WhenEvent[NMaximize[Derivative[1, 0][h][x, t]] > 10, tEnd = t; 
  "StopIntegration"]}, h[x, t], {x, -1/2, 1/2}, {t, 0, 10}, 
  Method -> {"MethodOfLines", "DiscretizedMonitorVariables" -> True, 
  "SpatialDiscretization" -> {"TensorProductGrid", "MinPoints" -> 100000}}];

Any ideas what I'm doing wrong?

$\endgroup$
  • 1
    $\begingroup$ Why Max[ ] .. ? $\endgroup$ – Dr. belisarius Mar 4 '16 at 15:38
  • $\begingroup$ I want to know the time at when the gradient goes over 10 anywhere $\endgroup$ – M Horsley Mar 4 '16 at 15:43
  • $\begingroup$ Not sure why to use Table? $\endgroup$ – M Horsley Mar 4 '16 at 15:45
  • $\begingroup$ is the idea to reject the entire solution in the event it has a large gradient? I think you need to compute the solution then analyze it. $\endgroup$ – george2079 Mar 4 '16 at 15:48
  • $\begingroup$ Tried amending code to use NMaximise. It returns an error: "The function value is not a number at {x}..." $\endgroup$ – M Horsley Mar 4 '16 at 15:48
6
$\begingroup$

After playing with this a bit: it seems WhenEvent[ h[x,t] ] gets passed the list of discrete x-grid values at each time step. So we can numerically compute the gradient:

rr = Reap[
   sol = NDSolve[{...,
      WhenEvent[
       Sow[Max[Abs[Subtract @@@ Partition[ h[x, t] , 2, 1]]] 100000] > 10, 
       "StopIntegration"]},
     h[x, t], {x, -1/2, 1/2}, {t, 0, 10},
     Method -> {"MethodOfLines", 
       "DiscretizedMonitorVariables" -> True, 
       "SpatialDiscretization" -> {"TensorProductGrid", 
         "MinPoints" -> 100000}}]];
 f = h[x, t] /. First@sol;

note I assume uniform grid spacing which in this case is set by your specified minpoints (But that's not guaranteed). You can verify that by looking at the interpolation function grid.

The critical time:

tmax = f[[0]]["Domain"][[2, 2]]

0.000172062

Plot3D[f, {x, -1/2, 1/2}, {t, 0, tmax}]

enter image description here

here is the derivative and you can see it becoming unstable..

g = D[f, x];
Plot3D[g, {x, -1/2, 1/2}, {t, 0, tmax}, PlotRange -> {-10, 10}]

enter image description here

here is the max gradient value captured by Reap/Sow

ListPlot[rr[[2, 1]]]

enter image description here

$\endgroup$
  • $\begingroup$ Thank you very much for the detailed response. This works brilliantly! $\endgroup$ – M Horsley Mar 9 '16 at 12:21

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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