add Jacobain elements for du3^j_err_dq^j and drhoa_err_dq^j

[sajjadazimi]

Purpose

• Improve stability by helping Newton converge,
• clean up code for sgs advection terms in manual sparse Jacobian.

In this PR we add Jacobian elements for $\partial u^{3j}_err/\partial q^j$ and $\partial \rho a_err/\partial q^j$ through variations of $\rho^j$ with respect to condensate masses $q$.

$$\rho = \frac{p}{R_m T}$$

so

$$\frac{\partial \rho}{\partial q} = -\frac{p}{(R_m T)^2}\frac{\partial R_m T}{\partial q} = - \frac{\rho^2}{p}\frac{\partial R_m T}{\partial q}.$$

We know that

$$\begin{align*} R_m T =& R_m T_0 + \frac{R_m}{c_{p, m}}\left[mse - \phi - (q_t - q_l - q_i) LH_{v0} + q_i LH_{f0}\right]\\\ =& \left[(1-q_t)R_d + (q_t - q_l - q_i)R_v\right]T_0\\\ +& \frac{ \left[(1-q_t)R_d + (q_t - q_l - q_i)R_v\right]}{(1- q_t)c_{p, d} + (q_t - q_l - q_i)c_{p, v} + q_l c_{p, l} + q_i c_{p, i}}\left[mse - \phi - (q_t - q_l - q_i) LH_{v0} + q_i LH_{f0}\right].\\\ \end{align*}$$

The derivate of $R_m T$ with respect to the condensate mass $q$ is

$$\begin{align*} \frac{\partial R_m T}{\partial q} &= -R_v T_0 + \frac{R_m}{c_{p, m}} LH + \left[\frac{-R_v c_{p, m} - R_m \Delta c_p}{c_{p, m}^2}\right]c_{p, m} (T-T_0)\\\ &= -R_v T_0 + \frac{R_m}{c_{p, m}} LH -R_v (T-T_0) - \frac{R_m}{c_{p, m}} \Delta c_p (T-T_0)\\\ &=\frac{R_m}{c_{p, m}} \left[LH - \Delta c_p (T-T_0)\right] - R_v T\\\ &= \frac{\kappa_m}{\kappa_m + 1} \left[LH - \Delta c_p (T-T_0)\right] - R_v T, \end{align*}$$

where $\Delta c_p$ is the difference of $c_p$ of the condensate and water vapor, $LH$ is the latent heat of the condensate ($LH_{v0}$ for liquid and $LH_{v0} + LH_{f0}$ for ice) and $\kappa_m = R_m/c_{v, m}$. Therefore,

$$\begin{align*} \frac{\partial \rho}{\partial q} = - \frac{\rho^2}{p}\left[ \frac{\kappa_m}{\kappa_m + 1} \left[LH - \Delta c_p (T-T_0)\right] - R_v T\right]. \end{align*}$$

---

• I have read and checked the items on the review checklist.

[szy21]

Do we need this (I think in the loop we are doing =)?

[szy21]

This sets it to grid-mean if ρa is smaller than eps - is that what you want? If so could you update the comment?

[szy21]

Thanks! Just some minor comments.