Take the 2-minute tour ×
Mathematica Stack Exchange is a question and answer site for users of Mathematica. It's 100% free, no registration required.

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}]

plot

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!!

share|improve this question
add comment

1 Answer 1

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}]

Fermi

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

share|improve this answer
    
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 NDSolve`ProcessEquations. –  Jens Jan 5 '13 at 22:23
add comment

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

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