引言

近年来，随着我国天气气候数值模式分辨率的不断提高(许晨璐等，2017；宫宇等，2018；牛嫣静等，2018)，对于计算机计算能力要求也日渐增长。为了满足模式运行时效，需要加强气象高性能计算能力建设(魏敏，2015a；2015b)，但目前由于通用CPU在降低功耗、提高主频上遇到瓶颈，而作为加速部件的GPU则表现出出色的浮点计算能力、低能耗高带宽、可靠的并行架构等方面的特点(Yang et al，2014)，并具备CPU相同种类的并行性，如能够采用合适的优化方法，将数值模式程序以CPU为主体执行逻辑运算部分、GPU执行并行运算部分的方式进行改写，可大大提升数值模式运行效率(王明清等，2017；桂彦等，2016)。2018年6月全球超算top500榜单中，越来越多的超算中心采用了异构架构，其中绝大部分使用的是CPU+GPU异构体系，充分展现了异构架构在高性能计算领域强劲的发展势头。因此，通用CPU中运行的模式程序如何能够适应异构架构系统，做好异构系统下模式优化，提升运行效率，已成为天气气候数值模式未来发展的重要课题。

BCC_AGCM3.0是国家气候中心开发的第三代大气环流模式，模式框架是在美国通用大气模式CAM3.0基础上发展而来的(顾聪等，2015；方宝辉等，2015)，引入质量通量型积云对流参数化方案、整层位温守恒干绝热调整方案、积雪面积覆盖参数化方案，并调整洋面感热和潜热通量参数化方案(Wu et al，2008；2010；吴捷等，2017)，支持共享内存及分布式并行架构，内部通过物理-动力耦合实现并行策略，是气候应用研究及气候预测的重要参考工具(颉卫华和吴统文，2010；段春锋等，2017；卢国阳等，2017)。本文在BCC_AGCM3.0原有并行程序基础上，选取计算量占比较大的程序段展开优化研究，针对GPU架构调整程序结构，改写优化程序段部分代码，达到提升程序段运行效率的目的。

1 研究环境简介 1.1 试验环境

GPU出色的计算能力，来源于其芯片架构组成。GPU芯片中的绝大部分为算术逻辑单元ALU，逻辑控制单元仅占一小部分。同样芯片面积的CPU，不仅包括占用面积较大的缓存区，还有大量的逻辑控制单元，承担计算任务的ALU只占小部分面积，因此通用CPU能够较好地处理复杂逻辑控制的程序，但在处理大规模科学计算程序时受其硬件条件限制，无法发挥优势。如果GPU中运行的程序能够做到取消分支预测、乱序执行等逻辑判断，充分进行大批量的计算，则相比CPU能够获得较好的性能提升(蒋华等，2016)。

本文选取NVIDIA Tesla K20 m作为GPU加速设备。优化试验环境的设备、硬件型号及配置信息如表 1所示。

1.2 OpenACC

OpenACC是一种编译器导语，它可以向编译器提供如何进行编译、如何优化代码的相关信息(孙晨等，2018)，在不同的高性能计算架构中具有高度的简便性、高性能和良好的可移植性(顾文静等，2018)。OpenACC的核心功能是能够把计算部分和数据从主机设备卸载到加速部件上，CPU(主机设备)与GPU(加速部件)都拥有各自的内存空间(江霞等，2016；曾文权等，2013)，在此种情形下，编译器在OpenACC导语编译运行时会分析导语意图并对主机设备及加速部件之间的数据部分进行管理，加速部件获得计算所需初始数据后进行计算，计算完成后结果数据返回主机设备中，并准备进行下一次计算(何沧平，2016)。

GPU加速器架构具有多个并行层次，以NVIDIA GPU为例，其内部拥有多个流处理器，每个流处理器各自包含了相应的线程块，线程块则由32位SIMD(单指令多数据)线程组构成。OpenACC中设计了gang、worker、和vector指导语句，它们分别对应了上述三层并行结构，能够把程序代码中的循环迭代映射到对应层次上。gang语句指示循环迭代部分将分成多块并行执行，worker语句指示在每个gang内使用多个worker并行执行，而vector语句则对应线程块内的SIMD线程组指定向量宽度，一般SIMD线程组的宽度为32，可设置vector的典型向量宽度为32的整数倍(罗布·法伯，2017)。

1.3 优化程序段选取

为了能够体现加速优化后的程序段对提高气候模式整体运行效率的促进作用，一般采用模式中运行耗时占比较大的程序段进行相关工作(肖洒和魏敏，2018)。BCC_AGCM3.0大气环流模式结构中主要依赖动力框架及物理过程部分的计算，且动力框架部分为模式核心部分，更适合进行程序移植，因此本文集中对模式动力框架中的程序段进行优化改写。通过对模式积分3 d后，统计出模式部分相关程序段耗时占比情况如表 2所示。

基于上述分析，按照模式内部运行时间占比排序，选取动力框架部分的grcalcs、grcalca及quad进行GPU加速优化，三段动力框架程序段的功能分别为通过逆勒让德变换，使用对称与反对称方法，将谱空间转换为格点空间和通过高斯积分完成波数至谱系数的转换。三段程序段中，逻辑判断内容较少，与计算相关的多重循环数量较多，循环次数较大，适于进行GPU加速优化。

2 优化方法

气候模式的运行加速效果主要来自于其计算密集部分的GPU并行化，通过采用OpenACC语言向编译器提供编译指令信息，指示编译器对模式代码进行并行化、数据传输、指定运行先后次序等操作，使模式代码在运行时达到性能提升目的。

要保证GPU能够正确计算出模式结果，需要将模式计算中所需的相关数据及时更新至GPU上，而此过程需要耗费一定时间，所以在CPU与GPU协同运行情况下，如何处理主机和加速设备之间的数据传输是能否提升模式程序段运行效率的关键问题。从宏观上看，加速设备数据采用由上而下、计算采用从下往上的管理方式是最好的，即在模式程序段外层的高级别位置创建和删除设备数据，在内层低级别位置进行计算并行化，这需要熟悉模式中过程间的调用流程和计算位置，保证模式计算结果正确前提下修改代码。图 1为加速设备端与主机计算在模式程序段内的生命周期。

2.1 异步执行策略

在CPU与GPU协同运行模式下，模式内程序的执行次序默认为串行执行，即加速设备→主机→加速设备→主机。在此执行次序下，加速设备与主机均需等待对方计算完成后方能开始计算，此种执行方式在模式内程序存在数据前后依赖时能够确保计算结果的正确性。但对于相互之间独立的模式程序段，通过使用异步执行方式，使加速设备与主机端同时进行计算，既能保证计算结果的正确，又能提升程序运行效率。

以BCC_AGCM大气环流模式动力框架部分grcalcs程序段为例，其内部包含4个主要循环段，其中循环段1与循环段2无数据依赖，循环段3、4之间有数据依赖，且循环段3依赖循环段1和2的计算结果。若按照任务串行方式执行，则为循环段1→循环段2→循环段3→循环段4流程运行。而将无数据依赖的循环段1和2设置为异步执行方式后，执行方式则如图 2所示，循环段1和2同时执行，相比前一种方式减少了运行时间。

2.2 循环内移

在优化过程中，需要将模式程序中计算量较大的部分放到GPU上进行并行化展开，利用GPU强大的浮点计算能力将计算耗时降至最低，同时尽量减少主机端与加速设备端的数据传输次数，避免不必要的时间损耗，程序的加速效果方能最佳。但分析动力框架格点空间转换部分grcalcs和grcalca函数段，存在外层函数多次循环调用情况，当对grcalcs和grcalca函数段内部计算循环段进行GPU并行化展开后，测试显示，由于外层循环调用次数较多，在每次调用时均会对GPU进行内核启动、数据同步等一系列操作，而这些操作耗时严重。

为了缩减这部分耗时，扩大计算并行度，采用循环内移方法，将外层函数调用的循环移至函数调用内部，并针对函数内程序内容，分别进行数组扩维、修改传递参数、精简数据传递、优化循环计算方式等操作，修改循环层及计算相关部分数组(图 3)。

2.3 OpenACC并行化

在循环内移后，对各个计算密集部分的多层循环段进行OpenACC加速试验。如果直接对各个循环段并行化，运行时会导致模式程序段性能下降，这是由于模式程序段中缺少数据管理导致。在每个OpenACC并行化计算开始时，编译器会默认将计算时所需数组变量从CPU拷贝至GPU中，在计算结束后将相关变量从GPU拷贝回CPU，这导致大量时间损耗在传输过程。为减少传输过程的时间损耗，在循环段计算开始前，设置加速数据区域，将循环计算所需要的数组通过!$acc data create语句在GPU加速器上初始化并分配空间，之后传输相应数据至GPU中，为计算做好准备。计算完成后将所需要的结果数据从GPU传回CPU中，便于CPU得到结果进行下一步运算，确保数据的完整性，避免数据的重复传输(图 4)。

对计算密集部分循环段进行加速，使用! ＄acc kernels/end kernels (莫德林等，2014)添加并行区域指示编译器并行化指定区域代码，为了减少CPU与GPU之间不必要的数据传输，使用present子句指定相应数组已存在于GPU中，无需再次重复传输。对内部每层循环进行向量化参数配置，循环段中gang、worker、vector根据实际试验调整为线程块数与每个线程块内线程数的最优组合，其中vector参数根据循环层中变量plev×8的值来确定，GPU在每次运行时向量长度包含循环迭代次数且尽可能减少向量浪费，运行效率为最佳，此处设置vector为256综合考虑了GPU中向量总量及运行效率情况(图 5)。

3 试验结果及分析

将优化后的模式程序段编译运行，测试使用16tasks，积分1 d，得到三段程序段的运行时间，与优化前相对比，效率分别提升3.3、5.8、3.1倍(图 6)。对模式进行结果正确性验证，测试运行1月，选取1985年1月7日的2 m温度做比较，quad程序段与标准结果完全一致，grcalcs与grcalca程序段存在细微误差(图 7)。使用模式输出中的全球涡度均方根数值进行误差对比，设相对误差=(测试试验结果-控制试验结果)/控制试验结果×100%(魏敏，2015a；2015b)，计算后的相对误差如表 3所示，三段程序段结果误差的绝对值均小于0.5%。

4 结论与讨论

(1) 结合数值模式适应未来异构计算系统的需求，本文对气候模式BCC_AGCM3.0中动力框架部分耗时占比较大的三段程序段进行GPU加速优化试验，在优化过程中，采取CPU与GPU协同运行的优化策略，选取程序段中计算密集部分进行GPU加速优化，通过采用循环内移、OpenACC并行化展开、数据管理、异步设置及向量参数化配置等方法改写，使得模式程序段能够获得一定加速效果。

(2) 运行结果显示，BCC_AGCM动力框架部分程序段优化后计算结果与CPU版本存在一定误差，其中包含众多因素导致，CPU与GPU的计算架构不同是一方面，程序优化过程中出现的数组计算精度损失亦是一方面。由于气候模式是中长期积分，对精度要求较高，因此下一步工作计划在对模式深入优化的同时考虑保证尽量减少积分过程中的计算精度损失，尽可能提高气候模式预报结果的准确性。

(3) 在模式优化中，由于试验环境限制，未能进行大规模任务的模式运行试验，未来中国气象局建成的“派”曙光高性能计算机集群包含了配置GPU与KNL(众核CPU)的试验子系统，可为异构架构下模式移植与优化提供计算性能更佳、运行更稳定的试验环境。本项研究将依托新的GPU集群，进一步深入气候模式内部算法，使用多线程、OpenMP调用GPU等方式(洪向共等，2016)，改进优化策略，提升气候模式的预报时效性。