There is one 3D-diffusion IBVP for sphere (Dirichlet problem with zero conditions on the surface of sphere).
The equation has the following view:
$$\frac{\partial u}{\partial t}=\operatorname{div}\left(A(r,\varphi,\psi)\nabla u(r,\varphi,\psi,t)\right)+f(r,\varphi,\psi,t) ,$$ where the matrix $A$ is skew with variable coefficients:
$$\left(\matrix{1 \quad a \quad c\\-a \quad 1 \quad b\\-c\quad -b \quad 1 }\right) $$
I've written a numerical scheme for this problem, but the program works strangely: if all non-diagonal coefficients are zeros, everything is cool, but when at least one coefficient is non-zero, my scheme begins to show great error in neighbourhood of points $\varphi=0$ and $\varphi=\pi$. I've checked program many times, moreover, I've shown it to expert --he said, that program code is correct. He also said that this problem with big errors is connected with singularities in coefficients of terms $\frac{\partial u}{\partial r},\frac{\partial u}{\partial \varphi},\frac{\partial u}{\partial \psi}$. Has anyone faced with similar problems? How can I avoid those singularities?
Here is one example of error. The surfaces' equation is $u(r,\varphi,\psi,t)=\cos(t)(1-r^2)$
[upd] Here is my numerical scheme: $$\frac{u^{n+1}_{i,j,k}-u^{n}_{i,j,k}}{\Delta t}=AA_{i,j,k}\frac{u^{n+1}_{i+1,j,k}-u^{n+1}_{i-1,j,k}}{2\Delta r}+BB_{i,j,k}\frac{u^{n+1}_{i,j+1,k}-u^{n+1}_{i,j-1,k}}{2\Delta \varphi}+CC_{i,j,k}\frac{u^{n+1}_{i,j,k+1}-u^{n+1}_{i,k,k-1}}{2\Delta \psi}+\frac{u^{n+1}_{i+1,j,k}-2u^{n+1}_{i,j,k}+u^{n+1}_{i-1,j,k}}{\Delta r^2}+\frac{u^{n+1}_{i,j+1,k}-2u^{n+1}_{i,j,k}+u^{n+1}_{i,j-1,k}}{r_{i}^{2}\Delta \varphi^2}+\frac{u^{n+1}_{i,j,k+1}-2u^{n+1}_{i,j,k}+u^{n+1}_{i,j,k-1}}{r_{i}^{2}\sin^{2}\varphi_{j}\Delta \psi^2}+f_{i,j,k} $$ Coefficients $AA,BB,CC$ are the following: $$AA_{i,j,k}=\frac{2}{r_{i}}+a_{i,j,k}\frac{\cot \varphi_{j}}{r_{i}}+\frac{1}{r_{i}}\left(\frac{\partial a}{\partial \varphi}\right)_{i,j,k}-\frac{1}{r_{i}\sin\varphi_{j}}\left(\frac{\partial c}{\partial \psi}\right)_{i,j,k} $$ $$BB_{i,j,k}=-\frac{a_{i,j,k}}{r_{i}^2}+\frac{\cot\varphi_{j}}{r_{i}^2}-\frac{1}{r_{i}}\left(\frac{\partial a}{\partial r}\right)_{i,j,k}+\frac{1}{r_{i}^2\sin\varphi_{j}}\left(\frac{\partial b}{\partial \psi}\right)_{i,j,k} $$ $$CC_{i,j,k}=\frac{c_{i,j,k}}{r_{i}^2\sin\varphi_{j}}+\frac{1}{r_{i}\sin\varphi_{j}}\left(\frac{\partial c}{\partial r}\right)_{i,j,k}-\frac{1}{r_{i}^2\sin\varphi_{j}}\left(\frac{\partial b}{\partial \varphi}\right)_{i,j,k} $$