efficient IMEX

All we need to modify are the Jacobian terms for ϕ\phi. From the IMEX preprint

Gm,j(gm,jϕ)=gm,j−Em,jϕ−gΔtA^j,jEmw+(gΔtA^j,j)2(1−μm,j) G_{m,j}(g_{m,j}^\phi) = g_{m,j} - E_{m,j}^\phi - g \Delta t \hat{A}_{j,j} E_m^w + (g \Delta t \hat{A}_{j,j})^2 (1-\mu_{m,j})

which gets translated into (dicarding any terms that disappear under differentiation),

C1−∑(δϕ)−C2−C3−(gΔtAj,j)2μm,j C_1 - \sum(\delta \phi) - C_2 - C_3 - (g\Delta t A_{j,j})^2\mu_{m,j} . The tricky thing is that we are calculating this quantity at model interfaces, which means and so we find μk=pk+1−pk12(dp3dk+1+dp3dk)=(dp3di)−1p0[(p0−1Rdθv,k+1dp3dk+1ϕk+1−ϕk)11−κ−(p0−1Rdθv,kdp3dkϕk−ϕk−1)11−κ]=(dp3di)−1p0[(p0−1Rdθv,k+1dp3dk+1)11−κ(ϕk+1−ϕk)11−κ−(p0−1Rdθv,kdp3dk)11−κ(ϕk−ϕk−1)11−κ] \begin{align*} \mu_k &= \frac{p_{k+1} - p_k}{\frac{1}{2}(\textrm{dp3d}_{k+1} + \textrm{dp3d}_k)}\\ &= (\textrm{dp3d}_i)^{-1} p_0 \left[ \left( p_0^{-1} \frac{R_d \theta_{v, k+1} \textrm{dp3d}_{k+1}}{\phi_{k+1}-\phi_{k}} \right)^{\frac{1}{1-\kappa}} - \left( p_0^{-1} \frac{R_d \theta_{v, k} \textrm{dp3d}_{k}}{\phi_{k}-\phi_{k-1}} \right)^{\frac{1}{1-\kappa}} \right]\\ &= (\textrm{dp3d}_i)^{-1} p_0 \left[ \left( p_0^{-1} R_d \theta_{v, k+1} \textrm{dp3d}_{k+1} \right)^{\frac{1}{1-\kappa}} \left( \phi_{k+1}-\phi_{k} \right)^{\frac{1}{1-\kappa}} - \left( p_0^{-1} R_d \theta_{v, k} \textrm{dp3d}_{k} \right)^{\frac{1}{1-\kappa}} \left(\phi_{k}-\phi_{k-1} \right)^{\frac{1}{1-\kappa}}\right] \end{align*}

and so to make everything as easy to read as possible

∂ϕk+1[μk]=(dp3di)−1p0[(p0−1Rdθv,k+1dp3dk+1)11−κ(ϕk+1−ϕk)11−κ1(1−κ)(ϕk+1−ϕk)]=−pk+1dp3di,k+1(ϕk+1−ϕk) \begin{align*} \partial_{\phi_{k+1}} [ \mu_k] &= (\textrm{dp3d}_i)^{-1} p_0 \left[ \left( p_0^{-1} R_d \theta_{v, k+1} \textrm{dp3d}_{k+1} \right)^{\frac{1}{1-\kappa}} \left( \phi_{k+1}-\phi_{k} \right)^{\frac{1}{1-\kappa}}\frac{1}{(1-\kappa)\left( \phi_{k+1}-\phi_{k} \right)} \right]\\ &= -\frac{p_{k+1}}{\textrm{dp3d}_{i, k+1}(\phi_{k+1} - \phi_k)} \end{align*} which gives us the template

a=(Δtg(ϕ1))2(1−κ)b=adp3dint,1c1=p1Δϕ1JL=bk+1ckJU,k=1=2bkckJU,k≠1=bkckJD,k=1=1−JU,1JD,k≠1=1−JL,k−1−JU,k \begin{align*} a &= \frac{(\Delta t g(\phi_1))^2}{(1-\kappa)}\\ b &= \frac{a}{\textrm{dp3d}_{\textrm{int}, 1}}\\ c_1 &= \frac{p_1}{\Delta \phi_1} \\ J_L &= b_{k+1}c_k \\ J_{U,k=1} &= 2b_kc_k \\ J_{U,k\neq 1} &= b_k c_k \\ J_{D, k=1} &= 1- J_{U, 1} \\ J_{D, k \neq 1} &= 1-J_{L,k-1} - J_{U, k} \end{align*}

where we can expand bmcn=(Δtg(ϕm))2(1−κ)⋅dp3dint,m⋅pnΔϕn=(Δtg(ϕm))2(1−κ)⋅dp3dint,m⋅pnΔϕn \begin{align*} b_m c_n &= \frac{(\Delta t g(\phi_m))^2}{(1-\kappa) \cdot \textrm{dp3d}_{\textrm{int}, m}} \cdot \frac{p_n}{\Delta \phi_n}\\ &= \frac{(\Delta t g(\phi_m))^2}{(1-\kappa) \cdot \textrm{dp3d}_{\textrm{int}, m}} \cdot \frac{p_n}{\Delta \phi_n}\\ \end{align*}

The trick comes from differentiating μ\mu with respect to ϕ\phi μ=pk−pk+1dp3dint,k=RdTkdp3dkΔϕk−RdTk+1dp3dk+1Δϕk+1dp3dint,k+1 \begin{align*} \mu &= \frac{p_k - p_{k+1}}{\textrm{dp3d}_{\textrm{int}, k}} \\ &= \frac{ \frac{R_d T_k \textrm{dp3d}_k }{\Delta \phi_k} - \frac{R_d T_{k+1} \textrm{dp3d}_{k+1} }{\Delta \phi_{k+1}}}{\textrm{dp3d}_{\textrm{int}, k+1}} \end{align*} and thus ∂Δϕkμ=−pkdp3dint,k+1[Δϕk]−1 \begin{align*} \partial_{\Delta \phi_k} \mu &= -\frac{p_k}{\textrm{dp3d}_{\textrm{int}, k+1}} [\Delta \phi_{k}]^{-1} \end{align*}

Parent post: