潮流与泥沙数值模拟研究进展

潮流与泥沙数值模拟研究进展

趋势与沉积物值模拟始于20世纪60年代，是一种新兴的综合模拟技术。根据数学的概念，它融合了具体的工程技术，如数学、沉积物动力学和计算。然后使用数值计算模型模拟特定的物理现象，以模拟自然或人工过程。与物理模型相比,数值模拟具有经济、快速、修改方便等优点,并且不受比尺的限制,已广泛应用于解决各种海岸、河口动力问题的研究,如河口整治工程、环境保护工程等。此外,还可以用来揭示许多海岸、河口物理现象,例如:通过潮流、泥沙数值模拟,可了解水质点在任意时刻的流速、流向和水量的输送、湍流的强度以及涡旋、环流状况和泥沙输移以及底床变化等。随着计算机技术的迅速发展和水动力学及泥沙动力学理论的不断完善,数学模型在理论上和计算技术上也发展甚快。从总体情况而言,20世纪70年代以河口一维潮流计算为主,进入80年代后,除继续应用推广一维潮流计算外,大多已采用二维潮流数学模型,并根据实际工程需要,配合以泥沙、温度、盐度和污染物等物质输移模型。90年代以来,随着工程项目的需要,三维模型的研究和应用日趋广泛,并取得了一些颇有价值的研究成果。但由于其复杂性,目前往往采用简单的模型来代替复杂的数学模型,即用在一个或两个方向上的积分平均,将三维的问题降为多个二维或一维问题以简化三维流动结构。有关潮流、泥沙数值模拟的研究工作,国内外学者对不同方面做了若干总结性的评述。李孟国、曹祖德、严世强、熊德琪分别就海岸河口地区潮流场的数值模拟方法进行了较为系统的归纳总结和评述,并对潮流数值模拟的发展趋势进行了分析。曹祖德、孔令双等对海口、河口水动力数值模拟的发展动向做了多方面的预测。Abbott和Davies等对近年来三维潮流数学模型的进展做了较为全面的概括。本文回顾总结了近年来国内外潮流、泥沙数值模拟的研究,评价分析了数学模型的多种数值解法,特别是有关二维、三维潮流数值模拟以及水沙数学模型近年来的进展,并根据目前潮流、泥沙数值模拟技术中存在的不足,讨论了进一步研究的方向和发展趋势。1非恒速模型及边界条件自然界中,河道或河网的水流或洪水激流运动、峡口或潮汐通道的潮波运动、三角洲网河口的潮波顶托、电站日调节的水流泄放、闸门启闭后的水波运动等,这类水利水电工程和海岸工程中经常遇到的问题都可用数值模拟加以解决。由于计算水域狭窄,沿流向方向上的尺度远大于其它两方向上的尺度,以上问题都可以采用一维数学模型。1871年,圣维南(Saint-Venant)根据Boussinesq提出的缓变流定义而建立的非恒定流方程(即圣维南水力方程组),至今仍被工程界采用,成为一维非恒定流数值模拟的基础。一维河道潮流数值模拟分单支和河网2种情况,在模拟单支明渠一维非恒定流时,通常用的计算格式有:蛙跳格式,Lax格式及Preissmann格式等。其中,蛙跳格式和Lax均采用显式差分,差分方程较为简单,使用方便,但稳定性和精度不够高。Preissmann采用空间平均的4点隐式解法,并将非线性代数方程组线性化以提高计算速度,在国内外均得到广泛应用。然而,Preissmann的空间平均4点隐式格式在理论分析和应用实践上都表明对于具有小扰动波速的河床变形计算,在小柯朗数的情况下常出现数值扰动和虚拟振荡。对此,人们通过引入空间权重因子,将格式改造成4点偏心隐格式(国际上常称GeneralPreissmannScheme)以达到削减数值振荡的目的,这实质上是通过引进数值权重系数增加人工阻尼效应来抑制数值失稳。在模拟河网非恒定流时,经常采用的格式的格式有:Preissmann格式、Cao格式及罗肇森格式等。应用Preissmann格式时,假设各汊道在汇流处水头相等,先沿汊道起始点到终点求出各断面的循环系数,然后结合边界条件,求出各断面的水位和流量。Cao格式假设在汇流处各汊道水位相等,流量连续,其计算步骤与Preissmann格式类似。罗肇森格式采用三级联解法,即将参加计算的方程分成微段、河段、汊点三级,采用逐级处理再联合运算的方法,求得河网中各计算断面的水位、流量等值,解决了分汊河道非恒定流时分流口和汇流口的水位和汊河的分流量的计算。一维数学模型的定解条件包括边界条件和初始条件。通常情况下,潮流河道一维模拟中,下游边界用水深或水位控制,而上游边界条件则由流量或流速控制。初始条件通常根据实测资料差补给出初始面处的水位及流量。由于边界的控制和阻力项的调节作用,使初值的误差随计算时段增长而逐步消失。因此,初值数据可用比较简单的办法近似地给出,而采用适当增长计算时段的办法消除初值误差的影响。2维数值模拟技术在海岸、河口、湖泊、大型水库等广阔水域地区,水平尺度远大于垂向尺度,水力参数(如流速、水深等)在垂向方向上变化要小于水平方向上的变化,其流态可用沿水深的平均流动量来表示,因此可采用平面二维水动力数值模拟技术。而在另外一些水域,如窄深潮汐通道、窄深河口地区,有关参量(如流速、温度、含盐量、含沙量等)的垂向变化要比水平横向的变化大,这时可采用垂向二维数值模拟技术。在解决二维潮流数值模拟的过程中,有多种数值解法可供选择。这些数值解法就划分标准的不同,可以大体分类如下:从离散方法上分,有差分法、有限元法和有限体积法;从适应物理域的复杂几何形状上分,有贴体坐标变换及σ坐标变换;从时间积分上分,有显式、隐式、半隐格式;从求解方法上分,有ADI法、迭代法、多重网格法以及并行计算技术;从干湿、露滩动边界的处理上分,有固定网格和动态网格技术;从悬沙输运重力沉降项的处理上分,有源项化和对流化的做法。上述方法在实践中已有尝试,下面给出的是几种常用且有效的计算方法,并分析其优缺点。2.1显式法和隐式法有限差分法是将微分方程中的各项微分离散成在微小网格上各临近格点的差商形式,得到一个以各节点上函数值为未知变量的代数方程,称差分方程,由此解得各网格点的函数值。用有限差分法离散求解方程可分为显式法、隐式法、半隐半显式法3类。显式差分格式简单,容易理解,但一般说来收敛性和稳定性较差,结果往往不令人满意。隐式差分克服了显式差分的稳定性差、时间步长受到限制、精度低等的缺点,相对于显式法,隐式法稳定性高,因而时间步长可以取得很大。半隐半显式法兼备了显式法和隐式法各自的一些优点,在实际应用中应用较广,其中具有代表性的是ADI法。相比于显式法,ADI法具有计算稳定、精度高等优点,又较隐式法减少了许多计算量,是目前进行二维潮流数值模拟比较经济和行之有效的一种方法,但由于使用正方形或矩形网格,ADI法处理固边界不灵活。采用矩形网格时,差分方程简单,但稳定性差,且不易拟合复杂固岸边界,因而较少使用。相比之下,三角形网格、四边形网格以及多边形网格等形式具有适合复杂地形和固岸边界的拟合、布设随意、稳定性及收敛性好、精度较高等优点,因而应用较广。2.2有限元法的应用差分法是通过用差分项代替微分项,直接对方程进行离散化的方法,而有限元法是将给定的连续介质划分成许多具有合适形状的微小单元,在单元内选定一些适当的节点,将微分方程中的变量改写成由各变量或其导数的节点值与所选用的插值函数组成的线性表达式,借助于变分原理或加权余量法,将控制微分方程转换成控制所有孤立单元的有限元方程。然后,将所有局部单元上的方程汇集成总体的微分方程组或代数方程组,再换上应有的边界条件和初始条件,便成一个完备的代数方程组,解此方程组,可以求得各节点的函数值,从而求得微分方程在整个计算域上的数值解。有限元法解题能力强,可以适应复杂的几何形状,统一处理各种边界条件。有限元法虽然在流体方面应用较晚,但目前已被广泛地应用于各种流体力学问题中。用有限元法求解二维浅水环流方程始于1973年。时至今日,国内外学者已提出了不少基于有限元法的数学模型,其主要思想都是旨在对隐式有限元繁杂求解加以改进,如吕玉麟等用分裂时间法作数值积分以取代有限元方程组求解所用的多次迭代,谭维炎等引入单元影响系数以便于计算有限元方程中的积分项等。有限元法具有网格划分灵活、拟合复杂岸边界和地形容易等优点,但是数学推导繁琐,公式众多,通常要求解庞大的代数方程组,不易编程,这些问题都限制了有限元法的广泛应用。2.3拟线性双曲型偏微分方程特征线法又称特征理论法或特征差分法,是解双曲型偏微分方程的最精确的数值方法,较多用于定常二维、平面和轴对称的、无旋超音速流动问题,目前已经有了完整的数值算法来确定超音速流场中内点上的流动参数。其基本思想在于潮流控制方程属一阶拟线性双曲型偏微分方程,利用二维空间的特征理论,可导出两族特征曲面和相应的特征关系式,对待征关系式进行离散求解可得到变量(速度、水深等)的数值解。对特征关系式的离散可采用特征偏心格式等多种方式进行。特征线法物理概念明确,数学分析严谨,计算精度较高,但其时间步长和空间步长的比值要受到一定稳定条件的限制,导致其时间步长只能取很小的值。如果分析过程较长,就需消耗很长的计算时间。特征线法适合应用于水位、流速变化剧烈的强潮河口、海湾潮流计算,尤其是在涌潮地区,这是其它各种方法所不能比拟的。因为在涌潮地区,潮头破裂,沿x方向发生不连续现象,这时用一般差分法不能求解,但它的特征线仍然存在,因此,仍可用特征线法求解。然而,和隐式差分一样,这种方法在多维问题中显得很繁琐而难以应用。2.4破开算子法的计算步骤破开算子法是苏联学者Yancnko等于20世纪50年代中期提出。到80年代初我国学者开始使用破开算子法解决实际潮流问题,并在破开算子法的基础上提出了各种不同的分步算法,例如:分步有限差分法、分步杂交法、准解析法等。分步法有可能导致较大误差,胡庆云等对此进行了研究,提出了分步误差和比较方程等概念。分步误差指由于算子分裂、分步计算引起的额外的截断误差项。以往研究者认为它是高阶微量,可以忽略,可是对于具体实际问题,固定的时间步长,分步误差不一定是微量。破开算子法计算模型的截断误差,在计算上的困难是含有中间变量,中间变量是联系分步计算格式的过渡值,没有具体的物理意义,在算子非线性时,中间变量不能关于未知解函数显式地求解。文献将破开格式的截断误差处理为不分步计算模型的截断误差与分步误差之和,这两种截断误差,由于其产生的途径不一,对近似解的影响也就不同。分步法是求解二维潮流的一种简便、易行而有效的数值计算技术,在处理上具有相当的灵活性。具体分成多少步,如何分,在理论上并没有严格的限制,而且在分步之后对各个分步采用什么样的格式进行计算,也是相当灵活的。原则上只要在分步之后,寻找出一个尽量正确、合理的数值离散方式,就能取得良好的计算成果。然而,破开算子法并不是对所有问题都是行之有效的,当计算的时间步长较大时,沿坐标轴方向分步将产生轴化效应。此外,分步法在处理边界条件时也不太灵活。2.5子域解析解的有限方程离散化有限分析法是由Cheng(1981)提出的一种方法,它克服了有限差分法在求解不可压粘性流体时在大雷诺数下的各种困难,避免了差分近似中的各种数值效应,是求解大雷诺数的各种流体力学问题行之有效的方法。有限分析法把求解域划分为多个有限小的子域,在每个子域内求解基本方程的解析解。例如不可压纳维-斯托克斯方程是非线性的,求解析解时需将方程在子域范围内局部线性化。子域内的解析解把该子域的一个内点上的参变量与该域边界上各点的参变量联系起来,以代数方程表示。由于每个子域边界上的点即为其邻域内的点,而每个内点有一个代数方程,因此全场所有内点的代数方程式组成一个完备的联立代数方程组,解此代数方程组,即可得原方程的数值解。求解域边界的参变量为求解联立代数方程组时嵌入的边界条件。有限分析法是九点格式,它不仅与一般的5点中心差分格式中的边界点有关,而且与4个角点也有关。子域边界各点对中点上的参变量的影响也随来流的不同而有所改变,这更真实地反映了流体运动的物理现象,正是有限分析法的独到之处。因此,用有限分析法可以在较大的雷诺数范围内得到精度较高的数值解。在潮汐河口水力计算中,除了上述几种方法外,还有很多种方法,这里不一一介绍。从上述介绍中,可以看出,各种计算方法都有其优缺点,一种计算方法对某一问题是有效的,而对另一问题则可能是不可取的。因此,在实践过程中,应就某一具体问题,进行具体分析,选择最有效的计算方法。2.6拟合坐标变换边界拟合坐标法是一种坐标变换法,它和三角形、四边形网格差分法、有限元法、单元积分法等一样,是为了克服经典的矩形网格差分法不能精确拟合天然水域不规则的水陆边界的困难而发展起来的。该方法的基本思想是使用泊松方程作为坐标转换方程,将笛卡尔坐标系(x,y)上具有复杂边界的计算域转换到拟合坐标系(ξ,η)(亦称曲线坐标系)上具有规则边界的矩形域,将(x,y)计算域上的疏密程度不同的曲线网格变换成(ξ,η)域上的矩形网格,将(x,y)上的方程组变换成(ξ,η)上的方程组并求解,利用两坐标系坐标之间和流速分量之间的关系得到(x,y)上计算域的解。从而完成方程组的求解,刘玉玲、张宗孝将拟合坐标变换技术应用于三角形网格生成中,形成一种方便、有效的三角形网格生成方法。该方法直接利用边界网格节点的分布信息来控制区域内部网格节点的分布,在数值模拟河道、湖泊等复杂边界的流速场时,可快速获得一种适应物理量场变化情形的三角形网格系统,使该网格系统的流场计算精度明显提高。边界拟台坐标法可以准确贴合边界,在计算域内布置灵活,在二、三维潮流计算中已得到广泛应用。尤其在计算域边界较为复杂的情况下,更能体现出此方法的优点。3数值方法在数值模拟中的应用20世纪90年代以前,由于受计算机计算能力和计算方法的限制,海岸河口潮流数学模型一般采用二维模型,三维模式还较少用来解决实际生产问题。随着计算机和数值模拟技术的发展,三维模型的应用也越来越广泛,因此,很有必要发展和完善海岸河口三维潮流数学模型。河口与海岸环境中的流动具有明显的空间三维特性,从理论上讲,可以从RANS方程(ReynoldsaveragedN-Sequations)出发建立三维水动力学数学模型。考虑到近海潮流的水平空间尺度远大于垂直空间尺度,因此可以看作一种准平行流动,这时水质点运动的垂直加速度比重力加速度小得多,垂直方向的流体运动方程退化为静压方程,即所谓的静压假定。这样,我们可以用水位变量ζ(x,y,t)代替压力变量,得到三维浅水方程,使得需求解的数学问题得到简化。基于三维浅水方程,Leendertse的工作具有开创性,他在垂直方向采用固定分层法,即将计算水域划分为若干固定层,在每层中沿水深积分使之成为二维问题,并用ADI格式进行数值离散。Kim等应用固定分层方法建立了海湾三维潮流、盐度模型。为了更好地模拟河床地形的变化,Burchard等将Philips提出的σ坐标变换应用到河口与海岸三维模型中。为了较严格地确定涡粘性系数和扩散系数,湍流模式理论也在潮流模型中得到了应用。以Princeton大学Mellor教授为首的海洋动力环境数值模拟小组从20世纪80年代开始,致力于三维数值模型的开发与应用研究,其代表性软件为POM(PrincetonOceanModel)。该模型在垂直方程采用σ坐标系,在水平方向采用正交曲线坐标系,但数值求解过程仍在物理平面上进行。为了提高计算效率,POM采用了模态分裂法,将自由面的三维流动问题分成表面波的传播问题(外模态)和内波的传播问题(内模态),对外模态仍用显式差分,Wang将三维水动力学模型的控制方程用普遍张量的形式表示,并复演了美国Galveston海湾的三维流动和盐度分布。为了适应我国大型河口研究和重大工程建设的需要,近年来国内学者广泛展开对河口三维水动力学数学模型的研究,针对三维模型,提出或改进了多种计算方法。赵士清采用与Leendertse类似的固定分层方法,建立了较简单的数值模式,对长江口南槽和口外海域的三维潮流进行了数值模拟。韩国其等采用σ坐标系和算子分裂技术建立了三维潮流数值模型。窦正兴等采用σ坐标系和模态分裂法对渤海湾的三维潮流作了数值模拟。宋志尧等基于模态分裂法和ADI格式建立了三维潮流的计算模式,并应用于海岸辐射沙洲的潮流场分析。刘桦等建立了涉及河口密度分层效应的三维潮流、盐度数学模型,垂直方向为σ坐标系,水平方向为直角笛卡儿坐标系,经模态分裂后,外模态用改进的ADI法求解。该模型首次复演了完整的长江口三维潮流场,并对三维盐度场进行了初步模拟。闫菊等运用一个三维正压湍流封闭数值模式,成功模拟了胶州湾M2分潮的潮汐与潮流分布。张越美、孙英兰基于河口,陆架和海洋模式(ECOM模型),引入干湿网格法模拟潮滩涨落的改进,建立了渤海湾三维变动边界潮流模型。杨陇慧等应用三维高分辨率非正交曲线网格河口海洋模式,模拟了以长江河口、杭州湾及领近海区作为整体的4个主要分潮。4确定条件和边界处理4.1有条件时的情形下,主要有三种方式显然,潮流模型的控制方程是非定常的,所以在求基本方程的定解时不仅要有边界条件,而且要有初始条件。初始条件可用2种方式给出:一种是由已有的实测资料用平面内插得到整个计算区域上初始时刻的各函数值,由此给出的初值比较精确,但计算准备工作量较大,另一种方式是选定某一时刻(通常是憩流时刻),将函数初始值近似地认为是常数。由于水流初始条件的误差在正确的边界条件控制下会很快消失,因此,通常我们选定后一种方式,取初始条件为常数,这样既简便,又能达到计算要求。4.2闭边界上的流速计算结果的正确与否与边界条件是分不开的,它不象初始条件那样,误差在计算中可以逐步消失,边界条件的误差会直接影响整个计算的精度。在潮汐河口的计算中有两种边界,即开边界(或称水边界)和闭边界(或称陆边界),它们有不同的边界条件。在闭边界上,根据流体在固壁上不可穿越的原理,在不考虑渗透的情况下,我们可以认为闭边界上的法向流速为0;另外认为流体是粘性的,水体与固体边界是无滑动的,即流体在闭边界上的切向流速也应为0。但是在实际计算中,由于我们不可能将边界节点取在实际边界上(即水深为0处),所以在闭边界上,特别在水深变化剧烈的地方,如狭窄水道,应用水流无滑动条件往往得不到令人满意的结果。一般情况下,我们认为闭边界的切向流速不为0,在某些特殊情况下,可以采用流速边壁函数作为闭边界条件。在开边界上,一般均采用实测水文资料作为开边界条件,但这里必须指出,水流计算中,可作为边界条件的有潮位、流速和流量。在有选择的情况下,通常应选水位作为边界条件较为适宜,因为水位沿横断面的变化往往较小,即所谓的横比降较小,例如,在长江口北槽,横比降只有1‰,而流速沿横断面变化剧烈,且与河床形态有关,不是一种单一的曲线分布,很难确定测点以外各点的流速过程线。4.3开边界与闭边界选择计算区域往往根据工程的要求、河道的自然条件和实测水文资料等加以综合考虑。计算区域的选取关键是开边界的确定。选取开边界时,要保证开边界必须离工程区足够远,使工程影响不到;开边界应选在地形比较平缓和流态为缓流的地方,不可选在地形和流态较复杂的地方,这样会使计算不稳定;开边界最好选在有流速、潮位资料的断面上;开边界断面应与流线保持正交。但在实际应用中,开边界不一定都能满足上述条件,这要根据实际情况酌情考虑。闭边界即水和陆地的交界,一般以海岸线为闭边界。通常对复杂海岸,局部加以光滑处理。4.4动边界在潮汐河口,常存在着大量的沙洲和边滩,这些滩地大都滩面宽广,坡度平缓,随着潮位的涨落,水边线在滩地上上下移动。如以1‰坡度的边滩为例,对于一般的寻常潮(潮差为3m左右),水边线移动距离可达数公里。这里滩地不仅起蓄水作用,而且会影响主槽的流场。在涨潮时水流上滩,分散了涨潮水流,使涨潮流速减少,涨潮历时缩短,而在落潮时滩水归槽,增加了落潮流速,延长落潮历时。因此,在滩地宽广的河口,这种现象需要有动边界处理。动边界处理方法有很多,通常采用的有:干湿网格法、水边线步进法、窄缝法等。目前,国内外最为广泛采用的动边界处理方法为干湿网格法。所谓干湿网格法即运用一定的干湿点判断规则对网格点进行判断,其计算域在每一个时间步长都是改变的。该方法原理简单,适用性强,可以处理复杂的边界条件,而且在二维与三维模型中都有所应用,但是该方法判断较为烦琐,计算量大。史峰岩等曾经提出了一种岸界适应干湿网格变边界模型,导出了极坐标下的连续移动边界模型,从而改善了对岸界弯曲较大海域的风暴潮漫滩计算,并对黄河三角洲海域风暴潮的漫滩进行了数值模拟。水边线步进法基本原理是:在水位变动过程中,以浅滩水边线最接近的网格节点为边界点。因为水边线是运动的,因此,作为边界的网格节点也是变化的。水边线步进法根据水面涨落的情况决定边界位置,没有考虑边界处的流速,而且边界是活动的,需要经常判定固边界的位置。窄缝法是处理动边界的一种较好的方法,它原理明确,计算简便、灵活。窄缝法的特点是将岸滩前的水域引入到岸滩上,通过引入化引水深的概念使原来没有水的区域也可参与到计算中来,从而可以运用固定的计算边界来进行计算。其最大优点在于它可以不必像其它动边界方法那样对每一个时间步长,每一个节点都进行判断,而是直接将动边界问题转化为固定边界计算。何少苓等将此方法运用到了杭州湾的潮流计算中,取得了较好的效果。王泽良等将此方法运用到了渤海湾中,较好地模拟了具有广阔潮间带海域污染物的迁移扩散。1998年以前,窄缝法的应用还仅限于二维模型中,在三维模型中的应用则始于近几年。江毓武等将窄缝法引入三维数学模型,对厦门前埔围海工程进行潮流数值模拟,取得了良好的效果。许金电等应用窄缝法对变边界进行模拟,由此建立三维σ坐标潮流模型,并建立三维σ坐标污染物扩散模型。此模型应用于厦门海域,应用结果表明,这种三维可变边界污染物扩散数值模型能够较好地对具有较大面积浅海滩涂的海区进行水质模拟和预测。5泥沙输运试验及模型在海岸地区,波浪起着掀沙的作用,是引起泥沙运动的主要动力因素,而泥沙输移则主要靠水流的挟带作用。水流可以是近岸地区的任何一种水流,包括潮流,风吹流,波浪生成的沿岸流、裂流等。在近岸区,潮流一般是相对较弱的,特别是在沙质海岸上,它不足以引起泥沙的起动,但在细颗粒泥沙的海岸上(如粉沙质海岸),单靠潮流的作用也可使泥沙起动。因此,对这种海岸泥沙运动的研究必须考虑波、流的共同作用。由于波浪水流是周期性的往复振动流,其作用下的泥沙运动也是往复运动,运动机理要比河流中单向恒定作用下的泥沙运动机理复杂的多。因此,研究近岸区泥沙的输移及岸滩演变具有十分重要的意义,这也是海岸研究中最困难的课题之一。用数学模型来模拟泥沙的输移是近几十年发展起来的,由于数学模型试验不存在比尺的限制,因此目前在解决海岸、河口工程泥沙方面己成为一种被广泛采用的试验方法。泥沙的输移依赖许多因素,如:波浪、水流(潮流、河川径流)、底部剪切力、紊动强度、水下地形、泥沙类型、底床的侵蚀及泥沙的沉积等,又由于现场观测困难,许多基本问题还远没有解决,至今还没有一套成熟的输沙公式,一般还是套用河流泥沙的研究成果,再加一些适当的修正。水沙数值模拟的成败与否,往往不完全由数学模型和数值模拟本身所决定,有时还经常由边界条件,初始条件以及选用参数等因素的处理是否合适所控制。泥沙理论的研究始于Brahms,他提出了泥沙的起动流速。Duboys首次提出了推移质运动的拖曳力理论,对悬移质的研究则始于Schmidt。此后众多学者对泥沙理论的研究做出了重要贡献。Vanonif主编的《泥沙工程》总结了20世纪70年代前国外学者关于泥沙输运理论和实践的工作,钱宁、万兆惠编撰的《泥沙运动力学》是对20世纪80年代前国内外泥沙理论和试验进展的较完整论述。随着对泥沙运动研究的深入,波浪作用下的泥沙运动也得到了深入研究。Bagnold根据试验得到了泥沙起动公式,Kalkanis模拟明渠水流中爱因斯坦床沙质函数,并导出了一个适用于波浪作用下的推移质输沙率公式。Bagnold从能量的角度导出波浪作用下的全沙输沙率。近年来,波、流作用共同作用下的泥沙运动研究有了较大发展,练继建对波、流作用下的边界层及床面剪切力进行了研究,Quick、窦国仁等给出了波、流共同作用下的输沙率公式。泥沙数学模型的研究,也经历了一个从一维、二维到三维,从非耦合到耦合的过程。VanRijn提出了只适用于渐变流的三维悬沙数学模型并探讨了模型中各参数对结果的影响。目前,二维泥沙模型广泛应用于悬沙及底沙输移和底床的演变研究中。在三维泥沙模拟中,由于参数确定比较困难,三维泥沙模型还很难适用于近岸区大范围的泥沙研究。近年来,商业化的软件包也逐步出现,如美国陆军工程兵团水道实验站的SWM、荷兰Delft水利研究所的TRISULA和DELWAQ。值得指出的是,Falconer曾对沿水深积分的平面二维流动与水质数学模型的细化与完善作了全面的研究,开发了DIVAST。国内在河口与海岸数学模型开发与应用方面也积累了许多经验,开发了一些颇具特色的模型,但在数学模型的集成化程度及通用软件包的研制方面尚落后于国外