Numerical Investigation of two Component Reaction-Diffusion System with Perturbative Advection

I am given two coupled differential equations and try to analyze them numerically because as far as I know, there is no analytical solution known(only the steady-state on an infinite domain is known, which I use to define initial conditions).

These two equations basically describe a mass conserving reaction-diffusion system to which I may add perturbations that are geometrically induced advection terms(still mass conserving).

When I activate the perturbation I have small advection terms that shift the pattern but do not deform the pattern in a relevant way. That's also what I observe numerically.

The advection equation(without diffusion, etc) reads du/dt + v du/dx = 0

v > 0 -> pattern moves to the right.

v < 0 -> pattern moves to the left.

In my equations, the dominant advection terms are proportional to dR/dx which I called Rfirstderiv(x). So v is proportional to Rfirstderiv(x) and also its sign is determined by Rfirstderiv(x).

Also, intuitively the pattern is only sensitive to advection perturbations by Rfirstderiv(x) where the interface is located because that's where c du/dx does not vanish. So you would suspect, that the direction of the pattern shift changes when Rfirstderiv(x) changes sign at the point of the interface.

But in my numerics, the interface just moves across the extremum of R as if there wasn't any change in the sign of the advection term. This has nothing to do with the boundaries: The very same behavior also occurs when you apply periodic boundaries.

Also note that I have a two-component system where the corresponding interfaces are at the same position in x, so the argument above should be valid.

Also, the reasons for that behavior aren't the error warnings. For example, if I choose BDF-method instead the warnings disappear but still the interface moves across the extremum of R.

So my question is: Can it be correct? Personally, I don't think it can be. The interface shouldn't move across the extremum of R. But then I don't understand what is going on with the numerics. I can't spot an error. Can someone help, please? I would be so glad if someone could.

Best, Luca

Here I add the corresponding notebook as a link to the wolfram cloud:

https://www.wolframcloud.com/obj/2668fe35-f942-42c8-aa09-129789854d0a

-The first cell defines the parameters, the reaction kinetics, and the stationary pattern without perturbation(coarsening neglected)

-the second cell defines the geometric perturbation and includes plots where you can see yourself, that the advection term undergoes a change in sign when crossing the extremum.

-third cell defines initial conditions

-fourth cell defines the dynamic equations and boundary conditions

-fifth cell just executes the simulation

-sixth cell shows the two components as a function of time. One of the plots also shows the geometry(radius) where the perturbation is multiplied by 1000 for visibility.

-seventh line checks mass-conservation

• It is better to post the code here so that this QA is self contained. Aug 18 '21 at 13:37