地下水渗流基本方程及数学模型.pptx

,安徽理工大学 地球与环境学院 水资源与规划系,Ch2,地下水渗流微分方程及数学模型,地 下 水 动 力 学,2-1,渗流连续方程,一、,含水层状态方程,含水层状态方程主要包括地下水状态方程和多孔介质状态方程。,1,、地下水状态方程,Hooke,定律：,式中：,E,体积弹性系数（体积弹性模量），,20,时，,E,=2.110,5,N/cm,2,。其倒数为压缩系数。,等温条件下，水压缩系数,(coef.of compressibility),为,第1页,第1页,积分（,pp,0,，,VV,0,）改写得：,体积：,密度：,按,Taylor,级数展开，得到近似方程：,和,因 （质量守恒），故有,第2页,第2页,2,、多孔介质状态方程,多孔介质压缩系数（,Coefficient of compressibility,）表示多孔介质在压强改变时压缩性指标，用,表示。,多孔介质压缩系数,表示式为：,式中，,Vb,Vs+Vv,多孔介质中所取单元体总体积；,Vs,单元体中固体骨架（,solid matrix,）体积；,Vv,为其中孔隙（,voids,）体积。,介质表面压强；,第3页,第3页,Vv,nVb,；,Vs,(1-n)Vb,式中,多孔介质固体颗粒压缩系数，表示多孔介质中固体颗粒本身压缩性指标，,s,p,；,多孔介质中孔隙压缩系数,（,Compressibility of the pores of a porous medium,），表示多孔介质中孔隙压缩性指标。,n,多孔介质孔隙度。,因 ，故 。,第4页,第4页,3,、贮水率和贮水系数,考虑承压含水层受力情况，取一水平横截面,AB,，按,Terzaghi,（,18831963,）观点：,式中,上覆荷重引起总应力（,total stress,）；,作用在固体颗粒上粒间应力,(intergranular stress),；,横截面面积中颗粒与颗粒接触面积所占水平面积比；,p,水压强。,Terzaghi,令,=,称为有效应力（,effective stress,）。很小，,(1-)p p,，因此有：,第5页,第5页,图,11,一个可压缩承压含水层（,J.Bear,）,第6页,第6页,在水位下降为,H,时，有 。,即作用于固体骨架上力增长了,H,。,作用于骨架上力增长会引起含水层压缩，而水压力减少将造成水膨胀。,含水层本来就充斥了水，骨架压缩和水膨胀都会引起水从含水层中释出，前者就象用手挤压充斥了水海绵会挤出水,样。,第7页,第7页,思考题,1、在涣散层覆盖岩溶地域，建立高楼群时，出现地表塌陷，怎样解释这一现象？,2、一列火车经过一个涣散含水层附近观测孔，该孔水位怎样改变？,3、我国不同地域因开地下水引发地面沉降怎样了解？,第8页,第8页,因,Vs=constant,，故 只在垂直方向上有压缩，,故,上两式表示垂直厚度改变、孔隙度改变与水压强改变关,系。,水头减少时含水层释出水特性，取面积为,1m,2,、厚度为,l m(,即体积为,l m,3,),含水层，考察当水头下降,1m,时释放水量。此时，有效应力增长了,H,g1=,g,。,介质压缩体积减少所释放出水量（,dV,b,）为,与水体积膨胀所释放出水量（,dV,）之和,第9页,第9页,上述两者之和所释放出水量为,或,式中,s,贮水率,释水率,（,specific storativity,），量纲,L-1,，为弹性释水,贮水,；,式中,M,含水层厚度,(m),；,*,贮水系数（,storativity,）。,*,=s,M,贮水系数,*,和贮水率,s,都是表示含水层弹性释水能力参数，在地下水动力学计算中含有主要意义。,第10页,第10页,贮水率,表示含水层水头改变一个单位时，从单位体积含水层中，因水体积膨胀（压缩）以及骨架压缩（或伸长）而释放（或储存）弹性水量。单位,1/L,。,贮水系数,又称,释水系数或储水系数,，为含水层水头改变一个单位时，从底面积为一个单位，高度等于含水层厚度柱体中所释放（或贮存）水量；指面积为一个单位、厚度为含水层全厚度,M,含水层柱体中，当水头改变一个单位时弹性释放或贮存水量，无量纲。既适合用于承压含水层，也适合用于潜水含水层。,贮水率,是描述地下水三维非稳定流或剖面二维流中水文地质参数，既适合用于承压水也适合用于潜水。对于平面二维非稳定流地下水运动，当研究整个含水层厚度上释水情况时，用贮水系数来表达。,第11页,第11页,第12页,第12页,第13页,第13页,*范围值：n10-3 n10-5；范围值：0.05 0.30。实际测出值往往小于理论值。,上述两参数之间不同，还在于潜水含水层存在滞后疏干现象。,弹性释水与重力给水:对于含水层而言，因为受埋藏条件限制，抽水时，水给出存在着不同。,潜水含水层在抽水过程中，大部分水在重力作用下排出，疏干作用于水位变动带（饱水带）和包气带两部分，因为包气带存在，使得饱水带中水释放存在延滞和滞后现象。,当水头下降时，可引发二部分水排出。在上部潜水面下降部位引发重力排水，用给水度表示重力排水能力；在下部饱水部分则引发弹性释水，用贮水率*表示这一部分释水能力。,必须区分二者之间不同，潜水含水层还存在滞后疏干现象。,第14页,第14页,承压含水层抽水时，水释放是由于压力减少造成，这一过程是瞬时完毕。只要水头下降不低到隔水顶板下列，水头减少只引起含水层弹性释水，可用贮水系数,*,表示这种释水能力。,4,、导压系数,描述含水层水头改变传导速度参数，其数值等于含水层导水系数与贮水系数之比或渗入系数与贮水率之比。,第15页,第15页,二、渗流连续方程,连续性方程就是质量守恒方程，也称为水均衡方程,水均衡基本思想：,对某一研究对象，流入 流出,V,研究对象能够是大区域，也能够是微分单元体。,大区域水均衡计算经惯用于区域水资源评价。,本课程基于微分单元体做水均衡，推导渗流连续性方程。,为反应含水层地下水运动普遍规律，我们选定在各向异性多孔介质中建立地下三维不稳定流动连续性方程。,第16页,第16页,因为渗流场中各点渗流速度大小、方向都不同，为了反应液体运动质量守恒关系，需要在三维空间中建立微分方程形式表示连续性方程。,在渗流场中任意取一点P(x,y,z)，以P为中心沿直角坐标轴取一微小六面体，体积为 ，称为特性单元体，设单元体无限小，但确保单元体穿过介质骨架和空隙。,设,v,x,，,v,y,，,v,z,分别为该点在,X,、,Y,、,Z,方向上渗流速度。,Abcd,面中点 。,沿,X,轴方向流入：,第17页,第17页,流出：,利用,Taylor,级数展开，略去二阶导数以上高次项，有：,单元体本身水质量在,t,时间内改变量,为液体密度。,由质量守恒定律,得到,渗流连续性方程,:,=,同理,第18页,第18页,或,上式即为非稳定流渗流连续方程，表明渗流场中任意体,积含水层流入、流出该体积含水层中水质量之差永远恒等于,该体积中水质量改变量。它表示了渗流区内任何一个“局,部”所必须满足质量守恒定律。,若把含水层看作刚体，,=constant,，,n,不变，即水和介质,没有弹性变形或渗流为稳定流，则渗流连续性方程为,第19页,第19页,上式表明，在同一时间内流入单元体水体积等于流出水体积，即,体积守恒。,连续性方程,是研究地下水运动基本方程，各种研究地下水运动微分方程都是依据连续性方程和反应质量守恒定律方程建立起来。,第20页,第20页,2-2,渗流基本微分方程,一、,承压水运动基本微分方程,基本假设：,（,1,）单元体体积无限小，为承压含水层；,（,2,）含水层侧向受到限制，,x,、,y,为常量,z,为变量,存在垂向压缩，水密度,、孔隙度,n,和随压力,p,而改变；,（,3,）由,引起改变,远小于单元体内液体质量改变量（含,）,，可忽略不计；,（,4,）水流服从,Darcys Law,；,第21页,第21页,（,5,）,K,不因 改变而改变；,（,6,）,s,和,K,也不受,n,改变（由于骨架变形）影响。,流体质量：,由于含水层侧向受到限制，可假设,x,、,y,为常量，只考虑垂向压缩。于是，只有水密度,.,孔隙度,n,和单元体高度,z,三个量随压力而改变，于是有：,由含水层状态方程，,=,=,第22页,第22页,由于 因此有 ，,Z,为定值，则,则可得到：,于是,连续性方程,变为：,又,则,第23页,第23页,令 则,依据连续性原理有：,则有：,即：,将 代入整理得：,第24页,第24页,因此有,上式为,三维流微分方程,，也可写成：,物理意义：,渗流空间内任一单位体积含水层在单位时间内流入与流出该体积含水层中弹性水量改变量，即单位体积含水层水量均衡方程。,第25页,第25页,基本微分方程,(Basic Differential Equation),是研究承压含水层中地下水运动基础。它反应了承压含水层中地下水运动质量守恒关系，表明单位时间内流入、流出单位体积含水层水量差（左端）等于同一时间内单位体积含水层弹件释放,(,或弹性贮存,),水量（右端）。它还通过应用,Darcy,定律反应了地下水运动中能量守恒与转化关系。可见，基本微分方程表示了渗流区中任何一个“局部”都必须满足质量守恒和能量守恒定律。,数学意义：,表示渗流空间内任一点任一时刻渗流规律。,在柱坐标系中，基本微分方程为：,第26页,第26页,或,由地下水流基本微分方程，在均质各向同性介质中，方程简化为：,对于各向异性介质，若把坐标轴方向和各向异性介质主方向定为一致，则有,在二维流情况下，基本微分方程可表示为：,第27页,第27页,上式即为,承压水平面二维流微分方程,，该方程是研究承压水含水层中地下水运动基础，反应了承压水含水层中地下水运动质量守恒关系，表明单位时间流入、流出单位体积含水层水量差等于同一时间内单位体积含水层弹性释放（或贮存）水量。,在实际渗流问题中若存在抽、注水及越流影响，只要在微分方程中左端中通过加、减,W,项，通常把该项称为源汇项。所谓源项表示在垂直方向上有水流入含水层，此时,W,为正；汇指在垂直方向上有水流出含水层，此时,W,为负。,此时微分方程变成：,第28页,第28页,二维流情况下：,在二维流情况下，令压力传导系数（导压系数），则,均质各向同性含水层,基本微分方程为：,非均质各向同性含水层中稳定流运动：,均质各向同性含水层中稳定流运动：,上式也称,Laplace,方程。稳定运动方程右端都等于零，意味着同一时间内流入单元体水量等于流出水量。这个结论不但适合用于承压含水层，也适合用于潜水含水层和越流含水层。,第29页,第29页,二、潜水运动基本微分方程,1,、,Dupuit,假设,在潜水面上任意取一点,P,，有：,图,1-16 Dupuit,假设,第30页,第30页,该点流速,v,方向与潜水面相切，则由达西定律有：,vs=-KJ=-Ksin,。,当,很小时，,tg=sin,。此时，,(1),潜水面比较平缓，等水头面呈铅直，水流基本水平，可忽略渗流速度垂直分量,v,Z,；,图,tg,q,、,sin,q,与角度,q,关系,第31页,第31页,（,2,）隔水底板水平，铅垂剖面上各点水头都相等，各点水力坡度和渗流速度都相等，能够近似地用代替，此即著名,Dupuit,假设。,渗流速度：，,H=H(x),通过宽度,B,铅直平面流量为,，,H=H(x),式中,Qxx,方向流量；,h,潜水含水层厚度；,h=H,（隔水层水平时）。,第32页,第32页,对于更普通情况，,H=H(x,，,y),有：,则得：,由于,Dupuit,假设引入，将垂直方向水流速度忽略，减少了,z,变量，简化了计算，但会产生一定误差，经验证实当 时，产生误差很小，误差表示式为：,第33页,第33页,Dupuit,假设无效地域：,(1),存在入渗潜水分水岭地段；,(2),渗出面附近。渗出面,(seepage surface),是在下游边界面上，潜水面下列、下游水面以上地段。渗出面上潜水面往往和边界面相切，有较大垂向分速度。,(3),垂直隔水边界附近。,图,1-17 Dupuit,假设无效地域,第34页,第34页,2,）,Boussinesq,方程,依据,Dupuit,假设，可建立相关潜水含水层中地下水流方程。,图,1-18,潜水非稳定运动,第35页,第35页,（,1,）潜水一维流方程（沿,x,方向运动）,在,t,时间内，上、下游流入、流出单元体水量差为：,在该段时间内，垂直方向补给量为,0,，故,t,时段水量总改变量为,由于水量改变引起潜水面升降，设其改变速率为 ，则在,t,时段，由于潜水面改变而引起小土体内水体积增量为 则有,第36页,第36页,将 （均质各向同性）代入上式，能够得到有入渗补给潜水含水层中地下水非稳定运动一维流方程，又称为,Boussinesq,方程：,式中,K,、,潜水含水层渗入系数、给水度；,W,含水层单位时间、单位面积上垂向补排量，补给为正，排泄为负。,第37页,第37页,(2)潜水二维流方程,均质各向同性含水层，Boussinesq方程为：,式中h=H-Z，Z=0时，h=H。,非均质含水层，Boussinesq方程为：,在推导潜水基本微分方程时应用了Dupuit假设，忽略了弹性储存，所选单元体是一个包含了整个含水层厚度在内土柱，这与承压水非稳定运动时选取无限小单元体不同。因此，应用潜水运动基本方程得到H(x，y，t)只能代表该点整个含水层厚度上平均水头近似值，不能用来计算同一垂直剖面上不同点水头改变。,第38页,第38页,(3),潜水三维流方程,若不用,Dupuit,假设，,Boussinesq,方程普通形式：,在上面潜水基本运动微分方程中右端项为,贮水率,而不是,给水度,，其原因在于，当不考虑,Dupuit,假设时，单元体位于渗流区内部，其贮存量改变只能是弹性释水而不是疏干排水，因此推导出潜水非稳定运动方程和承压水非稳定运动方程形式同样。在这种情况下，地下水非稳定运动特点由边界条件来反应。,第39页,第39页,对于,各向异性介质,，坐标轴方向同主方向，有：,假设固体骨架是不可压缩，同时假设忽略水压缩性,即,常数，有：,或：,第40页,第40页,（,4,）潜水稳定运动微分方程,没有入渗和蒸发时，潜水稳定运动方程式为：,非均质,或,均质,（,5,）地下水运动基本微分方程统一形式：,Z,含水层底板标高。,第41页,第41页,三、越流含水层中地下水运动基本微分方程,在自然界中，存在下列情况，承压含水层上、下岩层并不是绝对隔水，其中一个或两个也许是,弱透水层,(Aquitard),。,即使含水层会通过弱透水层和相邻含水层发生水力联系，但它还是处于承压状态，将其称为,半承压含水层,(Semi-confined aquifer),。,当该含水层和相邻含水层间存在水头差时，地下水就会从高水头含水层通过弱透水层流向低水头含水层。这种现象称为,越流,(Leakage),。,半承压含水层称为,越流含水层,(Leakage aquifer),。,第42页,第42页,假设：,主含水层渗入系数,K,远远不小于若透水层渗入系数,K1,；主含水层弹性释放水量、弱透水层越流量远远不小于弱透水层弹性释放水量。主含水层中水流近似地看作二维流问题，,对于均衡单元体，依据水均衡原理能够写出下列形式连续性方程,(continuty equation),：,第43页,第43页,第44页,第44页,式中，,v,1,，,v,2,分别为通过上部和下部弱透水层垂直越流速率或越流强度，即,其中，,H,1,(x,y,t),和,H,2,(x,y,t),分别为上含水层和下含水层中水头，如,T,表示主含水层导水系数，则得到不考虑弱透水层弹性释水条件下非均质各向同性越流含水层中非稳定运动基本微分方程：,对于均质各向同性介质来说，有：,第45页,第45页,式中,分别称为上、下两个弱透水层越流原因。,越流原因,B(Leakage Factor),量纲为,L,。弱透水层渗入性愈小，厚度越大，则,B,越大，越流量越小。在自然界中，越流原因值改变很大，能够从几米到若干公里。对于一个完全隔水覆盖层来说，,B,为无穷大。,另一个反应越流能力参数是越流系数,s,。其定义为：当主含水层和供应越流含水层间水头差为一个长度单位时，通过主含水层和弱透水层间单位面积界面上水流量。,因此，,K,1,、,m,1,分别为弱透水层渗入系数和厚度。,s,越大，相同水头差下越流量越多。,第46页,第46页,2-3 数学模型建立及求解,一、数学模型相关概念,同一形式偏微分方程代表了整个一大类地下水流运动规律，而对于不同边界性质、不同边界形状含水层，水头分布是不同。而且对于偏微分方程而言，方程本身并不包含反应特定渗流区条件全部信息，方程可能存在无数个解，如需要从大量可能解中求得与特定区域条件相对应唯一特解，就必须提供反应特定区域特性信息。,这些信息包含：,（1）微分方程中相关参数，当这些参数确定后，微分方程才干被确定下来。,第47页,第47页,（,2,）渗流区范围和形状，,当微分方程所相应区域被拟定之后才干对方程求解。,（,3,）边界条件,(boundary conditions),：,表示渗流区边界所处条件，用以表示水头,H,（或渗流量,q,）在渗流区边界上所应满足条件，也就是渗流区内水流与其周围环境互相制约关系。,（,4,）初始条件,(initial conditions),：,表示渗流区初始状态，某一选定初始时刻,(t=0),渗流区内水头,H,分布情况。,将边界条件和初始条件并称为,定解条件,(definite solution condition),，微分方程和定解条件一起构成渗流场数学模型。,第48页,第48页,1,、数学模型：,描述某一研究区地下水流运动数学方程与其定解条件共同构成表示某一实际问题数学结构。,亦即从物理模型出发，用简练数学语言，即一组数学关系式来刻画它数量关系和空间形式，从而反应所研究地质体地质、水文地质条件和地下水运动基本特性，达到复制或再现一个实际水流系统基本状态目的一个数学结构。,其中微分方程表示地下水流动规律，定解条件表明研究对象所处特定环境条件，即所研究地下水流真实状态。,定解问题,是给定了方程（或方程组）和相应定解条件数学物理问题。,建立模型,是指建立数学模型过程。,第49页,第49页,2,、定解条件,1,）定解条件,定解条件指水头、流量等渗流运动要素在流场边界上已知改变规律，这种改变规律是由流场外部条件引起，但它不断地影响流场内部渗流过程并在整个期间始终起作用。,定解条件包括边界条件和初始条件。,2,）边界条件,边界条件是渗流区边界所处条件，用以表示水头,H,（或渗流量,q,）在渗流区边界上所应满足条件，也就是渗流区内水流与其周围环境互相制约关系。,第50页,第50页,(1),第一类边界条件,(Dirichlet,条件,),：,假如在某一部分边界,(,设为,Sl,或,1),上，各点在每一时刻水头都是已知，则这部分边界就称为第一类边界或给定水头边界，表示为,:,或,给定水头边界不一定就是定水头边界。,能够作为第一类边界条件,来处理情况：,河流或湖泊切割含水层，,两者有直接水力联系时，这部分边界就能够作为第一类边界处理。此时，水头,是一个由河湖水位统计资料得到关于,t,函数。但要注意，一些河、湖底部及两侧沉积有一些粉砂、亚粘土和粘土，使地下水和地表水直接水力联系受阻，就不能作为第一类边界条件来处理。,第51页,第51页,在自然界，这种情况很少见。就是附近有河流、湖泊，也不一定能处理为定水头边界，还要视河流、湖泊与地下水水力联系情况，以及这些地表水体本身径流特性而定。在没有充足依据情况下，不要随意把某段边界拟定为定水头边界，以免造成很大误差。,区域内部抽水井、注水井或疏干巷道也能够作为给定水头内边界来处理。此时，水头通常是按某种要求事先给定，比如给定抽水井允许降深等。上面简介都只是给定水头边界。注意，给定水头边界不一定是定水头边界。,排泄地下水溢出带、冲沟或排水渠边界也可近似看作给定水头边界。,第52页,第52页,(2),第二类边界条件,(Neumam,条件,),：,当知道某一部分边界,(,设为,S,2,或,2,),单位面积,(,二维空间为单位宽度,),上流入,(,流出时用负值,),流量,q,时，称为第二类边界或给定流量边界。相应边界条件表示为,:,或,式中，,n,为边界,S,2,或,2,外法线方向。,q,1,和,q,2,则为已知函数，分别表示,S,2,上单位面积和,2,上单位宽度侧向补给量。常见这类边界条件：,隔水边界（流线、分水岭）：,第53页,第53页,抽水井或注水井：,补给或排泄地下水河渠边界上，如已知补给量。,（,3,）第三类边界条件：,某边界上,H,和 线性组合是已知，即有：,又称混合边界条件，,，,为已知函数。,边界为弱透水层（渗入系数为,K,1,，厚度或宽度为,m,1,），,第54页,第54页,第55页,第55页,在,s,3,上,在 上,浸润曲线边界条件：,当浸润曲线下降时，从浸润曲线边界流入渗流区单位面积流量,q,为：,式中，,为给水度，,为浸润曲线外法线与铅垂线间夹角。,第56页,第56页,3),初始条件,初始条件,:,某一选定初始时刻,(t=0),渗流区内水头,H,分布情况。,或,其中，,H,0,为,D,上已知函数。,第57页,第57页,3,、渗流数学模型解法,(1),解析法（,analytic method,）,用参数分析及积分变换等办法直接求解数学模型解办法。其解为准确解，使用简朴，但该办法存在一定不足，只适合用于含水层几何形状规则、方程式简朴、边界条件单一情况。,解析解（,Analytic solution,）又称准确解,，是用解析办法求解数学问题所得到解析表示式。,(2),数值法,(numerical method),用数值办法（离散化办法）求解数学模型办法，其解为近似解。该办法是求解大型地下水流问题主要办法。它把整个渗流区别割成若干个形状规则小单元，每个小单元近似处理成均质，然后建立每个单元地下水流动关系式。把形状不规则、非均质问题转化为形状规则均责问题。依据研究需要，拟定单元划分数量，对于非稳定流还要对时段进行划分。最后，把局部整合起来，加上定解条件。,数值解,(Numerical solution),是 用,数值办法求得数值解，是一个近似解。,(3),模拟法：,利用物理现象与水流相同性，在试验室内采用模拟办法求解。,第58页,第58页,1-8,如何建立,数学模型,水文地质概念模型,描述地下水运动方程,边界条件,初始条件,数学模型,第59页,第59页,例,-1,概念模型,数学模型,潜水含水层数学模型,第60页,第60页,下图为某供水水源地平面图和水文地质剖而图，图中,BC,、,ED,为流线，,BE,为天然隔水边界，,CD,为河流。场地有两个含水层，上部含水层是潜水含水层，下部含水层为承压含水层，中间为弱透水层。己知潜水开采强度为,1,，承压水开采强度为,2,，大气降水补给强度为,W,，承压含水层通过弱透水层补给上部潜水，部分参数见下图,，,试依据图示写出开采过程中地下水非稳定流（二维流）数学模型。,例,-2,越流下潜水,承压含水层数学模型,第61页,第61页,第62页,第62页,潜水含水层,承压含水层,第63页,第63页,下图为某供水水源地平面图和水文地质剖而图，图中,BC,、,DE,垂直等水位线，,CD,为湖边界，,BE,为河流，场地有两个含水层，上部含水层为潜水，下部含水层为承压水，中间隔水层不透水。己知潜水开采强度为,1,，承压水开采强度为,2,，大气降水补给强度为,W,，试依据图示写出开采过程中地下水非稳定流（二维流）数学模型。,独立潜水,承压含水层数学模型,例,-3,第64页,第64页,第65页,第65页,承压含水层,潜水含水层,第66页,第66页,多含水层系统数学模型,例,-4,第67页,第67页,第68页,第68页,第69页,第69页,第70页,第70页,数学模型,微分方程,定解条件,2,数 学 模 型,第71页,第71页,