Imex and operator splitting methods — Examples

With $a(y)=-\kappa y$ and $b=0$, [F1] gives $Y=(1-h\kappa )y$, so after $N$ steps the numerical value is $(1-h\kappa {)}^N$. If $h\kappa >2$, then $1-h\kappa <-1$ and $\mathrm{\mid }1-h\kappa {\mathrm{\mid }}^N>1$ for every $N\ge 1$. By [F2], $y(t)=\exp (-\kappa t)$ is the exact solution, and [F3] gives $\mathrm{\mid }y(Nh)\mathrm{\mid }=\exp (-\kappa Nh)<1$.

With $a=0$ and $b(y)=-\kappa y$, the equation in [F1] becomes $(1+h\kappa )Y=y$. Since $0<1<1+h\kappa$, [F4] gives $0<(1+h\kappa {)}^{-1}<1$. Thus the implicit treatment damps every such scalar step, while the explicit treatment grows whenever $h\kappa >2$.

At the boundary $h\kappa =2$, the explicit amplification is $-1$, whose magnitude is exactly one. Thus the strict instability condition is $h\kappa >2$; equality is neutral for this scalar test equation.

By [F4], the linear fields have global exact flows. Direct multiplication gives $BA=\left(\begin{array}{cc}0& 0\\ 0& 1\end{array}\right),AB=\left(\begin{array}{cc}1& 0\\ 0& 0\end{array}\right),C:=BA-AB=\left(\begin{array}{cc}-1& 0\\ 0& 1\end{array}\right)\mathrm{.}$ By [F1], $[a,b](x)=Cx$.

Multiplying once more gives $CA-AC=\left(\begin{array}{cc}0& -2\\ 0& 0\end{array}\right),CB-BC=\left(\begin{array}{cc}0& 0\\ 2& 0\end{array}\right)\mathrm{.}$ These are the matrices of $[a,[a,b]]$ and $[b,[a,b]]$ under the same convention.

At $e_1$, $C{e}_1=-e_1\mathrm{\ne }0$, so [F2] gives Lie local defect $-\frac{1}{2}{h}^2{e}_1+O(h^3)$, whose leading coefficient is nonzero. Thus its local error has exact order two.

The Strang leading matrix from step 2.1 is $-\frac{1}{24}(CA-AC)-\frac{1}{12}(CB-BC)=\left(\begin{array}{cc}0& 1\mathrm{/}12\\ -1\mathrm{/}6& 0\end{array}\right)\mathrm{.}$ It sends $e_1$ to $(0,-1\mathrm{/}6{)}^T\mathrm{\ne }0$, so [F2] gives a nonzero cubic leading local defect and exact local order three.

Let $H$ be any of $H_K$, $V$, or $H_K+V$. Its exact flow exists for all real time by [F7], and every solution obeys ${\psi }^{\mathrm{\prime }}=-iH\psi$. By [F3]–[F6], differentiating the real function ${\psi }^{\ast }\psi$ gives $\frac{d}{dt}({\psi }^{\ast }\psi )=({\psi }^{\mathrm{\prime }}{)}^{\ast }\psi +{\psi }^{\ast }{\psi }^{\mathrm{\prime }}=i{\psi }^{\ast }{H}^{\ast }\psi -i{\psi }^{\ast }H\psi =0,$ where [A1] supplies $H^{\ast }=H$. Therefore [F6] makes $\mathrm{\parallel }\psi (t){\mathrm{\parallel }}_2^2$ constant on every admitted interval, and each exact subflow and the total flow preserve norm.

Direct multiplication gives $C:=BA-AB=H_{K}V-V{H}_K=\left(\begin{array}{cc}0& -1\\ 1& 0\end{array}\right)\mathrm{\ne }0.$ This is the linear matrix of $[a,b]$ by [F1]. Also $CA-AC=i\left(\begin{array}{cc}0& 1\\ 1& 0\end{array}\right),CB-BC=2i\left(\begin{array}{cc}1& 0\\ 0& -1\end{array}\right)\mathrm{.}$

By [F1], a Lie or Strang step is a composition of exact component flows. Each factor preserves norm by step 1.1, so the composition preserves norm by repeated application of the same equality.

For $\psi =e_1$, the Lie coefficient $\frac{1}{2}C\psi$ is nonzero. The Strang coefficient $-\frac{1}{24}(CA-AC)\psi -\frac{1}{12}(CB-BC)\psi$ is also nonzero, since its first coordinate is $-i\mathrm{/}6$. Applying [F2] shows the exact subflows retain norm while the Lie and Strang maps have their respective nonzero local splitting defects.

For each mode, solving $u_k^{\mathrm{\prime }}=-k^2{u}_k+g_k$, $u_k(0)=0$, gives $u_k(t)=\frac{g_k}{k^2}(1-\exp (-k^{2}t))\mathrm{.}$ Differentiating verifies the equation and initial value. Also $\mathrm{\mid }{u}_k(t)\mathrm{\mid }\le \mathrm{\mid }{g}_k\mathrm{\mid }\mathrm{/}{k}^2\le k^{-3}$, so [F4] shows the reference state belongs to ${\mathrm{\ell }}^2$ for every $t\ge 0$.

Put $x=k^{2}h$. By [F1], each Strang coordinate obeys ${\tilde{u}}_{n+1,k}=\exp (-x){\tilde{u}}_{n,k}+h\exp (-x\mathrm{/}2)g_k\mathrm{.}$ Starting from zero and summing this scalar recurrence gives ${\tilde{u}}_{N,k}=h\exp (-x\mathrm{/}2)\frac{1-\exp (-Nx)}{1-\exp (-x)}{g}_k,$ where the finite geometric identity follows by multiplying the sum by $1-\exp (-x)$; the denominator is positive because $x>0$ and the exponential is strictly increasing.

Define $\sigma (x):=\exp (x\mathrm{/}2)-\exp (-x\mathrm{/}2)$ for real $x$, and $q(x):=1-x\mathrm{/}\sigma (x)$ for $x>0$. For $x>0$, ${\sigma }^{\mathrm{\prime }}(x)=\frac{1}{2}(\exp (x\mathrm{/}2)+\exp (-x\mathrm{/}2))>1$, because $\exp (x\mathrm{/}2)+\exp (-x\mathrm{/}2)-2=(\exp (x\mathrm{/}4)-\exp (-x\mathrm{/}4){)}^2>0.$ Thus [F6] gives $\sigma (x)>x$, so $0<q(x)<1$. Direct differentiation gives $\sigma (0)=0$, ${\sigma }^{\mathrm{\prime }}(0)=1$, and ${\sigma }^{\mathrm{\prime }\mathrm{\prime }}(0)=0$; [F7] therefore gives $\sigma (x)=x+O(x^3)$ on $0<x\le 1$. Hence $q(x)=(\sigma (x)-x)\mathrm{/}\sigma (x)\le C{x}^2$ there.

By [F2] and [F3], this constant Dirichlet forcing is outside the generator domain and has nonzero boundary value while the sine modes vanish at both endpoints. The modal comparison below therefore concerns the boundary-incompatible data identified in the Statement.

Since $Nx=k^{2}T$, the identity $\frac{h\exp (-x\mathrm{/}2)}{1-\exp (-x)}=\frac{1}{k^2}\frac{x}{\exp (x\mathrm{/}2)-\exp (-x\mathrm{/}2)}$ turns the numerical value into the exact modal value times $x\mathrm{/}\sigma (x)$. Therefore for odd $k$, $\mathrm{\mid }{\tilde{u}}_{N,k}-u_k(T)\mathrm{\mid }=k^{-3}(1-\exp (-k^{2}T))q(k^{2}h),$ and the error is zero for even $k$.

The quotient $\sigma (x)\mathrm{/}x$ is increasing for $x>0$: its derivative has numerator $p(x)=x{\sigma }^{\mathrm{\prime }}(x)-\sigma (x)$, with $p(0)=0$ and $p^{\mathrm{\prime }}(x)=x{\sigma }^{\mathrm{\prime }\mathrm{\prime }}(x)=x\sigma (x)\mathrm{/}4>0$. Therefore $x\mathrm{/}\sigma (x)$ decreases and $q$ increases. In particular $q(x)\ge q(1)>0$ on $1\le x\le 4$.

For sufficiently small $h\le 1$, set $H:=1\mathrm{/}\sqrt{h}$; [F8] gives $H\ge 1$. For odd $k\le H$, steps 2.2 and 1.3 give $\mathrm{\mid }{\tilde{u}}_{N,k}-u_k(T)\mathrm{\mid }\le C{h}^{2}k$. The sum of their squared errors is at most $C{h}^4{H}^3=C{h}^{5\mathrm{/}2}$. For $k>H$, use $q<1$ and $1-\exp (-k^{2}T)<1$; dyadic blocks $[2^{j}H,2^{j+1}H)$ contain at most $2^{j+1}H$ integers and contribute at most $2^{j+1}H(2^{j}H{)}^{-6}$. Summing the geometric series gives a tail at most $C{H}^{-5}=C{h}^{5\mathrm{/}2}$. Thus the squared ${\mathrm{\ell }}^2$ error is at most $C_T{h}^{5\mathrm{/}2}$.

Choose $h$ sufficiently small that $h\le 1\mathrm{/}64$ and $h\le T$. Then $\sqrt{h}\le 1\mathrm{/}8$ by [F8], so $H\ge 8$. The interval $[H,2H]$ contains at least $H-1$ integers, and among any $m$ consecutive integers at least $(m-1)\mathrm{/}2$ are odd; hence it contains at least $H\mathrm{/}4$ odd integers. For each such $k$, $1\le k^{2}h\le 4$, $k^{-6}\ge (2H{)}^{-6}$, $q(k^{2}h)\ge q(1)$ by step 2.3, and $1-\exp (-k^{2}T)\ge 1-\exp (-T\mathrm{/}h)>1\mathrm{/}2$ since $T\mathrm{/}h\ge 1$ and $e=1+1+{\sum }_{n\ge 2}1\mathrm{/}n!>2$. Therefore the sum of squared errors over these modes is at least $\frac{H}{4}(2H{)}^{-6}\frac{q(1{)}^2}{4}=c{H}^{-5}=c{h}^{5\mathrm{/}2}$ for a constant $c>0$. Taking square roots of the upper bound in step 3.1 and this lower bound, then using [F9], proves the claimed two-sided $h^{5\mathrm{/}4}$ error rate.

Direct multiplication gives $L{v}_0=0$ for $v_0=(1,1{)}^T$ and $L{v}_1=-2v_1$ for $v_1=(1,-1{)}^T$. These vectors are linearly independent, so they form a basis of ${\mathbb{R}}^2$ and diagonalize both $a(y)=y$ and $b(y)=10Ly$.

On each span $\mathbb{R}{v}_j$, every stage equation in [F1] remains in that span. The implicit denominators are $1$ on $v_0$ and $1+20\gamma h>0$ on $v_1$, so the stage solves are unique there. Hence [F2] applies with $z_a=h$ and $z_b=0$ or $z_b=-20h$, respectively.

For $v_0$ and $h=1\mathrm{/}10$, $z_b=0$, so [F2] gives $R(1\mathrm{/}10,0)=1+\frac{1}{10}+\frac{1}{200}=\frac{221}{200}>1.$ This mode grows under the positive reaction term.

For $v_1$, $z_b=-2$. Substituting [F3] into [F2] gives the displayed fraction, with denominator $(1+2\gamma {)}^2=(3-\sqrt{2}{)}^2=11-6\sqrt{2}$. By [F5], $7\mathrm{/}5<\sqrt{2}<3\mathrm{/}2$: squaring the positive endpoints gives $49\mathrm{/}25<2<9\mathrm{/}4$. Thus the numerator is greater than $1\mathrm{/}200$ and less than $9\mathrm{/}40$, while the denominator is greater than $2$. Consequently $0<R(1\mathrm{/}10,-2)<9\mathrm{/}80<1$.