divergence of static stability identical parameters

Question

All, The static stability of the atmosphere is defined as:
$\sigma =-\dfrac{T}{\theta} \dfrac{d\theta}{dp}= \dfrac{dT}{dp}-\dfrac{R}{c_p} \dfrac{T}{p}$

where: $T$ is temperature, $P$ is pressure, $R$ is ideal gas constant, $c_p$ is heat capacity at constant pressure, $\theta$ is potential temperature, which is defined as $\theta=T \left(\dfrac{p_o}{p}\right)^{R/c_p}$, $p_0$ is reference pressure.

The calculated statics stability using $\dfrac{dT}{dp}-\dfrac{R}{c_p} \dfrac{T}{p}$ or $-\dfrac{T}{\theta} \dfrac{d\theta}{dp}$ should be identical, yet they actually differ sometimes as shown below

(plz see attached Python code used to calculated statics stability )

I am wondering if there is analytical or numerical way to determine when both identical definition of static stability might diverge, though? Sometimes, the way that terms is calculated numerically can significantly change the results, for example, the physical properties of the Jacobian of the 2D advection equation was found to differ based on the numerical technique used in the calculation see Arakawa (1966)

Answer 1

From the code you provided, you used finite difference approximation to the differential expressions. That is the root of the difference. When you use finite difference schemes (here centered differencing), you have implicitly neglected the higher-order terms in Taylor series.

Generally, a differential term can be sometimes written in different forms. These different forms give somewhat different results when they are discretized. A well-known example is: $$u\frac{du}{dx}=\frac{1}{2}\frac{du^2}{dx}$$ However, if you using centered differencing, you cannot prove: $$u_i\frac{u_{i+1}-u_{i-1}}{2\Delta x}=\frac{1}{2}\frac{u_{i+1}^2-u_{i-1}^2}{2\Delta x}$$

because higher-order terms are dropped out. So you have to properly choose your expression for discretization. For diagnostic purpose, such difference can generally be ignored, as it is due to your finite difference scheme.