• 1 Example
  • Problem
  • Solving method
  • Solving problem
  • Solution
  • Checking the Solution
  • Solution in Cartesian coordinates
• 2 Example
  • Problem
  • Solving problem
  • Solution
  • Checking the Solution
  • Note
• 3 Example
  • Problem
  • Solving problem
  • Solution
  • Checking the Solution
• 4 Example
  • Problem
  • Solving problem
  • Checking the Solution
• 5 Example
  • Problem
  • Solving problem
  • Checking the Solution
• 6 Example
  • problem
  • theory
  • procedures
  • examples
• 7 Example
  • problem
  • solving problem
  • Checking the Solution
• 8 Example
  • problem
  • solving problem
  • Checking the Solution

Boundary value problems for three-dimensional elliptic equations

by Aleksas Domarkas

Vilnius University, Faculty of Mathematics and Informatics,

Naugarduko 24, Vilnius, Lithuania

aleksas@ieva.mif.vu.lt

NOTE: In this session we find  solutions of boundary value problems for three-dimensional

Laplace's , Poisson or Helmholtz equations in regions r<R,  r>R and R1<r<R2.

1 Example

Problem

Solve   boundary-value problem  

                             ,  r1< r < r2,

                                 u| r=r 1 = u1, u| r=r 2  = u2

here  ;     -- constants.

Solving method

We consider problem   

 ,  r1< r < r2,

u| r=r1 =0, u| r=r2 =0

and problem  

,  r1< r < r2,

  w| r=r1 =u1, w| r=r2 =u2.

Then u+w is solution of our problem.  Solutions of these problems are sought in the form

u = v(r)*cos(phi)*sin(theta)*cos(theta),

Solving problem

>

restart;

>	with(linalg,laplacian):

>	l:=expand(laplacian(u(r,theta,phi),[r,theta,phi], coords=spherical))-x*z;

>	tr:={x=r*cos(phi)*sin(theta),y=r*sin(phi)*sin(theta),z=r*cos(theta)};

>	l := subs(tr,l);

>	subs(u(r,theta,phi)=v(r)*cos(phi)*sin(theta)*cos(theta),l):

>	l1 := simplify(%/(cos(phi)*sin(theta)*cos(theta)));

>	spr := dsolve(l1,{v(r)});

>	solve({subs(r=r1,rhs(spr)),subs(r=r2,rhs(spr))},{_C1,_C2});

>	v(r) := subs(%,rhs(spr));

>	subs(u(r,theta,phi)=w(r),laplacian(u(r,theta,phi),[r,theta,phi], coords=spherical));

>	dsolve(%,w(r));

>	solve({subs(r=r1,rhs(%))=u1,subs(r=r2,rhs(%))=u2},{_C1,_C2});

>	w(r) := subs(%,rhs(%%));

>	sol:=v(r)*cos(phi)*sin(theta)*cos(theta)+w(r):

Solution

>

u=sol;

Checking the Solution

>	simplify(laplacian(sol,[r,theta,phi], coords=spherical))-subs(tr,x*z);

>	simplify(subs(r=r1,sol));

>	simplify(subs(r=r2,sol));

Solution in Cartesian coordinates

>	solve(tr,{r,theta,phi}):

>	allvalues(%):

>	subs(%[1],sol):

>	sol1:=simplify(%):

This is solution in Cartesian coordinates.

Checking the Solution:

>	simplify(laplacian(sol1,[x,y,z])-x*z);

>	subs(tr,r=r1,sol1):simplify(%,radical,symbolic):combine(%,trig);

>	subs(tr,r=r2,sol1):simplify(%,radical,symbolic):combine(%,trig);

>

2 Example

Problem

>	restart;with(linalg,laplacian):with(orthopoly,P):

Solve  three-dimensional interior boundary-value problem

,  r < R,

u| r=R  = g .

Here r=sqrt(x^2+y^2+z^2), u=u(x,y,z) or u=u(r,phi,theta).

The function g may be polynomial  with variables z and x^2+y^2 or with cos(theta) . 

>	g:=randpoly([z,x^2+y^2],degree=rand(2..5)(),terms=3);

or

>	#g:=randpoly([cos(theta)]);

>

Solving problem

>	n:=degree(g);

>	tr:={x=r*cos(phi)*sin(theta),y=r*sin(phi)*sin(theta),z=r*cos(theta)};

>	subs(tr,r=R,g);

>	f:=simplify(expand(%));

>	us:=sum(A[k]*r^k*P(k,cos(theta)),k=0..n);

>	k:='k': ug:=sum(A[k]*r^k*P[k],k=0..n);

>	koef:=solve(identity(subs(r=R,us)=f,theta),{seq(A[k],k=0..n)});

>	subs(koef,us):sol:=expand(%);

>	sol1:=subs(koef,ug);

In  Cartesian coordinates:

>	solve(tr,{r,theta,phi}):allvalues(%): subs(%[1],sol): solc:=simplify(%);

Solution

Problem:

>	print(laplacian(u(x,y,z),[x,y,z])=0,` if  r<R`);u=g,`  if  r`=R;

Solution:

>

sol;

or

>

sol1;

where  is the nth Legendre (spherical) polynomial.

In  Cartesian coordinates:

>

solc;

>

Checking the Solution

>	simplify(laplacian(sol,[r,theta,phi],coords=spherical));

>	simplify(laplacian(solc,[x,y,z]));

>	simplify(subs(tr,solc-g),{x^2+y^2+z^2=R^2,r=R}):

>	combine(expand(%));

>

sol1;

>	subs(seq(P[k]=P(k,cos(theta)),k=0..n),sol1):

>	simplify(%-sol);

Then sol=sol1 .

>

Note

In similar way can solve boundary-value problems in domains r>R and  R1<r<R2.

3 Example

Problem

>	restart;with(linalg,laplacian):

Solve  three-dimensional interior boundary-value problem

,  r < R,

   | r=R  = g .

Here ,  or .

>	g:=cos(theta);

Solving problem

>	with(linalg,laplacian);

>	assume(a>0);

>	eq:=laplacian(u(r,theta,phi),[r,theta,phi],coords=spherical)+ a^2*u(r,theta,phi)=0:

>	subs(u(r,theta,phi)=v(r)*cos(theta),eq):

>	l:=combine(%/cos(theta));

>	dsolve({l,v(0)<infinity,D(v)(R)=1},v(r));

>	sol:=subs(%,v(r)*cos(theta));

Solution

Problem:

>	Delta*u+a^2*u=0,` if  r<R`;u=g,`  if  r`=R;

Solution:

>

u=sol;

>

Checking the Solution

>	simplify(subs(u(r,theta,phi)=sol,eq));

>	diff(sol,r):

>	simplify(subs(r=R,%));

>

4 Example

Problem

u| r=1 =f .

We use package orbitals  from share libary  R5 or R4.

>	restart;with(share):readshare(orbitals,science): with(orbitals):

>	assume(theta,real);assume(phi,real);

See ?share and ?share,contents for information about the share library

>	#f:=randpoly([x,y,z],terms=3);

>

f:=x*y*z;

Solving problem

>	tr:={x=r*cos(phi)*sin(theta),y=r*sin(phi)*sin(theta),z=r*cos(theta)}:

>	fs:=subs(tr,r=1,f);

>	d:=degree(f);R:=1:

We search solution in the form:

>	u=Sum(Sum('a(i,j)*(r/R)^i*Y(i,j)','j'=-'i'..'i'),'i'=0..d);

Here     is    from orbitals .

>	value(rhs(%));

>	Y:=(l,m)->ComplexSurfaceHarmonic(l,m,theta,phi):

We compute coefficients:

>	a:=(i,j)->Int(Int(fs*conjugate(Y(i,j))*sin(theta),theta=0..Pi),phi=0..2*Pi):

>	convert(%%%,trig):value(%):

>

evalc(%):

>	expand(%,trig):

>

%;

Solution:

>	sol:=map(simplify,collect(%,r));

>

Checking the Solution

>	linalg[laplacian](sol,[r,theta,phi],coords=spherical):combine(%,trig);

>	simplify(fs-subs(r=R,sol)):combine(%,trig);

>

5 Example

Problem

=0, r<1,

u| r=1 =  .

Solving problem

>

restart;

>

Order:=2;

 We search solution in the form:

>	Sum(r^n*(A[0,n]*LP(n,0)+ Sum((A[k,n]*cos(k*phi)+B[k,n]*sin(k*phi))*LP(n,k),k=1..n)),n=0..Order);

Here LP(n,k)=P(n, k, cos(theta)),  P(n, k, x)  is associated  Legendre function. Function LP we define later.

>	so:=value(%);

>	fc:=combine(subs(r=1,so));

>	f:=cos(2*phi+Pi/3)*sin(theta)^2;

>	applyop(expand,1,f);

>	applyop(combine,1,%);

>	collect(%,[cos(2*phi),sin(2*phi)]);

>	dp:=collect(%,[sin(phi),cos(phi)]);

>	sk:=collect(fc-dp,[sin(phi),cos(phi),cos(2*phi),sin(2*phi)]);

>	LP:=proc(n::nonnegint,k::nonnegint) local i; orthopoly[P](n,t); for i to k do diff(%,t) end do; sin(theta)^k*subs(t=cos(theta),%); RETURN(%); end;

>

sk;

Next we compute coefficients:

>

K:={}:

>	coeff(sk,sin(phi));

>	indets(%,indexed);

>	solve(identity(%%,theta),%);K:=K union %;

>	coeff(sk,cos(phi));

>	indets(%,indexed);

>	solve(identity(%%,theta),%);K:=K union %;

>	coeff(sk,sin(2*phi));

>	indets(%,indexed);

>	solve(identity(%%,theta),%);K:=K union %;

>	coeff(sk,cos(2*phi));

>	indets(%,indexed);

>	solve(identity(%%,theta),%);K:=K union %;

>	subs(sin(phi)=0,cos(phi)=0,cos(2*phi)=0,sin(2*phi)=0,sk);

>	indets(%,indexed);

>	solve(identity(%%,theta),%);K:=K union %;

>	solution:=subs(K,so);

Other form of the solution:

>	sol:=r^2*cos(2*phi+Pi/3)*sin(theta)^2;

>	expand(solution-sol);

>

Checking the Solution

>	combine(linalg[laplacian](sol,[r,theta,phi],coords=spherical));

>	f-subs(r=1,sol);

>

6 Example

problem

| r=R =

theory

where

,     is associated  Legendre function.

procedures

>

restart;

Procedure for computing LP(n,k)=P(n, k, cos(theta)), here  P(n, k, x)  is associated  Legendre function.

>	LP:=proc(n::nonnegint,k::nonnegint) local i; orthopoly[P](n,t); for i to k do diff(%,t) end do; sin(theta)^k*subs(t=cos(theta),%); RETURN(%); end:

Procedure for solving problem   ,    | r=R =  :

>

lap3di:=proc(f,R) local A,B,so,lp,N; global LP;N:=degree(f); Sum(r^n*(A(0,n)*lp(n,0)/2+ Sum((A(k,n)*cos(k*phi)+B(k,n)*sin(k*phi))*lp(n,k),k=1..n)),n=0..N); so:=value(%): lp:=LP; A:=(k,m)->(2*m+1)*(m-k)!/2/(m+k)!/Pi/R^m* int(int(f*LP(m,k)*cos(k*phi)*sin(theta),phi=-Pi..Pi),theta=0..Pi); B:=(k,m)->(2*m+1)*(m-k)!/2/(m+k)!/Pi/R^m* int(int(f*LP(m,k)*sin(k*phi)*sin(theta),phi=-Pi..Pi),theta=0..Pi); RETURN(map(factor,value(so))); end:

examples

Example 1 ([vlad],16.20.3.)

>	f:=sin(theta)*(sin(phi)+sin(theta));R:=1;

>	sol:=lap3di(f,R);

Checking the Solution:

>	simplify(linalg[laplacian](sol,[r,theta,phi],coords=spherical));

>	simplify(expand(f-subs(r=R,sol)),trig);

Example 2 ([vlad],16.21.2.)

>	f:=sin(3*phi+Pi/4)*sin(theta)^3;R:='R':

>	sol:=lap3di(f,R);

Checking the Solution:

>	simplify(linalg[laplacian](sol,[r,theta,phi],coords=spherical));

>	simplify(expand(f-subs(r=R,sol)),trig);

Example 3

>	f:=randpoly([x,y,z],degree=2);R:=3;

>	tr:={x=r*cos(phi)*sin(theta),y=r*sin(phi)*sin(theta),z=r*cos(theta)}:

>	f:=subs(tr,r=R,f);

>	sol:=lap3di(f,R);

Checking the Solution:

>	simplify(linalg[laplacian](sol,[r,theta,phi],coords=spherical));

>	simplify(expand(f-subs(r=R,sol)),trig);

>

7 Example

problem

[vlad] 16.21.4.

( + )| r=R =

solving problem

>

restart;

>	f:=sin(theta)^2*(sqrt(2)*cos(2*phi+Pi/4)+2*cos(phi)^2);

We express   in Fourier series.

>	select(has,f,phi);

>	Order:=degree(%);

>	fseries:=proc(f)

>	local k,S;

>	S:=int(f,phi=-Pi..Pi)/2/Pi;

>	for k from 1 to Order do S:=S,

>	simplify(int(f*cos(k*phi),phi=-Pi..Pi)/Pi)*cos(k*phi),

>	simplify(int(f*sin(k*phi),phi=-Pi..Pi)/Pi)*sin(k*phi);od;

>	RETURN(sum(S[m],m=1..nops([S])));

>

end:

>

>	ff:=fseries(f);

>	simplify(expand(ff-f),trig);

Then ff=f.

 We search solution in the form:

>	Sum(r^n*(A[0,n]*LP(n,0)+ Sum((A[k,n]*cos(k*phi)+B[k,n]*sin(k*phi))*LP(n,k),k=1..n)),n=0..Order);

>	so:=value(%);

We define  LP(n,k)=P(n, k, cos(theta)), here  P(n, k, x)  is associated  Legendre function.

>	LP:=proc(n::nonnegint,k::nonnegint) local i; orthopoly[P](n,t); for i to k do diff(%,t) end do; sin(theta)^k*subs(t=cos(theta),%); RETURN(%); end:

>	subs(r=R,so+diff(so,r))-ff:

>	T:=collect(%,[seq(cos(k*phi),k=1..Order),seq(sin(k*phi),k=1..Order)]);

We need from identity T=0 compute coefficients.

>

K:={}:

>	for k to Order do

>	coeff(T,sin(k*phi));

>	indets(%,indexed);

>	solve(identity(%%=0,theta),%);

>	K:=K union %; od:

>	for k  to Order do

>	coeff(T,cos(k*phi));

>	indets(%,indexed);

>	solve(identity(%%=0,theta),%);

>	K:=K union %; od:

>	subs(seq(cos(k*phi)=0,k=1..Order),seq(sin(k*phi)=0,k=1..Order),T):

>	indets(%,indexed):

>	solve(identity(%%=0,theta),%):

>	K:=K union %;

Solution:

>	sol:=subs(K,so);

>

Checking the Solution

>	simplify(linalg[laplacian](sol,[r,theta,phi],coords=spherical));

>	simplify(expand(f-subs(r=R,sol+diff(sol,r))),trig);

>

>

8 Example

problem

[vlad] 16.24

| r=R1 = , | r=R2 =

solving problem

>

restart;

>	R1:=1;R2:=2;

>

#f1:=sin(theta)*sin(phi);f2:=0;#1 #f1:=3*sin(theta)^2*sin(2*phi);f2:=3*cos(theta);#2 #f1:=7*sin(theta)*cos(phi);f2:=7*cos(theta);#3 #f1:=sin(theta)^2*(3-sin(2*phi));f2:=4*f1;#4 #f1:=12*sin(theta)*cos(theta/2)^2*cos(phi);f2:=0;#5 #f1:=sin(2*phi)*sin(theta)^2;f2:=cos(2*phi)*sin(theta)^2;#6 #f1:=cos(phi)*sin(2*theta);f2:=sin(phi)*sin(2*theta);#7 #f1:=31*sin(2*theta)*sin(phi);f2:=31*sin(theta)^2*cos(2*phi);#8 f1:=cos(theta);f2:=cos(phi)*(12*sin(theta)-15*sin(theta)^3);#9

We express  1 and f2 in Fourier series.

>	Order:=max(degree(expand(f1)),degree(expand(f2)));

>	fseries:=proc(f)

>	local k,S;

>	S:=int(f,phi=-Pi..Pi)/2/Pi;

>	for k from 1 to Order do S:=S,

>	simplify(int(f*cos(k*phi),phi=-Pi..Pi)/Pi)*cos(k*phi),

>	simplify(int(f*sin(k*phi),phi=-Pi..Pi)/Pi)*sin(k*phi); end do;

>	RETURN(sum(S[m],m=1..nops([S])));

>

end proc:

>	ff1:=fseries(f1);

>	ff2:=fseries(f2);

 We search solution in the form:

>	Sum(r^n*(A[0,n]*LP(n,0)+ Sum((A[k,n]*cos(k*phi)+B[k,n]*sin(k*phi))*LP(n,k),k=1..n)),n=0..Order)+ Sum(r^(-n-1)*(C[0,n]*LP(n,0)+ Sum((C[k,n]*cos(k*phi)+D[k,n]*sin(k*phi))*LP(n,k),k=1..n)),n=0..Order);

>	so:=value(%):

We define  LP(n,k)=P(n, k, cos(theta)), here  P(n, k, x)  is associated  Legendre function.

>	LP:=proc(n::nonnegint,k::nonnegint) local i; orthopoly[P](n,t); for i to k do diff(%,t) end do; sin(theta)^k*subs(t=cos(theta),%); RETURN(%); end:

>	subs(r=R1,so-ff1):

>	T1:=collect(%,[seq(cos(k*phi),k=1..Order),seq(sin(k*phi),k=1..Order)]):

>	subs(r=R2,so-ff2):

>	T2:=collect(%,[seq(cos(k*phi),k=1..Order),seq(sin(k*phi),k=1..Order)]):

Next we compute coefficients  from identities  T1=0 and T2=0  .

>

K:={}:

>	for k to Order do

>	coeff(T1,sin(k*phi));

>	indets(%,indexed);

>	solve(identity(%%=0,theta),%);

>	K:=K union %; od:

>	for k  to Order do

>	coeff(T1,cos(k*phi));

>	indets(%,indexed);

>	solve(identity(%%=0,theta),%);

>	K:=K union %; od:

>	subs(seq(cos(k*phi)=0,k=1..Order),seq(sin(k*phi)=0,k=1..Order),T1):

>	indets(%,indexed):

>	solve(identity(%%=0,theta),%):

>	K1:=K union %:

>	T:=subs(K1,T2):

>

K:={}:

>	for k to Order do

>	coeff(T,sin(k*phi));

>	indets(%,indexed);

>	solve(identity(%%=0,theta),%);

>	K:=K union %; od:

>	for k  to Order do

>	coeff(T,cos(k*phi));

>	indets(%,indexed);

>	solve(identity(%%=0,theta),%);

>	K:=K union %; od:

>	subs(seq(cos(k*phi)=0,k=1..Order),seq(sin(k*phi)=0,k=1..Order),T):

>	indets(%,indexed):

>	solve(identity(%%=0,theta),%):

>	K:=K union %;

Solution:

>	sol:=map(normal,subs(K1,K,so));

>

Checking the Solution

>	simplify(linalg[laplacian](sol,[r,theta,phi],coords=spherical));

>	combine(expand(f1-subs(r=R1,sol)),trig);

>	combine(expand(f2-subs(r=R2,sol)),trig);

>

While every effort has been made to validate the solutions in this worksheet, Waterloo Maple Inc. and the contributors are not responsible for any errors contained and are not liable for any damages resulting from the use of this material.

Back to contents