Heat Equation on a surface sphere using NDSolve?

[georges...@gmail.com]

Hello Folks,

I am new to Mathematica, but it seems to be the most suitable tool for my issue.

I have to solve a equation which is similar to heat equation on a sphere surface. The exact expression in a spherical coordinates system is (in LaTeX):

\frac{\partial{h}}{\partial t} = \frac{1}{n} \left[ \frac{1}{r^2 \sin^2{\phi}} \frac{\partial}{\partial \theta} \left( K \frac{\partial h}{\partial \theta} \right) + \frac{1}{r^2 \sin \phi} \frac{\partial}{\partial \phi} \left( K \sin \phi \frac{\partial h}{\partial \phi}\right) + s(\theta,\phi,t)\right]

K is kind of a diffusion coeff.
s a source term (varying in time and space)
and n a constant (although it might depends on \phi and \theta).
\theta and \phi are defined here in a mathematical way.
r the radius (constant as we deal with a perfect sphere).
(no variation regarding the radius, we only look at the surface of the sphere).

I found that in Mathematica, and in a Cartesian frame, such equation (somewhat not exactly the same but similar) would be:

N = 1
K = 34
myfunc = NDSolve[{D[h[x, y, t], t] ==
1/N*(D[h[x, y, t], x, x] + D[h[x, y, t], y, y]), h[x, y, 0] == 0,
h[0, y, t] == K*t, h[9, y, t] == K*t, h[x, 0, t] == K*t,
h[x, 9, t] == K*t}, h, {x, 0, 9}, {y, 0, 9}, {t, 0, 6}];

I'm kind of stuck right now, because I do not know how can I switch to the spherical frame. I would appreciate any ideas/suggestions.

Thanks a lot and sorry for bothering you guys with this.

G.

[George]

I found something close to what I need, but in another system:

PDE := n*(diff(h(theta, phi, t), t)) = (diff(K*(diff(h(theta, phi, t),
theta)), theta))/(r^2*sin(phi)^2)+(diff(K*sin(phi)*(diff(h(theta, phi,
t), phi)), phi))/(r^2*sin(phi)) ;

sol := pdsolve(PDE, h(theta, phi, t));

but I really want to do it with Mathematica, as my other system's license will expire soon and only Mathematica will be available in my place.

Any ideas how I can convert this into Mathematica?

thanks.

[Roland Franzius]

Looks a bit confusing.

The Laplacian on the unit sphere with metrics

ds^2 = dtheta^2 + sin^2 theta dphi^2

is with s -> sin theta

Laplacian = s^-1 d_theta s d_theta + s^-2 d_phi d_phi

and the general heat equation with production and destruction terms is

d_t F(t,theta,phi) = diffconst Laplacian( F(t,theta,phi) )
- F(t,theta,phi)/lifetime + productionrate

to avoid the sine switch to cylinder coordinates

u = Cos theta, du = -Sqrt(1-u^2) dtheta, d_theta = - Sqrt(-1-u^2) d_u

with

Laplacian = d_u (1-u^2) d_u = d_u d_u - 2 u d_u + 1/(1-u^2) d_phi d_phi


equns=Simplify@
{
D[F[t,u,p],t] == K D[(1-u^2)D[F[t,u,p],u]],u] +
K/(1-u^2) Derivative[0,0,2][F][t,u,phi]
- a F[t,u,p] + b,
F[t,u,2 Pi]==F[t,u,0],
Derivative[F][t,-1.0,p]==0,
Derivative[F][t,1.0,p]==0,
F[0,u,p]==F0[u,p]}

NDSolve[equns,{t, 0.0, tmax}, {u,-1.0, 1.0}, { p, 0.0, 2.0 N[Pi]}]


The polar points u=+-1 are replaced by reflecting boundary conditions.
This guarantees smothness of the solution at the poles and conservation
of probability.

In case the start distribution F0 is rotational symmetric around some
axis you can choose this direction as north and discard the variable phi
completely.

--

Roland Franzius