CN116738806A - 含丁坝水流交汇区微塑料输运规律的模拟方法 - Google Patents
含丁坝水流交汇区微塑料输运规律的模拟方法 Download PDFInfo
- Publication number
- CN116738806A CN116738806A CN202310410747.3A CN202310410747A CN116738806A CN 116738806 A CN116738806 A CN 116738806A CN 202310410747 A CN202310410747 A CN 202310410747A CN 116738806 A CN116738806 A CN 116738806A
- Authority
- CN
- China
- Prior art keywords
- micro plastic
- particles
- micro
- plastic particles
- flow
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Pending
Links
- 238000000034 method Methods 0.000 title claims abstract description 81
- 238000004088 simulation Methods 0.000 title claims abstract description 28
- 229920003023 plastic Polymers 0.000 title claims abstract description 13
- 239000004033 plastic Substances 0.000 title claims abstract description 13
- 239000002245 particle Substances 0.000 claims abstract description 218
- 229920000426 Microplastic Polymers 0.000 claims abstract description 177
- XLYOFNOQVPJJNP-UHFFFAOYSA-N water Substances O XLYOFNOQVPJJNP-UHFFFAOYSA-N 0.000 claims abstract description 77
- 239000012530 fluid Substances 0.000 claims abstract description 49
- 230000033001 locomotion Effects 0.000 claims abstract description 46
- 238000004364 calculation method Methods 0.000 claims abstract description 43
- 230000008878 coupling Effects 0.000 claims abstract description 21
- 238000010168 coupling process Methods 0.000 claims abstract description 21
- 238000005859 coupling reaction Methods 0.000 claims abstract description 21
- 238000011160 research Methods 0.000 claims abstract description 20
- 230000008569 process Effects 0.000 claims abstract description 12
- 230000000704 physical effect Effects 0.000 claims abstract description 7
- 230000001133 acceleration Effects 0.000 claims description 20
- 125000000484 butyl group Chemical group [H]C([*])([H])C([H])([H])C([H])([H])C([H])([H])[H] 0.000 claims description 18
- 230000000694 effects Effects 0.000 claims description 12
- 238000002347 injection Methods 0.000 claims description 9
- 239000007924 injection Substances 0.000 claims description 9
- 230000005484 gravity Effects 0.000 claims description 8
- 230000010354 integration Effects 0.000 claims description 8
- 230000003993 interaction Effects 0.000 claims description 8
- 230000003068 static effect Effects 0.000 claims description 6
- 238000004422 calculation algorithm Methods 0.000 claims description 4
- 238000011438 discrete method Methods 0.000 claims description 4
- 230000004907 flux Effects 0.000 claims description 4
- 238000004519 manufacturing process Methods 0.000 claims description 4
- 238000013178 mathematical model Methods 0.000 claims description 4
- 230000021715 photosynthesis, light harvesting Effects 0.000 claims description 4
- 230000001052 transient effect Effects 0.000 claims description 4
- 230000005514 two-phase flow Effects 0.000 claims description 4
- 230000005540 biological transmission Effects 0.000 claims description 3
- 239000008187 granular material Substances 0.000 claims description 3
- 239000000178 monomer Substances 0.000 claims description 3
- 239000012798 spherical particle Substances 0.000 claims description 3
- 239000004698 Polyethylene Substances 0.000 description 10
- 229920000573 polyethylene Polymers 0.000 description 10
- 239000000463 material Substances 0.000 description 9
- 239000004952 Polyamide Substances 0.000 description 7
- 238000010586 diagram Methods 0.000 description 7
- 239000003344 environmental pollutant Substances 0.000 description 7
- 229920002647 polyamide Polymers 0.000 description 7
- 239000004793 Polystyrene Substances 0.000 description 6
- 238000005187 foaming Methods 0.000 description 6
- 231100000719 pollutant Toxicity 0.000 description 5
- -1 Polyethylene Polymers 0.000 description 4
- 229920002223 polystyrene Polymers 0.000 description 4
- 238000011084 recovery Methods 0.000 description 3
- 239000004743 Polypropylene Substances 0.000 description 2
- 230000009286 beneficial effect Effects 0.000 description 2
- 230000008859 change Effects 0.000 description 2
- 238000006073 displacement reaction Methods 0.000 description 2
- 229920000139 polyethylene terephthalate Polymers 0.000 description 2
- 239000005020 polyethylene terephthalate Substances 0.000 description 2
- 229920001155 polypropylene Polymers 0.000 description 2
- 238000004062 sedimentation Methods 0.000 description 2
- 239000000126 substance Substances 0.000 description 2
- 239000000725 suspension Substances 0.000 description 2
- 230000009471 action Effects 0.000 description 1
- 238000004458 analytical method Methods 0.000 description 1
- 230000007547 defect Effects 0.000 description 1
- 238000005516 engineering process Methods 0.000 description 1
- 239000000835 fiber Substances 0.000 description 1
- 230000014759 maintenance of location Effects 0.000 description 1
- 230000005012 migration Effects 0.000 description 1
- 238000013508 migration Methods 0.000 description 1
- 239000004417 polycarbonate Substances 0.000 description 1
- 229920000515 polycarbonate Polymers 0.000 description 1
- 238000012545 processing Methods 0.000 description 1
- 230000001737 promoting effect Effects 0.000 description 1
- 238000000926 separation method Methods 0.000 description 1
- 239000007787 solid Substances 0.000 description 1
- 238000012876 topography Methods 0.000 description 1
- 238000012546 transfer Methods 0.000 description 1
- 238000012795 verification Methods 0.000 description 1
Landscapes
- Management, Administration, Business Operations System, And Electronic Commerce (AREA)
Abstract
本发明公开了一种含丁坝水流交汇区微塑料输运规律的模拟方法，首先选择一组微塑料颗粒的成分进行记录；查询微塑料颗粒物理性质，随后建立出研究区域所含不同形状的微塑料颗粒；然后使用有限体积法建立微塑料颗粒相网格和流体域网格，对模型边界和水气交界面以及丁坝附近区域的网格进行加密，基于Hertz‑Mindlin无滑移接触模型，利用离散元DEM方法计算微塑料颗粒运动；利用计算流体力学CFD方法计算水流运动；耦合步骤3和步骤4的计算过程，得出微塑料颗粒随水流运动的数据。本发明解决了现有技术中存在的无法准确探究含丁坝水流交汇区中的微塑料颗粒输运规律的问题。
Description
技术领域
本发明属于水流与颗粒耦合运动模拟技术领域，具体涉及一种含丁坝水流交汇区微塑料输运规律的模拟方法。
背景技术
微塑料被认为是粒径或等效粒径小于或等于5mm的微塑料颗粒[1]，化学成分多样，其成分有聚乙烯(PE)、聚丙烯(PP)，对苯二甲酸聚乙烯(PET)、聚苯乙烯(PS)、聚酰胺(PA)等；形态较为丰富，主要有碎片状、纤维状、发泡状和薄膜状。微塑料颗粒密度的差异，使其在水中运动时伴随有悬浮和沉降现象，在研究时比较难以追踪。内陆河流是水环境中微塑料自诞生之后不断扩散污染的主要载体，而交汇区是河流中常见的结构之一，该区域水体流态复杂、掺混剧烈，河流中携带的微塑料颗粒流经交汇区时，微塑料颗粒的运动状态与平缓河流中的相比有很大差异，同时会伴随产生微塑料污染的差异；水流交汇区另一个特点是交汇水流对下游边壁的冲刷，为了减缓、控制这一现象，常常采取在交汇区布设丁坝的整治措施，将部分高速水流导向对岸，使原来被冲刷岸边的水流流速变小，从而达到保护的目的，但丁坝的存在会进一步扰乱交汇区水流，使水中微塑料颗粒的运动规律更加复杂。
目前关于微塑料颗粒在内陆河流中运移研究的模拟方法，有学者[2]将粒子跟踪模块与流动模块耦合，研究了河流中的微塑料受降雨影响的运动规律；也有学者使用MIKE 21流动模型来模拟计算水流运动，使用颗粒追踪模块确定每个微塑料颗粒的悬浮与沉降。这些研究中存在一些局限：此类方法的研究区域并未考虑支流汇入的影响，未考虑微塑料颗粒会受水动力条件复杂变化的影响从而改变运动规律的这种情况，而现有针对交汇区污染物的研究，有学者[3]采用了三维水动力-污染物耦合输运数学模型的方法对90°等宽明渠水流中的污染物混合规律进行探究，该方法对交汇区污染研究有促进作用，但这种方法研究的污染物属于均匀连续分布于水体中，和微塑料颗粒个体性的不均匀分布有很大差异，不能适用于微塑料颗粒在交汇区输运规律的研究。此外，在交汇区布设丁坝时，使用三维水动力-污染物耦合模型的方法也尚未考虑到丁坝对微塑料输运规律的影响。
综上，对于微塑料颗粒在含丁坝水流交汇区输运规律的研究和方法不够充分，大多针对于单一种类进行研究，现有的方法无法充分结合微塑料的理化性质、水流流态和河道地形条件，难以具体反映含丁坝水流交汇区中微塑料颗粒输运规律。本发明为了解决这些问题，提出了一种研究含丁坝水流交汇区中微塑料颗粒输运规律的模拟方法，真实得考虑微塑料颗粒的理化性质，耦合含丁坝交汇区水流作用特性，模拟研究微塑料颗粒的输运，具有创新性及价值。
[1]Cole M，Galloway T S.Ingestion of nanoplastics and microplastics byPacific oyster larvae[J].Environmental Science&Technology,2015，49(24)：14625—14632.
[2]Lu X,Wang X,Liu X,et al.Dispersal and transport of microplasticparticles under different flow conditions in riverine ecosystem[J].Journal ofHazardous Materials,2023,442:130033[2023-02-21].
[3]毛颂平,朱海,王玲玲,华祖林,陆圣杰,薛晓鹏.90°等宽明渠交汇水流污染物混合规律数值模拟[J].水利水电科技进展,2022,42(03):64-69+96。
发明内容
本发明的目的是提供一种含丁坝水流交汇区微塑料输运规律的模拟方法，解决了现有技术中存在的无法准确探究含丁坝水流交汇区中的微塑料颗粒输运规律的问题。
本发明所采用的技术方案是，含丁坝水流交汇区微塑料输运规律的模拟方法，具体按照以下步骤实施：
步骤1、选择一组微塑料颗粒的成分进行记录；查询微塑料颗粒物理性质，随后建立出研究区域所含不同形状的微塑料颗粒；
步骤2、使用有限体积法建立微塑料颗粒相网格和流体域网格，对模型边界和水气交界面以及丁坝附近区域的网格进行加密，重命名颗粒相和流体域网格的入口平面、边壁平面、顶面和底面；
步骤3、设置微塑料颗粒属性和网格各面属性，基于Hertz-Mindlin无滑移接触模型，利用离散元DEM方法计算微塑料颗粒运动；
步骤4、设置初始条件和边界条件，建立控制方程，基于三维紊动流体数学模型，利用计算流体力学CFD方法计算水流运动；
步骤5、耦合步骤3和步骤4的计算过程，得出微塑料颗粒随水流运动的数据。
本发明的特点还在于，
步骤1使用0.1mm～1mm为半径的单体球状颗粒进行密集拼接组成微塑料颗粒形状，组合体的尺寸为0.1mm～5mm。
步骤2具体按照以下步骤实施：
建立三维河流或三维水槽模型，在建立的模型中确保丁坝在模型中的形状、长度、位置和夹角与实际研究区域保持一致，生成网格时，设置最小单元尺寸，使用有限体积法依次连接模型各边节点完成非结构网格的绘制；进行网格加密时，将模型边界和水气交界面以及丁坝附近区域网格尺寸减小，各边节点之间的距离随之减小，依次连接加密区域各边节点完成加密区域非结构网格的绘制，达到对模型边界和水气交界面以及丁坝附近区域网格加密的效果；
步骤2中网格命名时，将干流和支流入口平面分别命名为inlet1和inlet2，将两河交汇后的出口平面命名为outlet；把接近空气的顶面命名为top；其他侧面和底面分别命名为side和bot。
步骤3中微塑料颗粒属性包括：微塑料颗粒大小、泊松比、弹性模量、密度、剪切模量、相互作用的恢复系数、静摩擦系数、动摩擦系数和几何环境与不同微塑料颗粒之间相互作用的参数，网格各面属性包括：虚体或实体；虚体允许微塑料颗粒通过，实体则不允许微塑料颗粒通过，因此需要在设置好的虚体平面上添加动态颗粒工厂，选择步骤1中建立的微塑料颗粒注入的模式即无限颗粒注入或一定数量颗粒注入，确定每秒从虚体平面射入的微塑料颗粒数量，颗粒入射速度为河流实测的平均流速，在设置时考虑重力影响，方向根据坐标系具体判断，大小设置为9.81m/s2。
步骤3中选择无限动态颗粒按干流和支流实测的平均流速入射。
步骤3具体如下：
在微塑料颗粒计算过程中，每个微塑料颗粒会产生平移运动和旋转运动，以牛顿第二定律为原理计算颗粒的平移和旋转加速度，用经过离散化的时间步长进行数值积分从而在网格之间更新颗粒速度，具体控制方程如下：
旋转运动：
式中，I是转动惯性矩，单位为kg·m2；ω是角速度，单位为rad/s；M是作用在质量上的总扭矩，单位为N·m；t是时间，单位为s；
平移运动：
式中，v是颗粒的平移速度，单位为m/s；m是颗粒的重量，单位为kg；Fg是作用在颗粒上的联合重力，单位为N；Fc和Fnc是颗粒与周围颗粒或壁面之间的接触和非接触力，单位为N；
数值积分：
式中，v(t)是速度，单位为m/s；x(t)是位置，单位为m；a(t)是颗粒在给定步长的加速度，单位为m/s2；△t是时间步长，单位为s；
最后采用离散方法求解微塑料颗粒运动的控制方程，完成对微塑料颗粒的模拟；
步骤4具体如下：
对河流进行模拟时，需要考虑空气对流动的影响，因此使用VOF多相流模型模拟水-气两相流；进行边界条件的设置时，将干流和支流入口平面设置为速度进口，即可设置水流和空气流入的速度；将两河交汇后的出口平面设置为压力出口；把接近空气的顶面设置为压力进口；其他侧面和底面均设置为无滑移的壁面；同时，河道或水槽的糙率也可以根据研究区域的实测数据进行设置；进行初始条件的设置时，将干流和支流入口的流速设置为研究区域水流的实测平均流速，顶面压力进口压强设置为一个大气压，两河交汇后压力出口的压强设置为一个大气压。
水流计算使用三维瞬态方法，基本控制方程有质量守恒方程、动量守恒方程，具体如下：
质量守恒方程：
式(4)中，V表示控制体，控制体为流场单个非结构网格的体积；A表示控制面；ρ为流体密度；v为流体流速；n为控制面数；等式左边第一项表示控制体V内部质量的增量，第二项表示通过控制表面流入控制体的净通量；
动量守恒方程：
式(5)中：Fbx，Fby，Fbz分别是单位质量流体上的质量力在三个方向上的分量；pxy，pyz，pxz等是该流体内应力张量的分量；
为保证赋予的重力加速度方向与实际情况保持一致，重力加速度方向依据网格本身坐标系位置确定，湍流模型使用RNG k-ε：
k方程：
ε方程：
式(6)与式(7)中，Gk表示由于平均速度变化梯度引起的湍动能产生；Gb是由于浮力影响引起的湍动能产生；YM为可压缩湍流脉动膨胀对总耗散率的影响；μ为湍流粘性系数。
以步骤2流体域网格的每个非结构网格为单元列出控制方程组，将流体域所有非结构网格单元的方程组联立，采用有限体积法建立控制方程组的离散代数方程，离散格式使用二阶迎风格式，将步骤4边界条件代入，将流体域全局初始压力、流速、紊动能、紊动能耗散率设置为0，完成数据初始化；使用压力耦合方程组的半隐式方法，即SIMPLE算法，采用计算机进行数值计算，求解离散化的控制方程。
步骤5具体按照以下步骤实施：
耦合计算步骤3和步骤4的过程，在一个时间步长内，将步骤4计算结果中水流对微塑料颗粒的拖曳力和虚拟质量力传递给步骤3中的微塑料颗粒，通过求解微塑料颗粒运动的控制方程，完成一个步长内微塑料颗粒运动的计算，在下一个步长开始时，以上一个步长中的结果作为这一步的已知条件，在多个时间步长内进行迭代计算，从而实现质量、动量、能量从计算流体力学到离散元中的传递，实现了微塑料颗粒状态的更新，在计算完毕时，即得微塑料颗粒模型中的分布位置。
本发明的有益效果是，含丁坝水流交汇区微塑料输运规律的模拟方法，综合考虑含丁坝交汇区中微塑料的输移，着重对水流流态最复杂的区域进行分析，在河道地形条件方面具有更接近真实地形的优越性。本发明可精确设置微塑料颗粒的多种性能参数，其中包括微塑料颗粒的密度、泊松比、弹性模量、剪切模量、恢复系数、静摩擦系数和动摩擦系数，使本发明的微塑料颗粒不同于其他方法所研究的污染物，更加清晰得模拟了微塑料颗粒在水流中的相互作用和运动规律，综合评估了在修建丁坝保护下游边壁的同时对微塑料颗粒运动的影响。
附图说明
图1是本发明方法的流程图；
图2是微塑料颗粒与水流耦合作用的流程图；
图3是本发明中的概化模型图；
图4是颗粒形状示意图；
图5是颗粒相和流体域的网格图；
图6(a)是无丁坝微塑料颗粒速度分布图；
图6(b)是有丁坝微塑料颗粒速度分布图；
图7(a)是无丁坝微塑料颗粒驻留时间分布图；
图7(b)是有丁坝微塑料颗粒驻留时间分布图；
图8(a)是工况3中微塑料颗粒速度分布图；
图8(b)是工况3中微塑料颗粒驻留时间分布图；
图9(a)是PE材料的微塑料颗粒速度分布图；
图9(b)是PS材料的微塑料颗粒速度分布图；
图9(c)是PA材料的微塑料颗粒速度分布图。
具体实施方式
下面结合附图和具体实施方式对本发明进行详细说明。
本发明含丁坝水流交汇区微塑料输运规律的模拟方法，具体按照以下步骤实施：
步骤1、调查内陆河流中微塑料颗粒的成分分布，选择具有代表性的一组微塑料颗粒的成分进行记录；查询微塑料颗粒物理性质，所述物理性质包括密度、泊松比、弹性模量、剪切模量、恢复系数、静摩擦系数和动摩擦系数，建立不同形状的微塑料颗粒；
步骤1应对目标区域的微塑料颗粒赋存情况进行资料查询，得到本区域微塑料的化学成分、物理性质、尺寸形状等数据。选取丰度大、密度分布适中的微塑料颗粒，进行模拟；
在内陆河流中，微塑料颗粒形状丰富，若对微塑料形状有关的规律进行研究，则需建立不同形状的微塑料颗粒，具体实施方法为：依据目标形状的最大尺寸，使用0.1mm～1mm为半径的单体球状颗粒进行密集拼接组成如图4所示的微塑料颗粒形状，组合体的尺寸在0.1mm～5mm范围内。
步骤2、使用有限体积法建立微塑料颗粒相网格和流体域网格，对模型边界和水气交界面以及丁坝附近区域的网格进行加密，重命名颗粒相和流体域网格的入口平面、边壁平面、顶面和底面；
步骤2具体按照以下步骤实施：
建立三维河流或三维水槽模型，在建立的模型中确保丁坝在模型中的形状、长度、位置和夹角与实际研究区域保持一致，生成网格时，设置最小单元尺寸，使用有限体积法依次连接模型各边节点完成非结构网格的绘制；进行网格加密时，将模型边界和水气交界面以及丁坝附近区域网格尺寸减小，各边节点之间的距离随之减小，依次连接加密区域各边节点完成加密区域非结构网格的绘制，达到对模型边界和水气交界面以及丁坝附近区域网格加密的效果。
步骤2中网格命名时，将干流和支流入口平面分别命名为inlet1和inlet2，将两河交汇后的出口平面命名为outlet；把接近空气的顶面命名为top；其他侧面和底面分别命名为side和bot。
步骤3、设置微塑料颗粒属性和网格各面属性，基于Hertz-Mindlin(无滑移)接触模型，利用离散元(DEM)方法计算微塑料颗粒运动；微塑料颗粒属性包括：微塑料颗粒大小、泊松比、弹性模量、密度、剪切模量、相互作用的恢复系数、静摩擦系数、动摩擦系数和几何环境与不同微塑料颗粒之间相互作用的参数，网格各面属性包括实体或虚体；
步骤3中虚体允许微塑料颗粒通过，实体则不允许微塑料颗粒通过，因此需要在设置好的虚体平面上添加动态颗粒工厂，选择步骤1中建立的微塑料颗粒注入的模式即无限颗粒注入或一定数量颗粒注入，确定每秒从虚体平面射入的微塑料颗粒数量，颗粒入射速度为河流实测的平均流速，在设置时考虑重力影响，方向根据坐标系具体判断，大小设置为9.81m/s2。
步骤3微塑料颗粒计算参数在颗粒相网格中的计算方程具体如下：
在微塑料颗粒计算过程中，每个微塑料颗粒会产生平移运动和旋转运动，以牛顿第二定律为原理计算颗粒的平移和旋转加速度，用经过离散化的时间步长进行数值积分从而在网格之间更新颗粒速度，具体控制方程如下：
旋转运动：
式中，I是转动惯性矩，单位为kg·m2；ω是角速度，单位为rad/s；M是作用在质量上的总扭矩，单位为N·m；t是时间，单位为s；
平移运动：
式中，v是颗粒的平移速度，单位为m/s；m是颗粒的重量，单位为kg；Fg是作用在颗粒上的联合重力，单位为N；Fc和Fnc是颗粒与周围颗粒或壁面之间的接触和非接触力，单位为N；
数值积分：
式中，v(t)是速度，单位为m/s；x(t)是位置，单位为m；a(t)是颗粒在给定步长的加速度，单位为m/s2；△t是时间步长，单位为s；
最后采用离散方法求解微塑料颗粒运动的控制方程，完成对微塑料颗粒的模拟；
步骤4、设置初始条件和边界条件，建立控制方程，基于三维紊动流体数学模型，利用计算流体力学(CFD)方法计算水流运动；
步骤4中对河流进行模拟时，需要考虑空气对流动的影响，因此使用VOF多相流模型模拟水-气两相流；进行边界条件的设置时，将干流和支流入口平面设置为速度进口，即可设置水流和空气流入的速度；将两河交汇后的出口平面设置为压力出口；把接近空气的顶面设置为压力进口；其他侧面和底面均设置为无滑移的壁面；同时，河道或水槽的糙率也可以根据研究区域的实测数据进行设置；进行初始条件的设置时，将干流和支流入口的流速设置为研究区域水流的实测平均流速，顶面压力进口压强设置为一个大气压，两河交汇后压力出口的压强设置为一个大气压。
步骤4中水流计算使用三维瞬态方法，基本控制方程有质量守恒方程、动量守恒方程，具体如下：
质量守恒方程：
式(4)中，V表示控制体，控制体为流场单个非结构网格的体积；A表示控制面；ρ为流体密度；v为流体流速；n为控制面数；等式左边第一项表示控制体V内部质量的增量，第二项表示通过控制表面流入控制体的净通量；
动量守恒方程：
式(5)中：Fbx，Fby，Fbz分别是单位质量流体上的质量力在三个方向上的分量；pxy，pyz，pxz等是该流体内应力张量的分量；
为保证赋予的重力加速度方向与实际情况保持一致，重力加速度方向依据网格本身坐标系位置确定，湍流模型使用RNG k-ε：
k方程：
ε方程：
式(6)与式(7)中，Gk表示由于平均速度变化梯度引起的湍动能产生；Gb是由于浮力影响引起的湍动能产生；YM为可压缩湍流脉动膨胀对总耗散率的影响；μ为湍流粘性系数。
步骤4以步骤2流体域网格的每个非结构网格为单元列出控制方程组，将流体域所有非结构网格单元的方程组联立，采用有限体积法建立控制方程组的离散代数方程，离散格式使用二阶迎风格式。将步骤4边界条件代入，将流体域全局初始压力、流速、紊动能、紊动能耗散率设置为0，完成数据初始化；使用压力耦合方程组的半隐式方法，即SIMPLE算法，采用计算机进行数值计算，求解离散化的控制方程。
步骤5、耦合步骤3和步骤4的计算过程，进行耦合计算；
步骤5具体按照以下步骤实施：
图2为水流和微塑料颗粒的耦合计算原理图，通过此耦合原理，耦合计算步骤3和步骤4的过程，在一个时间步长内，将步骤4计算结果中水流对微塑料颗粒的拖曳力和虚拟质量力传递给步骤3中的微塑料颗粒，通过求解微塑料颗粒运动的控制方程，完成一个步长内微塑料颗粒运动的计算，在下一个步长开始时，以上一个步长中的结果作为这一步的已知条件，在多个时间步长内进行迭代计算，从而实现质量、动量、能量从计算流体力学到离散元中的传递，实现了微塑料颗粒状态的更新，在计算完毕时，即得微塑料颗粒模型中的分布位置。
实施例
本发明含丁坝水流交汇区微塑料输运规律的模拟方法，现通过以下实例进行说明：
步骤1，查得目标区域具有代表性的三种微塑料颗粒为聚乙烯(PE)、聚苯乙烯(PS)、聚碳酸酯(PA)，材料物理性质如下表，同时建立不同形状的微塑料颗粒如图4所示，为提高计算效率，本例计算工况选择图4a中的发泡状微塑料颗粒进行验证；
材料属性表
步骤2，根据模拟区域的几何形状，对丁坝所处交汇区的具体位置进行了绘制，建立三维水槽模型，具体形状和尺寸如图3，丁坝截面为矩形，长0.8m，宽0.2m，丁坝布设角度为90°；网格采用非结构网格进行划分，颗粒相和流体域使用相同的网格，网格单元数量为263428，网格密度适中，网格质量较好，在流场变化剧烈的区域即模型边界和水气交界面以及丁坝附近区域的网格进行了加密处理。颗粒相和流体域使用的网格示意图见图5；最后重命名网格各个面，命名干流入口所在的面为inlet1，支流入口所在的面为inlet2，水流出口所在的面为outlet，底面、侧面所在的面分别为bot和side；
步骤3，根据步骤1所得微塑料颗粒成分，对PE、PS和PA的微塑料颗粒进行名称、泊松比、密度、剪切或弹性模量、相互作用的恢复系数、静摩擦系数和动摩擦系数的定义，并导入步骤1建立的发泡状颗粒；为了使微塑料颗粒能够顺利注入，将inlet1和inlet2平面设置为虚体属性，其他不允许水流通过的壁面设置为实体属性；在本实例中，如图3采用动态颗粒工厂在固定的inlet1和inlet2平面上以河流实测平均流速0.5m/s随机注入，重力方向与z轴正方向相反，大小取9.81m/s2；
在微塑料颗粒计算过程中，每个微塑料颗粒会产生平移运动和旋转运动，以牛顿第二定律为原理计算颗粒的平移和旋转加速度，用经过离散化的时间步长进行数值积分从而在网格之间更新颗粒速度，具体控制方程如下：
旋转运动：
式中，I(kg·m2)是转动惯性矩；ω(rad/s)是角速度；M(N·m)是作用在质量上的总扭矩；t(s)是时间；
平移运动：
式中，v(m/s)是颗粒的平移速度；m(kg)是颗粒的重量；Fg(N)是作用在颗粒上的联合重力；Fc(N)和Fnc(N)是颗粒与周围颗粒或壁面之间的接触和非接触力；
数值积分：
式中，v(t)(m/s)是速度；x(t)(m)是位置；a(t)(m/s2)是颗粒在给定步长的加速度；△t(s)是时间步长；
最后采用离散方法求解微塑料颗粒运动的控制方程，完成对微塑料颗粒的模拟；
步骤4，使用VOF多相流模型模拟水-气两相流；进行初始条件和边界条件设置时，将inlet1和inlet2平面设置为速度进口，速度设置为0.5m/s，将outlet平面设置为压力出口，压力大小为一个大气压，把top平面设置为压力进口，压力大小为一个大气压，将bot平面和side平面均设置为无滑移的壁面；同时，将bot平面和side平面糙率设置为0.005m；
水流计算使用三维瞬态方法，基本控制方程有质量守恒方程、动量守恒方程，具体如下：
质量守恒方程：
式(4)中，V表示控制体(控制体为流场单个非结构网格的体积)；A表示控制面；ρ为流体密度；v为流体流速；n为控制面数；等式左边第一项表示控制体V内部质量的增量，第二项表示通过控制表面流入控制体的净通量；
动量守恒方程：
式(5)中：Fbx，Fby，Fbz分别是单位质量流体上的质量力在三个方向上的分量；pxy，pyz，pxz等是该流体内应力张量的分量；
要保证赋予的重力加速度方向与实际情况保持一致，重力加速度方向依据网格本身坐标系位置确定，湍流模型使用RNG k-ε：
k方程：
ε方程：
式(6)与式(7)中，Gk表示由于平均速度变化梯度引起的湍动能产生；Gb是由于浮力影响引起的湍动能产生；YM为可压缩湍流脉动膨胀对总耗散率的影响；μ为湍流粘性系数；
以步骤2流体域网格的每个非结构网格为单元列出控制方程组，将流体域所有非结构网格单元的方程组联立，采用有限体积法建立控制方程组的离散代数方程，离散格式使用二阶迎风格式。将步骤4边界条件代入，将流体域全局初始压力、流速、紊动能、紊动能耗散率设置为0，完成数据初始化；使用压力耦合方程组的半隐式方法，即SIMPLE算法，采用计算机进行数值计算，求解离散化的控制方程。
步骤5，图2为水流和微塑料颗粒的耦合计算原理图，通过此耦合原理，耦合计算步骤3和步骤4的过程，在一个时间步长内，将步骤4计算结果中水流对微塑料颗粒的拖曳力和虚拟质量力传递给步骤3中的微塑料颗粒，通过求解微塑料颗粒运动的控制方程，完成一个步长内微塑料颗粒运动的计算，在下一个步长开始时，以上一个步长中的结果作为这一步的已知条件，在多个时间步长内进行迭代计算，从而实现质量、动量、能量从计算流体力学到离散元中的传递，实现了微塑料颗粒状态的更新，在计算完毕时，即得微塑料颗粒模型中的分布位置。
同等水动力条件下，三种计算工况：
计算工况表
| 微塑料形状 | 尺寸 | 有无丁坝 | 微塑料材料 | |
| 工况1 | 发泡状 | 1mm、2.5mm、5mm | 无 | PE |
| 工况2 | 发泡状 | 1mm、2.5mm、5mm | 有 | PE |
| 工况3 | 发泡状 | 1mm、2.5mm、5mm | 有 | PE、PS、PA |
为了验证本方法的可行性和有益效果，采用上表三种工况对该方法进行说明。在不含丁坝时，图6(a)中微塑料颗粒在最大流速区内产生了较高的速度，分离区微塑料的速度则较小，停滞区微塑料颗粒速度几乎为零；在含丁坝时，微塑料颗粒速度如图6(b)所示，由于丁坝的影响微塑料颗粒在支流侧拥有较大速度，丁坝前靠近左岸的区域有少量微塑料颗粒速度为零；图7(a)、图7(b)两图分别为工况1、2的微塑料颗粒驻留时间图，两图相比，驻留时间在丁坝附近存在明显差异，丁坝前长驻留时间的微塑料颗粒明显多于未设置丁坝的工况1；综合图6(a)、图6(b)、图7(a)和图7(b)可以发现：在相同时间内，有丁坝时微塑料颗粒的水平位移小于无丁坝时的水平位移，显然丁坝阻挡了微塑料颗粒在水流中的输运，造成一些微塑料颗粒聚集于丁坝前，这对微塑料颗粒的水平输运产生了很大的影响；而丁坝附近的水流状态对微塑料的垂直位移也有影响，因此在研究其规律时不能使用无丁坝的研究方法。
为了验证本发明方法也适用于多种微塑料颗粒材料，设置了工况表中的工况3，图8(a)、图8(b)为该工况下的微塑料颗粒速度图和驻留时间图，表明微塑料颗粒在水体中分布更加均匀，且在丁坝对岸侧底部发生了微塑料颗粒的聚集并呈带状分布；为了观察不同微塑料颗粒种类在含丁坝交汇区的输运规律，绘制了图9(a)、图9(b)、图9(c)的微塑料颗粒分布速度图，从微塑料颗粒水平运输的距离来看，PE最远，PS次之，PA最小，有随材料密度增大而减小的趋势；从微塑料颗粒垂直运输来看，不同材料密度的微塑料颗粒其垂直运移也不同。
上述工况结果表明，本发明提出的方法可以模拟接近真实情况下、含丁坝水流交汇区中微塑料颗粒的输运规律，避免了污染物过于笼统、研究区域未考虑交汇影响的弊端，对微塑料输运规律可以进行直观和全面的分析，使模拟结果更加符合实际情况。
Claims (9)
1.含丁坝水流交汇区微塑料输运规律的模拟方法，其特征在于，具体按照以下步骤实施：
步骤1、选择一组微塑料颗粒的成分进行记录；查询微塑料颗粒物理性质，随后建立出研究区域所含不同形状的微塑料颗粒；
步骤2、使用有限体积法建立微塑料颗粒相网格和流体域网格，对模型边界和水气交界面以及丁坝附近区域的网格进行加密，重命名颗粒相和流体域网格的入口平面、边壁平面、顶面和底面；
步骤3、设置微塑料颗粒属性和网格各面属性，基于Hertz-Mindlin无滑移接触模型，利用离散元DEM方法计算微塑料颗粒运动；
步骤4、设置初始条件和边界条件，建立控制方程，基于三维紊动流体数学模型，利用计算流体力学CFD方法计算水流运动；
步骤5、耦合步骤3和步骤4的计算过程，得出微塑料颗粒随水流运动的数据。
2.根据权利要求1所述的含丁坝水流交汇区微塑料输运规律的模拟方法，其特征在于，所述步骤1使用0.1mm～1mm为半径的单体球状颗粒进行密集拼接组成微塑料颗粒形状，组合体的尺寸为0.1mm～5mm。
3.根据权利要求2所述的含丁坝水流交汇区微塑料输运规律的模拟方法，其特征在于，所述步骤2具体按照以下步骤实施：
建立三维河流或三维水槽模型，在建立的模型中确保丁坝在模型中的形状、长度、位置和夹角与实际研究区域保持一致，生成网格时，设置最小单元尺寸，使用有限体积法依次连接模型各边节点完成非结构网格的绘制；进行网格加密时，将模型边界和水气交界面以及丁坝附近区域网格尺寸减小，各边节点之间的距离随之减小，依次连接加密区域各边节点完成加密区域非结构网格的绘制，达到对模型边界和水气交界面以及丁坝附近区域网格加密的效果。
4.根据权利要求3所述的含丁坝水流交汇区微塑料输运规律的模拟方法，其特征在于，所述步骤2中网格命名时，将干流和支流入口平面分别命名为inlet1和inlet2，将两河交汇后的出口平面命名为outlet；把接近空气的顶面命名为top；其他侧面和底面分别命名为side和bot。
5.根据权利要求4所述的含丁坝水流交汇区微塑料输运规律的模拟方法，其特征在于，所述步骤3中微塑料颗粒属性包括：微塑料颗粒大小、泊松比、弹性模量、密度、剪切模量、相互作用的恢复系数、静摩擦系数、动摩擦系数和几何环境与不同微塑料颗粒之间相互作用的参数，网格各面属性包括：虚体和实体；虚体允许微塑料颗粒通过，实体则不允许微塑料颗粒通过，因此需要在设置好的虚体平面上添加动态颗粒工厂，选择步骤1中建立的微塑料颗粒注入的模式即无限颗粒注入或一定数量颗粒注入，确定每秒从虚体平面射入的微塑料颗粒数量，颗粒入射速度为河流实测的平均流速，在设置时考虑重力影响，方向根据坐标系具体判断，大小设置为9.81m/s2。
6.根据权利要求5所述的含丁坝水流交汇区微塑料输运规律的模拟方法，其特征在于，所述步骤3中选择无限动态颗粒按干流和支流实测的平均流速入射。
7.根据权利要求5所述的含丁坝水流交汇区微塑料输运规律的模拟方法，其特征在于，所述步骤3具体如下：
在微塑料颗粒计算过程中，每个微塑料颗粒会产生平移运动和旋转运动，以牛顿第二定律为原理计算颗粒的平移和旋转加速度，用经过离散化的时间步长进行数值积分从而在网格之间更新颗粒速度，具体控制方程如下：
旋转运动：
式中，I是转动惯性矩，单位为kg·m2；ω是角速度，单位为rad/s；M是作用在质量上的总扭矩，单位为N·m；t是时间，单位为s；
平移运动：
式中，v是颗粒的平移速度，单位为m/s；m是颗粒的重量，单位为kg；Fg是作用在颗粒上的联合重力，单位为N；Fc和Fnc是颗粒与周围颗粒或壁面之间的接触和非接触力，单位为N；
数值积分：
式中，v(t)是速度，单位为m/s；x(t)是位置，单位为m；a(t)是颗粒在给定步长的加速度，单位为m/s2；Δt是时间步长，单位为s；
最后采用离散方法求解微塑料颗粒运动的控制方程，完成对微塑料颗粒的模拟。
8.根据权利要求7所述的含丁坝水流交汇区微塑料输运规律的模拟方法，其特征在于，所述步骤4具体如下：
对河流进行模拟时，需要考虑空气对流动的影响，因此使用VOF多相流模型模拟水-气两相流；进行边界条件的设置时，将干流和支流入口平面设置为速度进口，即可设置水流和空气流入的速度；将两河交汇后的出口平面设置为压力出口；把接近空气的顶面设置为压力进口；其他侧面和底面均设置为无滑移的壁面；同时，河道或水槽的糙率也可以根据研究区域的实测数据进行设置；进行初始条件的设置时，将干流和支流入口的流速设置为研究区域水流的实测平均流速，顶面压力进口压强设置为一个大气压，两河交汇后压力出口的压强设置为一个大气压；
水流计算使用三维瞬态方法，基本控制方程有质量守恒方程、动量守恒方程，具体如下：
质量守恒方程：
式(4)中，V表示控制体，控制体为流场单个非结构网格的体积；A表示控制面；ρ为流体密度；v为流体流速；n为控制面数；等式左边第一项表示控制体V内部质量的增量，第二项表示通过控制表面流入控制体的净通量；
动量守恒方程：
式(5)中：Fbx，Fby，Fbz分别是单位质量流体上的质量力在三个方向上的分量；pxy，pyz，pxz等是该流体内应力张量的分量；
为保证赋予的重力加速度方向与实际情况保持一致，重力加速度方向依据网格本身坐标系位置确定，湍流模型使用RNG k-ε：
k方程：
ε方程：
式(6)与式(7)中，Gk表示由于平均速度变化梯度引起的湍动能产生；Gb是由于浮力影响引起的湍动能产生；YM为可压缩湍流脉动膨胀对总耗散率的影响；μ为湍流粘性系数；
以步骤2流体域网格的每个非结构网格为单元列出控制方程组，将流体域所有非结构网格单元的方程组联立，采用有限体积法建立控制方程组的离散代数方程，离散格式使用二阶迎风格式，将步骤4边界条件代入，将流体域全局初始压力、流速、紊动能、紊动能耗散率设置为0，完成数据初始化；使用压力耦合方程组的半隐式方法，即SIMPLE算法，采用计算机进行数值计算，求解离散化的控制方程。
9.根据权利要求8所述的含丁坝水流交汇区微塑料输运规律的模拟方法，其特征在于，所述步骤5具体按照以下步骤实施：
耦合计算步骤3和步骤4的过程，在一个时间步长内，将步骤4计算结果中水流对微塑料颗粒的拖曳力和虚拟质量力传递给步骤3中的微塑料颗粒，通过求解微塑料颗粒运动的控制方程，完成一个步长内微塑料颗粒运动的计算，在下一个步长开始时，以上一个步长中的结果作为这一步的已知条件，在多个时间步长内进行迭代计算，从而实现质量、动量、能量从计算流体力学到离散元中的传递，实现了微塑料颗粒状态的更新，在计算完毕时，即得微塑料颗粒模型中的分布位置。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202310410747.3A CN116738806A (zh) | 2023-04-18 | 2023-04-18 | 含丁坝水流交汇区微塑料输运规律的模拟方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202310410747.3A CN116738806A (zh) | 2023-04-18 | 2023-04-18 | 含丁坝水流交汇区微塑料输运规律的模拟方法 |
Publications (1)
| Publication Number | Publication Date |
|---|---|
| CN116738806A true CN116738806A (zh) | 2023-09-12 |
Family
ID=87905116
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202310410747.3A Pending CN116738806A (zh) | 2023-04-18 | 2023-04-18 | 含丁坝水流交汇区微塑料输运规律的模拟方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN116738806A (zh) |
Cited By (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN117494596A (zh) * | 2023-10-26 | 2024-02-02 | 中国船舶集团有限公司第七一九研究所 | 船舶核动力二回路流体与运行姿态的联合仿真方法及系统 |
Citations (5)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US20070165225A1 (en) * | 2004-03-06 | 2007-07-19 | Michael Trainer | Methods and apparatus for determining the size and shape of particles |
| US20150213163A1 (en) * | 2012-12-20 | 2015-07-30 | Institute Of Modern Physics, Chinese Academy Of Sciences | Particle flow simulation system and method |
| CN109033605A (zh) * | 2018-07-18 | 2018-12-18 | 华中科技大学 | 一种基于多阶段划分和多单位线选择的流域汇流模拟方法 |
| CN115392156A (zh) * | 2022-09-20 | 2022-11-25 | 大连理工大学 | 一种基于dem-cfd耦合计算的滑坡-堵江-涌浪灾害链模拟方法 |
| CN115691709A (zh) * | 2022-09-19 | 2023-02-03 | 西北农林科技大学 | 一种研究水流交汇区微塑料颗粒输运规律的模拟方法 |
-
2023
- 2023-04-18 CN CN202310410747.3A patent/CN116738806A/zh active Pending
Patent Citations (5)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US20070165225A1 (en) * | 2004-03-06 | 2007-07-19 | Michael Trainer | Methods and apparatus for determining the size and shape of particles |
| US20150213163A1 (en) * | 2012-12-20 | 2015-07-30 | Institute Of Modern Physics, Chinese Academy Of Sciences | Particle flow simulation system and method |
| CN109033605A (zh) * | 2018-07-18 | 2018-12-18 | 华中科技大学 | 一种基于多阶段划分和多单位线选择的流域汇流模拟方法 |
| CN115691709A (zh) * | 2022-09-19 | 2023-02-03 | 西北农林科技大学 | 一种研究水流交汇区微塑料颗粒输运规律的模拟方法 |
| CN115392156A (zh) * | 2022-09-20 | 2022-11-25 | 大连理工大学 | 一种基于dem-cfd耦合计算的滑坡-堵江-涌浪灾害链模拟方法 |
Non-Patent Citations (6)
| Title |
|---|
| XIA SHEN ET AL: "Characteristics of secondary flow and separation zone with different junction angle and flow ratio at river confluences", JOURNAL OF HYDROLOGY, pages 1 - 18 * |
| 中国科学院水利部成都山地灾害与环境研究所: "山洪灾害及减灾技术", 31 August 2013, 四川科学技术出版社, pages: 182 - 183 * |
| 袁亚雄 等: "高温高压多相流体动力学理论与应用", 30 April 2016, 北京理工大学出版社, pages: 45 - 46 * |
| 贺博: "丁坝对水流交汇区影响的数值模拟研究", 中国优秀硕士电子期刊网, pages 012 - 108 * |
| 魏文礼 等: "60°弯道丁坝水流水力特性数值模拟研究", 水力发电学报, vol. 36, no. 09, pages 91 - 98 * |
| 黄思: "流体机械数值仿真研究及应用", 31 October 2015, 华南理工大学出版社, pages: 102 - 104 * |
Cited By (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN117494596A (zh) * | 2023-10-26 | 2024-02-02 | 中国船舶集团有限公司第七一九研究所 | 船舶核动力二回路流体与运行姿态的联合仿真方法及系统 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| Neidhold et al. | Interactive physically based fluid and erosion simulation | |
| CN115691709B (zh) | 一种研究水流交汇区微塑料颗粒输运规律的模拟方法 | |
| CN116738806A (zh) | 含丁坝水流交汇区微塑料输运规律的模拟方法 | |
| CN103136436A (zh) | 考虑泥沙吸附的核素输移扩散分析方法 | |
| Zhang et al. | A new capillary force model implemented in lattice Boltzmann method for gas–liquid–solid three-phase flows | |
| Matsumura et al. | Numerical simulation of fluid flow through random packs of polydisperse cylinders | |
| Nemati et al. | Convective heat transfer from two rotating circular cylinders in tandem arrangement by using lattice Boltzmann method | |
| KR20110072551A (ko) | 텐서형 와점성계수를 가진 2차원 하천흐름모형을 이용하여 천수흐름을 해석하는 방법 | |
| Zounemat-Kermani et al. | Coupling of two-and three-dimensional hydrodynamic numerical models for simulating wind-induced currents in deep basins | |
| Chen et al. | Three-dimensional modeling of pollutants transportation in the flow field around a spur dyke | |
| CN108256266A (zh) | 一种一维水动力模型和二维水动力模型耦合方法及系统 | |
| Ramón et al. | Simulation of turbulent flows in river confluences and meandering channels with a Cartesian 3D free surface hydrodynamic model | |
| CN116451535B (zh) | 水流交汇区底泥影响下微塑料输运模拟方法 | |
| Borthwick et al. | Shallow flow modelling using curvilinear depth‐averaged stream function and vorticity transport equations | |
| Chen | A comparison of hydrostatic and nonhydrostatic pressure components in seiche oscillations | |
| Ushijima et al. | Prediction method for local scour by warmed cooling-water jets | |
| Dimitrieva | Stratified flow structure near the horizontal wedge | |
| Üneş | Investigation of density flow in dam reservoirs using a three-dimensional mathematical model including Coriolis effect | |
| Wang et al. | Research on 3D Hydrodynamic Model Based on Adaptive Mesh | |
| Li et al. | 3D numerical simulation of deposition patterns due to sand disposal in flowing water | |
| Goudarzizadeh et al. | Experiment-Based CFD Simulation of Flow Pattern at Open-Channel Junctions | |
| Malik et al. | Mechanism of Scouring Around Group of Bridge Piers in Tandem Arrangement | |
| Gil et al. | Parallel multigrid detached eddy simulation algorithm for three-dimensional unsteady incompressible flows on unstructured meshes | |
| Chaplina et al. | Experimental and Analytical Study of Submerged Jet | |
| Ng et al. | Investigation of the Configuration of Small Hydropower using a Novel Smoothed Particle Hydrodynamics Method |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination |