# NDSolve for PDE with discontinuous initial/terminal condition

I have an issue with NDSolve for the case of a PDE with discontinuous initial/terminal condition. Consider the PDE solution

Z=z/.First[NDSolve[{
D[z[t,x],t] == z[t,x] D[z[t,x],x]-1/2 (1-t)^2 D[z[t,x],{x,2}],
z[1,x]==UnitStep[x],
z[t,-5000]==0,
z[t,5000]==1},
{z},{t,0,1},{x,-5000,5000}]]


The terminal condition is the discontinuous UnitStep function. When I plot the solution for t=1 now,

Plot[Z[1, x], {x, -10, 10}]


Mathematica has obviously smoothened the terminal condition (which was explicitly given by the UnitStep function).

Is it possible to disable this smoothing? Or at least, I want to have a smooth function that does not exceed 1 and does not go below 0.

Note: What I am not searching for is a solution like "increase the MaxSteps" or something similar. I am interested in options on how Mathematica processes the initial/terminal condition.

Thanks a lot for your help!!

-

## 1 Answer

For the given grid resolution, you're probably going to have to smooth the UnitStep function with a width equal to that grid spacing (at least). Here is one way to do that, using the Fermi function instead of UnitStep:

f[x_, d_] := 1/(E^(-(x/d)) + 1)
Z = z /. First[
With[{d = 5000./10000},
NDSolve[{D[z[t, x], t] ==
z[t, x]*D[z[t, x], x] - 1/2*(1 - t)^2*D[D[z[t, x], x], x],
z[1, x] == f[x, d], z[t, -5000] == f[-5000, d],
z[t, 5000] == f[5000, d]}, {z}, {t, 0, 1}, {x, -5000, 5000}]]]

Plot[Z[1, x], {x, -10, 10}]


Here, the parameter d is the smoothing of the vertical step.

-
Thanks for your quick reply. Smoothing the function in advance looks like a good workaround to me. However, I would like to know what kind of smoothing procedure Mathematica uses internally and if I can modify or disable it. –  Posch79 Jan 5 '13 at 19:12
@Posch79 Could just be numerical instability of whatever solution method is used, try tweaking the Method argument and see if/how it changes –  ssch Jan 5 '13 at 22:06
In general, though, I don't know the answer to how the interpolation can be influenced. No matter what you do, the grid limits how fast your spatial variations can be. And UnitStep will always be too steep. But one can probably construct a more customized solution using NDSolveProcessEquations`. –  Jens Jan 5 '13 at 22:23