I am trying to solve a particular form of the Euler / Navier-Stokes equations in 1D, with very strong and non-linear diffusion coefficients.
My system of equations is
\begin{cases} \frac{\partial n}{\partial t} + \frac{\partial}{\partial x}\left(nv \right ) = S_n\\ \frac{\partial \left( n v\right)}{\partial t} + \frac{\partial}{\partial x}\left(nv^2 + nT \right ) = S_M + \frac{\partial}{\partial x}\left(\nu(T)\frac{\partial v}{\partial x} \right)\\ \frac{\partial}{\partial t}\left(\frac{3}{2}nT + \frac{1}{2}n v^2 \right ) + \frac{\partial}{\partial x}\left(\frac{5}{2}nTv + \frac{1}{2}n v^3 \right ) - \frac{\partial}{\partial x}\left(D(T)\frac{\partial T}{\partial x} \right) - \frac{\partial}{\partial x}\left(v\nu(T)\frac{\partial v}{\partial x} \right) = S_E \end{cases} where $\nu(T)$ and $D(T)$ are strongly nonlinear coefficients. For instance, in some versions of my model, $\nu(T), D(T) \propto T^3$. This system is associated with a set of boundary conditions and initial conditions that I think are not essential to describe my issue. The temperature diffusion term ($D$) is by far (at least one order of magnitude) the fastest timescale in my system.
I've tried solving this using a purely time-explicit scheme, but this leads to very small time-steps (if the diffusion terms are removed from the system, I can use higher time-steps, so I guess they are the culprits). I also tried a fully implicit method (Crank-Nicholson or Backward-Euler, implemented with Petsc), but the non-linear Newton solve-step is quite slow. I would like to improve this.
I wondered if a scheme like the following could be more efficient:
- Solve the advection part of the equation
using $Q = \begin{pmatrix} n \\ nv \\ \frac{3}{2}nT + \frac{1}{2}nv^2 \end{pmatrix}$ and $F(Q) = \begin{pmatrix} nv \\ nv^2+nT \\ \frac{5}{2}nTv + \frac{1}{2}nv^3 \end{pmatrix}$,
we solve $\frac{\partial Q}{\partial t} + \frac{\partial F}{\partial x} = S$ (where $S$ is a vector of the source terms). This is a quite "classical" problem, that can be solved explicitly (I've tried a FD-WENO3 that I was happy with). This gives us $Q^\star$, from which we can compute $n^\star$, $v^\star$ and $T^\star$.
- We "push" $v$ implicitly, assuming other variables are frozen. We write
$n^\star \frac{\partial v}{\partial t} = \frac{\partial}{\partial x}\left(\nu(T^\star)\frac{\partial v}{\partial x} \right)$.
This is linear in $v$, so can be solved very effectively, and get $v^{\star\star}$
- We "push" $T$ implicitly, assuming other variables are frozen. $n^\star \frac{\partial T}{\partial t} = \frac{\partial}{\partial x}\left(D(T^\star)\frac{\partial T}{\partial x} \right) + \frac{\partial}{\partial x}\left(v^{\star\star}\nu(T^\star)\frac{\partial v^{\star\star}}{\partial x} \right) $.
which again should not be complicated to solve. We get $T^{\star\star}$.
- Solve again the explicit part using $n^\star, v^{\star\star}, T^{\star\star}$ as initial solutions.
This is then the end of the iteration, and we have as a solution at $t=t+\Delta t$, using for steps 1, 2, 3, 4 the time-steps $\Delta t /2 $, $\Delta t$,$\Delta t$,$\Delta t /2$, respectively.
As you can guess, I am not an expert in numerical methods, and I have a hard time determining if this idea makes any sense. Is it consistent? I have the feeling that it is plain wrong. If it makes sense, is the order of the steps correct ? Otherwise, would you have any ideas of methods I could explore to solve this problem ?