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

$\frac{du}{dt}=c\frac{d^2u}{dt^2}$

$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:

DSolve::litarg:

"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?

share|improve this question
    
@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... –  Albert Retey Dec 18 '12 at 19:32
add comment

2 Answers

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

share|improve this answer
    
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). –  Albert Retey Dec 18 '12 at 19:38
    
back in the old days (V4) I used to get around discontinuities or inconsistencies with an If in the conditions. –  Mike Honeychurch Dec 18 '12 at 21:06
    
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}] –  Martin Wijaya Dec 29 '12 at 3:04
add comment

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

Remove[sol]
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

share|improve this answer
    
Dear b.gatessucks, what does space boundaries means?what does newIC[z_] means?, I'm a bit confused with the new boundary condition. –  Martin Wijaya Dec 29 '12 at 3:06
    
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. –  b.gatessucks Dec 29 '12 at 9:42
    
Oh, I see..ok thx. –  Martin Wijaya Dec 30 '12 at 13:12
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.