Improving the time integration of implicit discretized PDE with a non-linear source term

Question

This might be a naive question, but when applying a implicit discretization to a PDE with a source term, should the source be averaged in time? For example if we take the diffusion equation with a non-linear source term,

$$ u_t = u_{xx}+s(x,t,u) $$

We can apply the following central difference implicit scheme to the differential term,

$$ \frac{u_j^{n+1} - u_j^n}{\Delta t} = \left[ (1-\theta) (u_{j-1} - 2u_{j} + u_{j+1}) + \theta (u_{j-1} - 2u_{j} + u_{j+1})\right] + s(x,t,u) $$

but how should $s(x,t,u)$ be treated? Should we simply take the value the $n$ time point (this is what I have always done in the past),

$$ s(x,t) = s_j^n $$

or averaged over time,

$$ s(x,t) = (1-\theta)s_j^{n+1} + \theta s_j^{n} $$

I am not sure it is possible to implement a time average in this way because, in general, the $n+1$ points in time are unknowns!

Is this a silly question? Or is there some way of improving the time integration of the above equation by taking averages in time?

Answer 1

If the right hand side were independent of $u$ then one would generally use the averaged form $$ (1-\theta)s_(x,t^{n+1}) + \theta s(x,t^n). $$ In the nonlinear case you can't do that easily, as you note, but you can at least use some kind of extrapolation, for example approximate $$ (1-\theta)s(x,t^{n+1},u^{n+1}) + \theta s(x,t^n,u^n) \approx (1-\theta)s\left(x,t^{n+1},u^{n}+\Delta t\frac{u^n-u^{n-1}}{\Delta t}\right) + \theta s(x,t^n,u^n) . $$ You will find more tricks like this in the ODE literature -- read up on semi-implicit methods for ODEs of the kind $\dot x(t) = A x(t) + f(x(t),t)$.

Answer 2

The most robust way is solving it implicitly, otherwise for stiff nonlinearity in the source function $s(u)$ you will have to use very small time step to maintain numerical stability.

Consider

$u_t = s(t,u)$

Then

$u^{n+1} = u^n + dt \; s(t+dt, u^{n+1})$

Here one can use average between n and n+1 in the RHS to increase the order, or use more temporal points for even higher order in time. But to make it simple let's assume first order implicit time-stepping.

Now the equation can be cast in the form

${F(u^{n+1})}=0$

which can be solved by the Newton method, and for the discretized system doing linear solves with the Krylov subspace is the best known approach, so this is the Newton-Krylov method.