laplace equation : spherical coordinates ' version 1.0'
Δ u = 0 , u = u(ρ,phi,z)
rem :Δ u (x,y,z) = 0 -> Δ u(ρ,φ,θ)= 0
peter.vlasschaert@gmail.com,21/12/2017 

ref:  schaum's outline series ,theory and problems
         of vector analysis and an introduction to tensor analysis
auth: Murray R.Spiegel 
rem: curvilinear coordinates. (ρ,φ,θ)
rem : Δu = grad(div(u))

range : θ : 0  =< φ =< pi
range : ρ : 0  =< r =<ρ
range : φ : 0  =< φ=< 2*pi
rem : (x,y,z) => (ρ,φ,θ)
rem : =>  e_r = vectorproduct(e_θ,e_φ)  , cyclic permutation alowed,right rot. 

restart
=====

kill(all)$

load 'packages'
=============

load("vector3d")$

menu:view :greek letters (make active)
=================================

assume(ρ>=0);

spherical coordinates : 3D 
===================

depends(u,[ρ,φ,θ]);

depends(x,[ρ,φ,θ]);
depends(y,[ρ,φ,θ]);
depends(z, [ρ,θ]);

p1a:x=ρ*cos(φ)*sin(θ);
p2a:y=ρ*sin(φ)*sin(θ);
p3a:z=ρ*cos(θ);

u ,depends on ρ,φ,θ
==================

depends(x,[ρ,φ,θ]);
depends(y,[ρ,φ,θ]);
depends(z, [ρ,θ]);

depends(u,[ρ,φ,θ]);

find :∂ρ/∂x

distance : ρ    ' point on the sphere'

p2:ρ=sqrt(part(p1a,1)^2+part(p2a,1)^2+part(p3a,1)^2);

p3:diff(part(p2,2),x);

rem : reduce 

p4:ev(p3,p1a,p2a,p3a);
p5:'diff(ρ,x) =trigsimp(p4);

find:∂θ/∂x   

find : tan(θ)

p6:trigsimp(sqrt(p1a^2+p2a^2)/p3a);

∂θ/∂x :

p7:diff(lhs(p6),x)=diff(rhs(p6),θ)*p8;
p9:ev(p7,p1a,p2a,p3a);
p10:trigsimp(solve(p9,p8))[1];

find:∂φ/∂x   

find :tan(φ)

p11:p2a/p1a;

∂φ/∂x :

p12:diff(lhs(p11),x)=diff(rhs(p11),φ)*p13;
p14:ev(p12,p1a,p2a,p3a);
p15:trigsimp(solve(p14,p13))[1];

find :∂ρ/∂y

p16:diff(part(p2,2),y);

rem : reduce 

p17:ev(p16,p1a,p2a,p3a);
p18:'diff(ρ,y) =trigsimp(p17);

find:∂θ/∂y   

p19:diff(lhs(p6),y)=diff(rhs(p6),θ)*p20;
p21:ev(p19,p1a,p2a,p3a);
p22:trigsimp(solve(p21,p20))[1];

find:∂φ/∂y  

p23:diff(lhs(p11),y)=diff(rhs(p11),φ)*p24;
p25:ev(p23,p1a,p2a,p3a);
p26:trigsimp(solve(p25,p24))[1];

find :∂ρ/∂z

p27:diff(part(p2,2),z);

rem : reduce 

p28:ev(p27,p1a,p2a,p3a);
p29:'diff(ρ,z) =trigsimp(p28);

∂θ/∂z :

p30:diff(lhs(p6),z)=diff(rhs(p6),θ)*p31;
p32:ev(p30,p1a,p2a,p3a);
p33:trigsimp(solve(p32,p31))[1];

gradient :  grad() =∂/∂x()*i+∂/∂y()*j+∂/∂z().k
i  =  unit vector  ' x ' = e_i
j  =  unit vector  ' y ' = e_j
k =  unit vector  ' z ' = e_k
|i| = |j| = |k| =1

how to find ,unit vector for spherical coordinate
i_ρ  =  unit vector  ' ρ '  = e_r
i_φ  =  unit vector  ' φ ' = e_φ
i_θ   =  unit vector  ' θ ' = e_θ
|i_ρ| = |i_φ| = |i_θ| =1

put all differential in a matrix : 'aa'
first row   = [ ∂ρ/∂x,∂θ/∂x,∂φ/∂x]
sec row    = [ ∂ρ/∂y,∂θ/∂y,∂φ/∂y]
thrid row  = [ ∂ρ/∂z,∂θ/∂z,∂φ/∂z]

aa: matrix(
 [rhs(p5),rhs(p10),rhs(p15)], 
 [rhs(p18),rhs(p22),rhs(p26)], 
 [rhs(p29),rhs(p33),0]
);

=>  last element  (3) in the first row (1).

aa[1,3];

relation : unit vector 'cartesian and spherical',found from the geometry of the figure above.

p34:e[x] = sin(θ)*cos(φ)*e[ρ]+cos(θ)*cos(φ)*e[θ]-sin(φ)*e[φ];
p35:e[y] = sin(θ)*sin(φ)*e[ρ]+cos(θ)*sin(φ)*e[θ]+cos(φ)*e[φ];
p36:e[z] = cos(θ)*e[ρ]-sin(θ)*e[θ]-0*e[φ];

find coef from above 'from e_x'

b11:diff(rhs(p34),e[ρ]);
b12:diff(rhs(p34),e[θ]);
b13:diff(rhs(p34),e[φ]);

find coef from above 'from e_y'

b21:diff(rhs(p35),e[ρ]);
b22:diff(rhs(p35),e[θ]);
b23:diff(rhs(p35),e[φ]);

find coef from above 'from e_z'

b31:diff(rhs(p36),e[ρ]);
b32:diff(rhs(p36),e[θ]);
b33:diff(rhs(p36),e[φ]);

matrix coef from
row 1 : e_x
row 2 : e_y
row 3  :e_z

cc: matrix(
 [b11,b12,b13], 
 [b21,b22,b23], 
 [b31,b32,b33]
);

find : gradient : grad() = e_ρ*1*∂/∂ρ+e_θ*1/ρ*∂/∂θ+e_φ*(1/ρsin(θ))∂/∂φ

dd: matrix(
 [aa[1,1]/b11,aa[1,2]/b12,aa[1,3]/b13], 
 [aa[2,1]/b21,aa[2,2]/b22,aa[2,3]/b23], 
 [aa[3,1]/b31,aa[3,2]/b32,0]
);

gradient(u) 'spherical coordinates'

p40:'grad(u) =[dd[1,1]*diff(u,ρ),dd[1,2]*diff(u,θ),dd[1,3]*diff(u,φ)];

laplace equation : ' spherical coordinates'

x - part : laplace equation.

∂(u)/∂x= ∂(u)/∂ρ*∂(ρ)/∂x+∂(u)/∂θ*∂(θ)/∂x+∂(u)/∂φ*∂(φ)/∂x

p41:'diff(u,x) =aa[1,1]*diff(u,ρ)+aa[1,2]*diff(u,θ)+aa[1,3]*diff(u,φ);

u  =>  replace    'diff(u,x,1)

find :∂(∂u/∂x) /∂x   ' x- part   laplace equation '
diff(u,x,2) =∂/∂x(∂u/∂x) =∂/∂ρ(∂u/∂x )∂ρ/∂x +∂/∂θ(∂u/∂x )∂θ/∂x+∂/∂φ(∂u/∂x )∂φ/∂x

p42:'diff(u,x,2) = ev(rhs(p41),u=rhs(p41),diff)$
p42a:trigsimp(ratsimp(p42));

y - part : laplace equation.

∂(u)/∂y= ∂(u)/∂ρ*∂(ρ)/∂y+∂(u)/∂θ*∂(θ)/∂y+∂(u)/∂φ*∂(φ)/∂y

p43:'diff(u,y) =aa[2,1]*diff(u,ρ)+aa[2,2]*diff(u,θ)+aa[2,3]*diff(u,φ);

u  =>  replace    'diff(u,y,1)

find :∂(∂u/∂y) /∂y   ' y- part   laplace equation '
diff(u,y,2) =∂/∂x(∂u/∂y) =∂/∂ρ(∂u/∂y )∂ρ/∂y +∂/∂θ(∂u/∂y )∂θ/∂y+∂/∂φ(∂u/∂y)∂φ/∂y

p44:'diff(u,y,2) = ev(rhs(p43),u=rhs(p43),diff)$
p44a: trigsimp(ratsimp(p44));

z - part : laplace equation.

∂(u)/∂z= ∂(u)/∂ρ*∂(ρ)/∂z+∂(u)/∂θ*∂(θ)/∂z+∂(u)/∂φ*∂(φ)/∂z

p45:'diff(u,z) =aa[3,1]*diff(u,ρ)+aa[3,2]*diff(u,θ)+aa[3,3]*diff(u,φ);

u  =>  replace    'diff(u,z,1)

find :∂(∂u/∂z) /∂z   ' z- part   laplace equation '
diff(u,z,2) =∂/∂z(∂u/∂z) =∂/∂ρ(∂u/∂z )∂ρ/∂z +∂/∂θ(∂u/∂z )∂θ/∂z+∂/∂φ(∂u/∂z)∂φ/∂z

p46:'diff(u,z,2) = ev(rhs(p45),u=rhs(p45),diff);
p46a: trigsimp(ratsimp(p46));

p48:ratsimp(p42a+p44a+p46a);
p49:trigsimp(p48);

'grad(div(u)) = lhs;

                                                 ===============================================

rem:

p50:rhs(p49);
p50a :part(p50,2);
p51:part(p50,1,1)/p50a;
p52:part(p50,1,2)/p50a;
p53:part(p50,1,3)/p50a;
p54:part(p50,1,4)/p50a;
p55:part(p50,1,5)/p50a;

 ============================================================================================================

find :'diff(u,ρ,1) =>'diff(u,ρ,2) ?,p51,p52
?

check  ?

p56:1/ρ^2*diff(ρ^2*diff(u,ρ),ρ)=p51+p52;
p57:ratsimp(p56);
p58:is(p57);

find :'diff(u,θ,1) =>'diff(u,θ,2) ?,p53,p54
?

check  ?

p59:1/(ρ^2*sin(θ))*diff(sin(θ)*diff(u,θ),θ)=p53+p54;
p60:ratsimp(p59);
p58:is(p60);

        ==============================================================================

laplace equation 3D ' spherical coordinates ' :    Δ u = 0   : Δ u (x,y,z)= 0  -> Δ u (ρ,φ,θ)= 0 ,see below

'Δu= p55+'1/(ρ^2*sin(θ))*diff(sin(θ)*diff(u,θ),θ)+'1/ρ^2*diff(ρ^2*diff(u,ρ),ρ);

Δu = 1/(ρ^2*sin(θ))*diff(sin(θ)*diff(u,θ),θ)+1/ρ^2*diff(ρ^2*diff(u,ρ),ρ)+'diff(u,φ,2)/(sin(θ)^2*ρ^2)=0

  ============================================================================================================
============================================================================================================