I have a PDE problem that I can't solve. and I wonder if you could help me with it.

I have a Terzaghi Consolidation equation (similar to the heat equation) that is written as follows:


$u$ is pore water pressure

$z$ is depth

$t$ is time

$c$ is a constant

pore water pressure is equal at all depth at $t = 0$, but it becomes $0$ at top and bottom directly after that.

u = u[z,t] ($z$ is depth and $t$ is time)

with boundary conditions as follows:

u[z, 0] = 100 ($u$ at all depth when $t = 0$)

u[2, t] = 0 ($u$ at top of the layer; in this case, $2$ is the height of the top layer)

u[0, t] = 0 ($u$ at bottom of the layer)

I have solved this kind of problem before using a finite difference technique in a spread sheet, I wanted to find out whether Mathematica can solve this problem directly, but unfortunately I failed.

This is what I entered into Mathematica:

NDSolve[{Derivative[0, 1][u][z ,t] == Derivative[2, 0][u][z, t],
 u[0, t > 0] == 0, u[2,t > 0] == 0, u[2 > z > 0, 0] == 100}, u, {z, 0, 2},{t, 0, 100}]

Unfortunately, evaluating this expression leads gives the error:


"To avoid possible ambiguity, the arguments of the dependent variable in u[0,t>0]==0; should literally match the independent variables"

I dont understand why I'm getting this message.

I found a heat problem similar to this. I tried to use the that equation but changed the boundary problem, but that did not work.

Can you help me detect which part of equation is wrong, and how I can fix the boundary values?

  • $\begingroup$ @Martin Wijaya: the boundary conditions as you have given them are not valid Mathematica syntax. The are also obviously not consistent with the initial conditions at u[0,0] and u[2,0]. See Nassers answer for the correct syntax... $\endgroup$ Dec 18, 2012 at 19:32

2 Answers 2


You can get rid of the warnings if you modify your initial condition close to the space boundaries :

limit = 3;
newIC[z_] = 100 (HeavisideTheta[z - 10^-limit] - HeavisideTheta[z - (2 - 10^-limit)])

sol[z_, t_] = NDSolve[{Derivative[0, 1][u][z, t] == Derivative[2, 0][u][z, t], 
  u[0, t] == 0, u[2, t] == 0, u[z, 0] == newIC[z]}, 
  u[z, t], {z, 0, 2}, {t, 0, 100}][[1, 1, 2]]

Plot[{newIC[z], sol[z, 0]}, {z, 0, 2}, PlotRange -> All, Exclusions -> None]

enter image description here

  • $\begingroup$ Dear b.gatessucks, what does space boundaries means?what does newIC[z_] means?, I'm a bit confused with the new boundary condition. $\endgroup$ Dec 29, 2012 at 3:06
  • $\begingroup$ By space boundaries I mean the points z=0 and z=2. The new initial condition I propose takes the value 100 everywhere except in the immediate vicinity of those two points and so there are no warnings. $\endgroup$ Dec 29, 2012 at 9:42
  • $\begingroup$ Oh, I see..ok thx. $\endgroup$ Dec 30, 2012 at 13:12

I need to look more into this, the warning about

NDSolve::ibcinc: Warning: boundary and initial conditions are inconsistent.

As I've seen this many times. Are you sure the IC's you are using are correct? (ie. consistent) ? Please double check the textbook or from the source you obtained this.

But just wanted to say that you still get a solution and can animate it just fine. A different method solver might be needed like MethodOfLines. But need to look more into this. Here is the solution

ClearAll[u, z, t]
eq = Derivative[2, 0][u][z, t] == Derivative[0, 1][u][z, t];
ic = {u[0, t] == 0, u[2, t] == 0, u[z, 0] == 10};

sol = u /. First@NDSolve[Flatten[{eq, ic}], u, {z, 0, 2}, {t, 0, 10}]

Plot3D[sol[z, t], {z, 0., 2.}, {t, 0., 10}, PlotRange->All, AxesLabel->{z, t, u[z, t]}]

Mathematica graphics

  • $\begingroup$ the boundary conditions can be made consistent with the initial conditions by using something like Exp[-a*t] with an a of users choice. From comparing with giving these bc explicitly I suspect that Mathematica automatically chooses these boundary conditions with a==1 (or something very close to that). $\endgroup$ Dec 18, 2012 at 19:38
  • $\begingroup$ back in the old days (V4) I used to get around discontinuities or inconsistencies with an If in the conditions. $\endgroup$ Dec 18, 2012 at 21:06
  • $\begingroup$ Dear Nasser, Why do I have to use u/.First@NDSolve and Flatten instead of directly typing it as Plot3D[NDSolve[{Derivative[0, 1][u][z, t] == Derivative[2, 0][u][z, t], u[0, t] == 0, u[2, t] == 0, u[z, 0] == 100}, u, {z, 0, 2}, {t, 0, 100}], {z, 0., 2.}, {t, 0., 100}] $\endgroup$ Dec 29, 2012 at 3:04

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.