海洋数值模型的理论及应用ppt课件.ppt

第一节　海洋数值模型的发第一节　海洋数值模型的发展概况展概况观测，试验观测，试验理论分析理论分析数值模型数值模型海洋的研究手段海洋的研究手段海洋数值模型发展的客观要求海洋数值模型发展的客观要求海水运动方程推导过程中的动力学假设，比如海海水运动方程推导过程中的动力学假设，比如海水静力学假设，都是建立在大部分空间尺度大于水静力学假设，都是建立在大部分空间尺度大于10km10km，时间尺度大于惯性周期（，时间尺度大于惯性周期（2 2ππ/f/f）的情况）的情况下的下的海水运动方程不存在精确解：（海水运动方程不存在精确解：（1 1）偏微分方程）偏微分方程本身是非线性的，（本身是非线性的，（2 2）将一个实际问题具体化）将一个实际问题具体化时所需要的环境场如海底地形、海岸线几何形状、时所需要的环境场如海底地形、海岸线几何形状、表面驱动力等都无法给出解析函数表达式，只能表面驱动力等都无法给出解析函数表达式，只能在离散化的空间和时间区域取得在离散化的空间和时间区域取得因此，求解海水运动方程只能通过近似求解的办因此，求解海水运动方程只能通过近似求解的办法来确定，而这些数值解又可能引入明显的近似法来确定，而这些数值解又可能引入明显的近似误差误差海洋数值模型发展的历史海洋数值模型发展的历史海洋模型到按其水平网格的离散方式以及所使用海洋模型到按其水平网格的离散方式以及所使用的垂向坐标系的不同大致经历了如下几个发展阶的垂向坐标系的不同大致经历了如下几个发展阶段段最早最早出现并且还在使用的海洋模型是出现并且还在使用的海洋模型是BryanBryan等人等人开发的基于原始方程的低阶精度的有限差分模型，开发的基于原始方程的低阶精度的有限差分模型，它在水深方向采用它在水深方向采用z z坐标系坐标系目前常见的目前常见的HOPS（（ Harvard Ocean Prediction System ）），，MOM（（ GFDL Modular Ocean Model）），，POP (Parallel Ocean Program)，NCOM（（NCAR Community Ocean Model））模型在某种程度上模型在某种程度上都可以认为是该模型改进版都可以认为是该模型改进版海洋数值模型发展的历史海洋数值模型发展的历史上个世纪上个世纪7070年代年代，，sigma坐标系开始应用于坐标系开始应用于海洋模型在水深方向，比如目前被广泛使用海洋模型在水深方向，比如目前被广泛使用的的POM（（Princeton Ocean Model））、、ECOM （（Estuarine Coastal and Ocean Model）） 、、ROMS（（Regional Ocean Modeling System））模型都属于这种类型的模型都属于这种类型的模型模型与传统的与传统的 z z 坐标系相比，坐标系相比，sigmasigma坐标可以更坐标可以更好地贴合海底地形变化好地贴合海底地形变化海洋数值模型发展的历史海洋数值模型发展的历史自上世纪自上世纪9090年代以来年代以来，非结构化网格技术开始，非结构化网格技术开始在海洋模型中得到应用，相应的也出现了基于在海洋模型中得到应用，相应的也出现了基于非结构化网格离散方式的有限元模型非结构化网格离散方式的有限元模型，，如如SEOM模型模型（（ Spectral Finite Element Ocean Model），和有限体积模型），和有限体积模型，，如如FVCOM虽然与传统的结构化网格相比虽然与传统的结构化网格相比, ,非结构化网格非结构化网格可以更好地拟合陆地边界，但是代码实现上的可以更好地拟合陆地边界，但是代码实现上的困难以及计算稳定性的问题使其迄今还没有得困难以及计算稳定性的问题使其迄今还没有得到非常广泛到非常广泛的的应用应用新一代的海洋数值模型新一代的海洋数值模型新一代的海洋模型广泛采用随地坐标系（新一代的海洋模型广泛采用随地坐标系（ terrain-following coordinates），进而促进了有），进而促进了有关时间步长，对流项和压力梯度项等数值算法关时间步长，对流项和压力梯度项等数值算法的改进的改进进入新世纪以来进入新世纪以来，下一代的海洋数值动力模型，下一代的海洋数值动力模型正在紧锣密鼓的研制中，代表性的是正在紧锣密鼓的研制中，代表性的是TOMS (Terrain-following Ocean Modeling System)，它融合了目前最先进的物理知识、数值方法，它融合了目前最先进的物理知识、数值方法和数据同化技术和数据同化技术 http://www.aos.princeton.edu/WWWPUBLIC/htdocs.pom/TOMS.htm 海洋模型的商业软件海洋模型的商业软件就海流的仿真建模研究而言，现阶段国外已形就海流的仿真建模研究而言，现阶段国外已形成了不少具有代表性的、具有强大前后处理功成了不少具有代表性的、具有强大前后处理功能及核心仿真建模技术的软件或系统能及核心仿真建模技术的软件或系统（商业软（商业软件）件）如加拿大如加拿大Dalhousie University研究研究的的DALCOAST河口海岸预报系统河口海岸预报系统；；丹麦的丹麦的DHI Water & Environment 机构开发的的机构开发的的MIKE系系列软件系统列软件系统；；荷兰的荷兰的WL Delft Hydraulics开开发的的发的的Delft－－3D软件软件第二节　第二节　基于基于POM模型的一般介绍模型的一般介绍水动力模型基本特点水动力模型基本特点n垂直方向的垂直方向的 ᵴ ᵴ 坐标坐标变换变换n水平网格采用正交曲线坐标和水平网格采用正交曲线坐标和ArakawaCArakawaC差分格式差分格式n水平时间差分采用显格式，而垂直差分为隐格式水平时间差分采用显格式，而垂直差分为隐格式n自由表面可以模拟水位变化自由表面可以模拟水位变化n垂直和水平方向的混合扩散分别采用垂直和水平方向的混合扩散分别采用2.52.5阶的阶的Mellor-YamadaMellor-Yamada湍流闭合模式和湍流闭合模式和SmagorinskiSmagorinski模式模式n内外内外模态分别处理速度较慢的内重力波和速度较模态分别处理速度较慢的内重力波和速度较快的外重力波以提高整个模式快的外重力波以提高整个模式计算效率计算效率n包含了海水的热动力过程包含了海水的热动力过程Sigma坐标变换：坐标变换：hs = -1z = 0z = H(x, y) s = 0 V (i, j+1) V (i, j) U (i+1, j) U (i, j) ᵳ (i, j)模型采用（模型采用（a a）垂向）垂向sigmasigma坐标系和（坐标系和（b b）水平）水平Arakawa CArakawa C网格网格(1)(1)控制方程控制方程连续方程连续方程海水运动方程海水运动方程(2)(2)(3)(3)(4)(4)温度方程温度方程湍流动能方程湍流动能方程盐度方程盐度方程湍流混合长度方程湍流混合长度方程(5)(5)(6)(6)(7)(7)(8)(8)其中，其中，ω是基于是基于sigma坐标系的垂向流速矢量，它坐标系的垂向流速矢量，它与笛卡儿坐标系下的垂向流速与笛卡儿坐标系下的垂向流速W之间的关系可表示之间的关系可表示为为另外，方程中另外，方程中另外，方程中另外，方程中 分别表示海水的实际密度，分别表示海水的实际密度，分别表示海水的实际密度，分别表示海水的实际密度，BoussinesqBoussinesqBoussinesqBoussinesq 近似密度以及密度扰动近似密度以及密度扰动近似密度以及密度扰动近似密度以及密度扰动 其中，其中， 是实际水深是实际水深(9)(9)M—Y湍流闭合模型（湍流闭合模型（profq））在湍流动能及混合长度方程中在湍流动能及混合长度方程中在湍流动能及混合长度方程中在湍流动能及混合长度方程中 分别分别表表示湍流混合系数示湍流混合系数，，热热扩散系数和湍流动能扩散系数和湍流动能的垂直的垂直扩扩散系数散系数；；墙墙近似近似函数函数 ，其中，其中 ，，κκ=0.4=0.4 是是vonvon KarmanKarman 常数常数湍流动能和混合长度方程组由下列等式关系闭合湍流动能和混合长度方程组由下列等式关系闭合 其中，其中， 是稳定函数是稳定函数(10)(10)在在2.52.5阶阶M M——Y Y湍流闭合模型中，湍流闭合模型中， 可表示为可表示为其中，其中， 。

在不稳定层次在不稳定层次 条件下，条件下， 的上限值为的上限值为0.0230.023，对应的，对应的 值分别为值分别为2.01452.0145和和2.44012.4401；在稳定层次；在稳定层次 条件下，它下限值为条件下，它下限值为-0.28-0.28，对应的，对应的 值分别为值分别为 0.04700.0470和和 0.04610.0461 是是2.52.5阶阶M-YM-Y湍流闭合模式参数湍流闭合模式参数，由实验测得，由实验测得(11)(11)M—Y湍流模型的适用条件湍流模型的适用条件M-YM-Y湍流模型已被广泛应用于浅海潮汐，风生表湍流模型已被广泛应用于浅海潮汐，风生表面混合层的模拟以及底边界层的研究；面混合层的模拟以及底边界层的研究；该模型在混合较弱的层化流体中对湍流混合系该模型在混合较弱的层化流体中对湍流混合系数计算效果不佳，应用于河口的可行性也有待数计算效果不佳，应用于河口的可行性也有待于探讨；于探讨；虽然该模型无论是在物理上还是在数学上都存虽然该模型无论是在物理上还是在数学上都存在着明显的缺陷，但是目前依然得到广泛应用。

在着明显的缺陷，但是目前依然得到广泛应用只有深入了解海水的微细结构以及海洋湍流结只有深入了解海水的微细结构以及海洋湍流结构和性质后，才有可能构建出更为完善的湍流构和性质后，才有可能构建出更为完善的湍流闭合模型闭合模型 海水运动方程的水平扩散海水运动方程的水平扩散项项•其中水平切应力其中水平切应力• 是水平涡粘性系数，由是水平涡粘性系数，由Smagorinsky模式模式计算计算(12)(12)(13)(13)Smagorinsky水平扩散水平扩散模模式式•其中参数其中参数C C 为无量纲值，其取值为为无量纲值，其取值为0.10.1或或0.20.2，如果计算网格划分得足够细的话，，如果计算网格划分得足够细的话，C C 可以可以取值为取值为0 0 •由由SmagorinskySmagorinsky模式计算的水平涡粘性系数随模式计算的水平涡粘性系数随着网格精度的改善和速度梯度的减小而减小着网格精度的改善和速度梯度的减小而减小(14)(14)温度盐度的水平扩散项温度盐度的水平扩散项•其中其中• 是水平热扩散系数，一般地，是水平热扩散系数，一般地， 是一个小是一个小值，取值，取0.10.1或或0.20.2，甚至在某些情况下可以取为，甚至在某些情况下可以取为0 0(15)(15)(16)(16)垂向边界条件垂向边界条件1. 1. 海面海面(17)(17)2. 2. 海底海底其中，其中， 分别是湍流闭合模型参数和摩擦速度；分别是湍流闭合模型参数和摩擦速度；分别是海面和海底的摩擦系数；分别是海面和海底的摩擦系数；(18)(18)模型的外模态模型的外模态•控制近岸环流的海水运动方程包含了运动速度较快的控制近岸环流的海水运动方程包含了运动速度较快的外重力波和速度较慢的内重力波，因此进行模式分裂外重力波和速度较慢的内重力波，因此进行模式分裂可以大大提高模型的计算效率可以大大提高模型的计算效率•模型的外模态就是计算自由水位和垂向平均的流速变模型的外模态就是计算自由水位和垂向平均的流速变化，需要较小的时间步长化，需要较小的时间步长•内模态主要模拟流速矢量、温度、盐度等三维结构变内模态主要模拟流速矢量、温度、盐度等三维结构变化，对时间步长的要求相对宽松化，对时间步长的要求相对宽松•采用模式分裂技巧可以有效地减少因外模态计算所消采用模式分裂技巧可以有效地减少因外模态计算所消耗的计算时间，从而提高计算效率耗的计算时间，从而提高计算效率垂向积分的海水运动方程垂向积分的海水运动方程其中：其中：连续方程连续方程海水运动方程海水运动方程(19)(19)(20)(20)(21)(21)水平扩散项水平扩散项水平数值耗散项水平数值耗散项(22)(22)(23)(23)模型代码中的一些符号模型代码中的一些符号•im, jm, kb, imm1, jmm1, kbm1, mode, isplit•dte, dti, days, umol•z, zz, dz, dzz•aam2d, art, aru, arv, dum, dvm, fsm•dx, dy, h, el, d, ua, va, ut, vt, et, cor•swrad•wusurf, wvsurf, wubot, wvbot, wtsurf, wssurf, rad•u, v, ω, t, s, rho, l, q2•km, kh, kq, aam, aah•rmean, tclim, sclim数数值值积积分分流流程程图图外模态和内模态的相互作用外模态和内模态的相互作用内模态的计算方法内模态的计算方法三维变量的计算（三维变量的计算（以以T为例为例）可以分解为垂直扩散项的计）可以分解为垂直扩散项的计算和水平对流扩散项的计算，以温度方程为例方程写为算和水平对流扩散项的计算，以温度方程为例方程写为求解过程分为两步，第一步计算求解过程分为两步，第一步计算水平对流扩散项水平对流扩散项，采用，采用中心差分中心差分格式，由格式，由advt子程序实现子程序实现(24)(24)(25)(25)水平对流扩散项水平对流扩散项垂直扩散项垂直扩散项第二步计算第二步计算垂向扩散项垂向扩散项，数值方法采用，数值方法采用“蛙跳蛙跳”格式格式，，由由proft子程序实现子程序实现由于由于“蛙跳蛙跳”格式在奇数时间步时的解与在偶数时间步时格式在奇数时间步时的解与在偶数时间步时的解会发生偏离，因此每一时间步的计算结束后还需对上的解会发生偏离，因此每一时间步的计算结束后还需对上述解进行平滑处理，即述解进行平滑处理，即α为一个小值，一般取为一个小值，一般取0.05。

处理后处理后，，(26)(26)(27)(27)计算网格的安排计算网格的安排二维外模态网格分布二维外模态网格分布三维内模态三维内模态网格分布网格分布水平对流项的计算水平对流项的计算■ ■ 虽然模型采用虽然模型采用有限差分有限差分计算格式，但是对于每个网格计算格式，但是对于每个网格对流项的计算都是按照对流项的计算都是按照有限容积有限容积的方法进行处理，即的方法进行处理，即温度对流过程可表示为温度对流过程可表示为■ ■ 速度的对流过程与温度相似，可写为速度的对流过程与温度相似，可写为 其中，其中， 是由是由z z坐标转化为坐标转化为σσ坐标后坐标后产生的弯曲项产生的弯曲项 (28)(28)(29)(29)垂直扩散项的计算垂直扩散项的计算■ ■ 垂直扩散项的计算公式（垂直扩散项的计算公式（第第k k层，层，1

公式（温度通量公式（3434）的解可以近似写为）的解可以近似写为■■ 这样，表层海水的温度只能根据式（这样，表层海水的温度只能根据式（3131）求解，即）求解，即■■ 上式可以进一步表示为上式可以进一步表示为(36)(36)(37)(37)■ ■ 式（式（3737）还可以进一步写为）还可以进一步写为■■ 再与温度的通解（再与温度的通解（3434）比较）比较■■ 可以得到可以得到(34)(34)(38)(38)(39)(39)■ ■ 当短波辐射通量当短波辐射通量 的计算公式的计算公式■■ r r，，ad1ad1，，ad2ad2 均为光辐射常数，根据不同水质取值如下均为光辐射常数，根据不同水质取值如下(40)(40)ntp 1 2 3 4 5Jerlov type I Ia Ib II IIIr 0.58 0.62 0.67 0.70 0.78ad1 (m) 0.35 0.60 1.0 1.5 1.4ad2 (m) 23.0 20.0 17.0 14.0 7.9来源：来源：Jerlov，，1976；；Paulson and Simpson，，1977■ ■ 当当k=kb-1k=kb-1时时，即底层海水，假定海底的热通量为，即底层海水，假定海底的热通量为0 0，，根据前面的推导方式同样可以得到底层海水温度的计根据前面的推导方式同样可以得到底层海水温度的计算公式算公式■■ 对于盐度方程垂向扩散项的计算与温度方程相同，对于盐度方程垂向扩散项的计算与温度方程相同，只是不考虑太阳短波辐射这一项只是不考虑太阳短波辐射这一项(41)(41)时间步长的时间步长的CFL限制条件限制条件● ● 由于模型水平方向采用显格式，因此时间步长的选取由于模型水平方向采用显格式，因此时间步长的选取必须要满足必须要满足CFLCFL稳定条件。

稳定条件● ● 对于外模态，时间步长的限制条件为对于外模态，时间步长的限制条件为 其中其中 ，， 可能预见到的最大流速可能预见到的最大流速● ● 对于内模态，时间步长的限制条件较外模态的情形对于内模态，时间步长的限制条件较外模态的情形宽松很多，主要是速度较快的外重力波已经在外模态宽松很多，主要是速度较快的外重力波已经在外模态中考虑了一般中考虑了一般 取值取值3030~5050即可即可 (42)(42)时间步长的其它限制条件时间步长的其它限制条件● ● 对于动量或标量还有其它一些时间限制条件对于动量或标量还有其它一些时间限制条件 其中其中 或或● ● 以及以及 其中其中 分别为地球角速率和地理纬度分别为地球角速率和地理纬度 (43)(43)(44)(44)侧开边界条件（侧开边界条件（bcond））陆地及岸线是由陆地及岸线是由dumdum、、dvmdvm、、fsmfsm控制的，在陆地控制的，在陆地上这些变量的值设为上这些变量的值设为0 0，有水的地方设为，有水的地方设为1 1子程序子程序bcondbcond（（idxidx），），idx=1idx=1对应的水位边界条对应的水位边界条件；件；idx=2idx=2对应垂向平均流速边界条件；对应垂向平均流速边界条件；idx=3idx=3对应三维水平流速边界条件；对应三维水平流速边界条件；idx=4idx=4对应温盐边对应温盐边界条件；界条件；idx=5idx=5对应垂向流速边界条件；对应垂向流速边界条件；idx=6idx=6对应湍流动能和湍流混合长度边界条件对应湍流动能和湍流混合长度边界条件 模型对开边界条件的要求很高，而开边界条件模型对开边界条件的要求很高，而开边界条件本身又具有很大的不确定性，因此有必要对模本身又具有很大的不确定性，因此有必要对模型的外模态和内模态的开边界分别进行处理型的外模态和内模态的开边界分别进行处理 数数值值积积分分流流程程图图FormulaBoundaryCode Inflow condition EASTuaf(im,j) = 2*bc(j)/(h(im,j)+elf(im,j) + h(imm1,j +elf(imm1,j))elf(im,j) = elf(imm1,j)vaf(im,j) =setWESTuaf(2,j) = 2*bc(j)/(h(1,j)+elf(1,j) + h(2,j)+elf(2,j))elf(1,j) = elf(2,j)vaf(1,j) = setNORTHvaf(i,jm) = 2*bc(i)/(h(i,jm)+elf(i,jm) + h(i,jmm1) + elf(i,jmm1))elf(i,jm) = elf(i,jmm1)uaf(i,jm) = setSOUTHvaf(i,2) = 2*bc(i)/(h(i,1)+elf(i,1) + h(i,2)+elf(i,2))elf(i,1) = elf(i,2)uaf(i,1) = set Elevation condition h = BC EASTelf(imm1,j) = bc(j)elf(im,j) = elf(imm1,j) cosmeticuaf(im,j) = uaf(imm1,j)vaf(im,j) = setWESTelf(2,j) = bc(j)uaf(2,j) = uaf(3,j)vaf(1,j) = setNORTHelf(i,jmm1) = bc(i)elf(i,jm) = elf(i,jmm1) cosmeticvaf(i,jm) = vaf(i,jmm1)uaf(i,jm) = setSOUTHelf(i,2) = bc(i)vaf(i,2) = vaf(i,3)uaf(i,1) = set外外模模态态开开边边界界条条件件 一一FormulaBoundaryCodeRadiationEASTuaf(im,j) = sqrt(grav/h(imm1,j))* el(imm1,j) + bc(j)elf(im,j) = elf(imm1,j)vaf(im,j) = setWESTuaf(2,j) = - sqrt(grav/h(2,j))* el(2,j)+bc(j)elf(1,j) = elf(2,j)vaf(1,j) = setNORTHvaf(i,jm) = sqrt(grav/h(i,jmm1))* el(i,jmm1) + bc(i)elf(i,jm) = elf(i,jmm1)uaf(i,jm) = setSOUTHvaf(i,2) = - sqrt(grav/h(i,2))* el(i,2)+bc(i)elf(i,1) = elf(i,2)uaf(i,1) = setRadiation EASTgae = dte*sqrt(grav*h(im,j))/dx(im,j)uaf(im,j) = gae*ua(imm1,j) + (1.-gae)*ua(im,j)elf(im,j) = elf(imm1,j)vaf(im,j) = setWESTgae = dte*sqrt(grav*h(2,j))/dx(2,j)uaf(2,j) = gae*ua(3,j) + (1.-gae)*ua(2,j)elf(1,j) = elf(2,j)vaf(1,j) = setNORTHgae = dte*sqrt(grav*h(i,jm))/dy(i,jm)vaf(i,jm) = gae*va(i,jmm1) + (1.-gae)*va(i,jm)elf(i,jm) = elf(i,jmm1)uaf(i,jm) = setSOUTHgae = dte*sqrt(grav*h(i,2))/dy(i,2)vaf(i,2) = gae*va(i,3) + (1.-gae)*va(i,2)elf(i,1) = elf(i,2)uaf(i,1) = set外外模模态态开开边边界界条条件件 二二流速水位流速水位流速流速FormulaBoundaryCodeRadiationEASTgae = dte*sqrt(grav*h(imm1,j))/dx(imm1,j)elf(imm1,j) = gae*el(imm2,j) + (1.-gae) *el(imm1,j)elf(im,j) = elf(imm1,j)uaf(im,j) = uaf(imm1,j)vaf(im,j) = setWESTgae = dte*sqrt(grav*h(2,j))/dx(1,j)elf(2,j) = gae*el(3,j) + (1.-gae)*el(2,j)uaf(2,j) = uaf(3,j)vaf(2,jm) = setNORTHgae = dte*sqrt(grav*h(i,jmm1))/dy(i,jmm1)elf(i,jmm1) = gae*el(i,jmm2) + (1.-gae) *el(i,jmm1)elf(i,jm) = elf(i,jmm1)vaf(i,jm) = vaf(i,jmm1)uaf(i,jm) = setSOUTHgae = dte*sqrt(grav*h(i,2))/dy(i,2)elf(i,2) = gae*el(i,3) + (1.-gae)*el(i,2)vaf(i,2) = vaf(i,3)，uaf(i,1) = set Cyclic EAST(I=IM)elf(im,j) = elf(3,j)uaf(im,j) = uaf(3,j)，vaf(im,j) = vaf(3,j)WEST(I=1)elf(1,j) = elf(imm2,j)，elf(2,j) = elf(imm1,j)uaf(2,j)=uaf(imm1,j)，vaf(2,j)=vaf(imm1,j)NORTH(J=JM)elf(i,jm) = elf(i,3)uaf(i,jm) = uaf(i,3)vaf(i,jm) = vaf(i,3)SOUTH(J=1)elf(i,1) = elf(i,jmm2)，elf(i,2) = elf(i,jmm1)uaf(i,2)=uaf(i,jmm1)vaf(i,2)=vaf(i,jmm1)外外模模态态开开边边界界条条件件 三三水位水位FormulaBoundaryCodeInflow conditionEASTuf(im,j,k) = bc(j,k)vf(im,j,k) = setU = BC WESTuf(2,j,k) = bc(j,k)vf(1,j,k) = set NORTHvf(i,jm,k) = bc(i,k)uf(i,jm,k) = setSOUTHvf(i,2,k) = bc(i,k)uf(i,1,k) = setRadiation:EASTgai = sqrt(h(im,j)/hmax) uf(im,j,k) = gai*u(imm1,j,k) + (1.-gai)*u(im,j,k)vf(im,j,k) = setWESTgai = sqrt(h(2,j)/hmax)uf(2,j,k) = gai*u(3,j,k) + (1.-gai)*u(2,j,k)vf(1,j,k) = setNORTHgai = sqrt(h(i,jm)/hmax)vf(i,jm,k) = gai*v(i,jmm1,k) + (1.-gai)*v(i,jm,k)uf(i,jm,k) = setSOUTHgai = sqrt(h(i,2)/hmax)vf(i,2,k) = gai*v(i,3,k) + (1.-gai)*v(i,2,k)uf(i,1,k) = set内模态开边界条件（一）内模态开边界条件（一）FormulaBoundaryCodeUpstream advection on T or SEASTuf(im,j,k) = t(im,j,k) -dti/(dx(im,j)+dx(imm1,j)) * ((u(im,j,k) + abs(u(im,j,k))) * (t(im,j,k)-t(imm1,j,k)) + (u(im,j,k) - abs(u(im,j,k))) * (tbe(j,k)-t(im,j,k)))WESTuf(1,j,k) = t(1,j,k) -dti/(dx(1,j)+dx(2,j)) * ((u(1,j,k) + abs(u(1,j,k))) * (t(1,j,k)-tbw(j,k)) + (u(1,j,k) - abs(u(1,j,k))) * (t(2,j,k)-t(1,j,k))) NORTHuf(i,jm,k) = t(i,jm,k) -dti/(dy(i,jm)+dy(i,jmm1)) * ((v(i,jm,k) + abs(v(i,jm,k))) * (t(i,jm,k)-t(i,jmm1,k)) + (v(i,jm,k) - abs(v(i,jm,k))) * (tbn(i,k)-t(i,jm,k)))SOUTHuf(i,1,k) = t(i,1,k) -dti/(dy(i,1)+dy(i,2)) * ((v(i,1,k) + abs(v(i,1,k))) * (t(i,1,k)-tbs(i,k)) + (v(i,1,k) - abs(v(i,1,k))) * (t(i,2,k)-t(i,1,k)))CyclicMuch the same as the external case except replace uaf with uf, etc. and t, s, q2 and q2l are handled similar to elf内模态开边界条件（二）内模态开边界条件（二）第三节　第三节　POM模型的实际应用模型的实际应用————渤海河流羽状流的数值模拟渤海河流羽状流的数值模拟渤海及入海河流渤海及入海河流河流径流输入模型河流径流输入模型•将河流视为点源，而点源的模拟采用的是一个水将河流视为点源，而点源的模拟采用的是一个水平网格的面积平网格的面积 •径流动量对动量方程的贡献可以忽略不计径流动量对动量方程的贡献可以忽略不计•受河流径流影响的二维连续方程受河流径流影响的二维连续方程 •三维连续方程可表示为三维连续方程可表示为 •内外模式相容内外模式相容30天时纯浮力羽流条件下渤海的盐度分布，等值线间隔为天时纯浮力羽流条件下渤海的盐度分布，等值线间隔为0.5 psu，，最外层等值线为最外层等值线为30 psu。

a）海水表层；（）海水表层；（b）海水底层）海水底层 30天时潮致羽流条件下渤海的盐度分布，等值线间隔为天时潮致羽流条件下渤海的盐度分布，等值线间隔为0.5 psu，，最外层等值线为最外层等值线为30 psua）海水表层；（）海水表层；（b）海水底层）海水底层 东风东风南风南风风场加入风场加入1010天后渤海天后渤海表层海水流场表层海水流场东风东风南风南风风场加入风场加入10天后渤海表层盐度的水平分布状况天后渤海表层盐度的水平分布状况 等值线间隔为线间隔为0.5 psu海河口附近站点（海河口附近站点（117117º4545´E E，，38 38 º5555´N N）盐度）盐度变变化化的的时间序列时间序列海水表层海水表层海水底层海水底层7 7月份渤海表层盐度的水平分布状况，等值线间隔为月份渤海表层盐度的水平分布状况，等值线间隔为1.0 1.0 psupsu计算结果结果( (左图）；实测结果（右图）左图）；实测结果（右图）学习要点学习要点了解基于了解基于sigmasigma坐标变换的海水运动方程坐标变换的海水运动方程了解了解POMPOM模型的主要特征模型的主要特征熟悉模型运算的流程熟悉模型运算的流程掌握掌握M-YM-Y湍流闭合模式的基本特点湍流闭合模式的基本特点针对不同研究对象学会选择针对不同研究对象学会选择POMPOM模式的开模式的开边界条件边界条件河流模型与河流模型与POMPOM模型的耦合模型的耦合。