2
$\begingroup$

The main motivation behind my next question is that I think I derived a higher order numerical scheme for linear advection equation that is unconditionally stable using Von Neumann stability analysis.

In the stability analysis of the numerical scheme I have to find the maximal value of so called amplification factor $|S|$ where $S$ depends on 5 parameters, it means $S=S(x,y,f,c,d)$ and $$ S= \left(2 \left(-c^2+c \cos (x) (c-d (2 f+1) \cos (y)+2 d f+d-1)+d \cos (y) (2 c f+c+d-1)-2 c d f-c d \sin (x) \sin (y)-c d-3 i c \sin (x)+c-d^2-3 i d \sin (y)+d+6\right)\right) \cdot \left(2 i c^2 \sin (x)-i c^2 \sin (2 x)+c^2-2 i c d \sin (x+y)+2 c d \cos (x+y)+2 i c d \sin (x)-2 c (c+d+5) \cos (x)+2 i c d \sin (y)-2 c d \cos (y)+2 c d+10 i c \sin (x)-2 i c \sin (2 x)+(c+2) c \cos (2 x)+8 c+2 i d^2 \sin (y)-i d^2 \sin (2 y)-2 d^2 \cos (y)+d^2 \cos (2 y)+d^2+10 i d \sin (y)-2 i d \sin (2 y)-10 d \cos (y)+2 d \cos (2 y)+8 d+12\right)^{-1} $$ Note that $i$ is imaginary number.

My hypothesis is that $|S| \le 1$ for $-\pi \le x \le \pi$, $-\pi \le y \le \pi$, $-1 \le f \le 0$ and $c \ge 0$ and $d \ge 0$.

Up to now I checked the values of $|S|$ for many tabulated values of input variables e.g. for $c \le 100$ and $d \le 100$. Moreover I used Mathematica to compute numerically the maximal values of $|S|$ with many starting values and several built-in optimization methods, I got always the maximal value of $|S|$ being $1$. Of course I used the numerical scheme to compute several examples of linear advetion equation with very large time steps and observed no instabilities.

I would appreciate any useful thoughts, helps or even your tries ;-) how to find the values of $(x,y,f,c,d)$ bounded by above constraints for which $|S|>1$. Or from the other view - what kind of methods and arguments would you use to support the hypothesis that $|S| \le 1$ ?

P.S. Here is the input for Mathematica

S=(2*(6 + c - c^2 + d - c*d - d^2 - 2*c*d*f + 
     d*(-1 + c + d + 2*c*f)*Cos[y] + 
     c*Cos[x]*(-1 + c + d + 2*d*f - d*(1 + 2*f)*Cos[y]) - 
     3*I*c*Sin[x] - 3*I*d*Sin[y] - c*d*Sin[x]*Sin[y]))/(12 + 8*c + c^2 + 8*d + 2*c*d + d^2 - 2*c*(5 + c + d)*Cos[x] + 
  c*(2 + c)*Cos[2*x] - 10*d*Cos[y] - 2*c*d*Cos[y] - 
  2*d^2*Cos[y] + 2*d*Cos[2*y] + d^2*Cos[2*y] + 
  2*c*d*Cos[x + y] + 10*I*c*Sin[x] + 2*I*c^2*Sin[x] + 
  2*I*c*d*Sin[x] - 2*I*c*Sin[2*x] - I*c^2*Sin[2*x] + 
  10*I*d*Sin[y] + 2*I*c*d*Sin[y] + 2*I*d^2*Sin[y] - 
  2*I*d*Sin[2*y] - I*d^2*Sin[2*y] - 2*I*c*d*Sin[x + y])

P.S.S. If for any reasons you need more insights, see my paper at http://dl.dropboxusercontent.com/u/386482/articles/frolkovic-algoritmy-corr.pdf where the scheme is given in (3.10) for $f=0$.

Improve this question
$\endgroup$
3
  • $\begingroup$ I know cylindrical algebraic decomposition sometimes works, but this problem looks like it's too large. $\endgroup$
    – Kirill
    Commented Oct 16, 2016 at 21:35
  • $\begingroup$ I tried it directly, without coming up with anything clever, and QEPCAD failed after an hour with "Prime list exhausted." $\endgroup$
    – Kirill
    Commented Oct 16, 2016 at 23:13
  • $\begingroup$ If you could clean up the equation by clustering function terms, it might be easier to come up with estimates. For instance I can't tell where the inversion begins. Presumably you can also collect the trigonometric terms. They can always be estimated by $\pm 1,0$, whichever gives you the largest or smallest value with the right sign. $\endgroup$
    – Deathbreath
    Commented Nov 1, 2016 at 16:06

2 Answers 2

Reset to default
1
$\begingroup$

In your expression for $S$, is $x = k_x \Delta x$ and $y = k_y \Delta y$ for some "wave number" vector $\mathbf{k} = [k_x \quad k_y]^T$ with $\Delta x $ and $\Delta y$ being grid spacings in $x$ and $y$ directions?

I assumed that $x = \pi$ and $y = \pi$ and obtained the following for $S$: $$ -{\frac {2\,cdp+{c}^{2}+{d}^{2}-c-d-3}{{c}^{2}+2\,cd+{d}^{2}+5\,c+5\,d +3}} $$ which satisfies $ |S| \le 1$. The quantity $p$ is equal to $2f +1$.

The values $x=\pi$ and $y=\pi$, of course, correspond to the minimum resolvable wavelengths in each direction.

I played with a few other values of $x$ and $y$ as well and in all those cases, $|S| \le 1$.

Improve this answer
$\endgroup$
0
$\begingroup$

Here is an example how I have tested the hypothesis that $|S| \le 1$ for given constraints. I have run the following code in Mathematica with many many values of variables nstart and nend and the method $M$:

(* amplification parameter *)
S=
(2*(6+c-c^2+d-c*d-d^2-2*c*d*f+d*(-1+c+d+2*c*f)*Cos[y]+c*Cos[x]*(-1+c+d+2*d*f-d*(1+2*f)*Cos[y])-3*I*c*Sin[x]-3*I*d*Sin[y]-c*d*Sin[x]*Sin[y]))
/
(12+8*c+c^2+8*d+2*c*d+d^2-2*c*(5+c+d)*Cos[x]+c*(2+c)*Cos[2*x]-10*d*Cos[y]-2*c*d*Cos[y]-2*d^2*Cos[y]+2*d*Cos[2*y]+d^2*Cos[2*y]+2*c*d*Cos[x+y]+10*I*c*Sin[x]+2*I*c^2*Sin[x]+2*I*c*d*Sin[x]-2*I*c*Sin[2*x]-I*c^2*Sin[2*x]+10*I*d*Sin[y]+2*I*c*d*Sin[y]+2*I*d^2*Sin[y]-2*I*d*Sin[2*y]-I*d^2*Sin[2*y]-2*I*c*d*Sin[x+y]);

(* numerical global maximum search *)
M="RandomSearch";(* RandomSearch NelderMead DifferentialEvolution SimulatedAnnealing*)
nstart=1;nend=10;

For[n=nstart,n<=nend,n++,
Print[NMaximize[{Abs[S],-Pi<=x<=Pi, -Pi<= y<= Pi, 0<= c <=n, 0<=d<=n, -1<=f<=0},
{x,y,c,d,f}, Method->{M,"SearchPoints"->20,RandomSeed->n}]]]

(* example of output *)
{1.,{x->3.87721*10^-7,y->-3.14159,c->0.996877,d->0.,f->-0.998816}}
{1.,{x->-3.14159,y->2.19618*10^-8,c->0.,d->1.71028,f->-0.941516}}
{1.,{x->-3.14159,y->-3.14159,c->0.,d->0.,f->-1.}}
{1.,{x->3.14159,y->3.14159,c->0.,d->0.,f->-1.}}
{1.,{x->4.71742*10^-7,y->3.14159,c->3.3804,d->0.,f->-0.997559}}
{1.,{x->-3.14159,y->-3.14159,c->0.,d->0.,f->-1.}}
{1.,{x->-3.14159,y->-2.4762*10^-7,c->0.,d->1.34533,f->-0.922121}}
{1.,{x->-0.000101745,y->-0.0000189795,c->3.48753,d->4.20828,f->-0.84976}
{1.,{x->3.14159,y->1.0698*10^-6,c->0.,d->8.98382,f->-0.63621}}
{1.,{x->0.000268716,y->-0.0000310003,c->0.593874,d->3.95438,f->-0.832499}}
Improve this answer
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.