interpolation combined with methods of characteristics can cause oscillations for the transport equation?

Question

I would like to know about the effect of using a higher order interpolator for the methods of characteristics. I am solving $$u_t+a(x,t)u_x=0$$ with some nonsmooth initial data $u_0(x)$ by the method of characteristics. It seems to be a great method, since I need not deal with any of the difficulties that arise from trying to solve this problem discretely by finite differences. The numerical solution will not have any oscillations at all!

However, that is for the case of "infinitely-refined grid" that is, if I know value of the function at any value of $x$ at the previous time step. Instead, consider a grid $x_i:i\in\{1,N\}$ and march forward in time. Since the characteristics don't cross the nodes perfectly I have to interpolate. I can choose linear interpolation but would like to take a higher order to increase the accuracy of interpolation.

Assume I choose a piecewise quadratic polynomial. Then the value at the interpolated point might be higher or smaller than the value at the end points of the interval and thus the extrema can get magnified from step to step. I can artificially restrict the interpolated value to be between $\min_{x\in I} u$ and $\max_{x\in I} u$ to prevent such growth (where $I$ is the interval the characteristic crosses at the previous known time).

There are two questions I would like to know the answer to:

1. Looks like if I take higher order than linear interpolator, methods of characteristics lose stability as the value can grow because the maximum or minimum of the parabola that I am fitting can go outside of the range of the endpoint values (Assume I have something like $\sin(x)$) as the initial profile? So it is not as effective in this case and if I keep the interpolation linear I can't do better than first order convergence.

2. Does the artificial restriction on the range of values affect the accuaracy of quadratic interpolation? If before the restriction the accuracy is second order, do I lower that by "squaeezing" the value to be between max and min value over the interval? I WOULD LIKE TO KNOW IF I CAN DO METHOD OF CHARACTERISTICS WITH SECOND ORDER ACCURACY?

Answer 1

This is an important and challenging issue.

1. Yes, using quadratic interpolation means that your solution values may lie outside the interval in which the initial data lie. This is not what we usually mean when we refer to numerical instability, but it is a potentially undesirable feature.

2. Yes, forcing the interpolated values to lie in an interval destroys the accuracy of your solution. Nevertheless, it is a commonly-used technique.

One may use slope-limiter methods or ENO/WENO interpolants to mitigate this problem, but with either of those approaches, you still have to give up either the formal 2nd-order accuracy or the strict maximum principle.

One very recent proposed approach that claims to give high order accurate, data-bounded interpolants is in this paper of Martin Berzins.

What you want to do can also be done if you advect cell averages instead of point values; see Xiangxiong Zhang's thesis.

A final remark, unrelated to the focus of your question: for general functions $a(x,t)$, you'll need to use a Newton iteration to implement the method of characteristics. This is a major drawback compared to traditional finite difference methods.