I have to find and plot a numerical solution for the following equation (I have to write a solver):

$$u_{t} = (u^2 u_x)_x$$ with the following conditions $u(0,t) =0, u(1,t) = \sqrt{\frac{2c-2}{t}}, u(x,0) = 0$, where $c$ is a speed of the wave

I try to use the following method:

$$\frac{u_{m}^{n+1} - u_{m}^n}{\tau} = \frac{F_{m +1/2}^n - F_{m-1/2}^n}{h}$$ where $$F_{m+1/2}^n = (u_{m+1/2}^n)^2(\frac{u_{m+1}^n - u_{m-1}^n}{h})$$ where $$u_{m+1/2}^n = \frac{u_{m+1}^n + u_{m}^n}{2}$$

tt = 0.1;

tnn = 10000;

xnn = 20;

tau = tt/tnn;

h = 1/xnn;

u0[x_] := x;

Do[v[m, 0] = u0[h*m], {m, 0, xnn}]

Do[v[0, n] = 0; v[xnn, n] = 0;
 Do[v[m, n + 1] = 

v[m, n] + tau/h*( ((v[m + 1, n] - v[m, n])/ h)(((v[m + 1, n])^2 + (v[m, n])^2)/ 2) - ((v[m, n] - v[m - 1, n])/ h)(((v[m, n])^2 + (v[m - 1, n])^2)/2) ), {m, 1, xnn - 1}], {n, 0, tnn - 1}]

u = ListInterpolation[
   Table[v[i, j], {i, 0, xnn}, {j, 0, tnn}], {{0, 1}, {0, tt}}];

Plot3D[u[x, t], {x, 0, 1}, {t, 0, tt}, PlotRange -> All]

I have changed the code and have got the correct solvation, but for the u(0,t)=u(1,t) = 0, u(x,0) = f(x). I don't actually understand how should I write the conditions like u(1,t) = f(t).

  • $\begingroup$ I doubt about your last expression. Would it be with a plus sign instead of a minus sign? $\endgroup$
    – nicoguaro
    Apr 4, 2017 at 21:17
  • $\begingroup$ yeap, sure. Misprint. $\endgroup$ Apr 4, 2017 at 21:30
  • $\begingroup$ What are the boundary and initial conditions for your problem? $\endgroup$
    – nicoguaro
    Apr 4, 2017 at 21:43
  • 1
    $\begingroup$ probably it's better to use implicit method instead of explicit for solving this equation. $\endgroup$ Apr 4, 2017 at 22:22
  • 3
    $\begingroup$ The nonlinear equation definitely needs a smaller time step. Think of $u^2$ as a conductivity coefficient -- the time step then needs to be proportional to it. You definitely want to use an implicit method. $\endgroup$ Apr 4, 2017 at 23:31

1 Answer 1


You got good advices in comments, but as it is only 1D problem, you can manage it using fully explicit method with very small time steps. I checked your code in Mathematica and there are several misunderstandings there. After correcting them, the code works for me.

Similarly, your description in the question contains, I suppose, some mistakes. As I guess that this is some student project, so let me only indicate what shall be corrected ;-).

Your method in the question can be viewed as finite volume method for computational cells $(x_{m-1/2},x_{m+1/2})$. The quantity $F_{m +1/2}^n$ approximates the flux at the right point, the approximation in your question is inappropriate, check it, it shall use so-called central difference at point $x_{m+1/2}$. Analogously the nonlinear conductivity coefficient $u_{m+1/2}^n$ at $x_{m+1/2}$ is approximated inproperly.

Your code anyway does not implement what you describe in your question. So, check the scheme and implement it exactly and I suppose your code can work.

Be careful about a correct input data for your model. In your implementation you use different boundary and initial conditions as you describe in the question.

I hope after these indications you manage to finish the task, good luck!

  • $\begingroup$ Peter, did you try to solve the equation with the conditions I've written in the question or in the code? Because I can't get any solvation with the conditions like $u = \sqrt{\frac{2c - 2}{t}}$. $\endgroup$ Apr 6, 2017 at 19:25
  • $\begingroup$ I tried the one in the code. The one you write in the question can cause some problems as it is unbounded for t=0, definitely for fully explicit scheme. $\endgroup$ Apr 7, 2017 at 7:12

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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