python - 使用numpy和metpy的湿通量发散不同

Answer 1

如果我正确理解了您的代码，我认为前者的错误是您在计算通量发散时包含了术语 d(qu)/dy 和 d(qv)/dx，我认为这是不正确的。基本上水分通量散度只有d(qu)/dx + d(qv)/dy，所以我认为你需要更改代码中的行来挑选qu梯度命令的x导数和qv的y导数。在我之前的回答中，我混淆了轴，我尝试了离线测试，如果我理解正确，你需要这个：

dqu_dx = np.gradient(qu,axis=1)  # take d/dx of qu 
dqv_dy = np.gradient(qv,axis=0)  # take d/dy of qv
divg = np.add(dqu_dx,dqv_dy)

或更简单地说：

divg=np.add(np.gradient(qu,axis=1),np.gradient(qv,axis=0))

如果我制作一个虚拟测试数组并计算梯度，那么按该顺序取轴似乎是正确的：

a=np.array([[1,2,3],[2,10,15],[1,2,3],[2,3,4]])
np.gradient(a,axis=0)

给

array([[ 1. ,  8. , 12. ],
   [ 0. ,  0. ,  0. ],
   [ 0. , -3.5, -5.5],
   [ 1. ,  1. ,  1. ]])

因此，如果我不感到困惑，我认为我的轴陈述是正确的……

希望这应该给出与metpy结果相似但可能不完全重现的东西，因为在metpy中计算梯度的方法可能不同（例如，它们可能使用上游差异）。

最后一点，这只是等压水平水分通量散度，要计算总水分通量散度，您也需要垂直分量......（或者您可以计算可降水的总柱散度，但这已经可用作为“单级”字段组中 ERA5 的直接输出）。由于对流层低层的湿度要高得多，您也可以尝试绘制与 ERA5 的总柱湿度差异，并将其用作 850hPa 计算的代理参考。