• {{article.mytitle}}

    • SAV方法

    1. 个人笔记(川大)

    笔记主页

    1.1. 复习回顾

    拉格朗日乘数法
    优势：

    • 二阶精度
    • 无条件稳定

    缺陷：

    • 只适用于双井模式
    • 要解变系数线性方程组

    1.2. IEQ方法

    ϕt=−GEϕ\phi_t = -\mathcal G E_\phi
    假设F(ϕ)F(\phi)有下界，取q=F(ϕ)+C0q = \sqrt{F(\phi)+C_0}，
    {ϕt=−Gμμ=Eϕ=−Δϕ+2qqϕqt=qϕϕt\begin{cases} \phi_t = - \mathcal G \mu \\ \mu = E_\phi = -\Delta\phi + 2q q_\phi\\ q_t = q_\phi \phi_t \end{cases}

    1.2.1. CN-2逼近

    {ϕn+1−ϕn=−ΔtGμn+1+μn2μn+1+μn2=−Δϕn+1+ϕn2+2qn+1+qn23qϕn+1−qϕn2qn+1−qn=3qϕn+1−qϕn2(ϕn+1−ϕn)\begin{cases} \phi^{n+1} -\phi^n= - \Delta t\mathcal G \frac{\mu^{n+1}+\mu^n}{2} \\[1em] \frac{\mu^{n+1}+\mu^n}{2}= -\Delta\frac{\phi^{n+1}+\phi^n}{2}+ 2\frac{q^{n+1}+q^n}{2} \frac{3q_{\phi}^{n+1}-q_{\phi}^n}{2}\\[1em] q^{n+1}-q^n = \frac{3q_{\phi}^{n+1}-q_{\phi}^n}{2} (\phi^{n+1} -\phi^n) \end{cases}

    上式分别与(μn+1+μn)/2,(ϕn+1−ϕn),2(qn+1+qn)(\mu^{n+1}+\mu^n)/2,( \phi^{n+1}-\phi^n), 2(q^{n+1}+q^n)做内积，分部积分即可证得修正能量的稳定性，其中修正能量定义为
    E~(ϕn+1)=12∣∣∇ϕn+1∣∣+∣∣qn+1∣∣2\tilde E(\phi^{n+1}) = \frac{1}{2}||\nabla \phi^{n+1}|| + ||q^{n+1}||^2
    该方法本质上是对拉格朗日乘子法的推广，因此它们的优点和缺点是一致的

    1.3. SAV方法

    假设∫Fdx\int Fdx有下界, 令r(t)=∫Fdx+C0r(t) = \sqrt{\int Fdx +C_0},于是
    E(ϕ)=12(ϕ,Gϕ)+(F,1)=12(ϕ,Lϕ)+r2−C0E(\phi) = \frac{1}{2}(\phi, \mathcal G\phi) +(F,1)= \frac{1}{2}(\phi, \mathcal L\phi) +r^2 - C_0
    方程修改为
    {ϕt=−Gμμ=Lϕ+2rF′2∫Fdx+C0rt=(F′,ϕt)2∫Fdx+C0\begin{cases} \phi_t = - \mathcal G \mu\\[1em] \mu = \mathcal L \phi + 2r \frac{F'}{2 \sqrt{\int Fdx +C_0}}\\[1em] r_t = \frac{(F',\phi_t)}{2 \sqrt{\int Fdx +C_0}}\\ \end{cases}
    同样可以证明上式能量是稳定的

    下面介绍离散格式

    1.3.1. CN-2

    {ϕn+1−ϕn=−ΔtGμn+1+μn2μn+1+μn2=−Lϕn+1+ϕn2+2rn+1+rn2F′(ϕ~n+12)2∫F(ϕ~n+12)dx+C0rn+1−rn=12∫F(ϕ~n+12)dx+C0(F′(ϕ~n+12),ϕn+1−ϕn)ϕ~n+12=32ϕn−12ϕn−1\begin{cases} \phi^{n+1} -\phi^n= - \Delta t\mathcal G \frac{\mu^{n+1}+\mu^n}{2} \\[1em] \frac{\mu^{n+1}+\mu^n}{2}= -\mathcal L\frac{\phi^{n+1}+\phi^n}{2}+ 2\frac{r^{n+1}+r^n}{2} \frac{F' (\tilde \phi^{n+\frac{1}{2}})}{2 \sqrt{\int F(\tilde \phi^{n+\frac{1}{2}})dx +C_0}}\\[1em] r^{n+1}-r^n =\frac{1}{2 \sqrt{\int F(\tilde \phi^{n+\frac{1}{2}})dx +C_0}} (F' (\tilde \phi^{n+\frac{1}{2}}),\phi^{n+1} -\phi^n)\\[1em] \tilde \phi^{n+\frac{1}{2}} = \frac{3}{2}\phi^n - \frac{1}{2}\phi^{n-1} \end{cases}
    上式中的第三行方程虽然含有ϕn+1\phi^{n+1}，但是因为它在积分项里面，因此可以事先把它当成一个与tt相关的常数Zn+1Z^{n+1}。从而整个方程组的求解只需要求解非线性项，因为常数Zn+1Z^{n+1}的存在，因此每一步求解需要两次，由于两次都是线性常系数方程，所以计算量非常小。

    矩阵的解释为
    (c1lG0−Lc2l∗∗0c3)(ϕn+1μn+1rn+1)=b‾n\left(\begin{array}{ccc}{c_{1} l} & {\mathcal{G}} & {0} \\ {-\mathcal{L}} & {c_{2} l} & {*} \\ {*} & {0} & {c_{3}}\end{array}\right)\left(\begin{array}{c}{\phi^{n+1}} \\ {\mu^{n+1}} \\ {r^{n+1}}\end{array}\right)=\overline{b}^{n}
    因此可以先求解前两行方程
    (c1lG−Lc2)(ϕμ)=b‾\left(\begin{array}{cc}{c_{1} l} & {\mathcal{G}} \\ {-\mathcal{L}} & {c_{2}}\end{array}\right)\left(\begin{array}{l}{\phi} \\ {\mu}\end{array}\right)=\overline{b}

    对于CH方程，需要求解如下线性方程
    {1Δtϕ+Δμ=fΔϕ+μ=g\begin{cases} \frac{1}{\Delta t}\phi + \Delta \mu = f\\[1em] \Delta \phi +\mu =g \end{cases}

    待定系数法，设
    ψ=aϕ+Δϕbψ−Δψ=f\psi = a\phi + \Delta \phi\\[1em] b\psi - \Delta \psi = f

    可以解得a+b=0,ab=1Δta+b =0, ab=\frac{1}{\Delta t}

    因此解耦需要加上稳定化因子

    1.4. 作业

    1.5. python代码

    点我下载SAV代码

    笔记主页