Short story: I am trying to specify the maximum step size for NDSolve function, but it appears to be working only if I do not specify boundary conditions using NeumannValue and DirichletCondition functions.
Long story: I have the following set of equations and boundary conditions
Gc0 = DirichletCondition[c[y] == cB, y == -1];
GmA = NeumannValue[- (Q[y] + (2 + c[y])/3), y == -1];
GmB = NeumannValue[(Q[y] + (2 + c[y])/3), y == 1];
GQA = NeumannValue[-(4/15) m[y], y == -1];
GQB = NeumannValue[(4/15) m[y], y == 1];
eqs1 = {
- c'[y] + m[y] == 0, Gc0,
- m''[y] + Q'[y] + c'[y]/3 + m[y] == GmA + GmB,
- Q''[y] + (4/15) m'[y] + Q[y] == GQA + GQB
}
So I set up the solver i the standard way, with the max step size being 2/100
nSol[eq_, par_] := NDSolve[eq /. par, {c, m, Q}, {y, -1, 1}, MaxStepSize -> 2/100]
res = nSol[eqs1, {cB -> 0.03}];
Plot[{c[y] /. res, m[y] /. res, Q[y] /. res}, {y, -1, 1}, PlotRange -> All]
The above works completely fine, but when I checked the grid of the interpolating function I find
c["Grid"] /. res
{{{-1.}, {-0.9}, {-0.8}, {-0.7}, {-0.6}, {-0.5}, {-0.4}, {-0.3},{-0.2}, {-0.1}, {0.}, {0.1}, {0.2}, {0.3}, {0.4}, {0.5}, {0.6},{0.7}, {0.8}, {0.9}, {1.}, {-0.95}, {-0.85}, {-0.75}, {-0.65}, {-0.55}, {-0.45}, {-0.35}, {-0.25}, {-0.15}, {-0.05}, {0.05}, {0.15}, {0.25}, {0.35}, {0.45}, {0.55}, {0.65}, {0.75}, {0.85}, {0.95}}}
Which does not look like spacing from -1 to 1 with the max step size being 2/100.
However, if I set my equations in the following way
eqs2 = {
- c'[y] + m[y] == 0,
- m''[y] + Q'[y] + c'[y]/3 + m[y] == 0,
- Q''[y] + (4/15)m'[y] + Q[y] == 0,
c[-1] - cB == 0,
m'[-1] - Q[-1] -(2+c[-1])/3 == 0,
m'[1] - Q[1] -(2+c[1])/3 == 0,
Q'[-1] - (4/15)m[-1] == 0,
Q'[1] - (4/15)m[1] == 0
}
and I repeat the calculation I get the same plots, but the MaxSetpSize now seem to be doing what I would expect
c["Grid"] /. res
{{{-1.}, {-0.999904}, {-0.999808}, {-0.995684}, {-0.99156}, {-0.987436}, {-0.974794}, {-0.962152}, {-0.94951}, {-0.936868}, {-0.916868}, {-0.896868}, {-0.876868}, {-0.856868}, {-0.836868}, {-0.816868}, {-0.796868}, {-0.776868}, {-0.756868}, {-0.736868}, {-0.716868}, {-0.696868}, {-0.676868}, {-0.656868}, {-0.636868}, {-0.616868}, {-0.596868}, {-0.576868}, {-0.556868}, {-0.536868}, {-0.516868}, {-0.496868}, {-0.476868}, {-0.456868}, {-0.436868}, {-0.416868}, {-0.396868}, {-0.376868}, {-0.356868}, {-0.336868}, {-0.316868}, {-0.296868}, {-0.276868}, {-0.256868}, {-0.236868}, {-0.216868}, {-0.196868}, {-0.176868}, {-0.156868}, {-0.136868}, {-0.116868}, {-0.0968684}, {-0.0768684}, {-0.0568684}, {-0.0368684}, {-0.0168684}, {0.00313162}, {0.0231316}, {0.0431316}, {0.0631316}, {0.0831316}, {0.103132}, {0.123132}, {0.143132}, {0.163132}, {0.183132}, {0.203132}, {0.223132}, {0.243132}, {0.263132}, {0.283132}, {0.303132}, {0.323132}, {0.343132}, {0.363132}, {0.383132}, {0.403132}, {0.423132}, {0.443132}, {0.463132}, {0.483132}, {0.503132}, {0.523132}, {0.543132}, {0.563132}, {0.583132}, {0.603132}, {0.623132}, {0.643132}, {0.663132}, {0.683132}, {0.703132}, {0.723132}, {0.743132}, {0.763132}, {0.783132}, {0.803132}, {0.823132}, {0.843132}, {0.863132}, {0.883132}, {0.903132}, {0.923132}, {0.943132}, {0.963132}, {0.981566}, {1.}}}
I also tried different values of MaxStepSize, but it seems to be doing something only if I define the boundary conditions as in eqs2 and not as in eqs1.
What I am missing here?
DirichletConditionandNeumannValuetrigger the use of the finite element method, which leads to a spatial discretization of your ode. Now,MaxStepSizeis an option for a temporal integrator inNDSolve. If you want to time integrate the ODE (this is probably what you want) then you should rewrite theDirichletConditionandNeumannValuein equational form. $\endgroup$