I am trying to solve the following Laplace PDE in MMA:

$$\frac{\partial^2 T}{\partial x^2}+\frac{\partial^2 T}{\partial z^2}=0$$

subject to the boundary conditions: $$\frac{\partial T}{\partial z}|_{z=0}=0,$$

and, at $z=1$ $$-\frac{\partial T}{\partial z}+A=BT,$$

This is my code:

DSolve[{D[T[x, z], {x, 2}] + D[T[x, z], {z, 2}] == 0, 
(D[T[x, z], z] /. {z -> 0}) == 0, 
-(D[T[x, z], z] /. {z -> 1}) + A == B T[x, 1]}, T[x, z], {x, z}]

But it just repeats my equation.

Thank you for kindly help in advance!

  • $\begingroup$ Can you check to make sure that you don't have any definitions for the symbols you use in your equations? $\endgroup$
    – MarcoB
    Sep 25, 2015 at 5:14
  • $\begingroup$ The syntax in D[T[x, 0],z]is not correct: you should evaluate in $z=0$ after taking the derivative. This yields 0==0 and so your error. Apart from this, are you sure about you boundary/initial conditions? Seems to me you are missing some. $\endgroup$
    – anderstood
    Sep 25, 2015 at 5:32
  • $\begingroup$ Hi, @anderstood, thanks. Please check my modified post. It seems that MMA has done nothing... $\endgroup$
    – Enter
    Sep 25, 2015 at 5:49
  • $\begingroup$ What if you solve only the pde (DSolve[{D[T[x, z], {x, 2}] + D[T[x, z], {z, 2}] == 0}, T[x, z], {x, z}]) and add the boundary conditions afterwards to identify the constants? $\endgroup$
    – anderstood
    Sep 25, 2015 at 5:55
  • $\begingroup$ Thanks, @anderstood. Your method looks feasible. But I have not realized it by now... $\endgroup$
    – Enter
    Sep 25, 2015 at 6:10

2 Answers 2


Another problem easy to do by separation of variables. However, BC's on x are definitely required to solve the problem. The assumption of exp[i a x] requires BC's that make sinusoidal solutions in x.

pde = D[T[x, z], x, x] + D[T[x, z], z, z] == 0;

With the BC's you give I add 2 more for x. The ones I chose will give sinusoidal solutions. You, of course can choose others

bc1 = T[0, z] == 0;
bc2 = T[1, z] == 0;
bc3 = (D[T[x, z], z] /. z -> 0) == 0;
bc4 = (D[T[x, z], z] /. z -> 1) == A - B T[x, 1];

Separate T in the form

T[x_, z_] = X[x] Z[z];

pde/T[x, z] // Apart

(* D[X[x],x,x]/X[x]+D[Z[z],z,z]/Z[z]==0 *)

One terms is dependent on x, the other z, so each must be equal to a constant.

xeq = D[X[x], x, x]/X[x] == -alpha^2

 X[x_] = (X[x] /. DSolve[xeq, X[x], x][[1]]) /. {C[1] -> c1, C[2] -> c2};

zeq = D[Z[z], z, z]/Z[z] == alpha^2

Z[z_] = (Z[z] /. DSolve[zeq, Z[z], z][[1]]) /. {C[1] -> c3,  C[2] -> c4}

(* c3 E^(alpha z)+c4 E^(-alpha z) *)

From the z BC's, it is obvious we want hyberbolics rather than exponentials.

(Z[z] // ExpToTrig) // Collect[#, {Sinh[alpha z], Cosh[alpha z]}] &

(* (c3-c4) Sinh[alpha z]+(c3+c4) Cosh[alpha z] *)

Z[z_] = % /. {c3 - c4 -> c3, c3 + c4 -> c4};
(* c3 Sinh[alpha z] + c4 Cosh[alpha z] *)

Apply the BC's


(* c1 (c3 Sinh[\[Alpha] z] + c4 Cosh[\[Alpha] z]) == 0 *)

From which

c1 = 0

T[x, z]
(* c2 Sin[alpha x] (c3 Sinh[alpha z] + c4 Cosh[alpha z])*)

Combine constants

c2 = 1;


(* Sin[alpha] (c3 Sinh[alpha z] + c4 Cosh[alpha z]) == 0 *)

From which

alpha = n Pi

with n required to be

$Assumptions = n \[Element] Integers && n > 0;

(* Pi c3 n Sin[Pi n x]==0 *)

From which

c3 = 0;

bc4 // Simplify
(* c4 Sin[Pi n x] (B Cosh[Pi n]+Pi n Sinh[Pi n])==A *)

Use orthogonality to solve for c4

Integrate[%[[1]] Sin[n Pi x], {x, 0, 1}] ==  Integrate[%[[2]] Sin[n Pi x], {x, 0, 1}];

c4 = c4 /. Solve[%, c4][[1]]

T[x, z]

(*-((2 A ((-1)^n - 1) Sin[Pi n x] Cosh[Pi n z])/(Pi n (B Cosh[Pi n] + Pi n Sinh[Pi n]))) *)

The actual T will be an infinite series in n, but obviously from the form of the solution, the even n terms are 0. In order to not compute the 0 terms we can do the following:

Tn[x_, z_] = (T[x, z] /. n -> 2 m - 1) /. m -> n // Simplify

(*(4 A Sin[Pi (2 n-1) x] Cosh[Pi (2 n-1) z])/(Pi (2 n-1) (B Cosh[Pi-2 Pi n]+Pi (2 n-1) Sinh[Pi (2 n-1)]))*)

The general solution is the above summed from n = 1 to infinity. To make a plot:

tn = Function[n, #] &@Tn[x, z];

Assign A and B some values.

A = 1;
B = 2;

Tterms[k_] := Compile[{x, z}, #] &@Total@tn@Range[k]

T = Tterms[50];

Plot3D[T[x, z], {x, 0, 1}, {z, 0, 1}, AxesLabel -> {"x", "y", "T"}]

enter image description here


There is a problem about your solution, it is that DSolve cannot handle it and, therefore, returns you the input unevaluated, as a signal. So, what you can do, if you do not want to solve it analytically (seems a simple task, does not it?) you might want to solve it numerically, that is, using NDSolve. This, of course, will reguire from you to somehow fix the values of A and B, otherwise, however, it can be solved. Try this:

   A = 1;
B = 2;

<< NDSolve`FEM`
mesh = ToElementMesh[Rectangle[{0, 0}, {1, 1}]];
nv1 = NeumannValue[A - B*T[x, z], z == 1];
nv2 = NeumannValue[0, z == 0];
dc = DirichletCondition[T[x, z] == 1, x == 0];
sl = NDSolve[{D[T[x, z], {x, 2}] + D[T[x, z], {z, 2}] == nv1 + nv2, 
   dc}, T[x, z], {x, z} \[Element] mesh]

(*  {{T[x, z] -> InterpolatingFunction[{{0., 1.}, {0., 1.}}, <>][x, z]}}  *)

To have a look at the solution evaluate this:

Plot3D[T[x, z] /. sl[[1, 1]], {x, z} \[Element] mesh, 
 ImageSize -> 300]

which should look as follows:

enter image description here

Take into account that I added a Dirichlet condition

 dc = DirichletCondition[T[x, z] == 1, x == 0];

otherwise there is nothing to look at in the solution, as it was.

Have fun!

  • $\begingroup$ @ Alexei Boulbitch, Thank you. But I have some troubles in debugging your code. I am using MMA 9. First, what is the meaning of << before your NDSolveFEM``. I guess you want to use Needs? Second, when I run it, MMA gives NDSolve::dvnoarg: "The function T[x,z] appears with no arguments. ". In addition, I think the original two boundary condition are sufficient to solve for the function T(x,z) even analytically. For example, assuming $T(x,z)=\theta(z)\exp[i\alpha x]$, we can obtain $T(x,z)=\frac{A \cosh \alpha z}{B\cosh \alpha + \alpha \sinh \alpha}$. $\endgroup$
    – Enter
    Sep 26, 2015 at 3:13
  • 1
    $\begingroup$ @ lxy Welcome. To address your questions. 1. << NDSolveFEM` is the short form of Needs... . You may use the latter instead. 2. I am not quite sure, if Mma 9 already has the NDSolveFEM package. If not (which is possible), it is the reason for the strange warnings like "The function T[x,z] appears with no arguments". 3. Indeed the original boundary conditions enable one to solve the equation. The solution, however, was a constant. That's why I added something to see that there is at least some variation, and in this way to check that there is, indeed, a solution. $\endgroup$ Sep 28, 2015 at 7:24
  • $\begingroup$ @ Alexei Boulbitch, Thanks, I confirmed that MMA9 does not have FEM package, and it can run on version 10 :) $\endgroup$
    – Enter
    Sep 28, 2015 at 9:07

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.