球谐函数展开快速算法及其并行算法研究.docx
《球谐函数展开快速算法及其并行算法研究.docx》由会员分享，可在线阅读，更多相关《球谐函数展开快速算法及其并行算法研究.docx（29页珍藏版）》请在冰豆网上搜索。
球谐函数展开快速算法及其并行算法研究
球谐函数展开快速算法及其并行算法研究
球
谐
函
数
展
开
快
速
算
法
及
其
并
行
算
法
研
究
国
防
学
技
术
大
学
研
究
生
院
分类号TP312学号
09060065
,,密级
公开
工学硕士学位论文
球谐函数展开快速算法及其并行算法研究
硕士生姓名王翔
学科专业计算机科学与技术
研究方向计算机应用技术
指导教师宋君强研究员
国防科学技术大学研究生院
二〇一一年十一月
ResearchonParallelAlgorithmsoftheFast
AlgorithmforSphericalHarmonic
Expansions
CandidateWangXiang
SupervisorProfSongJunqiang
Athesis
Submittedinpartialfulfillmentoftherequirements
forthedegreeofMasterofEngineering
inComputerScienceandTechnology
GraduateSchoolofNationalUniversityofDefenseTechnology
ChangshaHunancoma
November2011
国防科学技术大学研究生院工学硕士学
位论文
摘要i
ABSTRACTii
第一章引言1
11数值天气预报与谱模式1
com谱模式的发展与现状1
com球谐函数与Silberman方法2
com变换法4
com谱模式的优缺点5
12球谐函数6
com球谐函数的推导6
com截断问题6
com球谐函数展开算法的发展7
13GPU通用计算与数值天气预报9
14本文研究内容10
15论文结构11
第二章CPUGPU异构集群研究12
21GPU发展的历史与现状13
22CPUGPU异构集群14
comCPU与GPU的性能对比14
com典型CPUGPU异构集群16
23并行编程环境19
comMPI并行程序设计原理19
comCUDA架构与编程模型20
comMPICUDA混合编程模型研究22
第三章球谐函数展开快速算法研究24
31球谐函数展开简介24
32Rokhlin-Tygert球谐函数展开快速算法25
comSerSTEP1函数分解25
comSerSTEP2快速Fourier变换25
com1余弦变换25
com2正弦变换27
第I页
国防科学技术大学研究生院工学硕士学位论文
comSerSTEP3低阶连带Legendre多项式展开系数27
comSerSTEP4任意阶连带Legendre多项式展开系数29
330阶Legendre多项式展开系数的计算31
com直接法31
com变换法32
com两种方法比较34
34标准特征向量矩阵的计算35
35RT算法串行程序实现36
第四章RT算法的MPI并行程序设计38
41RT算法的并行策略38
42RT并行算法描述40
43RT并行算法的MPI实现42
com程序结构43
com通信机制44
com同步策略45
44RT算法的MPI并行程序实验及分析46
com实验方案I47
com实验方案II48
com实验方案III49
第五章RT算法的MPICUDA并行程序设计51
51GPU在RT算法中的加速策略51
comGPU在RT串行算法中的加速策略51
comGPU在RT并行算法中的加速策略51
52库函数52
comCuFFT52
comLAPACK53
comCULA55
comCulaBLAS55
53RT算法的异构混合并行实现56
comRT算法串行程序的CUDA加速56
comRT算法MPICUDA异构并行实现57
54RT算法的异构并行程序实验58
com实验方案I59
com实验方案II59
第II页
国防科学技术大学研究生院工学硕士学
位论文
com实验方案III60
com实验方案IV62
com实验方案V62
第六章结论与展望64
致谢65
参考文献67
作者在学期间取得的学术成果71
附录A公式推导72
第III页
国防科学技术大学研究生院工学硕士学位论
文
表目录
表31四个系数矩阵31
表32求0阶展开系数两种方法比较35
表33求标准特征向量矩阵四种方法比较36
表34串行程序在不同规模下各部分时间对比36
表41MPI版RT并行程序文件列表43
表42MPI版RT并行程序通信函数列表45
表43各进程在示例条件下所得行号46
表44使用48进程实验方案47
表45使用24进程实验方案48
表46三种进程数实验对比49
表51CuFFT函数参数说明53
表52CLAPACK函数参数说明54
表53culaSsyev参数说明55
表54culaSgemv参数说明56
表55CUDA加速RT串行程序文件列表57
表56MPICUDA异构并行RT算法文件列表58
表57SerSHE与CUDASerSHE实验对比59
表58MPIParaSHE在三种进程规模下的实验对比60
表59MPICUDAParaSHE程序64节点实验对比61
表510MPICUDAParaSHE程序512节点实验对比62
表511MPIParaSHE与MPICUDAParaSHE程序对比实验63
第IV页
国防科学技术大学研究生院工学硕士学位
论文
图目录
[1]
图11两种波数截断方式7
图12北极附近的风速向量8
[27]
图21CPU与GPU的结构对比12
[27]
图22CPU和GPU峰值浮点计算能力比较15
[27]
图23CPU和GPU的带宽比较16
[36]
图24DEGIMA的互连网络17
[55]
图25天河一号组织结构图18
图26CUDA的软件系统结构21
图27CUDA编程模型22
图28MPICUDA并行计算模型23
图31Rokhlin-Tygert球鞋函数展开快速算法26
图41RT并行算法流程42
图42MPI版RT并行程序中文件间的关系44
图43ParaSTEP4程序示例45
图44划分通信组程序示例46
图45两种方案中的对比48
图46求特征向量矩阵时间增长49
图51三种进程规模的并行效率60
图52MPICUDAParaSHE在两种进程规模下加速比变化63
第V页
国防科学技术大学研究生院工学硕士学位论文
第VI页
国防科学技术大学研究生院工学硕士学位论文
摘要
计算机技术的迅速发展促进了数值天气预报模式的不断推进并使偏微分方
程离散化相关的数值方法变得更加完善和精确从而减小了求偏微分方程数值解
的误差谱变换方法采用球谐函数的有限截断展开来逼近球面某垂直层上的状态
变量充分利用变量在物理空间和谱空间之间的适当变换有效降低偏微分方程
中以非线性项为主的计算量进而获得高精度的计算解其中物理空间由状态变
量在经纬度格点的直积组成谱空间由状态变量对应的截断球谐函数展开系数组
成谱模式具有精度高稳定性好等优点但其运算量和存储量较大
球谐函数展开作为一种数学方法被广泛应用于数值天气预报模式地球物
理化学物理偏微分方程数值解等众多领域球谐函数是球面纬向微分算子以
及球面Laplace算子的特征函数基于球谐函数的谱变换是谱模式计算中的核心
问题涉及Fourier变换和Legendre变换其中Fourier变换应用于经度方向
Legendre变换应用于纬度方向
本文基于Rokhlin-Tygert球谐函数展开快速算法深入研究了算法的并行化问
题主要内容包括
1分析数值天气预报谱模式以及球谐函数的特点研究了并行计算理论和
技术包括MPI并行CUDA并行MPICUDA异构并行
2深入研究并补充完善了Rokhlin-Tygert球谐函数展开快速算法增加了
0阶标准Legendre多项式展开系数的计算步骤
3针对球谐函数展开快速算法设计了RT并行算法
4实现了RT算法的MPI并行程序并分别在某五万亿次集群和天河一
号超级计算机上进行测试其结果均显示在三角截断波数M
1023时并行
效率可达87以上
5实现了RT算法的MPICUDA异构并行程序调用了CULA库函数
天河一号上的实验结果显示CUDA的加速效果在M
1023时逐渐显著单
节点双进程布局的并行效率随的增大而增长并超过单节点单进程在
M
M4095
时使用512节点1024进程可获得超过6700的加速比
主题词球谐函数展开并行算法GPU通用计算CULAMPICUDA
第i页
国防科学技术大学研究生院工学硕士学位论
文
ABSTRACT
ThecontinuousprogressofNumericalWeatherPridictionNWPmodelis
promotedbytherapicdevelopmentofcomputertechnologywhichmakesmore
consummateandaccuratethenumericalmethodofdiscretizationofcontinuouspartial
differentialequationsPDEandthereforedecreasestheerrorofnumericalsolutionsof
PDESpectraltransformmethodmeansusingsphericalharmonicexpansionsSHE
withlimitedtruncationstoapproximatesomevariableofaspecificverticallayerof
spheresurfacetakingadvantageofappropriatetransformationofthevariablebetween
physicalandspectralspacetogainaccuratesolutionsanddecreasethecomputation
whichmainlyattributestothenon-lineartermsofPDEThephysicalspaceisformedby
Cartesianproductsofthevariableonthepointsoflongitudinal-latitudinalgridwhile
thespectralspaceisformedbythecorrespondingcoefficientsoftruncatedSHEofthe
variableSpectralmodelisadvancedwithitshighaccuracyandstabilitywhilewiththe
drawbackoflargemassofcomputationandstorage
AsamathematicalmethodSHEisimplementedinmanydomains
suchasNWP
modelGeophysicsChemicalphysicsnumericalsolutionsofPDEetcSpherical
harmonicsaretheeigenfunctionsoflatitudinaldifferentialoperatorsandLaplace
operatorsonspheresurfaceBasedonsphericalharmonicsthespectraltransformisthe
kernelofthecomputationofspectralmodelconsistingofFourierTransformand
LegendreTransformwhichisappliedinlongitudinalandlatitudinaldirection
respectively
BasedonRokhlin-TygertsfastalgorithmsforSHERTAlgorithmthis
dissertationdeeplyinvestigatestheparallelizationofthealgorithmmainlyincluding
1AnalyzethecharacteristicsofthespectralmodelofNWPandspherical
harmonicsandstudythetheoriesandtechniquesofparallelcomputingincludingMPI
parallelismCUDAparallelismMPICUDAheterogeneousparallelism
2DeeplystudyandconsummateRokhlin-TygertsfastalgorithmsforSHEto
whichaddingaprocedureforcomputingthecoefficientsof0-ordernormalized
Legendrefunctions
3DesignRTparallelalgorithmaccordingtothefastalgorithmsforSHE
4ImplementMPIparallelprogramforRTAlgorithmwhichwastestedonthe
5-TerascaleClustersandTianhe-1Asupercomputerrespectivelytheresultsofwhich
showthatwhenthewavenumberoftriangulartruncationsM1023theparallel
efficiencyoftheprogramwouldreachover87
5ImplementMPICUDAheterogeneousparallelprogramforRTAlgorithm
usingCULAlibraryfunctionsTheresultofthetestonTianhe-1Ashowsthatfirstly
CUDAaccelerationwouldbeeffectivegraduallywhenM1023secondlythe
第ii页
国防科学技术大学研究生院工学硕士学位论文
parallelefficiencyof2-process-per-nodedistributionwouldincreaseasMgrowsand
eventuallytranscendthatof1-process-per-nodedistributionthirdlywhenM4095
thespeedupofusing1024processeson512nodescouldreachover6700
KeyWordsSphericalHarmonicExpansionsParallelAlgorithmGPU
GeneralComputingCULAMPICUDA
第iii页
国防科学技术大学研究生院工学硕士学位论文
第iv页
国防科学技术大学研究生院工学硕士学位论文
第一章引言
11数值天气预报与谱模式
com谱模式的发展与现状
1922年Richardson首次提出的关于数值天气预报Numerical
Weather
Prediction的设想在1950年由Charney等人在普林斯顿的高等研究所的ENIAC
[1]
计算机上实现他们应用正压过滤模式第一次成功地制作了24小时的数值预报
全球或半球范围的的大气数值模式的进展一方面表现在逐步减少模式方程组中
所取的各种近似改进模式中所包含的各种物理过程的参数化方案使之能够更
准确真实全面地反映实际大气中的各种动力和热力过程另一方面也表现在
模式的水平和垂直分辨率在不断地提高数值天气预报模式的进展与计算机技术
的迅速发展紧密地联系在一起而用于将连续偏微分方程离散化的数值方法也变
得更加完善和精确从而减小求偏微分方程数值解的误差数值天气预报和大气
环流数值模拟的实践已经证明有限差分法是求解流体动力学方程组的一种有效
方法
谱方法在数值天气预报模式中的应用可上溯至1942年当时前苏联学者
Blinova[2]提出利用球谐函数展开SphericalHarmonicExpansions的线性化模式来
制作长期天气预报随后Silberman[2]在1954年借助于球谐函数求解了无辐散的
正压涡度方程Lorenz[3]在1960年证明了正压无辐散流动的谱截断方程与其对应
的精确微分方程一样具有总动能和总涡度平方守恒的性质Platzman[4]于同年指出
这一性质可以自动地消除数值计算的非线性不稳定而这在当时是以有限差分法
为基础的模式网格点模式中存在的最主要的困难之一1965年Cooley[5]等人
提出的快速傅里叶变换FastFourierTransformFFT使得谱方法的计算效率大
大提高1966comrt[6]说明了如何将展开函数分解为较简单的形式以解决计
算机内存不足的问题1970年Eliasen[7]等人和Orszag[8]分别独立提出变换法
计算二次型的非线性项以取代相互作用系数法使得?
谱模式所需的计
算量与计算机内存容量大为降低?
更容易将各种物理过程直接引入谱模式中
变换法的提出为谱模式的发展和应用铺平了道路
[9]
我国对谱模式的研究始于20世纪70年代初郑庆林从1974年开始利用北
半球三层准地转谱模式进行了较多个例的预报实验其72小时内的预报效果是比
[10]
较好的他又于1976年做出了北半球七层初始方程谱模式的个例预报该预
报
[11]
模式考虑了地形和动量热量的水平扩散等物理过程1989年郑庆林将该模
第1页
国防科学技术大学研究生院工学硕士学位论文
式发展为全球七层大气环流谱模式并进行了30天的长期数值天气预报试验该模
式包含了较完整的物理过程方程组求解方案能有效地克服散度方程以及σ坐标
系中在大地形附近计算气压梯度力项时所存在的大量之间的小差问题其非线性
项的谱计算方法有优越性且在整个30天内模式的预报误差都比对应的持续性误
[12]
差小1991年韩卫清纪立人通过一个单向嵌套的有限区域浅水波方程谱模
式的具体方案预报试验验证了谱模式对于边界条件的敏感性并提出了一种有
[13]
助于谱模式中对陡峭地形处理的地形区计算方案1996年陈嘉滨舒静君提
出用热力学变量比容α代替温度变量T并将地面气压本身作为展谱预报变量解
决了模式大气总质量不守恒的问题且维持气压梯度项为二阶项不出现混淆误
[14]
差谷湘潜于1998年推导出大气自忆方程的谱形式采用T42L9谱模式为动
力
核建立了全球自忆T42模式对500hPa高度场的预报试验结果将月平均环
流距平
相关系数提高135月平均环流均方根误差降低294
com球谐函数与Silberman方法
在球坐标中相对涡度为
1v
ucos
11
acos
无辐散流场可用流函数表示为
1
u
a
12
1
v
acos
涡度表示为
2112
22cos13
acoscos
则正压涡度方程可写成[56]
21222
2214
tacosa
函数
YmPmsineim15
ll
是方程
22mm16
aYll1Y0
ll
的解其中Pmsin是连带Legendre多项式AssociatedLegendre
Functions
l
第2页
国防科学技术大学研究生院工学硕士学位论文
进而流函数可按球谐函数展开为
JmJ
2mm
17
tatY
ll
mJlm
从而涡度按球谐函数展开为
JmJ
2mm
18
ll1tY
ll
mJlm
将公式17和18代入14交换非线性平流项的求和附标ml和ml
1122
二次写出该项的表达式两式相加再除以2可得
JmJdmJmJi
lmmm
ll1Y2imYll
lll22
mJlmdtmJlm2cosm1l1m2l2
19
mm
dP2dP1
llimm
m1m2m12m2112
1ll1mPmPe
11l
升级会员 