NotebookDirectory[]

Equilibrium temperature distribution - 2D problem

The problem

The objective is to solve the PDE 

(((∂^2)u)/(∂(x^2)))(x,y)+(((∂^2)u)/(∂(y^2)))(x,y)=0  on  {(x,y) ∈  | 0≤x ≤K, 0 ≤y≤L},

subject to the conditions

u(x,0) = (Subscript[f, 1])(x), u(x,L)=(Subscript[f, 2])(x)   (call these BCx) 

u(0,y) = (Subscript[g, 1])(y), u(K,y)=(Subscript[g, 2])(y)  (call these BCy)

The trick is to split this problem into two problems

Problem 1

The objective is to solve the PDE 

(((∂^2)u)/(∂(x^2)))+(((∂^2)u)/(∂(y^2)))=0  on  0≤x ≤K, 0 ≤y≤L,

subject to the conditions

u(x,0) = (Subscript[f, 1])(x), u(x,L)=(Subscript[f, 2])(x)  (call these BCx) 

u(0,y) = 0, u(K,y)=0  (call these DBCy)

Step 1. First ignore BCx  and solve

(((∂^2)u)/(∂(x^2)))+(((∂^2)u)/(∂(y^2)))=0  on  0≤x ≤K, 0 ≤y≤L,

subject to the conditions

u(0,y) = 0, u(K,y)=0  (call these DBCy)

Using the Separation of Variables (SofV) method we find few solutions of this problem:

Sin[((n Pi)/K)x] ((Sinh[((n Pi)/K)y])/(Sinh[((n Pi)/K)L]))

Sin[((n Pi)/K)x] ((Sinh[((n Pi)/K)(L-y)])/(Sinh[((n Pi)/K)L]))

Test these solutions:

(D[#,{x,2}]+D[#,{y,2}])&[Sin[((n Pi)/K)x](( Sinh[((n Pi)/K)y])/( Sinh[((n Pi)/K)L]))]

FullSimplify[(Sin[((n Pi)/K)x] (( Sinh[((n Pi)/K)y])/( Sinh[((n Pi)/K)L])))/.{x->{0,K}}]

We need to tell Mathematica that n is an integer.

FullSimplify[(Sin[((n Pi)/K)x] (( Sinh[((n Pi)/K)y])/( Sinh[((n Pi)/K)L])))/.{x->{0,K}},n∈Integers]

(D[#,{x,2}]+D[#,{y,2}])&[Sin[((n Pi)/K)x] ((Sinh[((n Pi)/K)(L-y)])/( Sinh[((n Pi)/K)L]))]

FullSimplify[(Sin[((n Pi)/K)x](( Sinh[((n Pi)/K)(L-y)])/(Sinh[((n Pi)/K)L])))/.{x->{0,K}},n∈Integers]

Step 2. Now that we have few solutions we form many solutions.  This is the From-Few-Many (FFM) idea, which is commonly known as the superposition principle:

(∑_(n=1)^nn) (Subscript[a, n ])Sin[((n Pi)/K)x] ((Sinh[((n Pi)/K)y])/(Sinh[((n Pi)/K)L]))+(∑_(n=1)^nn) (Subscript[b, n ])Sin[((n Pi)/K)x](( Sinh[((n Pi)/K)(L-y)])/(Sinh[((n Pi)/K)L]))

Next we choose (Subscript[a, n]) and (Subscript[b, n]) such that the above function satisfies BCx conditions. First substitute y=0. This leads to the formula for (Subscript[b, n]).

(Subscript[f, 1])(x) =(∑_(n=1)^nn) (Subscript[b, n ])Sin[((n Pi)/K)x]

(Subscript[f, 1])(x)Sin[((j Pi)/K)x] =(∑_(n=1)^nn) (Subscript[b, n ])Sin[((n Pi)/K)x]Sin[((j Pi)/K)x]

Clear[K]; FullSimplify[Integrate[Sin[((n Pi)/K)x]Sin[((j Pi)/K)x],{x,0,K}],And[n∈Integers,j∈Integers,Or[j>n,j<n]]]

Clear[K]; FullSimplify[Integrate[Sin[((n Pi)/K)x]Sin[((n Pi)/K)x],{x,0,K}],And[n∈Integers]]

(∫_0^K)(Subscript[f, 1])(x)Sin[((j Pi)/K)x]dx =(∑_(n=1)^nn) (Subscript[b, n ])(∫_0^K)Sin[((n Pi)/K)x]Sin[((j Pi)/K)x]dx

(∫_0^K)(Subscript[f, 1])(x)Sin[((j Pi)/K)x]dx = (Subscript[b, n ])(∫_0^K)Sin[((j Pi)/K)x]Sin[((j Pi)/K)x]dx

(∫_0^K)(Subscript[f, 1])(x)Sin[((j Pi)/K)x]dx = (Subscript[b, j ])*(K/2)

Then substitute y=L. This leads to the formula for (Subscript[a, n]).

Problem 2

The objective is to solve the PDE 

(((∂^2)u)/(∂(x^2)))+(((∂^2)u)/(∂(y^2)))=0  on  0≤x ≤K, 0 ≤y≤L,

subject to the conditions

u(x,0) = 0, u(x,L)=0  (call these DBCx) 

u(0,y) = (Subscript[g, 1])(y), u(K,y)=(Subscript[g, 2])(y)  (call these BCy)

Step 1. First ignore BCy  and solve

(((∂^2)u)/(∂(x^2)))+(((∂^2)u)/(∂(y^2)))=0  on  0≤x ≤K, 0 ≤y≤L,

subject to the conditions

u(x,0) = 0, u(x,L)=0  (call these DBCx)

Using SofV method we find few solutions of this problem:

Sin[((n Pi)/L)y] ((Sinh[((n Pi)/L)x])/(Sinh[((n Pi)/L)K]))

Sin[((n Pi)/L)y] ((Sinh[((n Pi)/L)(K-x)])/(Sinh[((n Pi)/L)K]))

Test these solutions:

(D[#,{x,2}]+D[#,{y,2}])&[Sin[((n Pi)/L)y](( Sinh[((n Pi)/L)x])/(Sinh[((n Pi)/L)K]))]

FullSimplify[(Sin[((n Pi)/L)y] ((Sinh[((n Pi)/L)x])/(Sinh[((n Pi)/L)K])))/.{y->{0,L}},n∈Integers]

(D[#,{x,2}]+D[#,{y,2}])&[Sin[((n Pi)/L)y] ((Sinh[((n Pi)/L)(K-x)])/(Sinh[((n Pi)/L)K]))]

FullSimplify[(Sin[((n Pi)/L)y] ((Sinh[((n Pi)/L)(K-x)])/(Sinh[((n Pi)/L)K])))/.{y->{0,L}},n∈Integers]

Step 2. Now that we have few solutions we form many solutions.  This is the FFM idea, which is commonly known as the superposition principle:

(∑_(n=1)^nn) (Subscript[c, n ])Sin[((n Pi)/L)y] ((Sinh[((n Pi)/L)x])/(Sinh[((n Pi)/L)K]))+(∑_(n=1)^nn) (Subscript[d, n ])Sin[((n Pi)/L)y] ((Sinh[((n Pi)/L)(K-x)])/(Sinh[((n Pi)/L)K]))

Now we choose (Subscript[c, n]) and (Subscript[d, n]) such that the above function satisfies BCy conditions. First substitute x=0. This leads to the formula for (Subscript[d, n]). Then substitute x=K. This leads to the formula for (Subscript[c, n]).

A symbolic implementation

Here are the given quantities

Clear[lK1,lL1,f11,f21,g11,g21,nn1];

nn1=15; 

lK1 = 1; lL1 = 1;

f11[x_]=4(x^2)(1-x); 

g21[y_]=4y(((1-y))^2); 

f21[x_]=0; 

g11[y_]=0;

Clear[aa1];

aa1[n_]=FullSimplify[(2/(lK1 ))Integrate[f21[x]Sin[((n Pi)/lK1)x],{x,0,lK1}],And[n∈Integers,n>0]]

Clear[bb1];

bb1[n_]=FullSimplify[(2/(lK1 ))Integrate[f11[x]Sin[((n Pi)/lK1)x],{x,0,lK1}],And[n∈Integers,n>0]]

Clear[cc1];

cc1[n_]=FullSimplify[(2/(lL1 ))Integrate[g21[y]Sin[((n Pi)/lL1)y],{y,0,lL1}],And[n∈Integers,n>0]]

Clear[dd1];

dd1[n_]=FullSimplify[(2/(lL1 ))Integrate[g11[y]Sin[((n Pi)/lL1)y],{y,0,lL1}],And[n∈Integers,n>0]]

The solution is

Clear[uu1]; 

uu1[x_,y_]=(∑_(n=1)^nn1) aa1[n]Sin[((n Pi)/lK1)x] ((Sinh[((n Pi)/lK1)y])/(Sinh[((n Pi)/lK1)lL1]))+(∑_(n=1)^nn1) bb1[n]Sin[((n Pi)/lK1)x] ((Sinh[((n Pi)/lK1)(lL1-y)])/(Sinh[((n Pi)/lK1)lL1]))+ (∑_(n=1)^nn1) cc1[n]Sin[((n Pi)/lL1)y] ((Sinh[((n Pi)/lL1)x])/(Sinh[((n Pi)/lL1)lK1]))+(∑_(n=1)^nn1) dd1[n] Sin[((n Pi)/lL1)y] ((Sinh[((n Pi)/lL1)(lK1-x)])/(Sinh[((n Pi)/lL1)lK1]));

How good is our approximation for the function f1[x] in the boundary conditions? Here is a visual answer.

Plot[{f11[x],uu1[x,0]},{x,0,1},PlotStyle->{{Blue,Thickness[0.01]},{Red,Thickness[0.005]}}]

Plot[uu1[x,lL],{x,0,1}]

Plot[uu1[0,x],{x,0,1}]

Plot[{g21[y],uu1[lK,y]},{y,0,1},PlotStyle->{{Blue,Thickness[0.01]},{Red,Thickness[0.005]}}]

Plot3D[N[uu1[x,y]],{x,0,1},{y,0,1},Mesh->Automatic]

DensityPlot[N[uu1[x,y]],{x,0,1},{y,0,1},Frame->True,PlotRange->{{0,1},{0,1}},ColorFunction->(RGBColor[1,1-#,1-#]&)]

A numerical implementation

Here are the given quantities

Clear[lK2,lL2,f12,f22,g12,g22,nn2];

nn2=45; 

lK2 = Pi; lL2 = Pi;

f12[x_]=(((Sin[x]))^(1/10));g22[y_]=0; f22[x_]=0; g12[y_]=0;

Plot[f12[x],{x,0,lK2},PlotRange->{0,1.2}]

Clear[aa2];

aa2[n_]=(2/lK2)Integrate[f22[x]Sin[((n Pi)/lK2)x],{x,0,lK2}]

Clear[bbl2];

bbl2=Chop[Table[(2/lK2)NIntegrate[f12[x]Sin[((n Pi)/lK2)x],{x,0,lK2},
Method->{Automatic},MaxRecursion->200,AccuracyGoal->12,PrecisionGoal->16],{n,1,nn2}]]

Clear[cc2];

cc2[n_]=(2/(lL2 Sinh[((n Pi)/lL2)lK2]))Integrate[g22[y]Sin[((n Pi)/lL2)y],{y,0,lL2}]

Clear[dd2];

dd2[n_]=(2/(lL2 ))Integrate[g12[y]Sin[((n Pi)/lL2)y],{y,0,lL2}]

The solution is

Clear[uu2];

uu2[x_,y_]=(∑_(n=1)^nn2) aa2[n]Sin[((n Pi)/lK2)x] ((Sinh[((n Pi)/lK2)y])/(Sinh[((n Pi)/lK2)lL2]))+(∑_(n=1)^nn2) bbl2[[n]]Sin[((n Pi)/lK2)x] ((Sinh[((n Pi)/lK2)(lL2-y)])/(Sinh[((n Pi)/lK2)lL2]))+ (∑_(n=1)^nn2) cc2[n]Sin[((n Pi)/lL2)y] ((Sinh[((n Pi)/lL2)x])/(Sinh[((n Pi)/lL2)lK2]))+(∑_(n=1)^nn2) dd2[n] Sin[((n Pi)/lL2)y] ((Sinh[((n Pi)/lL2)(lK2-x)])/(Sinh[((n Pi)/lL2)lK2]));

Plot[{f12[x],uu2[x,0]},{x,0,lK2},PlotRange->{0,1.2}]

Plot[uu2[x,lL2],{x,0,lK2}]

Plot[uu2[0,x],{x,0,lK2}]

Plot[uu2[lK,x],{x,0,lL2}]

Plot3D[N[uu2[x,y]],{x,0,lK2},{y,0,lL2}]

DensityPlot[N[uu2[x,y]],{x,0,lK2},{y,0,lL2},Frame->True,PlotRange->{{0,lK2},{0,lL2}},ColorFunction->(RGBColor[1,1-#,1-#]&)]

A symbolic implementation with a problem

Here are the given quantities

Clear[lK3,lL3,f13,f23,g13,g23,nn3];

nn3=20; 

lK3 = 1; lL3 = 1;

f13[x_]=3(x^2)+1; 

g23[y_]=4-8y (y-1); 

f23[x_]=4+8 x (x-1); 

g13[y_]=1+3(y^2);

Clear[aa3];

aa3[n_]=FullSimplify[(2/(lK3 ))Integrate[f23[x]Sin[((n Pi)/lK3)x],{x,0,lK3}],And[n∈Integers,n>0]]

Clear[bb3];

bb3[n_]=FullSimplify[(2/(lK3 ))Integrate[f13[x]Sin[((n Pi)/lK3)x],{x,0,lK3}],And[n∈Integers,n>0]]

Clear[cc3];

cc3[n_]=FullSimplify[(2/(lL3 ))Integrate[g23[y]Sin[((n Pi)/lL3)y],{y,0,lL3}],And[n∈Integers,n>0]]

Clear[dd3];

dd3[n_]=(2/(lL3 ))Integrate[g13[y]Sin[((n Pi)/lL3)y],{y,0,lL3}]

The solution is

Clear[uu3]; 

uu3[x_,y_]=(∑_(n=1)^nn3) aa3[n]Sin[((n Pi)/lK3)x] ((Sinh[((n Pi)/lK3)y])/(Sinh[((n Pi)/lK3)lL3]))+(∑_(n=1)^nn3) bb3[n]Sin[((n Pi)/lK3)x] ((Sinh[((n Pi)/lK3)(lL3-y)])/(Sinh[((n Pi)/lK3)lL3]))+ (∑_(n=1)^nn3) cc3[n]Sin[((n Pi)/lL3)y] ((Sinh[((n Pi)/lL3)x])/(Sinh[((n Pi)/lL3)lK3]))+(∑_(n=1)^nn3) dd3[n] Sin[((n Pi)/lL3)y] ((Sinh[((n Pi)/lL3)(lK3-x)])/(Sinh[((n Pi)/lL3)lK3]));

How good are the approximations? Here are visual answers:

Plot[{f13[x],uu3[x,0]},{x,0,lK3},PlotStyle->{{Blue,Thickness[0.01]},{Red,Thickness[0.005]}},PlotRange->All]

Plot[{g23[y],uu3[lK3,y]},{y,0,lL3},PlotStyle->{{Blue,Thickness[0.01]},{Red,Thickness[0.005]}},PlotRange->All]

Plot[{f23[x],uu3[x,lL3]},{x,0,lK3},PlotStyle->{{Blue,Thickness[0.01]},{Red,Thickness[0.005]}},PlotRange->All]

Plot[{g13[y],uu3[0,y]},{y,0,lL3},PlotStyle->{{Blue,Thickness[0.01]},{Red,Thickness[0.005]}},PlotRange->All]

Plot3D[N[uu3[x,y]],{x,0,lK3},{y,0,lL3},Mesh->Automatic,PlotRange->{0,6.5}]

DensityPlot[N[uu3[x,y]],{x,0,lK3},{y,0,lL3},Frame->True,PlotRange->{{0,lK3},{0,lL3}},ColorFunction->(RGBColor[1,1-#,1-#]&)]