（计算机科学与技术专业论文）中国新一代全球数值天气预报模式切线性伴随模式技术.pdf

国防科学技术大学研究生院学位论文 摘要 四维变分同化技术是提高数值天气预报时效的重要方法，切线性模式和伴随模式是同 化系统的关键技术。本文针对我国自主研发的多尺度通用数值预报模式系统g r a p e s ，为 建立g r a p e s 全球四维变分同化系统的切线和伴随模式，开展了以下工作： 基于g r a p e s 全球模式绝热动力学框架和四维伴随变分资料同化理论，建立了 g r a p e s 全球模式的切线性模式( g 凡垤e st l ) 和伴随模式( g r a p e sa d ) ，切线性模 式精度高，伴随模式正确，满足业务要求。 设计了全球伴随模式动态数据存储库，采用断点存储技术，解决了g r a p e s 全球伴随 模式存储和恢复基态值的难点问题。 设计了全球半隐式一半拉格朗日中r i t c h e 算法在物理空间中连续的切线性伴随模式， 解决了r i t c h e 算法中切线性模式出现分母为零的问题；基于g c r 算法建立了h e l m h o l t z 方程的高精度切线性和伴随模式。 基于p e t s c 库建立了g r a p e s 模式中求解h e l m h o l t z 方程的k r y l o v 子空间方法和预 条件子算法库，设计了h e l m h o l t z 方程的切线性和伴随的p e s t c 框架。 关键词：g r a p e s ，切线性模式，伴随模式，代码检测，h e l m h o l t z 方程，断点存储， p e t s c 第i 页 国防科学技术火学研究生院学位论文 a b s t r a c t f o u rd i m e n s i o n a lv a r i a t i o n a ld a t aa s s i m i l a t i o nt e c h n o l o g yi s a ni m p o r t a n tm e t h o dt o e n h a n c et h ee f f e c to fn u m e r i c a lw h e t h e rf o r e c a s t ，t a n g e n tm o d ea n da d jo i n tm o d ea r et h ek e y t e c h n o l o g yo fa s s i m i l a t i o ns y s t e m i no r d e rt os e tu pt a n g e n ta n da d j o i n tm o d eo fg r a p e s g l o b a lf o u rd i m e n s i o n a lv a r i a t i o na s s i m i l a t i o ns y s t e m ，b e l o wr e s e a r c h e sw e r ed o n eb a s e do nt h e m u l t i s c a l eu s e dg r a p e s ( g l o b a la n dr e g i o n a la s s i m i l a t i o np r e d i c t i o ne n h a n c e ds y s t e m ) w h i c hw a sd e v e l o p e di n d e p e n d e n t l yb yc h i n a g r a p e s t la n dg r a p e s a do fg r a p e sg l o b a lm o d ew e r es e tu pb a s e do ng r a p e s g l o b a lm o d ea d i a b a t i cd y n a m i c sf r a m ea n df o u rd i m e n s i o n a la d j o i n ti n f o r m a t i o na s s i m i l m i o n t h e o r y g l o b a la d j o i n td y n a m i cd a t a b a s es t o r i n gl i b r a r yw a sd e s i g n e da n dc h e c k p o i n tt e c h n o l o g y w a sa d o p t e d ，w h i c hs o l v e dt h eo r i g i n a lv a l u e s s t o r i n ga n dr e s t o r i n gd i f f i c u l tp o i n t so fg r a p e s g l o b a la d j o i n tm o d e g l o b a ls e m i i m p l i c i ta n ds e m i - l a g r a n g i a na t m o s p h e r i cm o d e l sr i t c h ec o m p u t a t i o n s c o n t i n u o u st a n g e n ta n da d j o i n tm o d ei nt h ep h y s i c a ls p a c ew a sd e s i g nw h i c hs o l v e dt h ep r o b l e m t h a td e n o m i n a t o rc o m ef o r t ht ob ez e r oi nt a n g e n tm o d ei nr i t c h ec o m p u t a t i o n h e l m h o l t z e q u a t i o n sh i 曲p r e c i s i o nt a n g e n ta n da d jo i n tm o d e w a ss e tu pb a s e do ng c r d i f f e r e n tk r y l o vs u b s p a c em e t h o da n dp r e c o n d i t i o n a la l g o r i t h m sl i b r a r yo fh e l m h o l t z w a ss e tu pb a s e do np e t s c h e l m h o l t ze q u a t i o n st a n g e n ta n da d j o i n tf r a m ei sd e s i g n e di nt h e p e s t c l i b r a r y k e yw o r d s ：g r a p e s ，t a n g e n tm o d e ，a d j o i n tm o d e ，c o d ec h e c k i n g ，h e l m h o l t z e q u a t i o n ，c h e c k p o i n t ；p e t s c 第i i 页 国防科学技术大学研究生院学位论文 表目录 表4 1g r a p e s 全球模式切线模式整体精度。31 表5 1p e t s c 支持的k r y l o v 子空间方法及预条件子3 8 表5 2 矩阵a 大小为1 9 2 8 5 * 1 9 2 8 5 时运行时间4 1 表5 3 矩阵a 大小为3 9 6 7 2 * 3 9 6 7 2 时运行时间4 1 表5 4 矩阵a 大小为8 2 0 8 0 * 8 2 0 8 0 时运行时间。4 1 表5 5g m r e s j a c o b i 在不同数据规模下的s n 值4 2 表5 6c r s o r 在不同数据规模下的s n 值4 3 第1 1 1 页 国防科学技术大学研究生院学位论文 图目录 图2 1 四维变分同化示意图5 图2 2 开发切线性模式的两种方法1 1 图2 3 开发伴随模式的三种方法1 1 图3 1 新一代多尺度通用模式体系结构图1 3 图3 2g r a p e s 模式程序简化关系结构图2 0 图3 3g r a p e s 模式软件框架示意简图2 0 图3 4s o l v eg r a p e s 关键子程序流程示意简图2 0 图4 1g r a p e s 全球模式切线性伴随模式开发流程2 1 图4 2 切线性模式代码的构造过程2 2 图4 3 伴随模式代码的构造过程2 4 图4 4 伴随模式计算流程2 4 图4 5 伴随模式全计算策略2 5 图4 6 伴随模式全存储策略2 6 图4 7 伴随模式断点存储策略2 6 图4 8s o l v eh e l m h o l t z 程序流程简图2 9 图4 9s o l v eh e l m h o l t zb 程序流程简图3 0 图4 1 0t ls o l v eg r a p e s 关键子程序流程示意图3 l 图4 11 子程序a 流程简图。3 2 图4 1 2 嵌套断点存储3 2 图4 1 3a ds o l v eg r a p e s 简单流程图3 3 图5 1p e t s c 库结构层次图3 7 图5 2g m i 迮s j a c o b i 在不同数据规模下的s n 图4 2 图5 3g m i 也s s o r 在不同数据规模下的s n 图4 2 图5 4b i c g j a c o b i 在不同数据规模下的s n 图4 2 图5 5c r j a c o b i 在不同数据规模下的s n 图4 3 图5 6c r s o r 在不同数据规模下的s n 图4 3 图5 7h e l m h o l t z 方程p e t s c 库切线性模式框架4 4 图5 8h e l m h o l t z 方程p e t s c 库伴随模式框架4 4 第1 v 页 独创性声明 本人声明所呈交的学位论文是我本人在导师指导下进行的研究工作及取得 的研究成果尽我所知，除了文中特别加以标注和致谢的地方外，论文中不包含 其他人已经发表和撰写过的研究成果，也不包含为获得国防科学技术大学或其它 教育机构的学位或证书而使用过的材料。与我一同工作的同志对本研究所做的任 何贡献均已在论文中作了明确的说明并表示谢意。 学位论文题目：生国堑= 垡全壁熬焦丞氢亟塑搓惑塑垡丝佳随搓盛垫盛 学位论文作者签名：麴! 里歪日期：蝴年占月如日 ， 学位论文版权使用授权书 本人完全了解国防科学技术大学有关保留、使用学位论文的规定。本人授权 国防科学技术大学可以保留并向国家有关部门或机构送交论文的复印件和电子 文档，允许论文被查阅和借阅；可以将学位论文的全部或部分内容编入有关数据 库进行检索，可以采用影印、缩印或扫描等复制手段保存、汇编学位论文 ( 保密学位论文在解密后适用本授权书) 学位论文作者签名：勤旦至 作者指导教师签名：趣苎! ! 垒 日期：加易年占月孑。日 日期：加8 年占月弓口日 国防科学技术火学研究生院学位论文 第一章绪论 1 1 课题背景 气象数值预报是在当代大气科学理论成果、当代计算机技术与大气遥感技术基础上建 立起来的现代气象预报技术，是气象预报向定量、定时、定点方向发展，满足国家与社会 的各种需求的科学途径，对于国家经济建设与社会发展有重要意义l l j 。 西方发达国家如欧洲中期天气预报中心、美国、英国、法国、加拿大、日本等国家和 组织在进入新世纪以来，开始投入了比以往更多的人力、财力、物力，致力于加速发展数 值天气预报业务体系，并取得了显著的成效【2 j 。 准确的数值天气预报需要好的初始场。资料同化就是为数值天气预报提供一个良好的 初值。随着数值预报模式分辨率和准确率的逐步提高，资料同化技术也从简单的观测资料 插值、最优插值发展到现在的三维变分资料同化( 3 d v a r ) 和四维变分资料同化( 4 d v a r ) 。 在目前的业务化中，最常用的是o i 、3 d v a r ，在国外已实现4 d v a r 的业务化【3 j 。 在变分同化的理论研究方面，早在1 9 5 8 年，我国学者顾震潮就提出应在数值预报中 考虑历史演变信息的观点【4 】，1 9 7 4 年丑纪范详细论述了天气预报中使用过去资料的问题瞄j ， 并提出把初值问题改为历史演变问题的观点【5 j 。1 9 7 1 年m o r e l 提出了将全球模式的解向所 有的有效观测进行调整的方法【6 1 ，利用该方法对同化模式交替向前和向后积分，把观测不 断地加入到模式中，当观测资料足够精确时，这个解是协调的。1 9 8 1 年t a l a g r a n d 利用稳 定性理论确定了这种方法收敛的依据，解释了以前许多人所得到的大量实验结果。1 9 8 6 年 l ed i m e t 等在前人工作的基础上对四维变分同化进行了系统、严格的定义，引入了表示模 式方程的解和观测资料之间距离的目标函数概念，c o u r t i e r 等在p e n e k o 等人工作的基础上 系统地给出在资料同化中利用伴随方程的理论。在9 0 年代早期，四维变分同化被认为是 最有发展潜力的同化方法，人们大力发展新的基于变分方法的预报系统，但由于计算过程 复杂、计算量太大，仅能在一些数值模拟试验上实现变分同化。而随着高性能计算机性能 的不断提高，特别是1 9 9 4 年c o u r t i e r 提出增量形式的四维变分同化计算方法【2 9 1 ，大大降低 计算时间，使得四维变分同化业务化成为可能。 2 0 0 1 年，国家科技部启动了“十五”国家科技攻关项目“中国气象数值预报技术创新 研究”【l 】，该项目的实施标志着我国气象数值预报研究开发战略从引进国外数值预报模式 为主到自主研究开发具有自主知识产权的先进模式系统的重大转变。我国目前数值预报的 总体水平无论与发达国家相比或与业务及服务的需求相比都有很大差距，更不能满足当前 国家经济建设与社会发展对气象服务提出的更新、更高的要求，如2 0 0 8 奥运所要求的精 细预报、国家安全与现代国防所要求的全球天气事件预报和社会发展所要求的各种高影响 天气预报。通过数值预报系统的技术创新研究，建立起我国具有先进水平的数值预报系统， 第l 页 国防科学技术大学研究生院学位论文 已成为国家一项十分紧迫的工作。 中国新一代数值预报系统( g r a p e s ，g l o b a la n dr e g i o n a la s s i m i l a t i o np r e d i c t i o n e n h a n c e ds y s t e m ) 利用高性能巨型计算机和高时空分辨率气象观测资料，充分吸收国内外 数值预报与相关学科的最新研究成果，能够提供时效更长、更详细更精确的范围涵盖气候、 中期、短期的天气预报，以及满足中尺度天气预报、空气污染和火警预报等越来越“精细 和趋向于多样化的各种气象预报服务需求，从根本上缩短我国与发达国家在气象科学技术 上的差距。 1 2 研究的必要性、目的与难点 建立起我国首个能够直接同化卫星辐射观测资料的变分同化系统和我国首个在业务 运行环境下可以实施的四维变分同化系统，对解决我国数值预报长期面临的资料不能有效 利用问题具有重大意义【l 】。g r a p e s 全球模式切线性伴随模式是实现g r a p e s 全球模式四 维变分同化的关键技术。 本课题研究的主要目的正是解决g r a p e s 全球模式四维变分同化的关键技术问题，实 现g r a p e s 全球模式的切线性伴随模式，为四维变分同化系统的实施提供一个可靠的基 础。 g r a p e s 模式是我国第一个完全拥有自主知识产权的多尺度数值天气预报系统，正在 不断完善和发展阶段，四维变分同化系统的开发将是一个自主创新的过程，必需针对全球 模式的特点创新设计相应的切线性伴随模式。 1 3 课题主要研究工作 研究四维变分同化理论，针对g r a p e s 全球模式绝热动力学框架的特点提出相应的切 线性伴随模式框架模型。 设计实现g r a p e s 全球模式的切线性伴随模式，并使用实际数据场资料进行精度检 测。 设计全球半隐式一半拉格朗日中r i t c h e 算法【l 】的在物理空间中连续的切线和伴随模 式。基于g c r 算法建立h e l r n h o l t z 方程【i 】的高精度切线和伴随模式。 对切线性伴随模式实现中的海量数据存取问题，分析设计g r a p e s 全球模式切线性伴 随模式的数据存储方法。 基于p e t s c 库建立g r a p e s 模式中求解h e l m h o l t z 方程的k r y l o v 子空间方法和预条 件子算法库，在p e s t c 算法库中设计h e l m h o l t z 方程的切线和伴随框架。 1 4 论文结构 第2 页 国防科学技术人学研究生院学位论文 本文共分为六章： 第一章绪论，主要介绍论文研究背景和意义，以及论文的主要研究内容和所做的主 要工作。 第二章全球四维变分同化理论，主要介绍四维变分同化的发展过程及其基础理论， 指出切线性伴随模式在四维变分同化中的关键性技术地位，接下来讨论了切线性伴随模式 的相关理论基础。本章为论文的研究理论基础。 第三章中国新一代数值预报系统，主要介绍g r a p e s 全球模式的相关特点：介绍从 原始的连续非线性数值预报方程到离散的非线性模式方程的过程，介绍g r a p e s 全球模式 的软件框架，分析g r a p e s 全球模式正模式程序的核心计算程序结构。 第四章g r a p e s 全球模式切线性伴随技术，分成以下内容 设计全球半隐式一半拉格朗日中r i t c h e 算法的在物理空间中连续的切线和伴随模 式。 基于g c r 算法建立了h e l m h o l t z 方程的高精度切线和伴随模式。 设计并实现g r a p e s 全球模式切线性模式。 设计并实现g r a p e s 全球模式伴随模式。 最后对切线性和伴随模式进行实际数据检测。 本章是论文的核心部分，是第二章理论在第三章实际问题上的应用，是对理论在实际 应用中许多技术问题的分析和研究。 第五章基于p e t s c 库求解h e l m h o l t z 方程的算法库，主要设计h e l m h o l t z 方程在不同 k r y l o v 子空间方法和多种预条件子下的算法库，分析多种方法的计算性能，设计h e l m h o l t z 方程的切线性和伴随模式框架。 第六章结束语 第3 页 国防科学技术火学研究生院学位论文 第二章四维变分同化理论 2 1 引言 资料同化1 9 j 是将观测资料加入到数值预报模式中去，以得到更精确预报结果的一种复 杂的数据处理方法。资料同化技术是数值天气预报的核心技术之一【l ，它是做好数值预报 的重要基础。目前在我军和国家气象中心的数值天气预报业务系统中，应用的同化方法是 基于最优插值的问隙资料同化技术，它是一种分析、预报的循环流程，同化是通过对预报 模式的每次积分、在逐段时间内进行一系统的订正来实现的。在每次订正或分析时，把由 模式得到的背景场和新的观测结合起来。基于最优插值的间隙同化方案存在着根本性的缺 陷：第一，动力约束和资料约束是分别进行的，未能同时考虑：第二，不同时刻的资料是 分别考虑的，彼此分割开了；第三，模式只是近似地描述了实际过程，但被按绝对准确来 使用；第四，一些非模式变量的观测资料必须转换成模式变量才能使用，而这种转换常常 有较大的误差；第五，一些非模式变量的观测资料由于难以转换成模式变量而没有使用； 第六，不规则、不均匀分布的测站资料必须显式地转换成格点资料，误差较大；第七，有 时初始化导致的改变，超过了资料的误差。 为了避免上述同化方法的缺陷，改进数值天气预报的预报质量，提高预报的准确率， 延长预报时效，必须大力开展四维变分资料同化方法的研究。四维变分同化【1 4 j 将资料处理 问题描述为以动力模式为约束的极小化问题，它利用最优控制原理，通过调整控制变量( 通 常是初始场) ，使得在指定的时间窗口内由控制变量得到的模式预报结果与实际观测资料 之间的偏差达到最小【2 9 1 。四维变分同化，被普遍认为是当前最有发展潜力的同化方法。它 最主要的优势在于能将不同时刻、不同地区、不同性质的气象观测资料，包括以前客观分 析方法中很难应用的卫星、雷达等非常规观测资料作为一个整体同时进行考虑，从而得到 满足高分辨率“全物理 预报模式要求的，质量场和流场基本平衡的、理想化的初始场， 从而使得预报精确度和时效大大提高。它从根本上解决了一直困扰着数值预报发展的可用 观测资料少、误差难以处理等问题。另外，随着卫星资料的逐步增多，人们希望找到一个 很好的分析方法以便能充分利用卫星资料，这是发展四维变分同化的最初出发点【1 7 j 1 8 o 缺 乏有效处理卫星资料的方法是发展四维变分同化的最显著的动机。 2 2 四维变分资料同化基础理论 四维资料同化就是把不同时刻、不同地区、不同性质的气象资料，纳入统一的分析预 报系统中来，通过统计与动力关系( 包括预报模式) 使之在动力与热力上协调起来1 0 1 i t 3 ) 。 而变分方法能很自然地描述观测资料和模式要素之间的非线性和非直接的关系。因而，在 第4 页 国防科学技术大学研究生院学位论文 此基础上发展起来的四维变分资料同化方法，充分利用过去和将来的常规资料和非常规资 料，使各气象要素自然满足要求的协调条件，从而求得质量场和流场基本平衡的、理想化 的初始场；或者用来更新预报值，以得到质量更好的连续天气分析图及预报酬2 2 1 。 场值 初始场 校正场 梯度 迹一一 1 广、j 伴随模式反向积分轨迹 图2 1 四维变分同化示意图 随着计算机速度的加快和内存容量的增大，特别是1 9 9 4 年c o u r i e r 提出增量形式的四 维变分同化计算方法【2 9 1 ，使得全球尺度上的四维变分资料同化变得可行，并在e c m w f 的 i f s 系统中实现了业务化。它的优越性主要包含： ( 1 ) 充分利用各种气象( 不同时刻、不同地区、不同性质的气象资料) 资料i l 。四 维变分资料同化可用预报模式提供的信息来弥补海洋、高原等地带常规资料的不足，反演 现有观测场无法获取的某些场地( 如青藏高原、沙漠、海洋等地区) 信息。而且，可利用 不定时的卫星观测资料、雷达资料等非常规资料来弥补常规资料的不足。更重要的是它不 但可以通过包含背景场估计而利用过去的观测资料信息，而且充分利用在同化周期内不断 接收到的新资料( 或叫未来资料) 作为约束来提高分析质量。 ( 2 ) 将动力约束和资料约束纳入同一个方程。四维变分资料同化通过由动力模式和 观测资料构造的泛函极小化来实现。它使得动力模式和观测资料之间建立起一种可信又客 观的联系，使两者可以趋向于双向匹配：由动力模式同化观测资料，而由观测资料优化动 力模式的要素场和某些参数。因为它是用预报模式作为约束条件，所以同化效果更好。另 外，由于这种同化的结果是来自预报模式本身，因此，它可直接输入模式进行预报，而不 需要再进行初值化处理1 2 5 。 ( 3 ) 有效地控制多种误差。四维变分资料同化通过引入背景误差协方差矩阵、观测 资料误差协方差矩阵、预报模式协方差矩阵来有效地控制各种误差来源。背景误差协方差 矩阵主要包含用其它方法得到背景场时产生的误差信息，从而控制输入模式误差。观测资 料误差协方差矩阵主要包含观测资料中的误差信息，如各气象要素的观测精度都不一样， 常规资料精度高、而非常规资料精度要差些。预报模式协方差矩阵主要控制由于预报模式 第5 页 国防科学技术大学研究生院学位论文 模拟实际物理过程的不足而产生的误差。 2 2 1 最大似然估计方法 四维变分资料同化实际是个优化问题。求解优化分析问题主要有两种方法。一是基于 最小方差估计的方法，它是优化插值分析的基础。另一个是最大似然估计方法。它为变分 分析提供了理论框架。对于线性问题，两者可导出相同的结果。最大似然估计方法【4 7 1 如下： 我们希望由给定的观测资料y 寻求最大概率的模式场x ，也就是用b a y e s 定理发现最 大条件概率p ( x i y ) ： p ( x y ) 。p ( y x ) p ( x ) 其中p ( y x ) 是在给定模式场x 情况下得到观测资料y 的条件概率。 观测前x 的先验概率。如果我们假设误差是呈高斯状态分布的。则： 尸( y i x ) e x p - 1 2 ( y 日( x ) ) r 0 1 ( y 一日( x ) j 而且 ( 2 - 1 ) e ( x ) 是在使用所有 ( 2 - 2 ) 尸( x ) e x p - 1 2 ( x x b ) 1b 一1 ( x x b ) l ( 2 - 3 ) 其中b 为背景误差协方差矩阵，o 为观测误差协方差矩阵。另外，还应考虑到预报误 差协方差矩阵，o 应为o + f ，f 为预报误差协方差矩阵。但在实际中往往假设模式预报 不存在误差，所以一般( 2 2 ) 中仅出现观测误差协方差矩阵o 。 j ( x ) 是一个对所有观测资料和模式场x 的物理或动力学约束的标量衡量标准，表示模 式预报结果与实际观测资料之间的偏差，具体表达式如下： j ( x ) = 1 2 ( x x b ) 。b 一1 ( x x b ) + 1 2 ( y h ( x ) 10 1 ( y h ( x ) ) ( 2 4 ) 2 2 2 四维变分资料同化的目标函数 在推导四维变分资料同化的目标函数时，先引入两个重要的假设【2 9 】： ( 1 ) 预测模式可表示成各中间预测步的乘积 假定x o 为模式初始场，在四维变分资料同化中也称为控制变量。x f 为由模式初始场x o 演变到f 时刻的模式预报场。 x f = x ( t f ) = m ( t o ，7 f ) x o( 2 - 5 ) 第6 页 国防科学技术大学研究生院学位论文 其中，m ( t o , f ) 为, k t o 时刻到f 时刻给定模式预报算子。 预报模式可被表达成中间预测步的乘积，它实际上是自然界因果关系的反映。通常， 它是从o 时刻模式初始状态x o 开始的数值预报模式积分。如果以f 作时间序列指标，用m f 表示从f 一1 到f 时刻的积分步，即m 户m ( t i _ 1 , 以) ，则x f = m f x f l ，依此类推可得： x i = m i m i 一1 m l x o ( 2 - 6 ) ( 2 ) 切线性假设 设x g ( f 。) 是初始时刻的初猜场，x 6 ( ，。) 是初始时刻的背景场。x ? 是根据背景场完全非 线性模式积分得到的f 时刻的模式场。x 罗是根据初猜场完全非线性模式积分得到的i 时刻的 模式场。则 x 尹= x g ( f f ) = m ( o ，f ) x g ( o ) ( 2 7 ) x ? = x6 ( f f ) = m ( f o ，t i ) x 6 ( f o ) ( 2 - 8 ) x 尹一x ? m ( f o ，t f ) ( x 孑一x 8 ) ( 切线性假设) ( 2 - 9 ) 其中m ( o ，t f )m ( t o , f ) 的切线性算子。h f 是f 时刻的观测资料算子，它把模式场x f 转 化成与f 时刻观测资料y f 类型相同、位置相同的观测资料。假设观测资料算子是拟线性的， 则 h f ( x 笋) 一h f ( x ? ) h ；( x g ) ( x g x ? ) ( 2 1 0 ) 其中，h ：是观测资料算子h f 的切线性算子。 四维变分资料同化是三维变分资料同化在时间上的扩展。把( 2 - 4 ) 在时间维上进行扩 展，就可得到四维变分资料同化目标函数j ( x ) 的标准公式： j ( x ) ：l ( x - x b ) 丁b 一1 ( x - - x b ) r + 1m ( 2 一l1 ) + 去( h f ( x i ) 一y ，r _ 1 ( h f ( x i ) 一y f ) 于是四维变分资料同化就转化为寻求一个模式解x + ，使得j ( x + ) = m i n j ( x ) ，并满足动 力学条件e ( x ) = 0 。这表明四维变分资料同化实质上是个条件极值问题。原理上，通过求 解相应的欧拉一拉格朗日方程是能够得到最小约束问题的解【5 1 1 。然而在气象观测资料同化 国防科学技术大学研究生院学位论文 中，这些方程如此之大和复杂以致仅仅在极其近似的情况下才能直接求解，这就大大限制 了在实际情况下应用这种方法的可能性。但是，如果注意到模式的解是由给定时刻相应的 初值唯一确定，上述变分问题就能够用下述方式来表示：寻找初始状态，使得目标函数最 小。那么问题就转变成无约束最优化问题，因为没有特殊条件附加在初始条件上。并且， 现在仅在初始状态下而不是在模式整个时段实现最小化，从而使得数值维数大大减小。于 是四维变分资料同化就转化为： j ( x o ) = i 1 ( x o x 8 ) 丁b - 1 ( x o x 8 ) r + 1m ( 2 1 2 ) + 去( h f ( x i ) 一y f ) 7 r _ 1 ( h f ( x f ) 一y f ) i = 0 求x ：，使得j ( x ：) = m i n j ( x o )( 2 1 3 ) 其中x o 是模式初始状态和控制变量。x 8 是在同化周期初始时刻f o 时的背景模式( 或 预报) 状态。一般由初始时刻的预报求得。它是对初始条件的优化以及对过去信息的综合 归纳的一个先验估计。b 是描述预报质量x 8 的预报( 或背景) 误差协方差矩阵。y f 是在f f 时刻的观测。r f 为7 f 时刻的观测资料误差协方差矩阵。 由上式可看出，四维变分资料同化主要是在给定时间内，寻求与观测资料拟合最好的 模式演变轨迹。在完全模式假设下，模式演变轨迹完全由初始条件x 。决定。但在实际应用 中，由于计算量过大以及物理参数不可微导致( 2 1 2 ) 式很难实现的。尤其是物理参数不可 微妨碍人们采用有效的最小化算法计算目标函数梯度。为了解决这些问题，一种增量形式 的四维变分资料同化的目标函数被提出来了【2 9 】，它避免了问题的不可微，并通过引入可微 的切线性模式和伴随模式计算得到模式初始状态的校正量，从而使得四维变分资料同化计 算在实际应用中可行。另外，一些气象中心也在大力发展在增量形式下可微的物理参数化， 以便提高分析质量和预报质量。 假设x o = x 8 + 血o ，并把其代a ( 2 1 2 ) 式，利用式( 2 - 6 ) 、( 2 - 7 ) 进行简化就可得到四维变 分资料同化目标函数的的增量形式： j ( 斑o ) = 寺苏o r b - 1 盘o + 寺( h ：斑( ，f ) 一d f ) 7 r _ 1 ( h ：6 x ( t i ) 一d f ) ( 2 1 4 ) i = o 求次：，使得j ( 斑：) = m i n j ( 6 x o ) ( 2 1 5 ) 其中df _ y ，一hf ( xf ) 。 2 2 3 计算目标函数梯度 第8 页 国防科学技术大学研究生院学位论文 四维资料变分同化的优越性有目共睹，尽管早在七十年代就已经提出了此方法【6 】，但 仅到8 0 年代中期，把伴随理论引入到气象学中计算目标函数梯度【5 1 】，才使得四维资料变 分同化能在实际中应用，由此可见伴随理论在四维变分资料同化中的重要性。 v j ( 斑叫= b - 1 8 x 0 + f ：。m * ( t o , t i ) 吖r h ：m ( t o , t i ) 8 x o - d i ) ( 2 - 1 6 ) ：v j b + v j d 其中，宰表示伴随，h ；+ 是h 的伴随算子。m ( 7 0 ，7 f ) + 是h 的伴随算子m ( 7 0 ，t i ) m ( ，t ，) = ( m ：m ：一。m i ) = m ，- m o m ：。 v j 。= e m ( 珀，h ；+ f i i = 0 = m i + m 汨a ( 2 1 7 ) ：h b + f o + m i + 【h i + f l + m ! + 【h ! + f 2 + + m ：+ h ：f n 】 f f = r i l ( h i m ( ，t ，) s x o d ，) 由此，采用梯度法或共轭梯度法就可求解以上的最小化问题。对于非线性程度高的数 据，最小化一般采用两个循环实现。内循环主要针对于二次方程式( 增量形式的目标函数) ， 在低分辨率模式上，围绕背景场作线性化，进行逐次迭代得到最佳的扰动值。外循环主要 是用内循环得到的最佳扰动值更新背景场，并用高分辨率模式作预报，并和观测值进行比 较。 2 3 切线性伴随模式理论 考虑一个动态物理系统以及用一个模式来描述这个系统 3 2 】。假设d r 研( m n ) 是 一系列的观测资料，而模式能计算与这些观测资料对应的值y r 朋。我们的目的是如何使 得观测资料和对应的模式计算值能吻合得很好。 函数： ，：要( y d ，y d ) 2 、 为了量化两者之间的差别，定义一个目标 ( 2 - 1 8 ) ( ，) 为内积。目标函数越小，模式计算值与观测资料拟合地更好。 为了使用模式，定义一系列g n 参数x 。x 为控制变量( 或初始场) 。x 与y 之间 的关系为 第9 页 国防科学技术人学研究生院学位论文 f ：r 甩jr 肌 x 专y 因此，可被表示成 j ：r 玎一r ( 2 - 1 9 ) x 专知x ) - d , f ( x ) _ d ) q 乏 因此我们的问题就转化为求一个控制变量x ，使得目标函数，最小。 ，( x ) = j ( x o ) + ( v x j ( x o ) ，x x o ) + o ( i x - x 0 1 ) 上式可简化为 彤= ( v x j ( x o ) ，) 对于控制变量x o ，y 的扰动可近似为 砰= a ( x 0 ) o x ( 2 - 2 1 ) ( 2 2 2 ) ( 2 2 3 ) 其中a ( x o ) 为f 在x o 处的雅可比矩阵。根据内积的对称性可得 6 , = ( f ( x o ) 一d ，a ( x 0 ) o x )( 2 - 2 4 ) 定义伴随算子么 ( 1 ，a w ) = ( 么+ v ，们 则 ( 2 - 2 5 ) 甜= ( 彳+ ( x o ) ( f ( x o ) 一d ) ，6 x )( 2 2 6 ) 把( 2 - 2 6 ) 式和( 2 2 2 ) 式比较，可得到目标函数关于控制变量的梯度。 v x j ( x o ) = 彳+ ( x o ) ( f ( x o ) 一d ) 其中线性算子彳( x o ) 一般称为切线性模式，4 ( x o ) 为伴随模式。 2 4 切线性和伴随模式开发方法 ( 2 - 2 7 ) 开发非线性模式的伴随模式必须先开发相应的切线性模式，现在开发非线性模式的切 线性模式主要有两种方法，而开发非线性模式的伴随模式主要有三种方法【2 1 1 。 2 4 1 切线性模式开发方法 第1 0 页 国防科学技术大学研究生院学位论文 开发切线性模式的两种方法如图2 2 所示： r 一j l 连续的非线性l 方程 方法一方法二 图2 2 开发切线性模式的两种方法 ( 1 ) 方法一：先由连续的非线性方程，推导出其连续的切线性方程，再对切线性方 程进行离散、编程，得到切线性模式。 ( 2 ) 方法二：直接对连续的非线性方程对应的离散非线性模式代码进行转换得到切 线性模式代码。 2 4 2 伴随模式开发方法 现有的开发伴随模式的三种方法如图2 3 所示 r 一j l 连续的非线性 方程 l 连续的切线性连续的切线性离散的非线性 方程 方程模式 上 j rj r 离散的切线性 离散的切线性 连续的伴随方程 模式模式 00上 离散的伴随离散的伴随离散的伴随 模式模式模式 方法一方法二方法三 图2 3 开发伴随模式的三种方法 ( 1 ) 方法一是直接从连续方程出发，先推导出其切线性方程，然后再推导出伴随方 第l l 页 国防科学技术大学研究生院学位论文 程，最后对伴随方程进行离散，建立模式。这对于简单模式是较好的，但对于复杂模式， 是很难实现的。主要有以下两个原因：( a ) 对于一个原始方程模式来说，用分步积分方法 推导伴随方程是相当繁琐的，而且，对各种物理参数化的处理也是相当困难的。( b ) 即 使准确的推导出伴随方程，对伴随方程的离散仍然存在难题。离散过程中使用的方法可能 与原始模式不一致，这会限制了梯度求解的准确性。 ( 2 ) 方法二是从连续方程出发，先推导出其切线性方程，然后对切线性方程离散、 编程，形成切线性模式，最后直接对切线性模式的代码进行转换得到伴随模式代码【3 0 l 。由 于数值天气预报模式的方程非常复杂，即使推导出其切线性方程，切线性模式的生成仍存 在很大的难处。最困难的是，如果连续方程出现一些变化，那么又得重新推导切线性方程， 工作量很大。 ( 3 ) 方法三是直接从非线性预报模式代码进行转换得到切线性模式代码，然后再直 接对切线性模式代码进行转换得到伴随模式代码岭引。由于一个预报模式是长的或复杂的， 可把其分解成m 步，每步表示一个子程序或一个单循环或一些代码的集合，把其通称为一 个代码段或一个操作步。在欧几里德空间里，对于欧几里德内积，在给定的输入输出情况 下，每个代码段存在唯一的伴随算子( 或伴随代码段) 。所以我们可以独立并行地开发各 个代码段的伴随代码，最后再以相反的顺序连接起来，形成伴随模式。 开发伴随模式的第三种方法简化了构造一个复杂模式的伴随模式。依照此方法，伴随 模式可以逐渐地建立起来。发展切线性模式各代码段的伴随模式代码可并行进行，对模式 的各个部分可分开独立的开发其切线性和伴随模式代码，由此建立的伴随代码容易检查， 效率高，易于调试。 2 5 小结 由于四维变分资料同化计算过程复杂、计算量太大，把伴随方法引入变分分析后，它 也仅能在一些比较简单的模式上实现。随着超级并行计算机计算能力逐步提高，更重要的 是【1 4 】提出增量形式，逐步提高线性化的精度，大大降低计算时间，才使得四维变分资料同 化业务化成为可能。该方法主要是在简化的物理模式上，尽可能准确地在演变轨迹邻近区 域内线性化，从而求解一个近似的二次方程问题，这保证最小化解的唯一性，加快了收敛 速度。为了发展增量形式的四维变分同化系统，首先得开发非线性数值预报模式的切线性 模式和切线性模式的伴随模式。对于复杂的三维全球预报模式，切线性模式和伴随模式的 开发仍然是个难题。 第1 2 页 国防科学技术大学研究生院学位论文 第三章中国新一代数值预报系统 自英国科学家r i c h a r d s o n 时代以来，数值预报模式动力框架的发展经历了“粗”原始 方程模式、简化( 滤波、地转平衡) 模式、原始方程模式( 正压、斜压) 的发展过程；自 1 9 5 0 年数值预报实际应用首次获得成功以来，世界上许多气象中心的业务数值预报系统发 展也经历了单一预报系统、全球与有限区域两套预报系统、多套预报系统( 短期、中期、 中尺度、台风、环境预报等) 并存、多尺度统一系统的几个阶段。当前，许多国家的气象 业务中心( e c m w f 、加拿大英、德、法等) 都建立程度不同的多尺度统一模式( 中尺度、 短期和中期天气预报、月季年气候预测，统一动力模式框架、统一“c o d i n g ”) 。另 外，随着高性能计算机能力( 峰值速度、内存容量、网络通信速度等) 的迅速提高，使得 气象数值预报模式的分辨率越来越高，物理过程越来越复杂、越来越完善；在新一代高分 辨率的动力模式框架设计中，大都考虑采用“全可压的 ( “f u l l c o m p r e s s i b l e ”) 非静力 平衡假设，认为当模式水平格距小于5 k m ，非静力平衡动力过程作用将是显著的、不可忽 略的。总之，多尺度通用模式成为各国新一代数值预报模式的共同选择，高分辨率的、物 理过程完善的多尺度通用模式不仅是业务技术发展的需要，同时也是降低业务运行成本、 密切科研与业务之间的合作、加速研究成果业务的现实需要i l 6 。 我国的业务数值预报发展也经历了世界上其它主要气象中心走过的发展阶段，目前正 处于“多套预报系统并存 的阶段。