（环境工程专业论文）河口海岸垂向及平面二维水流数值模拟及应用.pdf

青岛理工大学工学硕士学位论文 摘要 实际工程中河口海岸、湖泊、河流区域边界几何形状的不规则以及水动力学 方程的非线性的性质，所以难以做出理论分析解，只能求助于数值解。河口海岸 地区的水流运动严格的说都是三维的，但由于三维运动本身的复杂性给三维数学 模型处理造成一系列困难，在实际工程中往往根据问题的性质和地形特征，将三 维问题简化为二维问题。 在窄深潮汐通道、窄深河口海岸地区，水流的垂向变化比水平横向的变化大， 其流场可用水宽的平均流动量来表示，可采用垂向二维水流数值模型。而在另外 一些水域，如河口海岸广阔水域地区，水平尺度远大于垂向尺度，水流在垂向变 化小于水平方向的变化，其流场可用沿水深的平均流动量来表示，可采用平面二 维水流数值模型。 本文主要工作： 1 】对于垂向二维非恒定水流数学模型，既处理了不规则域，又研究了其算法。 由于河口海岸水深随时间、空间不断发生变化，因此计算区域也相应发生变化， 这样就给计算造成不便。垂向做盯变换，将不规则计算域变为规则计算域，可精 确跟踪自由表面与非规则的地形边界。计算方法采用交错网格、交替方向、算子 分裂的有限差分法。 2 对于平面二维非恒定水流数学模型，微分方程的求解方法用a d i 法，它 兼备了显式和隐式格式的优点。因此可将二维问题在计算意义上分解为二个先后 求解的一维问题，使计算工作量大为减少，而且该方法还有相当的精度。 所建立的垂向二维水流数学模型和平面二维水流数学模型均在实际工程中得 到应用，其中： 1 以上海洋山港部分港口的潮汐通道的水流实测资料作为验证资料，运用 编写的程序对其进行水流模拟，并将结果与实测比较，模拟结果是合理的，表明 模型是实用有效的。 2 孤东及新滩油田地区则采用平面二维水流数值模拟。本文采用平面二维 潮流数学模型模拟该海区的潮流运动情况，在比较成功模拟原流场的基础上，设 青岛理工大学工学硕士学位论文 计了两个工程方案，并对方案后的流场进行了模拟计算。该工程措施具有一定的 护岸作用，达到一定的保滩促淤效果。 关键词：河口海岸；垂向二维；平面二维；数学模型；仃坐标变换；a d i 方法 i i 青岛理工大学工学硕士学位论文 a b s t r a c t i np r a c t i c a lp r o j e c t s ，t h eb o u n d a r ys h a p ei r r e g u l a r i t yo ft h ee s t u a r ya n dc o a s t ，l a k e ， f i v e ra r e aa n dt h en o n - l i n e a rc h a r a c t e ro fh y d r o d y n a m i ce q u a t i o nm a k ei td i f f i c u l tt ot h e t h e o r e t i c a ls o l u t i o no ft h ee q u m i o n s ot h en u m e r i c a ls o l u t i o nh a st ob et h eo n l ya c h o i c e t h ef l o wo ft h ew a t e ri nt h ee s t u a r ya n dc o a s tc a nb es t r i c t l ys a yt h r e e d i m e n s i o n i tm a k e sas e r i e sd i f f i c u l tt ot h et h r e ed i m e n s i o n a ln u m e r i c a lm o d e lb e c a u s e o ft h ew a t e rf l o wc o m p l e x i t yo fi t i na c t u a l l ye n g i n e e r i n g ，t h et h r e ed i m e n s i o n a l p r o b l e m sc o n v e r ti n t ot w od i m e n s i o n a lp r o b l e m su s u a l l yb a s e do np r o p e r t ya n dt h e g e o g r a p h yc h a r a c t e r i s t i c so ft h ep r o b l e m i nt h en a r r o w - d e e pt i d ec h a n n e la n dn a r r o w - - d e e pe s t u a r ya n dc o a s ta r e a , t h e v e r t i c a lv e l o c i t yv a r i e si sf a rg r e a t e rt h a nt h a ti nt h eh o r i z o n t a ld i r e c t i o n ，a n dt h ef l o w f i e l dc a nb ee x p r e s s e da st h ea v e r a g ef l o w i n gv a l u ea l o n gt h ew a t e rw i d t h ，s ot h e t w o - d i m e n s i o nv e r t i c a ln u m e r i c a ls i m u l a t i o nt e c h n i q u es h o u l db ea d o p t e di nt h e s ea r e a s b u t ，i ns o m eo t h e rw a t e ra r e a ss u c ha st h ee s t u a r ya n dc o a s t ，t h eh o r i z o n t a ls c a l ei sf a r l a r g e rt h a nt h ev e r t i c a ls c a l e ，t h ev e l o c i t yv a r i a t i o ni ss m a l l e ri nt h ev e r t i c a ld i r e c t i o n t h a ni nt h eh o r i z o n t a ld i r e c t i o n ，a n dt h ef l o wf i e l dc a nb ee x p r e s s e da st h ea v e r a g e f l o w i n gv a l u ea l o n gt h ew a t e rd e p t h ，s ot h et w o - d i m e n s i o np l a n en u m e r i c a ls i m u l a t i o n t e c h n i q u eo ft h ew a t e rf l o wm a yb eu s e di nt h e s ea r e a s t h em a i nc o n t e n t si nt h i sp a p e ra r ea sf o l l o w s ： s i m u l a t e di r r e g u l a rb o u n d a r ya n dr e s e a r c h e dv a r i e so fa l g o r i t h m sa n ds t u d i e d t h et w o d i m e n s i o nv e r t i c a lc o m p u t i n g t h ec a l c u l a t i o nd i s t r i c ta l s oc o r r e s p o n d st h e o c c u r r e n c ev a r i e t yb e c a u s et h ed e p t ho ft h ee s t u a r ya n dc o a s tt a k e sp l a c et h ev a r i e t y c o n t i n u o u s l ya ta n yt i m ea n ds p a c e ，g i v i n gc a l c u l a t i o nr e s u l ti nt h ei n c o n v e n i e n c et h u s i nt w o d i m e n s i o nv e r t i c a lm o d e l 仃c o o r d i n a t et r a n s f o r m a t i o nw e r eu s e dt of o l l o w i n g t h et r a c k so ff r e es u r f a c ea n di r r e g u l a rt o p o g r a p h yb o u n d a r y ，s ot h ei r r e g u l a rb o u n d a r y c a nb ec h a n g e di n t oi n t ot h er e g u l a ro n e i nt h i sm o d e l ，t h eo p e r a t o rs p l i r i n g t e c h n i q u e ，s t a g g e r e d 酊dt e c h n i q u ea n df i n i t ed i f f e r e n c em e t h o dw e r eu s e d 2 t h es q u a r ed i s t a n c eo fd i f f e r e n t i a lc a l c u l u ss o l v e st h em e t h o du s e dt h ea d i m e t h o df o r t h et w o - d i m e n s i o np l a n eu n s t e a d yc u r r e n tn u m e r i c a ls i m u l a t i o n ，i th a dt h e a d v a n t a g eo ft h ed i s p l a ya n dt h ei m p l yf o r m a t t h e r e f o r e i tc a nb er e s o l v e dt h e t w o d i m e n s i o np r o b l e m si nt h em e a n i n go fc o m p u t i n gf o ro n e - d i m e n s i o np r o b l e mt h a t s o l v es u c c e s s i v e l y ，d e c r e a s i n gt h ew o r k l o a dg r e a t l ya n dt h em e t h o ds t i l lh a sq u i t e i i i 青岛理工大学工学硕士学位论文 a c c u r a c y t h et w o - - d i m e n s i o nv e r t i c a lm o d e la n dt w o - d i m e n s i o np l a n em o d e ld e v e l o p e dh a v e b e e nu s e di np r a c t i c a lp r o j e c t s ，i nw h i c h ： 1 t h em o d e lh a sb e e nv e r i f i e db yt h ep r o g r a m ec o m p i l e d 而mt h et i d ef i e l dd a t a o b t a i n e df r o mt h ep a r th a r b o ro ft h es h a n g h a iy a n g s h a nh a r b o r t h es i m u l a t i o nr e s u l t s c o u l db ec o m p a r e dw i mt h ea c t u a lm e a s u r e s ot h en u m e r i c a ls i m u l a t i o nr e s u l t sa r e s a t i s f a c t o r y t h i ss h o w st h a tt h em o d e li sp r a c t i c a la n de f f e c t i v e 2 】t h ee a s to fg u d o n ga n dn e wb e a c ho i l f i e l dr e g i o n st h e na d o p tt h ec u r r e n t n u m e r i c a ls i m u l a t i o no f t h et w o d i m e n s i o np l a n e t h et w o - d i m e n s i o np l a n et i d a lc u r r e n t n u m e r i c a lm o d e lw a su s e dt os i m u l a t et h et i d a lm o v e m e n ti nt h es e aa r e ai nt h i sp a p e r ， o nt h eb a s i so fs i m u l a t et h ef o r m e r l yf l o wf i e l ds u c c e s s f u l l y ，d e s i g n e dt w oe n g i n e e r i n g p r o j e c t s ，a n dt h ef l o wf i e l dw a sc a r r i e do nt os i m u l a t et h ec a l c u l a t i o n t h ee n g i n e e r i n g s m e a s u r eh a sc e r t a i n l yp r o t e c tt h es h o r ef u n c t i o n ，a t t a i n i n gt h er e s u l to fc e r t a i n l yp r o t e c t t of i l lu p k e yw o r d ：e s t u a r ya n dc o a s t ；t w o d i m e n s i o n a lv e r t i c a l ；t w o - d i m e n s i o np l a n e ； n u m e r i c a lm o d e l ；仃c o - o r d i n a t et r a n s f o r m a t i o n ；a d im e t h o d i v 青岛理工大学工学硕士学位论文 1 1 研究的目的和意义 第1 章绪论 在流体力学及相关领域的科学研究中，随着基本理论、数值计算方法以及计 算机技术的进步，水流及物质输运的数学模型以其特有的优势得到越来越多的重 视，发展十分迅速。尽管尚不能完全取代物理模型实验，但已经在各项工作中占 有重要地位，发挥着巨大的作用，成为科学研究的主要方法之一。 模拟水流流场的数学模型有二维与三维之分。三维数学模型正处在逐步完善 的发展过程中，实际工程与科学研究的迫切需要是其发展的主要动力，计算机的 普及与进步为其提供了强有力的技术支持，只有三维数学模型才能真实完整的模 拟流体运动以及物质输运的过程。由于目前的流体力学理论尚不完善，因此还未 建立起普遍适用的统一的标准模型，不同研究者建立的三维数学模型多种多样， 数值求解方法各不相同。这些模型各具特色，在一定的条件下都可以取得满足要 求的结果，具有一定的理论价值与工程价值，但同时又都存在一定的局限性和适 用范围。但同时也应注意到，通过对已有模型的分析研究，不断加以完善，提高 其精度，扩大其应用范围也是一条可行的研究途径。 河口海岸地区的水流及水体中物质的运动严格的说都是三维的，但由于三维 运动本身的复杂性给三维数学模型处理造成一系列困难：如( 1 ) 水平尺度远大于 垂向尺度，因此离散后网格步长相差甚大，会产生较大的舍入误差，易出现计算 的不稳定；( 2 ) 复杂的边界地形处理，往往要在边界处加密网格，因而要消耗较 多的计算时间；( 3 ) 自由水面的非恒定运动，使得垂向网格每个步长都是变化的， 而且较难进行垂向积分计算；( 4 ) 粘性系数确定比较难，低阶紊流模式大多是经 验公式，而且公式中参数往往不同作者得到不同的实验结果，使得计算结果较粗。 高阶紊流模式中往往引入大量的非线性偏微分方程并出现新的未知量，其结果大 大增加计算时间；( 5 ) 盐度、温度特别是泥沙的底部边界条件很复杂，往往需要 实验确定。边界条件中的参数和涡粘性系数往往影响数值解的稳定性。由于这一 系列的问题，在实际工程中往往根据问题的性质和地形特征，将三维问题简化为 二维问题。与三维数学模型相比，二维模型输入数据少、计算效率高。有时用二 1 青岛理工大学工学硕士学位论文 维模型可以充分反映流动特征，与完全三维模型相比还可节省计算工作量。 在窄深河口海岸地区，流速的垂向变化比水平横向的变化大，可采用垂向二 维水流数值模型。而在另外一些水域，如河口海岸、湖泊、大型水库等广阔水域 地区，水平尺度远大于垂向尺度，流速在垂向变化小于水平方向的变化，其流场 可用沿水深的平均流动量来表示，可采用平面二维水流数值模型。平面二维模型 与垂向二维模型虽然有各种各样，也还是初步的，尚缺乏共识，有待进一步深入。 研究一些既计算简便，又能确切反映物理实际的方法与格式，开发一些有效的、 通用的算法与标准程序，是十分必要的，也将是很有益的工作。 河口海岸的防洪除工程措施以外，还要求建立健全非工程措施。所谓非工程 措施，是利用洪水预警预报系统，及时给出发生特大洪水时，堤坝溃决后洪水演 进情况，包括淹没范围、水深、流速，以便防汛指挥机关有计划、有步骤撤离人 员和财产至安全区域，减少洪灾所造成的损失。 河口海岸地区的自然资源丰富，人口众多，经济发达，同时它是大气圈、水 圈、岩石圈和生物圈之间各种过程相互作用最为活跃的地带，是各圈层之间进行 物质、能量、信息交换最为频繁和集中的区域，也是受人类活动影响突出的地区。 在这些相互作用的过程中，泥沙运动通常起着纽带作用，这使得河口海岸地区泥 沙运动的规律的研究成为相关学科共同关心的问题之一。然而对河口海岸的非恒 定水流的研究人又是泥沙运动的基础，也是开发利用河口海岸各种自然资源、从 事生产活动时所要面临的基础问题之一，通过对其研究指导相关的生产实践，具 有突出的应用价值和科学意义。 因此系统地研究河口海岸垂向二维及平面二维水流数值模型及其在实际工程 中的应用是有意义的。 1 2 国内外研究现状 水环境数学模型可以描述水环境中物质混合、输移和转化规律。它是在分析 水环境中发生的物理、生物和化学现象的基础上，依据质量、能量和动量守恒的 基本原理，应用数学方法建立起来的模型。通过水环境数学模型的求解计算可以 预报水文、水质在时间和空间上的变化，为水资源的管理和控制服务。由于实际 2 青岛理工大学工学硕士学位论文 工程中湖泊、水库、河流、河口区域边界几何形状的不规则以及水动力学方程的 非线性的性质，所以难以做出理论分析解，只能求助于数值解。数值模拟具有成 本低、速度快、资料完备、具有改变各种参数、模拟真实条件及理想条件的能力。 计算水力学中所求解是非线性偏微分( 积分) 方程组，其数值方法的现有数学理论 尚不够充分，严格的稳定性分析、误差估计和收敛性证明等理论工作的发展还跟 不上数值模拟方法的进展，所以需依靠对非线性数学方程的数值实验及数值解与 试验或典型算例结果的比较和物理特性分析，以验证数值计算结果，进而进一步 改进计算方法。 水利工程中各种复杂的流动现象的研究方法主要有：物理模型、理论分析和原 型观测等几种方法【1 叫。 物理模型十分直观，其相似准则可从各物理量表达式的对比关系中导出，往 往可为数值模拟提供一系列必须的计算参数。但物理模型投资大、周期长，精度 受比尺效应和观测仪器的影响。理论分析揭示的是普遍性的规律，如果能用数学 方法求解微分方程得到流场中物理量的解析表达式，就可以得到相应的初始和边 界条件下的全流场的准确信息，这比在某一特定条件下进行的物理模拟或数值模 拟所得到的局部近似值要全面、深刻得多。同时，为证明数值模拟是对某一水力 现象的准确模拟，应该用理论分析的方法对数学模型的相容性、收敛性和稳定性 进行推导证明。最后，数值模拟结果的整理和分析也要借助理论分析工具。理论 分析的缺点是受限于数学工具，只有少量特定条件下的问题，可根据求解问题的 特性对方程和边界条件做相应简化而得到其解析解，目前可解的问题很少。原型 观测是认识的基础和来源，即使只是局部的定性的认识，也是十分可贵的。它不 存在物模缩尺问题，也不用考虑数模计算参数是否适当的问题，它可以为物模及 数模提供设计和计算的参数及最后验证的数据。原型观测方法常受制于外部环境 和人力物力，耗时长，观测点少且不易进行。在计算流体力学中，原型观测、物 理模型、理论分析是研究水流水质的几种基本方法，它们的发展是相互依赖相互 促进相辅相成的【4 】。 相对于这几种研究方法，数值模拟方法有其突出的优点。数值模拟方法可称 为数学模型试验，所有实验条件都以数字形式给出，对流体无扰动，不存在缩尺 效应，不受实验场地和观测仪器影响，同时实验条件和参数容易进行控制和改变， 3 青岛理工大学工学硕士学位论文 能模拟多种因素在复杂条件下的物理过程，超高温、超高压、有毒等恶劣环境对 其并无妨碍。数值实验时间取决于计算机速度和计算方法，较物理模型要快得多。 数值模拟方法具有高效、经济、简便的优点。因此数值模拟作为一种行之有效的 方法愈来愈广泛地被应用于科学研究与工程设计的各个领域【5 堪】。数值模拟也有局 限性：一是依赖于基本方程的可靠性，结果只是有限个离散点上的数值解，并有 一定的计算误差；二是它不象物理模型实验一开始能给出流动现象并定性的描述， 却往往需要有原型观测或物理模型实验提供某些流动参数，并对建立的数学模型 加以验证。 在流体力学及相关领域的科学研究中，随着基本理论、数值计算方法以及计 算机技术的进步，水流及物质输运的数学模型以其特有的优势得到越来越多的重 视，发展十分迅速。尽管尚不能完全取代物理模型实验，但已经在各项工作中占 有重要地位，发挥着巨大的作用，成为科学研究的主要方法之一。由于国内外学 者的长期努力，水流及物质输运的数值模拟已取得了大量成果，并在环境保护、 水资源可持续利用、海岸工程、港口航道工程、水利工程、医疗科学、航空航天、 船舶制造、气候预测、灾害预报、汽车工业、交通运输、石油工业及日常生活等 诸多领域得到了广泛应用。 1 2 1 数值模型分类 从不同的角度出发，水环境数学模型进行不同的分类，一般来讲，可以从如 下几方面划分 9 川】： 1 、按解的过程可以分为确定性模型和随机性模型 对一组给定的输入条件，确定性模型只给出一组确定值，是使用最广泛的一 种数学模型。随机性模型的输入是随机的，其解不具有惟一性。 2 、 按时间的变化可分为稳态模型和动态模型 当水流运动要素和系统的输入不随时间变化时，这种模型称为稳态模型；水 流为非恒定流动，不论输入是否随时间而变，系统内的物质量都将随时间而变， 这种模型称为动态模型；当水流是恒定流动，由于随时间变化导致系统内的物质 量随时间变化，描述这种现象的模型称为准动念模型。 4 青岛理工大学工学硕士学位论文 3 、按空间变化主要分为零维模型、一维模型、二维模型【1 5 1 7 】和三维模型【1 8 2 0 】 ( 1 ) 零维模型：水体处于完个混合状态，x 、y 、z 三个方向的水动力、水质要 素都均匀分布，这种模型称为零维模型。 ( 2 ) 一维模型：系统内质点的水动力、水质要素只在一个方向有梯度存在， 在另外两个方向上均匀分布的模型称为一维模型。它包括垂向一维模型( 适用于温 度分层的湖泊) 和纵向一维模型( 适用于河流) 。 ( 3 ) 二维模型：系统内质点的水动力、水质要素在两个方向上有梯度存在， 另一个方向均匀分布，称为二维模型，二维模型分为沿水深平均的一维模型和沿 宽度平均的一维模型。对于河口、海岸、湖泊、水库等地区属于宽浅型区域，即 水平尺度大于垂向尺度，可将实际的流体运动的三维模型沿垂向积分而得到的平 面二维数学模型【2 1 2 5 l 得到了广泛的应用。而对于另外窄深区域，流速、水温或浓 度在深度方向上分布具有明显差异的情况，如垂向有回流旋涡、密度分层等，有 关参量( 如流速、温度、浓度) 的垂向变化要比水平横向的变化为大，用垂向二维 模型比较合适 2 锄7 1 。 ( 4 ) 三维模型：系统内质点的水动力、水质要素在三个方向上都有梯度存在 的模型称为三维模型。 此外，按水体类型可以分为湖泊及水库模型、河流模型、河口模型和海洋模 型；按污染物质在水体中的运动特性可分为对流模型、扩散模型、对流一扩散模型； 按反应动力学性质可分为纯转移模型、纯反应模型、转移及反应模型、生态模型 等。 1 2 2 常用离散方法 在水环境数学模型中常用的数值离散方法主要有有限差分法、有限元方法、 有限体积法等。 1 、 有限差分方法( f d m )有限差分方法【2 8 3 2 l 是计算机数值模拟最早采用的 方法，至今仍被广泛运用。该方法将求解域划分为差分网格，用有限个网格节点 代替连续的求解域，数学概念直观，表达简单，其解的存在性、收敛性和稳定性 5 青岛理工大学工学硕士学位论文 早已有较完善的研究成果，是比较成熟的数值方法，目前应用最广。由于实际应 用中采用的时间和空间差分形式不同，差分法又可以分为显式、隐式及显隐式交 错等方法。显式差分格式应用较早，相对简单，可避免试算，但为了保持其稳定 件，需严格遵守柯朗条件( c , x t 6 x 1 ，c 为小扰动波的波速，c = 办+ 材) ，时 间步长和空间步长受到限制。从理论上讲，隐式格式是无条件稳定的，但在实际 应用中，由于空间、时间步长为有限量，其时间步长也有一定的限制。隐式差分 格式的优点是时间步长可以取得较大，稳定性能好，但计算过程中需迭代，计算 量较大。交替方向隐格式法( a d i ) 3 3 1 是由d o u q l a c e 和r a c h f o r d 等1 9 5 5 年提出的， 后来被l e e n d e r t s e 结合交替网格建立起来并首次用于计算平面二维流场。a d i 方法 是一种显隐格式交替使用的有限差分格式，该方法同时具有显式和隐式两种差分 格式的优点，与完全隐格式相比较，它不必每一时间步骤都要求解一个大型代数 方程组，因而所需的内存少，计算量也相应减少。同时a d i 方法不像显格式那样， 在计算中易出现波动现象，因为显、隐格式在坐标轴上交替使用，使误差的增长 量相互抵消。因此a d i 方法有较好的计算稳定性和计算精度，目前己广泛应用。 然而，a d i 使用交替隐式格式，不能恰当的考虑不同方向流动之间的作用，数值 解存在流速向量向某坐标轴倾斜的趋势【3 4 】，因此就降低了其效果。此外发展有特 征线法【3 5 】。在2 0 世纪7 0 年代初，前苏联学者y a n e n k o 等提出破开算于法【3 6 】，按 维数或按方程的性质对控制方程进行算子破开，从而简化计算。但目前许多学者 对破开算子法提出了疑义，因而该方法应用不是很广。王船海、程文辉【37 】研究发 现破开算子法在实用上并非对任何情况均适用的，主要是由于破开算子法将方程 破开后引进了分步误差的缘故。二维模型计算中，由于使用t a y l o r 级数展开，f d m 一般只适用于矩形或正交曲线网格，在计算域概化和数值解精度方面，存在着根 本性的困难。 2 、 有限元法( f e m ) 有限元法【3 8 】的基础是极值原理和剖分插值，它吸收了 有限差分中离散处理思想，同时采用了变分计算中逼近函数及对任意形状( 三角形 或四边形) 的许多微小单元进行积分处理的合理方法，因而具有很广泛的适应性， 特别适合于几何、物理条件比较复杂的问题。该方法具有较强的适应性，计算精 度较高，但存在计算格式复杂、计算量大的问题。常见的有限元计算方法有直接 法、变分法、加权余量法及能量平衡法等。其中变分法类的里兹法( 1 9 0 9 年) 、加权 6 青岛理工大学工学硕士学位论文 余量法类的g a l e r k m 法和最小二乘法常用于河口数值模拟。采用常规的有限元方 法时，对于对流效应比较强的情况，常常由于有限元网格不恰当而造成数值解的 失真或震荡。在有限元方法中通常是附加人工粘性和采用迎风格式，但这样做使 方程失去了加权余量的数学意义，从而近似方程不满足相容性，即在保证稳定的 条件下，失去了精度。后来又提出了流线迎风p e l r o v g a l e r k i n 有限元法【3 9 1 ，简称 s u p g ，该方法在稳定性、收敛性及精确度等方面已有很大的提高。在多级数学模 型计算中，因有限元方法贮存量比较大，且大型系数矩阵求解较困难等，直接影 响着计算速度，因而在非恒定问题及其它对计算速度要求比较高的问题中应用不 是很多。 3 、 有限体积法( f v m ) 有限体积法1 1 4 0 - - 4 3 】又称控制体积法，它同有限元法 一样将计算区域划分为若干互相连接又不重复的单元或控制体，每个网格点周围 有一个控制体积；将待解的微分方程对每一个控制体积积分，得出一组以计算结 点上物理量为未知数的离散方程，其中的未知数是网格点上的因变量的数值。有 限体积法的基本思想易于理解并能得出直接的物理解释。 在数值模拟的过程中，各种方法均有其自身的优点和适应性，在实际计算时 选择什么数值方法应根据所研究问题的特点和计算精度要求，以及研究者的习惯 而定。 1 2 3 数值模拟两个常见问题 ( 1 ) 不规则边界的处理 目前比较常用的数值方法主要有有限差分法、有限分析法、有限元方法、有 限体积法以及特征线方法等。但对于边界形状比较复杂的区域，一些数值格式的 应用受到了限制。对于这一问题，很多学者作了许多工作。在2 0 世纪6 0 年代 w i n s l o w 最早提出用伯微分方程变换生成网格，随后不少学者对这种方法进行了改 进。在2 0 世纪7 0 年代t h o m p s o n 等提出了贴体坐标系计算网格的生成方法，它使 计算网格与计算域边界相重合，从而使边界条件能准确地引用到计算网格节点上， 并使得变换后的网格成为规则的矩形网格。但数值网格生成方法计算比较繁琐， 变换后的基本方程也相当复杂，因而至今未能在二维数值模型计算中推广应用。 7 青岛理工大学工学硕士学位论文 ( 2 ) 动边界问题 由于河床冲淤变化和水位的波动，在河岸及洲滩等区域都会遇到动边界问题。 潮流的运动受各种因素的影响，如天文因素、气象因素、海洋因素和地形因素等， 而这些因素都是随机的。因此在模型中应考虑输入的随机特性，进而对输出进行 概率统计分析。孔令双】曾在胶州湾流场的数值模拟中考虑了随机因素的影响。 但这只是一个初步的尝试，需要更加深入地研究。l e e n d e r t s e 在1 9 7 0 年首先提出 了动边界的方法。由于河床冲淤变化和水位的波动，在海岸区域会遇到动边界问 题。动边界数值模拟的因难主要在于沿动边界法向的流动不同于明渠均匀流，非 恒定和非均匀性强，常用的曼宁摩阻公式等在形式上难以套用，再者是因水深很 小，对离散格式的要求很高，要求数值解不产生假振，保证水深总是大于零。 对 于动边界的处理一般有两种方法：一种是追踪动边界的准确位置进行模拟，曹祖 德【4 5 】利用边滩上活动水边界的流速来确定边界的位置；另外一种方法是把整个计 算区域作为固定计算域，不管有水无水照常进行计算，在计算过程中采取一些方 法对露出水面部分进行处理，汪德灌【4 6 】采用附加人工粘性项的方法处理这个问题。 李狮、何少苓等h 7 】引入了不透介质法处理动水边界，自由水面边界和床底边界， 将动边界问题转化为定边界问题求解，对大尺度水体中的恒定与非恒定流动进行 了成功的模拟计算。另外，不少学者提出了不同的动水边界处理方法，如窄缝法、 水边线步进法、井点法等。 1 2 4 数值模拟的计算网格 流体数值计算是在离散的计算网格上进行的，因此采用的计算网格是否适合 对计算的成败起着关键的作用，好的计算网格应该能完整准确地反映固体边界的 情况，疏密适当，既能够较好的控制计算量，使计算有较好的收敛性，获得经济 性，又能使计算结果误差不致过大，能真实地反映流体的运动情况。随着计算手 段的不断改进，计算网格的生成技术也在不断的发展，以至于可以独立成为一项 网格生成技术5 4 1 。 传统的计算流体力学多采用有限差分法，有限差分的一个重要问题就是计算 网格的问题。早期的差分采用矩形网格，不需要对控制方程进行过多的转换，计 算简便，易于编程，但对于复杂边界，采用矩形网格无法紧密贴合，将会形成阶 r 青岛理工大学工学硕士学位论文 梯状边界，边界条件难以处理。由此产生了变步长网格、三角形网格和四边形网 格等方法。 边界拟合坐标【5 5 】的计算网格仍然是曲线坐标下的矩形或长方体网格，也称为 结构网格。无结构网格一般采用三角形或不规则四边形网格构成，早期用于有限 单元法，直至二十世纪八十年代与控制体积法结合并经多年研究，逐步改进才得 以广泛应用。无结构网格与边界结合较好，便于控制网格密度，网格生成较边界 拟合法容易。但其网格排列不规则，节点多，计算量大，同时也要对控制方程进 行处理。从实际应用来说，可以依照实际情况，灵活的将有结构网格和无结构网 格结合使用，即结构非结构杂交网格，合理的采用结构非结构杂交网格可以使 得两种网格各展所长，优势互补，得到更为理想的计算结果。 自适应矩形网格被用来计算几何形状复杂的流场，该方法可以一次性生成计 算所需的网格，简单省时，无需从物理空间到计算空间的转换，易于在流场计算 中实现，网格容易加密，可提高计算精度。但由于采用的仍是矩形网格，在描述 外形的精度上要较贴体网格和无结构网格低，不能做到与边界完全一致。 水流计算中垂向网格的划分常用的是盯坐标法【5 以5 9 1 。在计算区域的水深、自 由表面、水底地形随时空变化，不便于做差分计算时，为了使计算网格连续跟踪 自由表面运动，精确拟合自由表面和水底边界，垂向做仃坐标变换。其基本原理 是：通过垂向的坐标变换，将水面至河床( 海床) 床面的坐标z 变换为仃，使计算 域垂向变化值在0 一l 之间，再在垂向划分网格。该方法原理简单，易于理解，计 算方便，对控制方程变换的复杂程度较低，易于实现编程计算，在一般情况下可 以获得较为理想的结果。仃坐标变换的主要优势是自由水面的处理相对简单，但 也存在一定不足，如变换后的控制方程形式更为复杂，计算结果需加以转换，使 计算难度加大还可能引入新的误差，某些情况下容易发散等。因此对于起伏剧烈 的地形，流体随时间变化剧烈的情况，此计算方法不太适合。 总之，计算网格的生成对计算结果非常重要，每种方法都有其价值和实用性。 本论文在垂向二维水流数值模拟中采用了盯坐标法，通过这种坐标变换，将不规 则区域变为规则区域，使计算网格能连续的跟踪自由表面运动，精确的拟合自由 表面和水底边界。 9 青岛理工大学工学硕士学位论文 1 3 数值模拟发展趋势 从国内外进行的数值模拟的进展分析，数值模拟将向以下方面发展【6 嘶2 1 。 ( 1 ) 建立数值模拟与物理模型相结合的复合模拟 数值模拟与物理模型是目前最常用的两种模拟方法，它们各有其模拟方法和 特点，这两种方法均得到了很大的发展。自七十年代以来又形成了一种新的模拟 技术一复合模拟，它是吸收了数值模拟和物理模型各自的优点，将两者结合在一 起的模拟技术。目前，复合模拟己显现其特有的优越性和生命力。 ( 2 ) 发展多因素全过程的综合数值模拟系统 水环境数值模拟是极端复杂的，其中必然包含着各种动力因素和流体的运动 状况，另外还有具体的地形地貌。动力因素有风吹流、径流、潮流等；流体运动 有悬移质、推移质和浮泥流等：地貌演变有冲刷和淤积过程等。因此，建立一个 多因素全过程的综合数值模拟系统是今后发展的方向。 ( 3 ) 与“3 s ( 地理信息系统g i s 、遥感系统r s 、全球定位系统g p s ) 结合 建立显明精确的仿真图形动态显示系统 仿真图形动态显示系统克服了数值模拟不能直观显示的缺点，g i s 可把复杂多 变的自然、社会变化以及变化过程以图形、图像的方式进行数字化处理。在其空 间和属性库输入水域基本数据、水文及污染源数据，利用其空间数据库采集、管 理和分析能力，可使水质监测与评价更科学、准确。通过数值模拟计算，可得出 反映水质变化特性的断面位置，并以图像显示水域水流变化的空间特性、统计特 性和变化趋势等。目前r s 、g p s 已与g i s 相结合进入了能及时提供各种对地观测 的具有整体性的动态资料，并对这些资料进行分析与处理的新阶段。 1 0 青岛理工大学工学硕士学位论文 第2 章垂向二维水流数值模型 在海岸河口、高山峡谷河道型水库地区，其水较深，水面较窄，宽深l k e , j , ， 有关参量( 如流速、温度、含盐度、含沙量等) 的垂向变化要比水平横向的变化为 大，如潮汐通道、河口深渊段等，此时可采用垂向二维水流数值模型 2 】 2 弛7 】【4 3 】。 垂向二维水流数值模拟的基本方程是将n a v i e r s t o k e s 方程沿河宽积分，并取 横向平均而得。在推导垂向二维水流运动的基本方程时则有如下假定： ( 1 ) 水质点横向流速对时间空间的偏导数可以忽略不计； ( 2 ) 流动要素( 流速) 在控制断面上沿横向分布足够均匀； ( 3 ) 两岸测向摩阻力相等，无滑移； ( 4 ) 两岸边壁无滑动条件，即甜= 1 ，= w = 0 ； ( 5 ) b = 8 ( x ) ，忽略柯氏力和质量力。 2 1 控制方程 垂向二维数学模型的基本方程组是将n a v i e r s t o k e s 方程沿河宽积分，并将横 断面概化为矩形，假定河宽不随时间变化，运用含参变量的牛顿一莱布里兹微积分 公式，并忽略一些二阶小量而得到的，系统坐标系如图9 - 1 所示。 连续方程： _ a ( b u ) + b 娑：0 ( 2 1 ) c o z 沿z 方向积分连续方程可得垂向积分方程，即水位方程： 丝+ 土旦f b ( u d z ) ：0 ( 2 2 ) 0 tb 反、j - # 动量方程： 要+ 鼍墨学+ w 罢+ 土罢一1a ( b e ，罢) + 昙( b e ：罢) ( 2 3 )+ 一立二+ w + 一二= 一l + i l ( z 一3 ) a tba ) ca z p 瓠 ba ) ca ) c 。8 z 、 2 8 z 。 压力方程： 青岛理工大学工学硕士学位论文 p = p o 七g 毫p d z 式中，甜、w 分别为横向平均流速的x ，z 向分量； 占，、占：分别为x ，z 方向的紊动粘性系数； 孝水位( 水面至基准面的垂向距离) ； h 河底至基准面的距离； p 水密度； p 。水面压力； p 浑水压力； g 重力加速度； ( 2 4 ) b 河宽。 方程( 2 1 ) 一( 2 4 ) 即为垂向二维水流数值模拟的基本方程组。 2 2 定解条件 图2 - 1 系统坐标系 1 2 x f 甜1 青岛理工大学工学硕士学位论文 数值计算的定解条件有二：一为边界条件；二为初始条件。垂向二维数学模 型的边界条件有自由表面条件、水底边界条件及上下游边界条件。 2 2 1 边界条件 ( 1 ) 自由表面边界条件 有风时：生h 丝o z = 7 r w x ， f 。2 成s 1 w 无风时，f 。= 0 。 式中，f 。风对水面剪切力在x 方向的分量； “。风速在x 方向的分量； 成空气密度； w 风速。 ( 2 ) 水底边界条件 有滑移时：鲁老= 7 c b x f 如= 昭万矧c 2 c ：1 h 1 1 6 n 无滑移时： u = w = 0 式中，底部剪切应力； 石垂向平均流速； c 谢才系数； 卜水深； n 粗糙系数。 1 3 青岛理工大学工学硕士学位论文 ( 3 ) 上下游开边界条件 已知水位过程线：孝= 孝( x ，z ，f ) 已知流速过程线： “= u ( x ，z ，f ) 2 2 2 初始条件 初始条件虽然对稳定解的值没有什么影响，但却直接影响到稳定解所需计算 时间的长短，因此选择合适的初始条件对数值计算同样很必要。 2 3 坐标变换 由于水深随时间、空间不断发生变化，因此计算区域也相应发生变化，这样 就给计算造成不便。垂向可做仃变换，将不规则计算域变为规则计算域。仃坐标 法的基本原理是通过垂向的坐标z 变换为盯，使其值在o 一1 之间，再在垂向上划分 网格。该方法原理简单，计算方便，对控制方程变换的复杂程度较低，相对较易 实现编程计算，在一般情况下可获得较为理想的结果。坐标变换后，可使计算网 格连续的跟踪流体的自由表面运动，精确的拟合自由表面和水底边界，这样原来 的垂向二维方程组也要作相应的变化。经过仃坐标变换之后，方程增加了若干非 线性项，增加了方程求解的难度，而在差分离散范畴内，这种处理方法是比较适 合的，否则，原方程将无法通过差分法求解。垂向二维水流数学模型的原坐标系 与新坐标系如图2 2 所示，其中左边图为原坐标系，右边图为经过坐标变换后的 新坐