用粒子有限元方法解决.docx

用粒子有限元方法解决流体-土壤-结构相互作用问题的可能性摘要我们提出一些粒子有限元方法（PFEM）的发展以针对在涉及到流体-土壤-结构相互作用（FFSI）的复杂耦合问题的分析。该PFEM采用修正后的拉格朗日描述来模拟节点（粒子）同时在流体与土壤域（后者包括土壤/岩石和结构）中的运动。联结这些粒子（节点）的网定义了离散域，其中主要的方程的每个组成元素都能用标准有限元分析解决。处理不可压缩连续体的稳定性将通过有限微分（FIC）法介绍。增量迭代计划作为解决非线性瞬态耦合的FSSI问题的方法已经介绍过了。这种用来模拟摩擦接触条件和在流体-土壤及土壤-土壤界面的材料腐蚀的程序已经介绍过了。我们提出一些通过PPEM法来解决FSSI问题的一些例子，比如岩石受到水流作用的运动，临近河床的桥梁基础受到的腐蚀，防洪堤和海上结构的稳定性以及山体滑坡的研究。1 介绍对于包括流体、土壤/岩石和结构之间的相互作用的问题的分析在许多工程领域上是相关的。在对山体滑坡及其对水库和相邻结构的影响，大波浪对海上及港口结构，受洪水和海啸侵袭的结构，在溢流情况下填石坝的土壤侵蚀和稳定性等问题的研究上的例子具有广泛性。这些研究可以看作是所谓的流体-结构相互作用（FSI）问题的一种扩展46。在自由表面流使用FEM的欧拉或者ALE模拟的FSI分析中的典型困难包括对流项和在流体方程中不可压缩性约束的修正，建模以及对流体中自由表明的追踪，在流体和移动的土壤域间通过表明接触的信息的转移，对波浪飞溅的建模，在流体域中处理多实体的大运动的可能性，对于结构和流体的有限元网格的有效更新，等。通过ALE和时空FEM来对FSI问题进行3D分析的范例报告在以下4,6,26,27,31,34,40。如果用朗格朗日描述来模拟土壤和流体域的控制方程，上述大部分问题就得到了解决。在拉格朗日模拟中，单个粒子的运动可以被追踪，所以，在有限元网格中的节点可以看作移动的材料点（这里称作“粒子”）。因此，在网格离散总域（包括流体和土壤部分）中的运动在瞬态解决方案中是可以遵循的。一个针对FSI分析的成功的朗格朗日法就是所谓的Soboran Grid CIP技术，其已被成功地应用在解决不同种类的3D问题中44。研究人员已经在之前工作中成功地开发了特定种类的拉格朗日配方以解决涉及到（自由表明流体）与土壤间的复杂相互作用的问题。该方法被称为有限元法（PFEM，朗格朗日配方的一个优点是，对流项从流体方程中消失了11,48。然而，难题则转移到了如何充分（并且有效）地移动网格节点的问题上。我们用混合了不同形状的元素一个网格再生过程，这些形状使用一种特殊形状函数的扩展Delaunay tessellation17,19。该PFEM的理论和应用报告在如下文献中2,7,10,18,20,21,23,26,32,3439。（不可压缩的）流体流动问题的有限元解决方案意味着解决了动力和不可压缩性方程。这并不是一个简单的问题，因为不可压缩条件限制了FE近似值对于速度和压力的选择于众所周知的div-稳定条件11,48。在我们的作品中，我们使用基于有限元微积分（FIC）的稳定混合有限元方法，其允许对速度和压力变量取线性近似值15,2931,33,34。在其他的有类似特征的稳定有限元方法我们提到有PSPG法41，多尺度法3,6,8和CBS法8,48。本文的目的是介绍目前针对流体-土壤-结构相互作用（FSSI）问题的PFEM法的最新进展。这些问题在土木、水利、海洋与环境工程及其他许多领域中都是相关的。它表明，PFEM法提供了一个通用的分析方法，通过简单并且有效的方式来修正这些复杂的问题。本文的布局如下。在接下来的部分对PFEM的主要思想进行了概况。接下来，使用拉格朗日描述和有限元微分配方的对于可压缩/不可压缩连续体的基本方程通过示意图给出。然后将简要介绍针对瞬态解决方案的一种算法。对耦合FSSI问题及网格生成和自由面节点的识别方法的修正进行概述。对于修正流体、土壤和结构表面间的摩擦接触的程序已进行过说明。我们提出几个应用PFEM法来解决FSSI问题的例子，例如岩石受水流作用的运动，临近河床的桥梁基础的侵蚀，防护堤和受海浪结构的稳定性以及山体滑坡的研究。2 粒子有限元的基本方法让我们考虑一个既含有流体子域也含有固体子域（固体子域可能包括土壤/岩石材料和/或者结构元件）的域。移动的琉璃粒子与固体边界相互作用，从而导致固体的变形，其反过来影响流动运动，因此，该问题是完全耦合的。在PFEM中流体域与固体域都是通过更新的拉格朗日配方进行建模47。也就是说，所有的变量都假设为在t时刻的当前配置。在两个域中的变量的新设定在下一个或在t+t时刻更新的配置中寻找。有限元法（FEM）用于解决针对每个子域的连续介质力学的方程。因此，必须生成一个离散这些域的网格以解决在标准FEM领域中针对各个子域的控制方程问题。数值解的质量取决于标准FEM法中的离散值选择。自适应网格细化技术可用于改善立体或结构发生大运动区域中的解决方案。2.1 PFEM的基本步骤为清楚起见，我们将定义集合或节点云（C）为流体域及固体域，用体积（V）来定义流体和固体的分析域，并用网格（M）定义这两个域的离散。一个PFEM的的典型解决方案包括下列步骤。1. 在每个时间点的起始点是在流体域和固体域的点云。例如nC表示在时刻t=tn（图1）的云计算。2. 用流体域和固体域的确定界限定义流体与固体中的分析域为nV。这是重要的一步，因为一些边界（例如流体中的自由表面）可能会在解决过程中（包括分离和重新进入节点）被严重扭曲。阿尔法形状法12用于边界定义。3. 通过有限元网格nM离散流体域和固体域。在我们的作品中，我们使用的是基于扩展Delaunay tesselation的创新网格生成方案17,19,20。4. 解决流体域和固体域中耦合的运动的拉格朗日方程。计算在这两个结构域中在下一个（更新的）针对t + t的组合的状态变量：速度、压力以及在所述流体和位移下的粘性应力，还有固体中的应力和应变。5. 将网格节点移动到一个新的位置+1C，就时间增量大小而言，在该处n+1表示时刻tn+t。这个步骤是步骤4解决方案过程的典型结果。6. 回到步骤1，对下一个时间步长重复求解过程 得到n+2C (Figure 2)。3 对拉格朗日连续体进行有限元微分/有限元计算3.1 控制方程要解决的方程是连续介质力学中的标准方程，写在拉格朗日参照系中：动量压力-速度关系在上述公式中vi是沿着全球（笛卡尔）轴的第i个速度，p是压力（假定为可压缩） 和K分别代表密度和材料的体积模量，特别地，bi和ij是body forces和（柯西）应力。公式（1）和（2）完全符合本构constitutive关系：不可压缩连续体可压缩/准不可压缩连续体其中ij为应力张量的分量其中S是第二Piola-Kirchhoff应力张量，F是变形梯度张量并且J=detF22,47。参数 和 采用以下针对流体或固体材料的取值：其中v是泊松比，G为剪切模量并且t为时间增量。在方程（3）和（4）中，sij是偏应力，ij是变形率，是粘度并且ij是Kronecker delta。t()表示t时刻的值。在方程（1）（4）中指数的范围从ij=1.nd，其中nd是问题中空间维度的数量（即nd=2就是二维问题）。这些方程完全符合机械问题规定的速度和表面牵引力的标准边界条件11,36,47,48。3.2 方程离散在方程（1）（4）的数值解的一个关键问题就是对不可压缩情况下质量平衡条件的满足（即在方程（2）中K = ）。若干用来解决这些问题的程序已存在于有限元文献中11,48。在我们的方法中我们使用基于所谓的有限积分程序的配方法15,2931,33,34。这种方法的精华是一种改进的质量平衡方程的解决方案，表述如下其中q为权重函数，T是一个由34给定的稳定化参数。在以上的描述中，h是每个有限元的特征长度，|v|是速度矢量的模载体。在方程（5）中i是辅助压力头像变量选择以确保方程（5）中的第二项可以被解释为动量方程的残差的加权总和，因此可以在精确解中消去。控制方程的设置完全是通过加入下列约束方程32,36其中wi是任意的加权函数。积分方程的其余部分通过应用标准的加权残值计算到控制方程（1），（2），（3）和（5）中与其相应的边界条件获得的11,22,48。我们插入介绍下一个标准有限元领域的变量设置问题。在三维问题中，存在三个速度vi，压力p，温度T以及三个压力梯度投影i。在我们的作品中，我们对所有3节点三角网格（二维）和4节点四面体网格（三维）情况下所有变量使用同阶线性插值法。由此产生的离散方程的设置使用的是标准Galerkin技术，如下表所示在公式（8）（10）中， ()表示节点的变量，() =。不同的矩阵和向量由以下给出22,34,36。在使用公式（8）（10）时的解法可以使用任意时刻的积分方案的典型更新拉格朗日有限元法来进行36,47。一中基本算法在以下第二节的概念过程中进行描述，呈现在Box 1中。4 新的网格生成PFEM成功的一个关键点是在根据在该空间域中的节点位置每一个时间步长的网格快速重生。事实上，任何快速网格算法都可用于该目的。在我们的作品中，在每一时间步长的网格生成都使用了所谓的扩展Delaunay tesselation（EDT），见文献17,19。CPU时间要求网格随着节点数量线性增长。CPU时间用于解决方程超过所需的网格节点增长数量。这种情况已经在所有用PFEM解决的问题中都有出现。作为对每一步消耗总CPU大约15%的大型三维问题网格一般规则，而方程的总结果（通常3次迭代至达到在一个时间步长内的收敛）和该系统组装消耗接近70%并且每个时间步长分别消耗CPU的15%。这些数字指向由标准单处理器Pentium IV PC对所有的计算和证明获得结果，并且证明了网格的生成在PFEM中具有可以接受的成本。网格重划的成本类似于以下报道24。事实上，相当速度能够用平行计算技术来获得。5 边界表明的识别在PFEM法中的一个主要任务是对边界域进行正确的定义。边界节点有时会明确标出。在其他情况下，节点的总体设置是唯一的可用信息并且其算法必须先确认出边界节点。在我们的作品中，我们使用一个扩展的Delaunay分区来识别边界节点19。考虑到该节点遵循变量h（x）的分布，其中h（x）一般是两个节点间的最小距离。在一个半径比h大的空心球体上的所有节点，被认为是边界节点。在实践中是一个接近但大于1的参数。取值于1.3至1.5之间已被发现是在所有分析实例中的最佳取值。这个标准与阿尔法形状概念相符12。一旦边界上的节点被确定，该边界表面就由所有的多面体表面（或在二维情况下为多变形）具有边界上所有节点并且只属于一个多面体定义出来。该方法据描述还允许确定在主要流体域外的隔离流体粒子。这些粒子被视为外部边界处的部分，其压力值固定位大气压力值。我们回想一下，每个粒子是由其所属的固体域或流体域密度表征的一个质点。当边界元由于从主分析领域分离出发成为一个节点而被淘汰时，质量消失后会在该“飞”节点下落并且生成一个新的由阿尔法形状算法确定的边界元重新生成时恢复。边界识别法同样对检测流体域和一个固定边界间的接触条件以及不同固体间的彼此相互作用适用，具体将在下一节进行详述。我们强调的主要区别PFEM法和经典有限元法的仅仅是重新网格化技术和在每个时间步长上的域边界的鉴定。6 在PFEM法中对接触前提的修正6.1 流体和固定边界之间的接触在PFEM法中固定边界上的规定速度的条件被强读式施加在边界节点上。这些节点可能属于固定的外部边界，或者连接到相互作用固体上的移动边界。流体粒子和固定边界间的接触是由于不可压缩条件，其自然地防止流体节点渗透到固体边界32,36。6.2 固体-固体界面间的接触两个固体界面之间的接触简单地通过引进一个两相互作用的固体界面间的接触元件层来看待。这个层是在网格生成步骤期间通过设定两个固体边界间的最小距离（hc）自动生成的。如果该距离超过最低值（hc），则所生成的元被视作流体元。否则该元被视为接触元，其中切线和法向量之间的关系以及相应的位移已经介绍过了（图2）。该算法已被证明是非常有效的，并且它允许对两个或更高相互作用的在水中移动的物体之间复杂的摩擦接触条件通过非常简单的方式进行识别和建模。该算法也可有效地用于结构力学应用中的模拟刚性或弹性之间的摩擦接触条件7,36。7 河床冲刷的建模河床冲刷和明渠流中泥沙的移动的预测是河流与环境工程的许多领域中的重要任务。河床冲刷可以导致流域斜坡的不稳定。它也可以破坏桥梁桩基础，从而造成结构失效。河床冲刷的建模也与牵扯到溢流情况下的土坝的表面材料的发展的预测相关。河床冲刷是洪灾中环境破坏的主要原因之一。河床冲刷模型传统地基于侵蚀率和剪应力水平之间的关系25。在最近的工作中我们提出了一中PFEM扩展的方法来模拟河床冲刷35,36。侵蚀模型是基于河床表面源于流体剪切应力的摩擦效果。由此产生的侵蚀模型类似于典型应用与摩擦接触条件下的表面磨损建模的Archard法1。该算法用于模拟在流化床中土壤/岩石粒子的侵蚀，如下所示：1. 计算每一点的河床冲刷表面有流体运动产生的切向应力。在三维问题中 =(2s+ t)2其中s和t是由河床节点正常方向n定义的平面中的切向应力。在二维问题中的值可估计如下：其中vkt是切向速度在节点k的模，hk是一个沿正常河床点k的规定距离。通常，hk是相邻节点k的流体元件的最小数量级（图3）。2. 计算源自河床表面的切向应力的摩擦作用为公式（12）将时间综合为3. 当nwf超过临界域值wc时河床点的侵蚀开始出现。4. 如果在一个河床节点上nwf Wc，那么该节点将被从床区分离，并被允许移至流体流中。这样的话，该围绕床节点的河床元的片的质量则会消失并且转移到新的流体节点上。该质量随后转移到流体中。5. 沉积物的沉淀可通过一个先前步骤中描述的逆转程序来建模。因此，靠近河床表面的以低于阈值的速度暂停的节点附着到河床表面。图3示出了河床冲刺算法介绍的示意图。8 范例8.1 水流对岩石的拖动预测岩石被水流拖动的临近速度在液压、港口、土木和环境工程领域的许多问题上是十分重要的。PFEM法已经成功地应用于1Tn quasi-spherical岩石由水流引起的运动的研究中。该岩石符合岩石集合保持刚性的性质。所分析的岩石与其余岩石之间的摩擦条件已经被假定。图4a显示，1m/s速度的水流无法移动独立的岩石。水流速度增加到2m/s才导致岩石的运动，在图4b中示出。8.2 海浪对码头及防波堤的影响图5和图6显示的破碎波在两个不同方位对包含钢筋混凝土块的防波堤（44 MTS各一个）的影响的分析。这些数值与对位于西班牙A Coru na的Langosteira海港的PFEM法研究相对应。8.3 土壤侵蚀问题图7示出了PFEM法对河床中土壤侵蚀、泥沙输移以及材料沉积的建模能力。土壤粒子首先在喷流作用下从河床表面分离出来。然后其由流体输送，并最终由于重力作用下降到河床表面的下游点。图8示出了西班牙A Coru na的Langosteira海港防波堤斜坡的未保护部分的逐步侵蚀。无保护的上肩区在海浪的作用下被逐步侵蚀。图9显示了河床中临近桥墩底部由于受到水流作用（图中没有示出水流）受到逐步侵蚀并且土壤粒子受到拖曳。注意桥梁基础由于临近土壤颗粒因为侵蚀的移动所显示的信息。关于PFEM法对河床冲刷问题的其他应用可以在35,36中找到。8.4 海浪对道路边坡侵蚀导致的货车落海图10显示了临近海岸的大量土壤由于受到海浪作用而导致的逐步侵蚀并且随后导致用二维物体代表货车的部分落入海中的代表性实例。该对象作为刚性固体进行模拟。这个例子中，虽然还是很简单的示意图，但是显示出PFEM模拟涉及到土壤侵蚀、自由表面波和刚性/可变形结构的复杂FSSI问题的可能性。8.5 模拟滑坡该PFEM法特别适用于对滑坡及其对周边结构的影响进行建模和仿真。图11显示了一个二维在落在两个相邻建筑物之间滑坡的示意图。滑坡物质被模拟为粘性不可压缩流体。8.6 Lituya湾滑坡当山体滑坡出现在水库附近是一个很有趣的例子43。落入水库中的材料碎片通常会引起大波，可能高出大坝的水波并引起意想不到的洪水，进而导致对结构和下游地区的人群造成严重破坏。在这个例子中，我们提出了一些对发生正1958年7月9日的Lituya湾滑坡进行三维分析的结果。滑坡源自一场地震和90亿吨的岩石坠落到海湾进而引起斜坡对岸达到524米高的大浪。图13示出用PFEM法仿真的滑坡图像。该滑坡体被模拟成具有规定剪切模量的连续体。滑坡体和地下土壤间的无摩擦影响已被考虑。该分析也没有考虑到土壤材料由滑坡体在运动期间引起的侵蚀和拖曳。该PFEM结果已与临近水库的北部山的最高水位相比较。这座山通过PFEM法获得的最大水位是551米。这比524米的数值还高5%。观测实验由13,14。该最大高度位置与观测值相差300米13,14。在南坡的最大水位高度观测值是208米，而PFEM结果（此处未示出）是195米（6%的误差）。关于这个例子的PFEM法的更多信息可以在38,39中找到。9 结论粒子有限元法（PFEM）是一种针对求解涉及到流体-固体和固体-固体界面以及河床冲刷及其他复杂现象中流体和固体粒子、表面波、溅水、摩擦接触流体-土壤-结构相互作用（FSSI）问题的很有前途的数值技术。该PFEM的成功在于对不可压缩连续体的方程的准确有效的解法采用更新的拉格朗日配方和稳定的有限元法，使得其允许对所有变量的低阶元使用同阶插值法。其他基本的解决方案成分是高效的有限元网格重新生成，通过阿尔法形状技术确定边界节点以及用简单算法来修正摩擦接触，通过网格生成确定流体-固体及固体-固体界面的侵蚀/磨损。给出的例子显示了PFEM对解决工程中一大类实际存在的FSSI问题的潜力。致谢This research was partially supported by project SEDUREC of the Consolider Programme of the Ministerio de Educacin y Ciencia (MEC) of Spain and the projects SAFECON and REALTIME of the European Research Council of the European Commission (EC).Thanks are also given to the Spanish construction company Dragados for financial support for the study of harbour engineering problems with the PFEM .