# Which scheme for inhomogeneous convection-diffusion equation with highly variable coefficients?

I have a 1D convection-diffusion equation $$\sigma_t = a(x,t) \sigma_{xx}+b(x,t)\sigma_x+f(x,t)$$ defined on the unit interval, with nonzero Neumann boundary conditions at both ends. It should be noted that the coefficients are not strictly positive, which explains the negative Peclet number below.

From what I understand, choosing an appropriate finite difference scheme depends on the Peclet number, or the ratio between the convective and diffusive coefficients, but this ratio $$b(x,t)/a(x,t)$$ has a strong dependence on time and varies between -700 and 1100, with the jump between these two approximate values being discontinuous in time. Here is an example plot showing this ratio at $$x=0.1$$, plots at other locations in space have a similar shape:

My question is, what sort of numerical scheme would be best suited to solve this problem? The region where $$Pe=0$$ can be discarded; the physics of the problem imply that the solution is identically zero in this region. My initial guess is to use a different scheme for each regime; i.e. when $$Pe\approx -700$$, use one scheme, and when $$Pe\approx 1100$$, use a different one, and use the final result of one scheme as the initial data for the next, but I'm unsure about the stability of such an approach. I've never seen convection-diffusion equations in the context of a negative Peclet number, so I'm unsure what scheme to use for the first regime. I've read elsewhere that higher-order upwinding is useful for large Peclet numbers, which might be an possible candidate for the second regime. Thoughts?

• Do I understand correctly that the viscosity can be negative? Because sign of the velocity does not matter when computing the Peclet number. – Abdullah Ali Sivas Jul 14 '20 at 1:41
• Yes, the viscosity can be negative. – kieransquared Jul 14 '20 at 4:21
• Given that the viscosity is negative over part of the domain, do you know that the problem is well-posed? It's not clear to me that it would be (though I'm not immediately convinced it couldn't be; it just seems potentially problematic). – David Ketcheson Jul 14 '20 at 6:51

It's just a linear 1D problem, you can easily do implicit time stepping here so that numerical stability would not be a problem. The accuracy should not be an issue either since for such a simple problem you should be able to use any spatial resolution you need.

More specifically, let the equation be discretized in space by any scheme, e.g., low-order central difference, and let the time step be implicit Euler:

$${\sigma}^{n+1} = \sigma^{n} + \tau (\hat{M} \sigma^{n+1} + \hat{F}),$$

where $${\sigma}^{n+1}$$ is the predicted state vector, $${\sigma}^{n}$$ is the old state vector, $$\hat{M}$$ expresses the differential operators and $$\hat{F}$$ expresses the source term discretized on the grid. The operators $$\hat{M}$$ and $$\hat{F}$$ should be evaluated at the mid-point in time, $$t_{n+1/2}$$.

Now, to perform a time step, solve the linear system for the predicted state:

$$(\hat{I} - \tau \hat{M} ) \sigma^{n+1} = \sigma^{n} + \tau \hat{F},$$

This scheme is unconditionally stable for any advection velocity and for any positive-definite diffusion coefficient; but for accuracy of time integration the time step will have to be chosen sufficiently small.

• This may not be a good advice, depending on the velocity field and viscosity. "Negative Peclet number" the asker mentions, makes me think that maybe the viscosity can be negative. And if the velocity field is changing rapidly or if the initial condition is discontinuous, this scheme will suffer from undershoots and overshoots. Either upwinding or some sort of slope limiting will be necessary. – Abdullah Ali Sivas Jul 14 '20 at 3:23
• @Abdullah The coefficient $a$ must be positive definite, otherwise the problem is mathematically ill posed. The coefficient $b$ can be positive or negative, does not matter. Since one can use very high spatial resolution on a 1D problem like this, just pushing the resolution will likely be enough to overcome the accuracy problems such as overshoots etc. Using upwind etc. may be beneficial but it is problem specific, depending what is more important - overshoots or numerical diffusion. Either way, it is probably at the level of fine tuning. – Maxim Umansky Jul 14 '20 at 3:29
• I agree with you completely. But if viscosity is positive definite, then Peclet number cannot be negative. That is my problem with the question. – Abdullah Ali Sivas Jul 14 '20 at 3:32
• $U$ should not be the velocity but the characteristic speed. Otherwise, it is meaningless because the characteristics of a uniform flow in one direction is not different than the characteristics of a uniform flow in the opposite direction. Direction of the velocity does not matter, only the magnitude. – Abdullah Ali Sivas Jul 14 '20 at 3:40
• @Abdullah The author said he/she used for $Pe$ the ratio $b/a$, so $b$ must be negative, that's the most likely explanation - because if $a$ is negative then the problem is not solvable, just ill-posed. – Maxim Umansky Jul 14 '20 at 3:43