# Trouble with differential equation

I tried to solve this differential equation:

$$\epsilon y''(x)+xy'(x)=-\epsilon \pi^2 \cos(\pi x)-\pi x\sin(\pi x)$$

with boundary conditions: $y(-1)=-2, \space y(1)=0$. If we take $\epsilon=0.1$, Mathematica can solve it without any trouble

Block[{e = 0.1, min = -1, max = 1},
Plot[Evaluate[
y[x] /. NDSolve[{e y''[x] +
y'[x] x == -e Pi^2 Cos[Pi x] - Pi x Sin[Pi x],
y[min] == -2, y[max] == 0}, y, {x, min, max}]], {x, min, max}]
]


But if we want a smaller $\epsilon$, let say 0.01, Mathematica seems unable to handle it. Is there any options to invoke or methods to employ to get the desired result? Anyway, this is the solution for $\epsilon=0.0001$.

Thank you.

• Piece of advice: don't use Block to inject values into parameters. Use With instead. – m_goldberg Jun 25 '15 at 14:24
• The fundamental problem would appear to be that, in the limit of small e, the order of the equation drops from second to first, with the result that there is one too many boundary conditions. – bbgodfrey Jun 25 '15 at 14:30

DSolve can handle this.

Clear[y];
y[x_, e_] = y[x] /. DSolve[{
e y''[x] + y'[x] x == -e Pi^2 Cos[Pi x] - Pi x Sin[Pi x],
y[-1] == -2, y[1] == 0}, y[x], x][[1]] // Simplify


Manipulate[
Plot[y[x, e], {x, -1, 1}],
{{e, 0.01}, 0.0001, 0.1, Appearance -> "Labeled"}]


• Thank you for your answer, but I'm looking for numerical solution. – Deco Jun 25 '15 at 14:59
• @Deco Bob's symbolic solution can be evaluated numerically quickly, and with better accuracy than a mere numerical approximation. In fact DSolve can solve the problem in full generality: DSolve[{e y''[x] + y'[x] x == -e Pi^2 Cos[Pi x] - Pi x Sin[Pi x], y[min] == -2, y[max] == 0}, y, x]. Is there some reason you do not think this is a superior solution? – Michael E2 Jun 25 '15 at 15:55

As noted by @bbgodfrey, the "shooting" algorithms that Mathematica tends to use are not well-adapted to this particular equation. Better would be some kind of relaxation method, which is what Mathematica uses (I think) for solving PDEs on a mesh. And an ODE is just a PDE in one dimension, so let's try solving this equation on a one-dimensional mesh:

Needs["NDSolveFEM"]
truey[x_,e_] = Cos[\[Pi] x] + Erf[x/(Sqrt[2] Sqrt[e])]/Erf[1/(Sqrt[2] Sqrt[e])];

e = 0.0001; min = -1; max = 1;
mesh = ToElementMesh[Interval[{min, max}], MaxCellMeasure -> 0.05];
bcs = {DirichletCondition[y[x] == -2, x == min], DirichletCondition[y[x] == 0, x == max]}
soln = NDSolve[{e y''[x] + y'[x] x == -e Pi^2 Cos[Pi x] - Pi x Sin[Pi x], bcs}, y, Element[x, mesh]]
Plot[{Evaluate[y[x] /. soln], truey[x, e]}, {x, min, max}]
Clear[e, min, max]


I've plotted here the result from Mathematica's finite element solver (in blue) versus the true solution found by @BobHanlon (in yellow). I've actually used a coarser mesh than is ideal to show the difference between the two; if you set MaxCellMeasure -> 0.01 in the above code (instead of 0.05), the two curves are indistinguishable at this resolution.

To see why NDSolve has difficulty with this problem for very small e, consider that NDSolve solves this two-point boundary value problem by some form of shooting. In other words, it varies y'[min] until one is found that yields the desired y[max]. However, as e becomes very small, the sensitivity of y[max] to y'[min] becomes great, because the differential equation becomes singular in that limit.

y''[x] + y'[x] x/e == -Pi^2 Cos[Pi x] - Pi x Sin[Pi x]/e


This is also apparent from Bob Hanlon's symbolic solution. To illustrate the sensitivity of y[max] to y'[min], consider

slope = D[Cos[π x] + Erf[x/(Sqrt[2] Sqrt[e])]/Erf[1/(Sqrt[2] Sqrt[e])],  x] /. x -> -1;
LogLogPlot[N[slope], {e, 0.0001, .1}, PlotRange -> All]


You can work around this sensitivity by increasing WorkingPrecision, e. g. for e = .01,

With[{e = 1/100, min = -1, max = 1},
Plot[Evaluate[y[x] /. NDSolve[{e y''[x] + y'[x] x == -e Pi^2 Cos[Pi x] -
Pi x Sin[Pi x], y[min] == -2, y[max] == 0}, y, {x, min, max},
WorkingPrecision -> 50, MaxSteps -> 50000]], {x, min, max}]]


but doing so rapidly becomes prohibitively expensive as e is further reduced.