I am trying to solve the 2D advection-diffusion equation in cylindrical coordinates:
$$ \frac{\partial c}{\partial t} = D\left(\frac{\partial^2 c}{\partial r^2} + \frac{1}{r}\frac{\partial c}{\partial r} + \frac{\partial^2 c}{\partial z^2}\right) - u_z\frac{\partial c}{\partial z} $$
However, in my situation, both $D$ and $u_z$ vary piecewise in $r$, so applying $Z = D\frac{\partial^2}{\partial z^2} - u_z\frac{\partial}{\partial z}$ to $c$ becomes an operator $\mathcal{\hat{Z}}$ that not only applies the derivatives but also needs to split $c$ into the different regions based on $r$ that have the different values of $u_z$ and $D$. If $C = [c_{ij}]$, with $i$ and $j$ corresponding to $r$ and $z$, respectively, for a small system with three distinct zones:
$$ \mathcal{\hat{Z}}(C)=\begin{bmatrix}1&0&0\\0&0&0\\0&0&0\end{bmatrix}CZ_1 + \begin{bmatrix}0&0&0\\0&1&0\\0&0&0\end{bmatrix}CZ_2 + \begin{bmatrix}0&0&0\\0&0&0\\0&0&1\end{bmatrix}CZ_3 $$
where (my actual system is much bigger, so it's not just one row in each section). This is fine to treat explicitly but becomes a problem when I'm solving this implicitly, since $\mathcal{\hat{Z}}$ is now non-invertible. However, I'd like to keep using an implicit scheme, at least for the derivatives in $r$. Is it valid to set up a finite difference scheme that is fully implicit in $r$ and fully explicit in $z$? Something like:
$$ \frac{c_{ij}^{n+1}-c_{ij}^n}{\Delta t} = D\left(\frac{c_{i-1,j}^{n+1}-2c_{ij}^{n+1}+c_{i+1,j}^{n+1}}{\Delta r^2}+\frac{1}{r}\frac{c_{i+1,j}^{n+1}-c_{i-1,j}^{n+1}}{2\Delta r} + \frac{c_{i,j-1}^n-2c_{ij}^n+c_{i,j+1}^n}{\Delta z^2}\right) - u_z\left(\frac{c_{i,j+1}^n-c_{i,j-1}^n}{2\Delta z}\right) $$
which is something like one half-step of an ADI scheme (but of course taking into account $\frac{\partial{D}}{\partial{r}}$, etc, at the transition zones). I know each half-step is supposed to be conditionally stable, but is this generally a reasonable scheme to use without the other half-step that would treat $z$ implicitly and $r$ explicitly? Is there a better way I can be treating $z$ implicitly that is actually able to be implemented, or should I give up on using any implicit method and just do the whole thing explicitly?
