1.前言 数值求解海洋与大气方程组可以看做是初边值问题。即给定初始条件和边界条件，方程封闭，即可求解。因此，如果想要得到精准的结果，需要保证满足两点，其一是控制方程能很好的描述海洋大气的物理运动机制，其二是需要初边值条件精确。从前面讲过的有限差分方程离散的形式可以看到，如果初始条件和实际情况有误差，这个误差会随着时间迭代一步步的积累下去。这样是我们看天气预报，能发现短期预报通常都很准，而长时间预报经常不准。不仅仅是初始条件的误差，整个模式系统不能收敛的误差都会随着时间迭代逐渐累积下去，使得结果变得渐渐无法接受。
$ Q2 _; `5 h0 j4 n- K 资料同化技术，可以一定程度上解决这一问题。资料同化，就是利用观测资料来修正模式结果，使得模式结果不偏离观测的结果。比如，针对上述问题，如果时间迭代每进行一定步数，就用以掌握的观测资料去修正结果的偏差，就能很好的解决上述问题。不仅是对模式结果的修正，在实际预报中，如何给出准确合理的初始条件（初始场）也是一个需要解决的问题。资料同化技术对模式预报的结果有显著影响，是业务化预报中不可缺少的一环。6 e- H! e p9 {6 C+ N8 a8 B
如今，在业务化预报中数据同化技术发展成熟，基本是以3DVar（三维变分同化）、4DVar（四维变分同化）和集合卡尔曼滤波为主。要想理解3DVar和4DVar，就需要理解贝叶斯理论和极大似然估计，以及一些基本的变分学理论。若想理解集合卡尔曼滤波，就需要懂得最小方差估计等理论的思想。这些方法需要有充分的数学基础，才能更好的理解。因此，本文计划从同化的传统方法讲起，分几篇文章，试图用相对易懂的语言来讲解这些数据同化技术。在最后，从全局的角度分析一下这些方法的相同与不同。( U! I; `0 a! J5 Y1 n" y. H h
如果想深入了解气象资料同化，推荐《资料同化理论和应用（上册）》这本书，关于该书的详细信息见文末参考资料。本文部分内容参考该书，因此沿用了该书的符号表示。" x& @& h2 P G& {, k
2.最小方差估计 最小方差估计是本文接下来要讲到的逐步订正法和最优插值法的理论基础。之所以在这些方法的视角下，资料同化是统计问题，是因为我们不知道真值是什么，只能通过观测值和模式结果值去估计真值，从而把模式计算的结果趋向于真值方向修正。
4 s9 l" k( Q% V. w4 |: ` 假设一个真值 ，不妨视为某时刻某个地区的温度值，而对这个地区的温度有 个观测资料， 。首先需要假设这 个观测资料的平均误差是0，误差之间不相关。
3 W. z& J# z e5 S7 ]
8 U r/ O2 K z 即可得到线性估计的形式如下。即使用这 个观测资料线性组合出一个分析值 。# c2 b, J( F9 L' ]( u. V
J7 W8 T0 ~% I& z* ]
2 t+ ]/ d! l! f5 y S3 \) @( {% d, H2 N8 k- U
因为要满足无偏条件，因此系数的和为1。% L/ \- Y T+ R
/ G! z4 k" J3 f* O6 y' }* T" J& Q3 T8 h# @6 s8 x, c& K
紧接着就可得出分析值的误差方差和观测值的误差方差的关系。5 W* a* y) C% u. x+ H; O, ^7 O
" {, d' G$ v( \. v6 K7 _4 C z9 m
$ s$ Z" b$ P, f& g3 v 这样我们就得到了一个有约束的优化问题，形式如下。在满足无偏的条件下，使得分析值的误差方差尽可能小，这样分析结果就变得更可能精准。
8 F6 X9 w8 ^/ t0 w& A% e4 b( A5 O q9 h
. h$ D% ?0 H& K8 b
在高数中，面对求条件极值的问题，学过一种方法叫拉格朗日乘数法。把上述优化问题即可从有约束优化问题转化成无约束优化问题，求该拉格朗日函数的极小值即可。
; T# }. ~" l4 C
% o$ N( P p% j- P+ x
0 S- V7 |8 w- I* B* } 求极小值的方法就是使得其对各个量的偏导数同时为0。4 \2 ^$ s1 h6 ?1 C$ k L
- A! b4 }3 g: T' [
. O9 r2 X* L! i3 v) Z 这样就求出了各个观测值的权重值，即
# H h$ Q' ~7 t$ |, G$ q2 ]* ?& q5 \1 z& r: {- A+ O
3 ]7 ?; m- c3 P& R: L
如果想分析结果的优劣，可以推导出分析值的误差如下。
; j2 B2 C) Y, S, e8 b( y
$ O: {4 b3 Q7 r D& P; Y: W' V) A5 ^5 {( O2 |' s. w' ^* r
从上式可以看出，每增加一个观测点，就算其很不精确（即方差很大），分析值的误差方差也是朝着变小的方向发展，而且误差方差很大的观测值会被给予很小的权重。5 a" Z0 G2 M/ y$ U0 K% g% n
然而，理想很好，但现实却总是和理想有偏差。因为无法知道真值，进而无法计算背景场的误差方差，从而无法精确得到加权系数。因此，针对不同类型的问题，产生了经验权重函数。9 G; b2 s( W H/ ~" ^& V" r2 K$ j- j
3.逐步订正法 把上述思想应用到具体问题，就可以得到如下模型。) A/ Z4 F2 B' A* G: I! _4 Y
6 O( x6 V& f* W. K- G
/ H5 l3 ]( k. x5 H1 r2 o. ^ 该形式是由最小方差估计推导而来，其中 为分析场， 为背景场。不妨把背景场理解为待同化的场，比如用模式计算的某一步的结果，用这一步的结果结合 观测场去得到一个更准确的分析场 ，以此来修正计算结果。可以看出，该方法的思想是，使用背景场和周围的k个观测点的值与背景场之差加权求和得来。1 x2 k0 ?2 H; L+ `/ D! \
4 \% q. K9 H/ S 上面是从场的角度说明该同化方法的思路，为了更便于理解，从场中某一个点来描述该方法的思路。为了同化某一点的值，首先需要定义一个影响半径 ，意思是要使用在这个半径范围内的观测点来进行同化，比如要同化东方明珠处的气温值，影响半径可以设置为上海这么大，这样就只需要借助上海几个气象站的观测值进行同化即可，而太远的观测点的观测值则意义不大，而且即使使用了权重也会很小。将影响半径 内的总观测点数计为k。两种经验权重函数如下所示，这些权重都是随着观测站的距离变远而减小的。% ~. N& P R+ O$ `" `! x5 N
/ Y, _5 ?$ C J, M. |5 l# V
1 Y' b$ Y' o9 ?( u L) `
$ Z' H2 l7 m" \ K1 m 综上，该方法的计算公式如下所示。
" l9 |9 O+ h$ p- X# j, j- h4 r, n6 Q1 M
6 s* q/ E5 A- ]- M- J 总结一下，该方法为单点分析方案，即在同化不同点的过程是独立的，只与自己周围的k个点做计算。而且是背景场和观测值的增量做加权平均。
2 e/ \" Z, L1 W/ w! j4.最优插值法 该方法的同化模型和上述逐步订正法的形式一致，只是在计算权重时不再使用经验函数。写成矩阵形式则如下所示。还是延续上文所使用的变量记法， ; K9 L2 u2 G6 i, L/ B/ Y
+ L7 K+ W; k+ @6 Y y% h9 z" h: N- p6 ^% L7 M" b% L4 ?
6 z/ v# |0 i0 A5 U0 Q" c 背景场协方差矩阵如下- e$ I- K, x8 _; g, E
观测场协方差矩阵如下- _$ \6 d' o0 v3 E, U
第i个背景场协方差的向量如下) Z1 X* f1 t5 u3 l- y
使用最小方差估计的方法可推导得到如下所示的权重函数表达式，若对推导过程感兴趣可查阅文末参考资料。
6 Q6 a# C* S& o/ Q) U; G4 D! h3 f8 n3 w2 k7 ?- K
& i# X* g6 R0 p* s8 m
因不知道所谓的真值，因此无法准确快速的求解协方差矩阵。因此，最优插值法的主要难点在于协方差矩阵的估计。估计协方差矩阵因特定问题而异，本文只对其基本思想简单概述，想了解更多信息可查阅相关资料。
+ y* ?# x4 J* ^& |; d( ] 集合卡尔曼滤波和变分同化方法所涉及的内容较多，之后会单独更新几篇文章来对其进行介绍。7 n2 e$ G) |7 E$ m. O9 R) a* f% d
版权声明 本文创作的初衷是用于帮助数值模式的学习者。欢迎转载，转载请私信并注明作者和出处，请勿用于任何商业用途。
- o$ d: x, B# f. Q参考资料资料同化理论和应用（上册）. 邹晓蕾著. 气象出版社, 2009. }4 F8 o7 O& @) `' O6 D8 k
大气模式、资料同化和可预报性. Kalnay E, 蒲朝霞译. 气象出版社, 2005.2 R3 @0 u. B. x" t- L
e, W3 a L) D+ A1 j) |# Z
|
) K* u$ r3 `! `; c9 F* k: C
& O; K) x7 |4 [9 x
" W: c% C3 ^2 F0 L