数值天气预报实习报告.doc

正压原始方程模式预报 实 习 报 告 姓 名： 昊缪 学 号：20111336086一、 实习目的与要求：1、 实习目的通过正压原始方程模式预报实习，加深理解该章所学的基本原理，学会综合运用所学的知识设计求解正压原始方程模式的数值计算方案；提高编制大而复杂的程序的能力。2、 实习要求以1973年4月29日08:00（北京时）我国东北、华北地区500hPa等压面位势高度场作为初始场，采用地转风初值、固定的水平侧边界条件和二次守恒平流格式，应用正压原始方程模式制作未来24小时有限区域500hPa位势高度场和风场的预报，并写一份实习报告。实习报告的主要内容是：（1）说明正压原始方程模式的特点。（2）分析数值预报结果，并讨论正压原始方程模式的预报性能。二、 模式的性能介绍实际大气运动包含有两种最基本的动力学过程，即准地转演变过程和地转适应过程。正压原始方程模式既有缓慢移动的大气长波，又有快速移动的重力惯性波。因此，这种模式不但可以模拟准地转演变过程，而且还可以模拟地转适应过程。虽然原始方程模式比准地转模式能够更好地描述实际大气中的物理过程，但是在这种模式的数值积分中也出现了一些新的亟待解决的问题：1、原始方程模式中包含快波解，为保证计算稳定性，时间步长必须取得很短，而使得计算量过大；2、当时间步长取得很短时，非线性计算不稳定的问题尤为突出；3、原始方程模式中包含多种波动解，对资料的初始化要求较高；4、对边界条件很敏感，要求给出适当的边界条件。三、 模式的基本假定、基本方程、技术方法、计算流程等的介绍1、模式的基本假定：（1）大气为正压状态：水平风速不随气压变化；（2）大气为均匀不可压缩的流体：密度为一常数；（3）大气处于静力平衡：水平气压梯度不随高度变化；（4）大气上界为自由面：有水平辐合辐散；（5）不考虑摩擦及非绝热作用。2、模式的基本方程：3、技术方法（1）为了解决测的风场与高度场之间的不平衡，以及风场、高度场与模式之间的不协调，直接用观测的风场和高度场作为原始方程模式的初始值容易产生高频振荡，使数值积分变得不稳定这个问题就必须对资料进行初始化，们采用的是静力初始化，初始条件采用地转风初值;（2）对有限区域的预报，在边界上必须人为地给出水平侧边界条件，本次选择固定边界条件；（3）在时间积分的最初阶段，采用欧拉-后差格式，来阻尼由于初始资料不平衡而引起的高频振荡，积分一段时间后，再改用“三步法”起步的中央差格式进行积分；（4）为滤除短波扰动，抑制非线性计算不稳定，可穿插进行空间平滑，包括：边界内第一圈格点“边界点”平滑：9点平滑；内点平滑：5点平滑。4、计算流程（1） 输入初始位势高度场；（2） 计算地图投影放大系数和科式参数；（3） 计算初始地转风场（应用书中4.120）；（4） 计算风场和位势高度场的变化项；（5） 时间积分，最终得出预订时刻的位势高度场和风场的预报.流程图四、 试验结果及分析：1、 试验结果： （1）初始时刻500hPa等压面实况图： （2）预报时刻500hPa等压面图：2、 实验结果分析根据初始场和预报场的图形进行对比，可以看出有明显变化的是鄂霍次克海附近的低压中心发生移动，以及贝加尔湖附近的低压槽线发生东移。由预报场可知，低压中心向东北方向移动，原东北有一高压脊570 hPa线，在预报时刻已经减弱，并在其东部生成更强的575hPa线高压脊，我国东北地区有正变高，且出现高压脊；贝加尔湖附近，槽线从它西侧移动到东侧，可能会影响到东北。预报场跟实况场的误差不是很大，误差产生原因：1、模式考虑的是正压方程，而实际大气是斜压的，计算过程中都是考虑在理想条件下；2、实际计算中由于步长选取的原因，计算本身就存在误差5、 小结与讨论：第一次上数值预报实习课，我遇到了很多困难，但是经过老师的讲解和指导，总算圆满完成了。这次实习加强了我们对正压原始方程模式的了解，主要考察了我关于五点平滑，九点平滑和地转风的编程计算，加深了我对数值天气预报这门课的了解与掌握。此次接触到的是正压原始方程模式，考虑条件都是比较理想的，算是个比较简单的模式。在500hPa高度上，大气基本是符合正压，准地转，准静力平衡的，所以，经过正压方程预报所得的结果基本与实况场差距不大。但是实际大气并没有假想的那样理想，实际大气是斜压的，复杂的系统。在以后的学习中，斜压方程预报模型可能才是以后使用更多的，难度也更大。在Fortran程序编译上，我的能力还有所欠缺，有些计算过程需要我们仔细考虑，通过对比后，再去选择计算过程更简单的，思路更明确的编程方法，也就是所说的程序最优化。附录：（1） 计算地转风初值：subroutine cgw(ua,va,za,rm,f,d,m,n)dimension ua(m,n),va(m,n),za(m,n),rm(m,n),f(m,n),dp=m-1q=n-1do i=1,mua(i,1)=-rm(i,1)*9.8*(za(i,2)-za(i,1)/(f(i,1)*d) ua(i,n)=-rm(i,n)*9.8*(za(i,n)-za(i,n-1)/(f(i,n)*d)do j=2,qua(i,j)=-rm(i,j)*9.8*(za(i,j+1)-za(i,j-1)/(2.0*f(i,j)*d)end doend dodo j=1,nva(1,j)=rm(1,j)*9.8*(za(2,j)-za(1,j)/(f(1,j)*d)va(m,j)=rm(m,j)*9.8*(za(m,j)-za(m-1,j)/(f(m,j)*d)do i=2,p va(i,j)=rm(i,j)*9.8*(za(i+1,j)-za(i-1,j)/(2.0*f(i,j)*d) end doend doend subroutine cgw （2）区域内5点平滑subroutine ssip(a,w,s,m,n,k,l)dimension a(m,n),w(m,n)if(l=1)thendo i=2,m-1do j=2,n-1 w(i,j)=a(i,j)+s*(a(i-1,j)+a(i+1,j)+a(i,j-1)+a(i,j+1)-4*a(i,j)/4.0enddoenddodo i=2,m-1 do j=2,n-1a(i,j)=w(i,j) enddoenddoreturnelsedo i=2,m-1do j=2,n-1 w(i,j)=a(i,j)+s*(a(i-1,j)+a(i+1,j)+a(i,j-1)+a(i,j+1)-4*a(i,j)/4.0 enddoenddodo i=2,m-1 do j=2,n-1a(i,j)=w(i,j)enddoenddo do i=2,m-1do j=2,n-1 w(i,j)=a(i,j)+(-s)*(a(i-1,j)+a(i+1,j)+a(i,j-1)+a(i,j+1)-4*a(i,j)/4.0enddoenddodo i=2,m-1 do j=2,n-1 a(i,j)=w(i,j)enddoenddoendifreturn end （3）边界九点平滑 subroutine ssbp(a,w,s,m,n) dimension a(m,n),w(m,n) do i=1,m do j=1,nw(i,j)=a(i,j) end doend do j=2 do i=2,m-1 w(i,j)=a(i,j)+s*(1-s)/2*(a(i-1,j)+a(i+1,j)+a(i,j-1)+a(i,j+1)-4*a(i *,j)+s*s/4*(a(i-1,j-1)+a(i-1,j+1)+a(i+1,j-1)+a(i+1,j+1)-4*a(i,j) end do j=n-1 do i=2,m-1 w(i,j)=a(i,j)+s*(1-s)/2*(a(i-1,j)+a(i+1,j)+a(i,j-1)+a(i,j+1)-4*a(i *,j)+s*s/4*(a(i-1,j-1)+a(i-1,j+1)+a(i+1,j-1)+a(i+1,j+1)-4*a(i,j) end do i=2 do j=2,n-1 w(i,j)=a(i,j)+s*(1-s)/2*(a(i-1,j)+a(i+1,j)+a(i,j-1)+a(i,j+1)-4*a(i *,j)+s*s/4*(a(i-1,j-1)+a(i-1,j+1)+a(i+1,j-1)+a(i+1,j+1)-4*a(i,j) end do i=m-1 do j=2,n-1 w(i,j)