物理仿真方法、装置、计算机设备和存储介质与流程

1.本技术涉及仿真技术领域，特别是涉及一种物理仿真方法、装置、计算机设备和存储介质。


背景技术：

2.随着科学技术的发展，出现了仿真技术，仿真是指通过数学方程组描述某一过程，并通过数值计算方法在计算平台上求解方程组，从而预测该过程的结果的一种技术手段。
3.目前，在很多场景中都需要对多物理现象进行仿真。多物理仿真涉及多于一种物理现象，且不同物理场之间会相互产生影响的仿真问题。其中物理现象包括但不限于流体、气体、固体、弹性体、生物组织的力学、声学、电学、光学以及磁学效应。例如，需要对子弹打出后，在空气中的移动进行仿真，也需要确定子弹打入某个物体的过程中，物体的变化状态与物体周围的空气的变化状态进行仿真。
4.传统技术中，在对多物理现象进行仿真时，通常是按顺序依次逐个地仿真涉及的每一种物理现象，求解某一物理现象时，其他物理现象的影响会被设定为已知，存在仿真准确度低的问题。


技术实现要素：

5.基于此，有必要针对上述技术问题，提供一种物理仿真方法、装置、计算机设备和存储介质。
6.一种物理仿真方法，所述方法包括：在针对目标固体对象与目标气体对象的物理仿真过程中，获取所述目标固体对象在当前时间步的当前固体状态以及所述目标气体对象在当前时间步的当前气体状态；基于所述当前固体状态对所述目标固体对象进行位置移动仿真，得到中间固体状态；基于所述当前气体状态对所述目标气体对象进行对流仿真，得到中间气体状态；基于所述中间气体状态、所述中间固体状态以及交界约束条件进行耦合仿真，得到所述目标固体对象在下一时间步的目标固体状态以及所述目标气体对象在下一时间步的目标气体状态；其中所述交界约束条件为所述目标固体对象与所述目标气体对象之间的交界对应的约束条件，所述交界约束条件包括交界压强约束条件。
7.一种物理仿真装置，所述装置包括：当前状态获取模块，用于在针对目标固体对象与目标气体对象的物理仿真过程中，获取所述目标固体对象在当前时间步的当前固体状态以及所述目标气体对象在当前时间步的当前气体状态；位置移动仿真模块，用于基于所述当前固体状态对所述目标固体对象进行位置移动仿真，得到中间固体状态；对流仿真模块，用于基于所述当前气体状态对所述目标气体对象进行对流仿真，得到中间气体状态；耦合仿真模块，用于基于所述中间气体状态、所述中间固体状态以及交界约束条件进行耦合仿真，得到所述目标固体对象在下一时间步的目标固体状态以及所述目标气体对象在下一时间步的目标气体状态；其中所述交界约束条件为所述目标固体对象与所述目标气体对象之间的交界对应的约束条件，所述交界约束条件包括交界压强约束条件。
8.在一些实施例中，所述耦合仿真模块用于：将所述中间气体状态作为气体压强约束方程的已知变量，得到气体约束条件；其中，所述气体压强约束方程为中间气体状态与目标气体压强之间的约束方程；将所述中间固体状态作为固体压强约束方程的已知变量，得到固体约束条件；其中，所述固体压强约束方程为中间固体状态与目标固体压强之间的约束方程；根据所述气体约束条件、所述固体约束条件以及所述交界约束条件进行耦合仿真，得到所述目标固体对象在下一时间步的目标固体状态以及所述目标气体对象在下一时间步的目标气体状态。
9.在一些实施例中，所述气体压强约束方程为将气体对应的压强推进方程代入气体动量修复方程得到，所述固体压强约束方程为将固体对应的压强推进方程代入固体动量修复方程得到。
10.在一些实施例中，所述目标固体状态包括所述目标固体对象对应的目标固体速度，所述目标气体状态包括所述目标固体气体对应的目标气体压强；所述耦合仿真模块用于：根据所述气体约束条件、所述固体约束条件以及所述交界约束条件进行耦合仿真，得到所述目标固体对象对应的目标固体压强以及所述目标气体对象对应的目标气体压强；根据所述目标固体压强以及所述固体动量修复方程修复得到所述目标固体对象对应的目标固体速度；根据所述目标气体压强以及所述气体动量修复方程修复得到所述目标气体对象对应的目标气体速度。
11.在一些实施例中，所述气体压强约束方程为：所述固体压强约束方程为：其中，δt为时间步长，sf代表气体对应的刚性系数，代表下一时间步气体的密度，为拉普拉斯算子，代表下一时间步气体的目标气体压强，代表气体的中间压强状态；代表空间散度算子，代表气体的中间速度状态，ss代表固体对应的刚性系数，代表下一时间步固体的目标固体密度，代表下一时间步固体的目标固体压强，代表固体的中间压强状态；代表固体的中间速度状态。
12.在一些实施例中，所述交界约束条件包括交界压强约束条件以及交界速度约束条件，所述交界压强约束条件包括交界处固体压强与气体压强一致，所述交界速度约束条件包括交界接触面所述目标固体对象与所述目标气体对象法向速度一致。
13.在一些实施例中，所述位置移动仿真模块用于：在拉格朗日视角下，基于所述当前固体状态对所述目标固体对象进行位置号移动仿真，得到以拉格朗日视角记录的中间固体状态；所述对流仿真模块用于：将以拉格朗日视角记录的中间固体状态投射到所述目标气体对象对应的气体网格中，并转化为欧拉视角记录的中间固体状态；基于欧拉视角记录的中间固体状态，以及欧拉视角记录的离散到所述气体网格中的当前气体状态，对所述目标气体对象进行对流仿真，得到欧拉视角记录的中间气体状态。
14.一种计算机设备，包括存储器和处理器，所述存储器存储有计算机程序，所述处理器执行所述计算机程序时实现物理仿真方法的步骤。
15.一种计算机可读存储介质，其上存储有计算机程序，所述计算机程序被处理器执
行时实现物理仿真方法的步骤。
16.上述物理仿真方法、装置、计算机设备和存储介质，由于可以先分别对固体对象以及气体对象分别进行仿真，得到中间的状态信息，然后再基于中间气体状态、中间固体状态以及交界约束条件进行耦合仿真，即一起进行求解，以对中间的状态进行修复，因此能对多物理现象更为准确的描述，使得得到的目标固体状态以及目标气体状态更加准确，提高了物理仿真的准确性。
附图说明
17.图1为一些实施例中物理仿真方法的应用环境图；
18.图2为一些实施例中物理仿真方法的流程示意图；
19.图3为一些实施例中耦合步骤的流程示意图；
20.图4为一些实施例中物理仿真装置的结构框图；
21.图5为一些实施例中计算机设备的内部结构图；
22.图6为一些实施例中计算机设备的内部结构图。
具体实施方式
23.为了使本技术的目的、技术方案及优点更加清楚明白，以下结合附图及实施例，对本技术进行进一步详细说明。应当理解，此处描述的具体实施例仅仅用以解释本技术，并不用于限定本技术。
24.本技术实施例提供的物理仿真方法，可以应用于如图1所示的应用环境中。其中，终端102通过网络与服务器104进行通信。终端102可以向服务器104发送仿真指令，服务器104接收到仿真指令，执行本技术实施例提供的物理仿真方法，得到物理仿真结果，服务器104可以存储物理仿真结果，也可以向终端102返回物理仿真结果，终端102可以通过动画的形式展示物理仿真结果对应的仿真画面。例如，以动画的形式播放子弹打入某个物体的过程中，子弹的变化状态、物体的变化状态以及物体周围的空气的变化状态，终端102可以但不限于是各种个人计算机、笔记本电脑、智能手机、平板电脑和便携式可穿戴设备，服务器104可以用独立的服务器或者是多个服务器组成的服务器集群来实现。可以理解，本技术实施例提高的物理仿真方法，也可以是终端执行的，或者由终端与服务器协同执行。
25.本技术实施例涉及物理仿真，首先对物理仿真进行解释。物理仿真中，仿真的结果是对空间中的相关物理量储存为一个随时间的变化数据序列，仿真顺序与时间流逝顺序相同，例如将第n个时间步的结果作为输入，计算出第n 1时间步结果的输出。
26.一般地，用户会输入最开始的即第0个时间步的数据(称为初始条件)，然后指定计算到某一个时间步m为止，一个物理仿真算法，其是不断重复从上一时间步到下一时间步的过程，即是第0步到第1步、第1步到第2步、
……
、第m-1步到第m步不断重复的过程，那么0，1，2,
…
,m-1,m的结果序列全体构成输出结果，即物理仿真结果。
27.而上文的“相关物理量储存为一个随时间的变化数据序列”也叫“物理量在时空间中的离散”，真实世界中，物理量是连续分布于空间的，即其在空间中的分辨率是无限大的，其随时间的变化也是连续的，分辨率也是无限大的，计算机无法储存无限多的数据，因此，仿真算法会首先将相关物理量在时空间中离散为可数个值。举例而言，某个仿真问题涉及
到一个三维正方体内部的密度、温度、与速度在时间0～1秒内的变化，如果将一个正方体空间分为10x10x10＝1000个网格，每个网格中心各放置一个密度、温度以及速度值，则在任意时刻(时间步)，此时刻对应的1000个密度、1000个温度、与1000个速度即对物理场的近似，其中空间点坐标与数据储存的位置可以由指定的映射关系l(x,y,z)
→
索引位置给出，又如果0～1秒被离散成了了十个时间步，一个网格对应1个密度、1个温度以及1个速度，有1000个网格，则每一个时间步的结果是1000个密度、1000个温度、和1000个速度的数据集合，总共10个不同时间步的数据集合构成的序列，是物理仿真的结果。
28.对于申请实施例的物理仿真方法，气体侧涉及到的物理量可以包括：密度ρf、压强pf以及气体速度uf。气体速度可以包括三个方向的速度分量：x轴方向上的速度u
fx
、y轴方向上的速度u
fy
以及z轴方向上的速度u
fz
。气体可以被离散到一套规则的背景网格上，例如上文提到的10*10*10个网格，固体可以分解为一个个的粒子，即多个固体粒子，粒子具体大小根据需要可以设置，固体侧涉及到的物理量可以包括：粒子位置x、粒子质量ms、粒子速度、以及粒子形变梯度。其中在三维空间中，粒子速度可以包括三个方向的速度分量：x轴方向上的速度u
sx
、y轴方向上的速度u
sy
以及z轴方向上的速度u
sz
。
29.在一些实施例中，如图2所示，提供了一种物理仿真方法，以该方法应用于图1中的服务器为例进行说明，包括以下步骤：
30.步骤s202，在针对目标固体对象与目标气体对象的物理仿真过程中，获取目标固体对象在当前时间步的当前固体状态以及目标气体对象在当前时间步的当前气体状态。
31.其中，目标固体对象是指需要进行物理仿真的固体，具体根据需要设置，可以为一个或者多个，目标固体对象具有弹性，可以发生形变，例如可以是子弹。目标气体对象是指要进行物理仿真的气体，目标气体对象可以被压缩，例如可以是空气，可以为一个或者多个。目标固体对象与目标气体对象相互作用，因此需要确定目标固体对象与目标气体对象在相互作用的过程中，各个时间步的状态变化。
32.在仿真的过程中，一个时间步一个时间步的进行仿真，基于当前时间步中各个物理现象的状态得到下一个时间步各个物理现象的状态。当前时间步是指当前的时间步。当前固体状态是当前时间步固体的状态。当前气体状态是当前时间步气体的状态。当前固体状态以及当前气体状态可以是输入的初始的状态，也可以是仿真得到的。例如，在开始进行仿真时，将初始的固体状态作为当前时间步的当前固体状态，将初始的气体状态作为当前时间步的当前气体状态，第一个时间步为下一时间步，即需要进行仿真得到物理现象的状态的时间步。基于初始的固体状态以及初始的固体状态仿真得到第一个时间步的固体状态以及气体状态后，则当前时间步更新为第一个时间步，第二个时间步为下一时间步。
33.固体状态可以包括粒子位置、粒子质量以及粒子速度。气体状态可以包括密度、压强以及气体速度。其中，粒子质量以及气体密度可以认为是不变的。
34.时间步长指跨越时间的长度，仿真到同样的终止时间，如果时间步长越小，则需要仿真的步骤数也就越多，整个仿真所需要计算的时间也会越长。时间步的步长可以根据需要设置。也可以是计算目标固体对象与目标气体对象这两种物理现象各自能够容许的最大时间步长，取其中的较小值作为耦合系统的时间步长。计算步长的算法可以根据需要设置，在满足cfl条件下计算得到时间步长。cfl条件又称克朗条件，英文全称为courant-friedrichs-lewy condition，cfl条件是有关数值计算稳定性和收敛的重要条件。
35.具体地，在针对目标固体对象与目标气体对象的物理仿真过程中，当仿真得到一个时间步的固体状态以及气体状态后，可以将该时间步的固体状态作为当前固体状态以及将该时间步的气体状态作为当前气体状态，以继续针对下一个时间步的固体状态以及气体状态进行仿真。
36.步骤s204，基于当前固体状态对目标固体对象进行位置移动仿真，得到中间固体状态。
37.位置移动仿真是指对目标固体对象的各个固体粒子进行移动仿真，得到粒子的位置以及形变梯度。位置移动仿真是在不考虑压强下进行的，其中气体的影响可以忽略或者将其作为边界条件。位置移动仿真的算法可以根据需要设置，例如可以是前向欧拉积分法。中间固体状态代表的是固体的中间态，由于其是在不考虑压强的影响下进行的，因此可以认为是还需要进一步进行修正的下一时间步的固体状态，即通过对中间固体状态进行修改，可以得到目标固体对象在下一时间步的目标固体状态。
38.具体地，服务器可以获取当前固体状态，基于预设的固体的位置更新算法(位置移动仿真算法)进行仿真，得到中间固体状态。其中，预设的位置移动仿真算法可以使用固体状态对应的控制方程进行求解，由于位置移动仿真是在不考虑压强下进行的，因此该控制方程为去除压强有关的方程项的方程。
39.以下对位置移动仿真的原理进行解释：一般而言，弹性固体的控制方程为公式(1)以及公式(2)，其中脚标s代表固体，ρ、u、p、e、σ、g分别代指密度、速度、压强、焓、应力张量、重力加速度，是物理仿真中的常用物理量，其表达的字母选取与本发明无关，代表空间梯度算符，代表散度算符。为时间离散符，其中f指代任意物理状态量。
[0040][0041][0042]
而在进行固体的物理仿真时，可以去除其中的压强项，则步骤s204可以是指求解方程(3)以及(4)，即忽略了公式(2)中的压强项，因此求解得到的是一个中间状态，后续需要进行修复。
[0043][0044][0045]
步骤s206，基于当前气体状态对目标气体对象进行对流仿真，得到中间气体状态。
[0046]
对流仿真与压强无关，固体的影响可以忽略或者将其作为已知的边界条件，对流仿真考虑的是气体速度对于质量、动量以及能量的输运，因此称为对流步骤。
[0047]
具体地，服务器可以基于预先设置的算法进行对流仿真，得到中间气体状态，其中，气体可以离散至多个气体网格中，即可以将气体看成是通过多个气体网格的气体单元组成的，可以获取每个气体网格对应的气体状态。其中，这里的“网格”代表的是一个单元，其形状可以根据需要设置，例如是可以正方形也可以是圆形。气体网格可以是结构网格，结构网格在二维的形状为矩形或在三维的形状为长方体(三维)，且网格的边与坐标轴方向平
行。即通过步骤s206可以将气体状态从当前时间步推进至一个待修复的中间状态。预先设置的算法可以包括有限差分格式与eno(essentially nonoscillatory schemes，基本无振荡)/weno(加权基本无振荡，weighted essentially non-oscillatory schemes)重构方法。
[0048]
在一些实施例中，可配合“多时间步积分方法”例如龙格库塔方法来提升中间气体状态和中间固体状态的精度。
[0049]
在一些实施例中，由于对流仿真与压强无关，因此在进行对流仿真之后，服务器还可以对压强进行更新，得到中间压强状态，加入到中间气体状态中，其中，中间压强状态用本构方程更新，例如在步骤s206后，利用气体的已知的状态方程得到中间态的压强p
*
。
[0050]
以下对对流仿真的原理进行解释：一般而言，可压缩气体(用脚标f代指)的控制方程为公式(5)、(6)以及(7)。
[0051][0052][0053][0054]
而在进行气体的对流仿真时，可以去除其中的压强项，则步骤s206可以是指求解方程(8)到(10)，即忽略了公式(6)以及(7)中的压强项，因此求解得到的是一个中间状态，后续需要进行修复。
[0055][0056][0057][0058]
其中，公式(8)～(10)的另一形式可以表示如公式(11)～(13)。其中上标“*”代表中间态，“n”代表当前时间步。例如代表气体密度的中间状态。ρfuf称为动量，ρ
fef
称为能量，则(ρ
fef
)n表示当前时间步的能量。
[0059][0060][0061][0062]
在对流仿真求解结束后，可以基于得到气体的中间状态以及气体的状态方程得到中间态的压强p
*
。
[0063]
可以理解，固体和气体满足本身物理模型所给出的本构方程。本构方程取决于具体材料。上述气体的控制方程为守恒格式，根据具体的本构模型，其亦可转换成另一种等价的非守恒格式。
[0064]
步骤s208，基于中间气体状态、中间固体状态以及交界约束条件进行耦合仿真，得
到目标固体对象在下一时间步的目标固体状态以及目标气体对象在下一时间步的目标气体状态。
[0065]
其中，耦合仿真是指一起进行仿真，即目标固体状态以及目标气体状态是一起解出来的，并不是分离的。在进行耦合仿真时，考虑了压强的影响，因此考虑了步骤s204以及s206中不考虑压强的作用带来的错误，因此能够对中间气体状态以及中间固体状态进行修复。
[0066]
交界约束条件为目标固体对象与目标气体对象之间的交界对应的约束条件，交界约束条件包括交界压强约束条件。例如，在气体和固体的交界面处，气体由于不能穿透固体，故气体与固体在交界处的压强大小一致,即交界压强约束条件包括气体与固体在交界处的压强大小一致。交界约束条件还可以包括交界速度约束条件，交界速度约束条件包括交界接触面目标固体对象与目标气体对象法向速度一致。
[0067]
具体地，在进行耦合仿真时，中间气体状态、中间固体状态可以被看出是一种待修复的状态，服务器可以基于与压强有关的气体控制方程、与压强有关的固体控制方程、中间气体状态、中间固体状态以及交界约束条件一起求解，得到满足约束条件的下一时间步的目标固体状态以及目标气体对象在下一时间步的目标气体状态。可以理解，得到下一时间步的目标固体状态以及目标气体对象在下一时间步的目标气体状态后，则下一时间步更新为当前时间步，可以循环执行本技术步骤s202～s208，直至所有时间步都得到仿真得到的状态。
[0068]
在一些实施例中，基于当前固体状态对目标固体对象进行位置移动仿真，得到中间固体状态包括：在拉格朗日视角下，基于当前固体状态对目标固体对象进行位置号移动仿真，得到以拉格朗日视角记录的中间固体状态。基于当前气体状态对目标气体对象进行对流仿真，得到中间气体状态包括：将以拉格朗日视角记录的中间固体状态投射到目标气体对象对应的气体网格中，并转化为欧拉视角记录的中间固体状态；基于欧拉视角记录的中间固体状态，以及欧拉视角记录的离散到气体网格中的当前气体状态，对目标气体对象进行对流仿真，得到欧拉视角记录的中间气体状态。
[0069]
欧拉视角追踪空间中固定的采样点的状态随时间的变化，拉格朗日视角则随着代表物质的点运动，并追踪这些点携带的状态变化。气体会充满空间中尚未被占据的地方，因此适合使用欧拉视角，而固体可能会在空间中大范围运动，因此适合拉格朗日视角。而本技术实施例中，在进行位置移动仿真时，是在拉格朗日视角下进行的，从而使得对固体的物理仿真更加准确。
[0070]
投射也可以称为投影。投影是转换视角的过程，视角转化的过程可以理解为将记录在粒子信息所代表的离散物理场1，重新分配到记录在附近网格格点(也即附近的另一格点)上。转换后格点记录的信息代表了一个新的离散物理场2，这两个离散的物理场其实都是对同一个真实物理场的近似，因此他们之间越接近投影的效果就越好。记投影前的点为起始点ps，即将投影到的点为目标点pd，投影的任意物理量为f，则有公式(14)如下，这个公式的意思是目标点的物理量的值为附近所有起始点的值的加权w
s,d
求和σs。w权重被称为通道函数，起始点为记录离散的固体信息的粒子，目标点为网格格点。
[0071]
f(pd)＝σ
sws,d
f(ps)
ꢀꢀ
(14)
[0072]
通过在固体位置移动以后，将记录于拉格朗日视角下的相关物理信息，包括但不
限于：质量，速度，体积，压强，投射到与气体同一套的气体网格上，因此能够方便地在同一套网格上构建耦合系统。因此步骤s208可以是在欧拉视角下进行的。当耦合仿真结束后，可以将目标固体对象的目标固体状态从欧拉视角转换为拉格朗日视角。
[0073]
本技术实施例中，对于可压缩气体的气体仿真，采用欧拉视角进行仿真。而弹性固仿真采纳混合视角：即固体的位置移动使用拉格朗日视角，状态数据记录在代表固体质量点的粒子上，与气体的耦合步骤使用欧拉视角，并在气体网格上进行，固体的状态数据记录在和气体同一套网格，避免重新绘制与重构耦合的网格，两视角之间采用投射技术来保证准确性，且耦合步骤采用欧拉视角并且在一套气体网格上进行仿真，从而使得仿真数学系统的构建更加方便，且能够减少网格重构的时间。
[0074]
在一些实施例中，在移动固体位置后，即使不考虑气体的压强影响，其速度及其形变梯度分布通常不符合本构方程限制，因此可以采纳速度及其形变梯度分布通进行修复。修复算法可以包括但不限于牛顿法或者拟牛顿法等，其中该修复既可以在拉格朗日视角下修复，也可以在投射到欧拉视角后，即投射到气体网格后，在对应的网格上修复。当修复成功后，再执行步骤s204。
[0075]
上述物理仿真方法中，由于可以先分别对固体对象以及气体对象分别进行仿真，得到中间的状态信息，然后再基于中间气体状态、中间固体状态以及交界约束条件进行耦合仿真，即一起进行求解，以对中间的状态进行修复，因此能对多物理现象更为准确的描述，使得得到的目标固体状态以及目标气体状态更加准确，提高物理仿真的准确性。
[0076]
以下对耦合仿真的原理进行解释，根据步骤s204～s208可知，在对气体和固体相互作用的条件下进行物理仿真时，气体的物理仿真采用了时间分离技术，分为压强无关的步骤(步骤s206)与压强有关步骤(步骤s208)，固体的物理仿真采用了时间分离技术，分为压强无关的步骤(步骤s204)与压强有关步骤(步骤s208)。
[0077]
而在气体和固体的交界面处，气体由于不能穿透固体，即两者在接触面法向(法向向量用w代指)速度一致，而且在此处的压强大小也一致，因此可以添加两种交界约束条件，如公式(15)的速度约束条件以及公式(16)的压强约束条件，在步骤s208进行耦合仿真时，可以获取与压强有关的控制方程以及交界约束条件进行耦合仿真。
[0078]
(u
f-us)
·
w＝0
ꢀꢀ
(15)
[0079]
p
f-ps＝0
ꢀꢀ
(16)
[0080]
其中，气体的时间分离第二步为压强相关的步骤，描述如公式(17)以及(18)
[0081][0082][0083]
其中公式(18)能量方程，可以等价变化转换为一个压强方程如公式(19)。此转化步骤的是为了能够直接利用耦合的压强约束条件。其中sf为刚度放缩矩阵，取决于气体的本构方程。
[0084][0085]
而基于与气体相似的原理，弹性固体的时间分离第二步的控制方程(20)以及(21)
与上述固体的时间分离控制方程(17)以及(18)形式一致，因此亦可以利用耦合的压强约束条件。因此，公式(15)至(21)构成了完整的强耦合系统，可以进行联立求解。可以理解，此处的强耦合系统同时包含了压强与速度未知量，其离散系统在数学中一般被称作kkt(karush-kuhn-tucker)系统，可以直接进行求解。也可以将变换为只包含压强未知量的spd系统或者pd(正定)系统，先求解压强，基于压强求解速度，以提高求解效率。其中，kkt系统、spd系统以及pd系统为代数系统。
[0086][0087][0088]
通过上述分析，步骤s208采用了强耦合的方式进行求解。所谓“强耦合”是指将各个物理现象及它们间的联系包含在一个完整的数学方程组中，一起求解。强耦合是对多物理现象更为准确的描述。时间分离指的是将一个同时涉及多个相互影响项的时间步分解为，逐次求解单个或部分影响项的子步骤。时间分离与完全分离求解不同，完全分离求解的分解步骤是依据物理现象，每一次求解结束后输出的即为该物理现象在下一时间步的结果。而时间分离中，每一个子步骤的结果会作为下一个子步骤的输入，待全部子步骤解完后得到下一个时间的结果。单物理系统使用时间分离技术可以降低求解难度。由于可压缩气体流动过程中，传递的信息是流动与声波的叠加，表现在数学方程中则是压强无关项与压强项，如不将两者分离，则时间步长被限制为网格尺寸除以声速与流速之和，远小于不可压缩流动的时间步长。而通过时间分离技术，对可压缩气体部分采用时间分离，将流动与声波传递分离，能够使得时间步长不再受声速限制，能够提高时间步长。
[0089]
此外，对于连续的数学方程组，物理状态量以及用于其的操作算符是在空间中连续，因此可以先进行离散，然后再进行求解，再被进行求解。离散包括空间离散和时间离散。对于离散，物理状态量的离散取决于网格的种类，以及将这些物理量放置于何处，放置位置包括但不限于网格的格点，网格的中心等。空间离散的操作算符分别代表空间梯度和散度算符，其离散除了网格类型，还取决于离散系统的方法及方法的特殊格式。离散系统的方法包括但不限于有限差分法、有限体积法或者有限元法。特殊格式包括但不限于不同的阶数等。对于时间离散，时间离散作用于操作算符其中f指代上述提到的任意物理状态量，可以将操作算符其转化为离散的时间推进算符其中δf代指两个离散量之间的差值，对于时间分离第一步，物理量的差值为待修复的中间状态减去当前时间步；对于时间分离第二步，物理量的差值为本时间步结束后的值减去待修复的中间状态；δt则为时间步长。
[0090]
对于将kkt系统转化为spd系统来提高求解效率的问题，记上标n、*、n 1分别为当前时间步的、待修复的中间态、下一时间步；记g,d分别为离散梯度算符与离散散度算符，根据数学通识d＝-g
t
，即离散散度算符对应的矩阵为负的离散梯度算符的转置，则控制方程都变换为成如下离散形式。其中，下标f，s同时存在是指可以适用于气体和固体。
[0091]
[0092][0093]
将上述两方程中的公式(22)代入公式(23)，得到只包含压强未知量的公式(24，即转换为pd系统或者spd系统。
[0094][0095]
因为在公式(24)只包含压强变量，在求解结束后，可以根据压强来修复得到速度等变量。气体与固体的速度变量可依照公式(25)修复。其中,u代表动量，m为质量运算符，m-1
等于密度的倒数。下标f以及s同时存在表示可以适用于气体以及固体。
[0096][0097]
此后，气体的能量项可以按照公式(26)修复，以便下一时间步中的时间分离第一步。
[0098][0099]
以下结合图3对如何得到只包含压强未知量的公式以及如何得到最终的修复完成的目标固体状态以及目标气体状态进行解释。
[0100]
在一些实施例中，如图3，基于中间气体状态、中间固体状态以及交界约束条件进行耦合仿真，得到目标固体对象在下一时间步的目标固体状态以及目标气体对象在下一时间步的目标气体状态包括：
[0101]
步骤s302，将中间气体状态作为气体压强约束方程的已知变量，得到气体约束条件；其中，气体压强约束方程为中间气体状态与目标气体压强之间的约束方程。
[0102]
其中，气体压强约束方程是中间气体状态与目标气体压强之间的约束方程，约束了中间气体状态与目标气体压强之间的对应关系。可以根据气体的压强推进方程以及动量修复方程得到，例如可以是将气体对应的压强推进方程代入气体动量修复方程得到。
[0103]
气体压强约束方程的原理描述如下：由于在进行位置仿真以及对流仿真时，没有考虑压强，因此中间状态的物理量需要修复，需将对流步骤中没有考虑的压强引入其中，其中，修复物理量时，以气体为例，动量的数学描述可以描述为公式(27)，能量修复的数学描述可以描述为公式(28)
[0104][0105][0106]
利用压强修复气体的动量与能量的方程中，因为需要n 1步的压强，在求解前尚未得知，因此引入可压缩气体所满足的压强推进方程，表示如公式(27)。将气体对应的压强推进方程代入气体动量修复方程(29)，可以得到气体压强约束方程，如公式(30)，称为扩散方程。公式(30)中的符号为拉普拉斯算子，有刚性系数sf＝(ρfc2)n，ρ为密度，c为声速。即为对一个物理场先做梯度操作，再做散度操作。上文已经指出梯度和散度算子可以被离散化，因此拉普拉斯算子即为连续进行两个离散操作，也是离散操作。如果记离散梯度算子为g，离散散度算子为d，那么上述扩散方程变成离散形式，可以表示为公式(31)。
[0107][0108][0109][0110]
具体地，可以将中间气体状态代入公式(31)，例如公式(31)中的以及为中间压强状态以及中间速度状态。从而使得公式(31)中的未知变量为下一时间步(n 1时间步)的压强
[0111]
步骤s304，将中间固体状态作为固体压强约束方程的已知变量，得到固体约束条件；其中，固体压强约束方程为中间固体状态与目标固体压强之间的约束方程。
[0112]
与步骤s302原理相似，可以将记录在网格上的固体的中间速度状态记录为以及与中间压强则也可推导出与固体压强约束公式(32)以及离散形式的公式(33)，这里不再赘述，其中固体的刚性系数ss与固体的材料模型有关。
[0113][0114][0115]
具体地，可以将中间固体状态代入公式(33)，例如公式(33)中的以及为中间压强状态以及中间速度状态。从而使得公式(33)中的未知变量为下一时间步(n 1时间步)的压强
[0116]
在一些实施例中，对于固体约束条件，也可以转化为只包含压强未知量的公式(34)以及只包含压强未知量的公式(35)，其中公式(34)为交界速度约束条件，公式(35)为交界压强约束条件。其中“m”代表质量运算符。
[0117][0118]
p
f-ps＝0
ꢀꢀ
(35)
[0119]
步骤s306，根据气体约束条件、固体约束条件以及交界约束条件进行耦合仿真，得到目标固体对象在下一时间步的目标固体状态以及目标气体对象在下一时间步的目标气体状态。
[0120]
具体地，可以根据气体约束条件、固体约束条件以及交界约束条件进行联立求解，得到目标固体对象在下一时间步的目标固体状态以及目标气体对象在下一时间步的目标气体状态。例如，联合公式(31)、(33)、(34)以及(35)进行求解，得到目标固体对象对应的目标固体压强以及目标气体对象对应的目标气体压强。
[0121]
在一些实施例中，例如，可以根据气体约束条件、固体约束条件以及交界约束条件进行耦合仿真，得到目标固体对象对应的目标固体压强以及目标气体对象对应的目标气体
压强。根据目标固体压强以及固体动量修复方程修复得到目标固体对象对应的目标固体速度，根据目标气体压强以及气体动量修复方程修复得到目标气体对象对应的目标气体速度。
[0122]
例如，以气体为例，求解pd系统后，n 1时间步固体与气体的压强已经取得，因此对于气体，可以利用动量修复方程(36)来得到n 1时刻的动量与速度。
[0123][0124]
此后再利用能量修复方程(37)得到n 1时刻的能量。
[0125][0126]
对于固体，还可以根据压强对形变梯度进行修复，通过根据固体投影在网格上的物理量进行速度以及形变梯度的修复，能够修复固体本身的应力不平衡。
[0127]
本技术实施例中，通过只包含压强的pd系统，求解压强之后再进行速度的求解，相较于直接将气体侧的压强与固体的速度场耦合的传统办法所得的kkt系统，能够提高求解效率。
[0128]
在一些实施例中，在求解时，固体与气体各个待求解的速度，压强自由度可以被标记出来，标记是指即将空间中的位置编号通过一个索引转化为一个平铺的排号的过程，固体与气体的交界处也可以被标记以便施加离散的约束条件。上述界面的标记既可以沿着网格的边缘进行，也可以使用网格切割/合并技术来提高精度。网格切割/合并技术本质上是通过对结构化网格内部进行建模来提升精度，本质依然是结构网格。
[0129]
在一些实施例中，由于耦合仿真是在气体网格上进行的，在修复结束后，固体在网格上的速度需投影回粒子，即从欧拉视角回到拉格朗日视角，以便继续下一时间步时固体的位置移动仿真。此处的透射方法为将拉格朗日视角的固体状态投射到欧拉视角的固体状态的反过程。
[0130]
在应用本技术实施例提供的方法进行仿真时，可以在仿真软件中输入以下初始条件数据，即仿真系统的初始输入可以包括：气体与固体的初始状态条件、边界条件、气体与固体的本构模型、计算时间步长的限制条件，时间分离方法中气体的压强无关与压强相关项、固体的位置更新算法、固体的自身形变修复算法、强耦合系统的构造与离散算法、离散系统的求解器。仿真过程中重复每一个时间步的仿真步骤，直到到达设定的终止时间，得到仿真结果。在仿真的每一个时间步，数据的输入可以包括：气体当前的气体状态例如密度、速度以及压强，固体当前的固体状态例如速度以及形变梯度，仿真系统的输出可以包括：当前时间步长、气体下一个时间步的状态例如密度、速度以及压强，固体下一个时间步的状态例如速度以及形变梯度。
[0131]
以下对时间离散和空间离散进行一步解释，在已经将物理量离散为每个离散的时间步下，空间中限个数据点的前提下，控制方程中依然含有连续的：时间算子和空间梯度算子空间散度算子它们需要被转化为离散的，计算机可以处理的算子。对于任意的时间步，记时间步长为δt，则方程中的时间导数们转化为离散的差分算子记上一步结束时的物理量带有角标n，即将求解的物理量带有角标n 1，则具体而言有，中间态与目标
态的关系可以表示如公式(38)～(42)
[0132][0133][0134][0135][0136][0137]
方程左边第一项即为差分的时间导数，计算机可以直接处理。空间算子取决于空间是被什么样的网格所填充，也取决于离散算法本身，空间梯度算子如公式(43)，其中，记其作用的任意物理量为f。
[0138][0139]
空间散度算子表示如公式(44)，其中散度作用于矢量，即在三个空间维度上各有一个值的量，例如速度为矢量，记其作用的任意物理矢量为f＝(f
x
,fy,fz)
[0140][0141]
通过上述离散处理，记录在离散空间中的物理量可以被映射到计算机内存空间中，因此空间离散算子就变成了内存中数据的加减乘除基本操作。而且方程中的续算子都被转化为离散形式，物理量也记录在离散的时间步下的，离散的空间数据点上，因此计算机可以直接处理。
[0142]
应该理解的是，虽然上述流程图中的各个步骤按照箭头的指示依次显示，但是这些步骤并不是必然按照箭头指示的顺序依次执行。除非本文中有明确的说明，这些步骤的执行并没有严格的顺序限制，这些步骤可以以其它的顺序执行。而且，上述流程图中的至少一部分步骤可以包括多个步骤或者多个阶段，这些步骤或者阶段并不必然是在同一时刻执行完成，而是可以在不同的时刻执行，这些步骤或者阶段的执行顺序也不必然是依次进行，而是可以与其它步骤或者其它步骤中的步骤或者阶段的至少一部分轮流或者交替地执行。
[0143]
在一些实施例中，如图4所示，提供了一种物理仿真装置，包括：
[0144]
当前状态获取模块402，用于在针对目标固体对象与目标气体对象的物理仿真过程中，获取目标固体对象在当前时间步的当前固体状态以及目标气体对象在当前时间步的当前气体状态；
[0145]
位置移动仿真模块404，用于基于当前固体状态对目标固体对象进行位置移动仿真，得到中间固体状态；
[0146]
对流仿真模块406，用于基于当前气体状态对目标气体对象进行对流仿真，得到中间气体状态；
[0147]
耦合仿真模块408，用于基于中间气体状态、中间固体状态以及交界约束条件进行耦合仿真，得到目标固体对象在下一时间步的目标固体状态以及目标气体对象在下一时间
步的目标气体状态；
[0148]
其中交界约束条件为目标固体对象与目标气体对象之间的交界对应的约束条件，交界约束条件包括交界压强约束条件。
[0149]
在一些实施例中，耦合仿真模块用于：将中间气体状态作为气体压强约束方程的已知变量，得到气体约束条件；其中，气体压强约束方程为中间气体状态与目标气体压强之间的约束方程；将中间固体状态作为固体压强约束方程的已知变量，得到固体约束条件；其中，固体压强约束方程为中间固体状态与目标固体压强之间的约束方程；根据气体约束条件、固体约束条件以及交界约束条件进行耦合仿真，得到目标固体对象在下一时间步的目标固体状态以及目标气体对象在下一时间步的目标气体状态。
[0150]
在一些实施例中，气体压强约束方程为将气体对应的压强推进方程代入气体动量修复方程得到，固体压强约束方程为将固体对应的压强推进方程代入固体动量修复方程得到。
[0151]
在一些实施例中，目标固体状态包括目标固体对象对应的目标固体速度，目标气体状态包括目标固体气体对应的目标气体压强；耦合仿真模块用于：根据气体约束条件、固体约束条件以及交界约束条件进行耦合仿真，得到目标固体对象对应的目标固体压强以及目标气体对象对应的目标气体压强；根据目标固体压强以及固体动量修复方程修复得到目标固体对象对应的目标固体速度；根据目标气体压强以及气体动量修复方程修复得到目标气体对象对应的目标气体速度。
[0152]
在一些实施例中，气体压强约束方程为：固体压强约束方程为：其中，δt为时间步长，sf代表气体对应的刚性系数，代表下一时间步气体的密度，为拉普拉斯算子，代表下一时间步气体的目标气体压强，代表气体的中间压强状态；代表空间散度算子，代表气体的中间速度状态，ss代表固体对应的刚性系数，代表下一时间步固体的目标固体密度，代表下一时间步固体的目标固体压强，代表固体的中间压强状态；代表固体的中间速度状态。
[0153]
在一些实施例中，交界约束条件包括交界压强约束条件以及交界速度约束条件，交界压强约束条件包括交界处固体压强与气体压强一致，交界速度约束条件包括交界接触面目标固体对象与目标气体对象法向速度一致。
[0154]
在一些实施例中，位置移动仿真模块用于：在拉格朗日视角下，基于当前固体状态对目标固体对象进行位置号移动仿真，得到以拉格朗日视角记录的中间固体状态；对流仿真模块用于：将以拉格朗日视角记录的中间固体状态投射到目标气体对象对应的气体网格中，并转化为欧拉视角记录的中间固体状态；基于欧拉视角记录的中间固体状态，以及欧拉视角记录的离散到气体网格中的当前气体状态，对目标气体对象进行对流仿真，得到欧拉视角记录的中间气体状态。
[0155]
关于物理仿真装置的具体限定可以参见上文中对于物理仿真方法的限定，在此不再赘述。上述物理仿真装置中的各个模块可全部或部分通过软件、硬件及其组合来实现。上
述各模块可以硬件形式内嵌于或独立于计算机设备中的处理器中，也可以以软件形式存储于计算机设备中的存储器中，以便于处理器调用执行以上各个模块对应的操作。
[0156]
在一些实施例中，提供了一种计算机设备，该计算机设备可以是服务器，其内部结构图可以如图5所示。该计算机设备包括通过系统总线连接的处理器、存储器和网络接口。其中，该计算机设备的处理器用于提供计算和控制能力。该计算机设备的存储器包括非易失性存储介质、内存储器。该非易失性存储介质存储有操作系统、计算机程序和数据库。该内存储器为非易失性存储介质中的操作系统和计算机程序的运行提供环境。该计算机设备的数据库用于存储物理仿真结果。该计算机设备的网络接口用于与外部的终端通过网络连接通信。该计算机程序被处理器执行时以实现一种物理仿真方法。
[0157]
在一些实施例中，提供了一种计算机设备，该计算机设备可以是终端，其内部结构图可以如图6所示。该计算机设备包括通过系统总线连接的处理器、存储器、通信接口、显示屏和输入装置。其中，该计算机设备的处理器用于提供计算和控制能力。该计算机设备的存储器包括非易失性存储介质、内存储器。该非易失性存储介质存储有操作系统和计算机程序。该内存储器为非易失性存储介质中的操作系统和计算机程序的运行提供环境。该计算机设备的通信接口用于与外部的终端进行有线或无线方式的通信，无线方式可通过wifi、运营商网络、nfc(近场通信)或其他技术实现。该计算机程序被处理器执行时以实现一种物理仿真方法。该计算机设备的显示屏可以是液晶显示屏或者电子墨水显示屏，该计算机设备的输入装置可以是显示屏上覆盖的触摸层，也可以是计算机设备外壳上设置的按键、轨迹球或触控板，还可以是外接的键盘、触控板或鼠标等。
[0158]
本领域技术人员可以理解，图5以及6中示出的结构，仅仅是与本技术方案相关的部分结构的框图，并不构成对本技术方案所应用于其上的计算机设备的限定，具体的计算机设备可以包括比图中所示更多或更少的部件，或者组合某些部件，或者具有不同的部件布置。
[0159]
在一些实施例中，提供了一种计算机设备，包括存储器和处理器，存储器中存储有计算机程序，该处理器执行计算机程序时实现以下步骤：在针对目标固体对象与目标气体对象的物理仿真过程中，获取目标固体对象在当前时间步的当前固体状态以及目标气体对象在当前时间步的当前气体状态；基于当前固体状态对目标固体对象进行位置移动仿真，得到中间固体状态；基于当前气体状态对目标气体对象进行对流仿真，得到中间气体状态；基于中间气体状态、中间固体状态以及交界约束条件进行耦合仿真，得到目标固体对象在下一时间步的目标固体状态以及目标气体对象在下一时间步的目标气体状态；其中交界约束条件为目标固体对象与目标气体对象之间的交界对应的约束条件，交界约束条件包括交界压强约束条件。
[0160]
在一些实施例中，提供了一种计算机可读存储介质，其上存储有计算机程序，计算机程序被处理器执行时实现以下步骤：在针对目标固体对象与目标气体对象的物理仿真过程中，获取目标固体对象在当前时间步的当前固体状态以及目标气体对象在当前时间步的当前气体状态；基于当前固体状态对目标固体对象进行位置移动仿真，得到中间固体状态；基于当前气体状态对目标气体对象进行对流仿真，得到中间气体状态；基于中间气体状态、中间固体状态以及交界约束条件进行耦合仿真，得到目标固体对象在下一时间步的目标固体状态以及目标气体对象在下一时间步的目标气体状态；其中交界约束条件为目标固体对
象与目标气体对象之间的交界对应的约束条件，交界约束条件包括交界压强约束条件。
[0161]
在一些实施例中，提供了一种计算机程序产品或计算机程序，该计算机程序产品或计算机程序包括计算机指令，该计算机指令存储在计算机可读存储介质中。计算机设备的处理器从计算机可读存储介质读取该计算机指令，处理器执行该计算机指令，使得该计算机设备执行上述各方法实施例中的步骤。
[0162]
本领域普通技术人员可以理解实现上述实施例方法中的全部或部分流程，是可以通过计算机程序来指令相关的硬件来完成，所述的计算机程序可存储于一非易失性计算机可读取存储介质中，该计算机程序在执行时，可包括如上述各方法的实施例的流程。其中，本技术所提供的各实施例中所使用的对存储器、存储、数据库或其它介质的任何引用，均可包括非易失性和易失性存储器中的至少一种。非易失性存储器可包括只读存储器(read-only memory，rom)、磁带、软盘、闪存或光存储器等。易失性存储器可包括随机存取存储器(random access memory，ram)或外部高速缓冲存储器。作为说明而非局限，ram可以是多种形式，比如静态随机存取存储器(static random access memory，sram)或动态随机存取存储器(dynamic random access memory，dram)等。
[0163]
以上实施例的各技术特征可以进行任意的组合，为使描述简洁，未对上述实施例中的各个技术特征所有可能的组合都进行描述，然而，只要这些技术特征的组合不存在矛盾，都应当认为是本说明书记载的范围。
[0164]
以上所述实施例仅表达了本技术的几种实施方式，其描述较为具体和详细，但并不能因此而理解为对发明专利范围的限制。应当指出的是，对于本领域的普通技术人员来说，在不脱离本技术构思的前提下，还可以做出若干变形和改进，这些都属于本技术的保护范围。因此，本技术专利的保护范围应以所附权利要求为准。