CN102262689A

CN102262689A - 能够实现实体/sph耦合作用的混合单元

Info

Publication number: CN102262689A
Application number: CN2011101323942A
Authority: CN
Inventors: 陈皓; 许敬晓; 王季先
Original assignee: Livermore Software Technology LLC
Current assignee: Livermore Software Technology LLC
Priority date: 2010-05-26
Filing date: 2011-05-20
Publication date: 2011-11-30
Anticipated expiration: 2031-05-20
Also published as: EP2390801A1; JP5792995B2; JP2011248882A; CN102262689B

Abstract

本发明公开了一种能够实现SPH粒子和FEM实体之间的耦合作用(510)的混合单元(110，120，130，140，150，160)。根据本发明的一个方面，配置混合单元(110，120，130，140，150，160)来实现基于FEM的实体单元与基于SPH的一个或多个相应粒子之间的耦合作用(510)。混合单元(110，120，130，140，150，160)在计算机辅助工程网格模型(220)中被定义为位于SPH粒子和FEM实体之间的缓冲器或者接口。例如，网格模型(220)的一部分由于具有产生大形变的可能性而包含有SPH粒子，而模型的其余部分包括FEM实体单元。混合单元被置于实体和粒子之间。每个混合单元包括两层：实体层和粒子层。

Description

能够实现实体/SPH耦合作用的混合单元

技术领域

本发明一般涉及计算机辅助工程分析，特别涉及使用基于有限元法(finiteelement method，FEM)的实体单元与基于平滑粒子流体动力学(smoothedparticle hydrodynamics，SPH)的粒子相结合来对产生大形变的结构进行时间推进式仿真(例如，汽车碰撞或爆炸仿真)的方法和系统，其中在所述实体单元和SPH粒子之间建立至少一层混合单元以实现耦合作用。

背景技术

连续介质力学已被用于仿真诸如固体和流体(即液体和气体)等连续物质。微分方程被用于解决连续介质力学中的问题。大量的数值方法已经在使用。最常见的方法之一是有限元分析(finite element analysis，FEA)或有限元法(FEM)，这是一种广泛应用于工业上的计算机化的方法，用于对与诸如三维非线性结构设计和分析等复杂系统相关的工程问题进行建模和求解。有限元分析的名称源于对所关注的目标物体的几何特征进行描述的方式。随着现代数字计算机的出现，有限元分析已在有限元分析软件中实现。基本上，有限元分析软件提供有关于几何描述的基于网格的模型、以及该模型中每一点处的相关材料性能。在此模型中，所分析的系统的几何特征被表示为各种大小的实体、壳体和梁，这些被称为单元。各单元的顶点被称为节点。该模型是由有限数目的单元组成的，这些单元都被分配有一材料名称以便于与材料性能相关联。因此，该模型表示了被分析的目标物体所占的物理空间以及它的周围环境。有限元分析软件随后涉及一列出了每种材料性能(例如，应力-应变本构方程、杨氏模量、泊松比和热传导率)的表格。此外，指定了目标物体的边界条件(即负荷、物理限制等)。遵循此种方式，建立了目标物体及其环境的模型。

一旦定义了模型之后，可以使用有限元分析软件对特定负荷或初始条件下的物理性能进行数值仿真。有限元分析软件广泛应用在汽车行业中以仿真汽车正面和侧面碰撞、乘坐的假人与安全气囊的相互作用、以及由金属薄板到车体部件的成形过程。这种仿真为工程师提供了宝贵的见解，使得工程师能够提高汽车的安全性，并且更快速地为市场带来新的车型。该仿真是在时域中进行的，这意味着有限元分析是在从首个求解周期开始的多个求解周期中计算的，在每个后续的求解周期中，仿真时间以时间步长ΔT递增。这种仿真被称为时间推进式仿真。

最具挑战性的有限元分析的任务之一是对涉及到产生极大形变的结构的撞击事件进行仿真，例如，汽车碰撞或爆炸仿真。随着现代计算机的改进，工程师不仅要仿真遭到结构破坏的撞击事件中的性能，他们也想在撞击事件完全破坏之前对产生弯曲之后的结构性能进行仿真。但是，很难用有限元分析和实体单元来仿真这些现象。例如，表示保险杠的泡沫材料的实体单元可能被挤压或压缩，变得过于扭曲或压扁，从而成为零体积或负体积，这将导致在仿真过程中产生数值问题(例如，因为数字计算机中的无效数字而仿真失败)。

为了解决零体积或负体积的问题，采用基于平滑粒子流体动力学(SPH)的粒子来替换这些无效的实体单元。然而，FEM和SPH的数学公式不同。为了使粒子和实体单元并存在同一模型中，必须建立一些联系以便于连接粒子和实体单元。现有技术的方法一直采用刚性连接特定粒子与实体单元的连接接口。然而，这种方法因为连接接口的随意放置(即刚性连接)通常会导致仿真结果非常不真实。例如，粒子和实体单元最初可能合理地被连接在一起。但是，当他们以不可预测的形式产生形变时，随意放置的这些刚性连接可能会导致连接非常不符合实际情况。

因此，希望在计算机辅助工程分析模型中提供一种更加真实的接口以使SPH粒子和FEM实体单元可以共存从而避免现有技术方法的问题和缺陷。

发明内容

本发明公开了能够实现SPH粒子和FEM实体之间的耦合作用的混合单元。根据本发明的一个方面，配置混合单元来实现基于FEM的实体单元与基于SPH的一个或多个相应粒子之间的耦合作用。混合单元在计算机辅助工程(CAE)网格模型中被定义为位于SPH粒子和FEM实体之间的缓冲器或者接口。例如，网格模型的一部分由于具有产生持久大形变的可能性而包含有SPH粒子，而该模型的其余部分包括FEM实体单元。混合单元被置于实体和粒子之间。每个混合单元包括两层：实体层和粒子层。

首先，混合单元的耦合作用是通过以下过程实现的：基于FEM计算实体层中的节点加速度、速度和位移以及单元应力；将计算出的该节点量与单元应力状态一起映射到粒子层，其中单元应力状态包括单元的应力值和当前材料状态(例如，弹性、塑性、屈服、应变硬化等)；基于SPH计算粒子层中的内力，并将该内力转移回实体层以计算节点力用于下一求解周期。

在接下来求解周期中，实体层的节点位移得到更新，并被映射到相应SPH粒子。基于SPH计算出内力，然后将该内力转移回实体层以计算节点力用于下一求解周期。

根据本发明的另一方面，一种计算机辅助分析模型可以包括位于周边的实体有限单元以及位于其余位置的SPH粒子。位于周边或边缘的实体单元被配置用于满足边界条件。

根据另一方面，SPH粒子被用于置换已受到超出其屈服极限的破坏的实体单元。SPH粒子可以使用表示材料的应变硬化效应的较柔软的材料模型进行建模。

通过结合附图详阅接下来对实施例的详细描述，本发明的其它目标、特征和优点将是显而易见的。

附图说明

本发明的这些以及其它特征、方面和优点通过以下描述、权利要求和附图能得到更好的理解：

图1为根据本发明实施例的各种示例性混合单元的示意图；

图2A-2B为根据本发明实施例的可以使用混合单元进行数值仿真的受到大形变的示例性结构的示意图；

图3A-3D为根据本发明一个实施例中激活混合单元的耦合作用的示例性顺序的一些流程图；

图4为根据本发明一个实施例可用于数值仿真屈服后结构性能的应力-应变曲线图；

图5A-5C共同示出了根据本发明实施例使用混合单元来对响应冲击负荷的屈服后结构性能和大形变进行数值仿真的示例性处理过程的流程图；

图6为可用于实施本发明实施例的示例性计算设备的主要部件的功能图。

具体实施方式

首先参照图1，为根据本发明一个实施例的各种示例性混合单元的示意图。混合单元包括两部分：实体层和相应的粒子层。实体层包括基于FEM的实体单元，而相应的粒子层包括基于SPH的一个或多个粒子。实体单元包括但不限于六面体、楔形和四面体。混合单元110是含有一个相应粒子的六面体，单元120是含有一个粒子的楔形单元，而单元130是含有一个粒子的四面体。含有一个以上粒子的示例性混合单元如含有8个粒子的单元140、含有6个粒子的单元150以及含有4个粒子的单元160所示。其它数目的粒子也可在根据本发明的另一实施例中实现。例如，含有27个粒子的六面体(未示出)。

混合单元的耦合作用是通过将实体层与粒子层相关联来实现的。内力计算方法的细节在如图5C所示的流程图中示出。例如，实体层用作对相应粒子层中粒子的约束条件。在一个实施例中，实体层的体积被设定为包围粒子的区域。

图2A-2B示出了目标物体210(即刚性球形式的弹丸)以较高的速度(箭头215指出)接触一结构220(即部分显示为网格模型的板)产生一冲击负荷的连续过程。结构中受到冲击负荷的部分使用单元225建模(如虚线所示)。单元225最初可以使用FEM实体单元建模。当球210对板220产生冲击时，单元225可能会发生破坏或者弯曲(见图4以及相应地关于材料破坏和弯曲的定义的说明)。失效的单元采用SPH粒子来替代以便于继续仿真。因为SPH粒子和FEM实体单元采用不同的公式，所以在SPH粒子和FEM实体单元之间建立了以混合单元形式出现的接口以便于实现耦合作用。

为了进一步证明上述的实例，图3A-3D示出了结构(板220)的连续平面图。首先，在图3A中板220被表示为全部的实体单元。接下来，在图3B中，中心实体单元失效并采用SPH粒子来替换(表示为包含中心点的带阴影的圆圈)。这可能是由于弹丸/球210猛烈地触碰到板220所导致的。创建至少一层混合单元(表示为带阴影的单元)作为实现SPH粒子和实体单元之间的耦合作用的接口。随后，在图3C中，中心单元周围的更多单元失效，并被SPH粒子所替代。可以看出，混合单元的接口是动态调整的，以便于永远位于SPH粒子和实体单元之间。

此外，混合单元可以放置在CAE模型的边界上，而SPH粒子位于模型的其余部分。图3D中所示的设置证明本发明的该方面。

图4示出了根据本发明一个实施例可用于确定结构中的屈服后结构性能的应力-应变曲线。该曲线400具有表示应力的纵轴402和表示应变的横轴404。材料包括两个区域：弹性区域406和塑性区域408。塑性区域408被进一步分为三类：屈服区域424、应变硬化区域426和颈缩区域428。在应力-应变曲线400的弹性区域的顶端为屈服点414，屈服应力与之相对应。极限应力对应于极限强度点416，而断裂或破裂应力对应于破坏位置418。根据一个实施例，FEM实体单元用于对材料的弹性性能进行建模。只要材料超出屈服点，就生成SPH粒子来取代实体单元。使用较为柔软的材料模型来对替换的SPH粒子进行建模，从而能更加真实地对应变硬化效应进行仿真。

现在参考图5A，为根据本发明实施例使用混合单元来对响应冲击负荷的屈服后结构性能和大形变进行数值仿真的示例性过程500的流程图。过程500优选在软件中实现。

过程500从在步骤502中定义结构(例如，汽车、飞机)的计算机辅助分析网格模型(例如，FEM网格模型)开始。该网格模型包括一个或多个混合单元，用于表示结构中最可能受到大形变的部位，例如，碰撞仿真中的汽车保险杠。该网格模型被用于时间推进式仿真中。接下来，在步骤504中，在时间推进式仿真的初始阶段(即时间等于0，或者第一个求解周期)对所有的单元和参数进行初始化。随后，过程500在判定步骤506中检查混合单元的耦合作用是否已经激活。如果“否”，过程500移到步骤508中通过将混合单元看做为FEM实体单元进行仿真。步骤508的细节在图5B和相应的说明中进行了详细描述。换言之，在耦合作用没有激活前，时间推进式仿真是使用FEM来进行的。否则，如果“是”，过程500移到步骤510，使用具有耦合作用的混合单元进行时间推进式仿真。图5C和相关的描述是针对步骤510的。

过程500移到步骤516，递增下一求解周期的仿真时间。然后，在判定步骤518中，确定时间推进式仿真是否已经结束。例如，与预设总仿真时间相对照来检查仿真时间。如果不是，过程500移回判定步骤506重复执行其余步骤以进行下一求解周期，直至判定步骤518为真。过程500在此之后结束。

图5B示出了步骤508的进一步细节。在步骤522中，过程500获取包括混合单元在内的每个单元的节点加速度、速度和位移。在一个实施例中，该节点量是在FEM的显式求解器(例如，f＝mxa，其中“f”为节点力，“m”为节点质量，“a”为节点加速度)中获得的。接下来，在步骤524中，根据FEM的实体公式计算单元内力。最后，在步骤526中，可以计算出节点力以进行下一求解周期，该节点力包括单元内力的贡献作用。任何给定的节点可以接收所有相连单元的贡献。

图5c示出了步骤510的更多细节。在步骤532中，处理过程获取混合单元的实体层中的节点量(即节点加速度、速度和位移，以及单元应力状态)，大体与FEM下的步骤522相似。接下来，在步骤533中，在耦合的开始阶段将取得的节点量和单元状态映射到相应的粒子层。然后在后续的周期中，只将节点位移映射到相应的粒子层。实体层的位移得到更新并映射到对应的SPH粒子作为约束条件。换言之，SPH粒子受到基于FEM计算出的节点位移的限制。计算出粒子层中的内力，并将其转移到实体层用于基于FEM计算节点力以进行下一求解周期。换言之，使用SPH来获取内力，在这些交换期间实体层内的内力计算在混合单元中被暂停。单元应力状态至少包括计算出的该单元的应力值。此外，单元应力状态包括在单元应力状态历史变量中该单元的当前状态。该变量或者其他等效的方法被用于跟踪整个时间推进式仿真过程中单元的状态。换言之，单元的屈服后状态可以从该单元的历史变量中确定。

随后在步骤534中，基于SPH方程计算相应粒子层中粒子所在的混合单元的内力。接下来，在步骤535中，该内力被转移回实体层。换言之，实体层的单元内力被从相应粒子层中的粒子计算获得的内力所替换。最后，步骤536与步骤526大致相同，计算出用于下一求解周期的节点力，该节点力包括了单元内力的贡献。

根据一个方面，本发明涉及能够实施本文所述功能的一个或多个计算机系统。一个计算机系统600的实例如图6所示。计算机系统600包括一个或多个处理器，例如处理器604。该处理器604与计算机系统内部通讯总线602相连。根据该示例性计算机系统来对各种软件实施例进行实施。在阅读该说明后，如何使用其它计算机系统和/或计算机体系来实施本发明，对于本领域技术人员而言是显而易见的。

计算机系统600还包括主内存608，优选随机存取存储器(random accessmemory，RAM)，还可能包括辅助存储器610。该辅助存储器610可以包括，例如一个或多个硬盘驱动器612和/或一个或多个可移动存储驱动器614，典型的有软盘驱动器、磁带驱动器、光盘驱动器等等。可移动存储驱动器614以众所周知的方式从可移动存储单元618读取数据和/或向其写入数据。可移动存储单元618，典型的有软盘、磁带、光盘等，通过可移动存储驱动器614读取和写入数据。应当理解，可移动存储单元618包括一计算机可用存储介质，其中存储有计算机软件和/或数据。

在另一实施例中，辅助存储器610可包括允许将电脑程序或其它指令加载到计算机系统600中的其它类似装置。这些装置可能包括，例如，可移动存储单元622和接口620。这样的例子可能包括程序盒和盒式接口(如视频游戏设备中存在的)、可移动内存芯片(如可擦除可编写只读式存储器(EPROM)、通用串行总线架构(USB)快闪存储器或PROM)及相关插座，以及其它允许将软件和数据从可移动存储单元622传送到计算机系统600的可移动存储单元622和接口620。一般来说，计算机系统600是由操作系统软件来控制和协调的，它执行诸如流程调度、内存管理、网络和I/O服务的任务。

也可能存在连接到总线602的通信接口624。通信接口624允许在计算机系统600和外部设备之间传输软件和数据。通信接口624的实例可能包括调制解调器、网络接口(如以太网卡)、通信端口、符合个人计算机存储卡国际协会(Personal Computer Memory Card International Association，PCMCIA)标准的插槽和卡等。通过通信接口624传输软件和数据。计算机600通过根据一套特定的规则(即协议)，通过数据网络来与其它计算设备进行通信。常见的协议之一是经常在因特网上使用的TCP/IP(传输控制协议/互联网协议)。一般来说，通信接口624将数据文件打包成较小的数据包并通过数据网络传输该数据包，或者对接收的数据包重新组装来获得原始数据文件。此外，通信接口624处理每个数据包的地址部分，以便使其到达正确的目标端口，或者以拦截最终目标为计算机600的数据包。在本文中，术语“计算机程序介质”、“计算机可读介质”、“计算机记录介质”和“计算机可用的介质”被用来泛指诸如可移动存储驱动器614(例如，闪存驱动器)的介质，和/或安装在硬盘驱动器612中的硬盘。这些计算机程序产品是用于为计算机系统600提供软件的装置。本发明涉及这样的计算机程序产品。

计算机系统600可能还包括输入/输出(I/O)接口630，用于来计算机系统600提供显示器、键盘、鼠标、打印机、扫描仪和绘图仪等等。

计算机程序(也称为计算机控制逻辑)作为应用模块606存储在主存储器608和/或辅助存储器610中。计算机程序也可通过通信接口624来接收。这种计算机程序，在执行时，使计算机系统600能够执行本文所讨论的本发明的功能。特别是，计算机程序，在执行时，使处理器604能够执行本发明的功能。因此，这些计算机程序象征计算机系统600的控制器。

在通过软件实施本发明的实施例中，软件可以存储在计算机程序产品中，并可以使用可移动存储驱动器614、硬盘驱动器612或通信接口624加载到计算机系统600中。应用模块606，在被处理器604执行时，使处理器604执行本发明所述的功能。

为完成想要做的任务，主存储器608可以载入有一个或多个应用模块606，在无论有否通过输入/输出接口630输入使用者指令的情况下，此模块可以被一个或多个处理器604执行。在运行时，当至少有一个处理器604执行应用模块606其中一个时，计算出结果并存储在辅助存储器610(即硬盘驱动器612)中。分析结果根据用户的指令通过I/O接口630以文本或图形演示的方式报告给客户。

虽然本发明是根据特定实施例进行描述的，但这些实施例仅仅是对本发明进行说明，而非限制的。对具体公开的示例性实施例进行的改动或改变会对本领域的普通技术人员起到暗示和提醒作用。例如，尽管产生大形变的示例性结构被表示并描述为受到弹丸冲击的板，其它受到冲击负荷的结构也可以使用本发明声明的方法进行数值仿真，例如，汽车碰撞中的汽车保险杠。此外，尽管实体单元被表示并描述为六面体、楔形或四面体，也可以使用其他类型的实体单元来替代，例如，五面体。总之，不能用所述的和特指的典型实施例来限定本发明的保护范围，并且，所有可轻易地对本领域普通技术人员能起到暗示和提醒作用的所有改动仍应涵盖在本专利申请的精神和所附加的权利要求的范围之内。

Claims

1.一种在计算机系统(600)中执行的使用混合单元(110，120，130，140，150，160)来实现基于有限元法FEM的实体单元和基于平滑粒子流体动力学SPH的粒子之间的耦合作用的方法(500)，其特征在于，所述方法包括：

定义(502)代表一结构的网格模型(220)，所述网格模型包括多个FEM实体单元、至少一个SPH粒子以及位于所述FEM实体单元和所述至少一个SPH粒子之间的至少一层混合单元(110，120，130，140，150，160)，每个所述混合单元(110，120，130，140，150，160)包括实体层和粒子层；

使用所述网格模型(220)对处于负荷条件下的结构执行(504，506，508，510，512，516，518)时间推进式仿真，所述仿真包括通过在实体层中计算出的一组节点量和在粒子层中计算出的相应内力之间的交换(532，533，534，535，536)所实现的耦合作用(510)；以及

获取(506，536)每个求解周期的时间推进式仿真结果，其中所获得的结果被保存到存储设备内的文件中，并根据用户指令以图形方式显示到显示器上。

2.根据权利要求1所述的方法，其特征在于，所述在实体层中计算出的一组节点量和在粒子层中计算出的相应内力之间的交换(532，533，534，535，536)是通过以下步骤来完成的：

基于FEM计算(532)所述实体层中的一组节点量以及单元应力状态；

将所述一组节点量以及单元应力状态从实体层映射(533)到粒子层；

基于SPH计算(534)所述粒子层中的相应内力；

将所述相应内力转移(535)到实体层中以基于FEM计算出节点力用于时间推进仿真的下一求解周期。

3.根据权利要求2所述的方法，其特征在于，所述一组节点量(533)在耦合作用的初始阶段(506，510)包括节点位移、速度和加速度。

4.根据权利要求3所述的方法，其特征在于，所述节点量(533)在耦合作用的初始阶段之后(510)包括节点位移。

5.根据权利要求4所述的方法，其特征在于，所述节点位移(533)被映射到所述相应SPH粒子上作为约束条件。

6.根据权利要求4所述的方法，其特征在于，在所述交换(532，533，534，535，536)期间，所述计算相应内力被暂停。

7.根据权利要求2所述的方法，其特征在于，所述粒子层包括一个或多个SPH粒子。

8.根据权利要求2所述的方法，其特征在于，所述单元应力状态包括跟踪所述每个混合单元的材料状态(400)的应力状态历史变量。

9.根据权利要求8所述的方法，其特征在于，所述每个混合单元的材料状态(400)包括屈服后区域(424)、弹性区域(406)和塑性区域(408)。

10.一种使用混合单元(110，120，130，140，150，160)来实现基于有限元法FEM的实体单元和基于平滑粒子流体动力学SPH的粒子之间的耦合作用的系统(600)，其特征在于，所述系统包括：

存储用于一个或多个应用模块(606)的计算机可读代码的存储器(608)；

连接至所述存储器(608)的至少一个处理器(604)，所述至少一个处理器(604)执行存储器(608)中的计算机可读代码以使得所述一个或多个应用模块(606)执行以下操作：

11.根据权利要求10所述的方法，其特征在于，所述在实体层中计算出的一组节点量和在粒子层中计算出的相应内力之间的交换(532，533，534，535，536)是通过以下步骤来完成的：

基于SPH计算(534)所述粒子层中的相应内力；

12.一种含有指令的计算机可读介质，所述指令用于控制计算机系统(600)通过执行根据权利要求1或2所述的方法(500)使用混合单元(110，120，130，140，150，160)来实现基于有限元法FEM的实体单元和基于平滑粒子流体动力学SPH的粒子之间的耦合作用。