CN111488670B - A nonlinear simulation method of mass spring soft tissue deformation - Google Patents

Abstract

The invention discloses a nonlinear mass spring soft tissue deformation simulation method, comprising the steps of: collecting medical image data of soft tissue, segmenting organ images, and using Delaunay triangular mesh to obtain three-dimensional surface mesh data; On the basis of the subdivision, the mesh nodes are used as mass points, and the edges connected between the nodes are spring components composed of structural springs and damping springs in parallel to form a mass-spring soft tissue model. Whether there is a collision between the surgical instrument and the mass-spring soft tissue model , if there is a collision, calculate the force of each particle, update the position and velocity of the particle, calculate the deformation value of the surface particle of the particle spring soft tissue model, until the surgical instrument and the particle spring soft tissue model do not collide; then output the simulation results and carry out End the entire process after visual rendering. The invention can more accurately simulate the deformation of soft tissue in real time.

Description

A nonlinear simulation method of mass spring soft tissue deformation

technical field

The invention relates to the technical field of soft tissue deformation modeling, in particular to a nonlinear mass spring soft tissue deformation simulation method.

Background technique

Physical motion simulation of organs and tissues is a central issue in the fields of anatomy teaching, surgical simulation, and surgical training. Human tissues and organs generally have soft tissue characteristics, usually manifested in material properties such as inhomogeneity, anisotropy, quasi-incompressibility, nonlinearity, plasticity, and viscoelasticity. Accurate and real-time simulation of such soft tissue deformation processes has always been very important. Challenging research topic. The mass-spring model is one of the most widely used models in human soft tissue modeling. Most of the existing interactive soft-tissue deformation simulation models based on the mass-spring model simplify the stress-strain constitutive relationship of soft tissue into a linear relationship. Its computational complexity is low and it is easy to discretize time and then iteratively calculate, but its computational accuracy is not high.

Due to the need of real-time simulation calculation, the existing methods mostly use linear elasticity theory to describe the soft tissue deformation, that is, it is assumed that there is a linear hysteresis relationship between the stress component and the infinitesimal strain. The advantage of using a linear elastic model is that the global stiffness matrix is constant throughout the simulation and can be pre-computed, inverted, and stored, improving computational efficiency. In fact, biological soft tissues are nonlinear materials that can only be approximated as linear materials with small deformations. In addition, large deformation, large rotation, and large displacement in soft tissue deformation will produce geometric nonlinearity, which cannot be solved by classical linear elasticity theory.

In a word, the existing interactive simulation system based on the mass-spring model can improve the computational efficiency by simplifying the deformation of the soft tissue into a linear elastic model in order to ensure that the real-time solution models are all simplified. The assumption of using linear elasticity as the base model, while reducing the amount of computation at runtime, limits the accuracy with which physical materials can be modeled. More importantly, the linear elastic model can only approximate the soft tissue deformation under the condition of small displacement and small deformation. However, biological soft tissue has material nonlinearity and geometric nonlinearity in the case of large deformation.

SUMMARY OF THE INVENTION

The purpose of the present invention is to provide a nonlinear particle spring soft tissue deformation simulation method for the technical defects existing in the prior art. , that is, the Euler elastic energy is used to describe the spring deformation in the mass-spring model, and then the spring length can be calculated to obtain a nonlinear stress-strain relationship, which is then used to characterize the nonlinear constitutive relationship of soft tissue.

The technical scheme adopted for realizing the purpose of the present invention is:

A nonlinear mass spring soft tissue deformation simulation method, comprising the steps of:

Collect medical image data of soft tissue, segment the organ image, and obtain 3D surface mesh data by Delaunay triangular mesh division;

On the basis of 3D mesh division, the mesh nodes are used as mass points, and the edges connected between nodes are spring components composed of structural springs and damping springs in parallel to form a mass-spring soft tissue model;

Detect whether there is a collision between the surgical instrument and the mass spring soft tissue model. If there is a collision, calculate the force of each mass point, update the mass position and mass velocity, and calculate the deformation value of the mass point on the surface of the mass spring soft tissue model until the surgical instrument and the mass spring soft tissue model. until there is no collision between them;

The whole process ends after outputting the simulation results and rendering them visually.

The present invention takes into account the nonlinearity of biological soft tissue, and can simulate the deformation of soft tissue in real time more accurately. The nonlinear mass spring soft tissue deformation simulation method of the present invention is closer to the real soft tissue biomechanical characteristics, and can effectively ensure the simulation accuracy and computational complexity.

Description of drawings

Fig. 1 is a flow chart of a nonlinear mass spring soft tissue deformation simulation method;

Fig. 2 is a visualization result of liver triangulation meshing;

Figure 3 is a schematic diagram of a spring assembly between two mass points;

Fig. 4 is the schematic diagram of a neighborhood of a certain particle;

Figure 5 is a schematic diagram of the parametric angle for calculating the total area of all triangular patches in a certain particle-neighborhood;

Fig. 6 is the visualization result after visual rendering and the schematic diagram of motion simulation (the upper layer is in a static state, and the lower layer is in a motion state);

Figures 7-8 show the effects of sphere and liver parameter a on stress-strain nonlinearity:

Figures 9-10 are graphs of the effects of sphere and liver parameter b on stress-strain nonlinearity, respectively.

Detailed ways

The present invention will be further described in detail below with reference to the accompanying drawings and specific embodiments. It should be understood that the specific embodiments described herein are only used to explain the present invention, but not to limit the present invention.

In the present invention, Euler elasticity is used to describe the spring deformation in the particle spring model to characterize the nonlinearity of soft tissue, that is, the Euler elastic energy is used to describe the spring deformation variable in the particle spring model, and then the calculated spring length can be obtained The nonlinear stress-strain relationship is satisfied, and then it is used to characterize the nonlinear constitutive relationship of soft tissue.

Mathematician Euler proposed the nonlinear constitutive relation of elastic energy and elastic line according to the principle of minimum potential energy when he studied the steady-state problem of elastic soft rod under the action of external force. Since then, mathematician Mumford introduced Euler elastic energy to the field of computer vision. The core idea is to describe the deformation process of cork sticks with two geometric quantities of length and curvature at the same time, that is, consider the following energy:

E(l)=∫ _l (a+bκ ² (x))ds,

where a, b are two positive parameters, κ is the curvature of curve l at position x, and ds is the arc length.

As shown in Figure 1, the nonlinear mass spring soft tissue deformation simulation method of the present invention comprises the following steps:

1. 3D mesh division

Collect medical image data such as CT/MR of soft tissue, segment the organs to be studied, and use Delaunay triangular mesh to obtain 3D surface mesh data. Taking the liver as an example, Fig. 2 shows the visualization result of triangulation of the liver. The liver surface is formed by the combination of points, edges and triangular patches.

On the basis of 3D meshing, the mesh nodes are used as mass points, and the edges connected between nodes are spring components composed of structural springs and damping springs in parallel, as shown in Figure 3, to form the initial soft tissue deformation model. Figure 3 is a schematic diagram of the spring assembly, x _i , x _j are two connected mass points.

2. Collision detection

Detect whether there is a collision between the surgical instrument and the soft tissue. If there is a collision, it is considered that the soft tissue is deformed by the force and then executes the next step. Otherwise, the simulation result is output to end the whole process.

3. Traverse all the particles and calculate the force on the particles

Assuming that there are a total of n particles, for any particle x _i , i=1,...,n, its mass is mi ∈ R and the resultant force f _{i ∈} R ³ acting on this point. Particles on the surface of soft tissue satisfy the following kinetic equation, which is Newton's second law equation:

是位移关于时间的二阶导数，即质点的加速度。where M is a 3n × 3n diagonal mass matrix,

is the second derivative of the displacement with respect to time, the acceleration of the particle.

1. Calculate the external force on the particle

The external force on the particle includes the gravity of the particle and the external force exerted by man, namely

F _i ^ext =m _i g+ ^fu ,

Among them, g=9.8N/kg is the standard gravitational acceleration, and f ^u is the artificially applied external force.

2. Calculate the damping force on the particle

为：Since energy dissipation occurs during deformation, the viscous force is represented by the spring damping force, the damping force between the mass points x _i and x _j

for:

弹簧的阻尼常数，v_i表示质点i的速度。in,

The damping constant of the spring, v _i represents the velocity of particle i.

3. Calculate the spring force on the mass

弹簧的弹簧系数，l⁰和l′分别代表弹簧的原始长度和形变后的弹簧长度。in,

The spring coefficients of the spring, l ⁰ and l' represent the original length of the spring and the deformed spring length, respectively.

The present invention uses Euler elastic energy to propose a new formula for calculating the spring length l _(i,j) between two mass points x _i , x _j , as follows:

的平均值来估计，即where κ _(i,j) is the average curvature of the particles x _i , x _j

to estimate the average value of

Unlike the original Euler elastic energy, in the present invention, the curvature on the curved surface is used instead of the curvature of the curve for calculation, which simplifies the calculation complexity. There are many ways to estimate vertex normals and curvatures for triangular meshes that describe smooth surfaces.

The present invention uses the following method to calculate the particle normal vector and the mean curvature.

Estimate the particle normal vector:

The area composed of the triangular patches composed of the particle _xi and the adjacent particles is regarded as the 'one neighborhood' of _xi , as shown in Figure 4, the 'one neighborhood' with _xi as the vertex in the triangular mesh model.

表示质点

的法向量，

为质点x_i周围的三角面片集合，

为质点x_i周围的三角面片法向量的集合，|x_i|是x_i周围相邻质点的总数。质点x_i的法向量用为其周围三角面片法向量的加权和来估计：make

Represents a particle

the normal vector of ,

is the set of triangular patches around the particle _xi ,

is the set of triangular patch normal vectors around the particle _xi , and | _xi | is the total number of adjacent particles around _xi . The normal vector of a particle x _i is estimated as a weighted sum of the normals of its surrounding triangular patches:

是x_i到x_j的边矢量，

为叉积。in,

is the edge vector from x _i to x _j ,

is the cross product.

Estimate particle mean curvature:

In the present invention, the discretization method based on the mean curvature manifold proposed by Matllieu-Desbrun is used to estimate the mean curvature of the particle _xi

为任意质点x_i关于其坐标(x,y,z)的梯度算子，A为质点x_i一邻域内所有三角面片的总面积，α_j,β_j的定义如图5所示，α_j,β_j分别为两个夹角，为两个相邻三角面片的对角。in,

is the gradient operator of any particle x _i with respect to its coordinates (x, y, z), A is the total area of all triangular patches in the neighborhood of the particle x _i , α _j , β _j are defined as shown in Figure 5, α _j , β _j are two included angles respectively, which are the opposite angles of two adjacent triangular patches.

4. Calculate the deformation value of the soft tissue surface particles

According to the above 1, 2, 3, the resultant force of the surrounding adjacent particles on the particle is calculated, and the deformation value of the particle on the surface of the soft tissue is obtained by using the dynamic equation of the particle on the surface of the soft tissue.

Among them, the dynamic equation of the soft tissue surface particles is:

质点x_i所受合力为f_i＝f_i ^ext-f_i ^s-f_i ^d.in,

The resultant force on the particle x _i is f _i =f _i ^ext -f _i ^s -f _i ^d .

f _i ^d is the resultant damping force on the mass point _xi , and f _i ^s is the resultant spring force on the mass point _xi ;

The particle position and velocity are updated according to the following explicit framework Verlet integration, and the deformation state of the soft tissue surface particle is obtained, that is, the deformation value.

Assign x _i (t+Δt) to x _i (t), and v _i (t+Δt) to vi ₍ t) to determine whether there is a collision.

In the formula, x _i (t+Δt) represents the position of the particle at time t+Δt, v _i (t+Δt) represents the velocity of the particle at time t+Δt, and _xi (t) and vi ₍ t) represent respectively The position and velocity of the particle at time t.

5. Visual rendering

Use OpenGL to visually render the real-time updated data, and give the visualization results and motion simulation diagrams, taking the liver as an example, as shown in Figure 6.

6. System stress-strain analysis

Taking a sphere and a liver as examples, Figures 7-8 and 9-10 show the nonlinear relationship between force and displacement for different parameters. For different a and b, the deformation relationship between the spherical surface and the liver under external force was tested and compared. It can be seen from the results that as a and b increase, the nonlinearity of the system becomes stronger.

It should be noted that, for the nonlinearity of soft tissue, the spring coefficient in the deformation process can also be set according to the stress-strain relationship of the soft tissue, but the premise of this method is to obtain the real stress-strain relationship of the soft tissue through experimental tests. This is complex and costly, especially for the simulation of human soft tissue.

Compared with the prior art, the present invention has the following advantages:

1. The nonlinear strain of the spring is realized by introducing Euler elasticity, which improves the limitation of the existing mass-spring system in simulating soft tissue motion, and solves the linear relationship between the spring elastic force and the spring deformation variable in the traditional mass-spring model, that is, stress-strain. The problem that it cannot meet the requirements of the actual simulation of human tissues and organs.

2. Since the spring in the mass spring system is elastic, its deformation is not only related to the length of the spring, but also to the bending degree of the spring. The curvature is a geometric quantity that describes the bending degree of the spring. Therefore, the present invention proposes a method for describing the spring deformation by using the spring length and curvature, which can more accurately simulate the soft tissue deformation process.

In the present invention, Euler elastic energy is used to describe the spring deformation, that is, the deformation of the spring is controlled by the spring length and the curvature at the same time. Since the curvature is a nonlinear function of the length, the elastic force calculated by Euler elasticity and the spring length change amount has a nonlinear relationship.

The invention estimates the average curvature on the curved surface in a neighborhood of the surface grid, which can effectively reduce the calculation cost, is suitable for parallel computing, and can ensure real-time simulation in real simulation.

The mass-spring model with nonlinear stress-strain relationship proposed in the present invention can be used to simulate the deformation of soft tissues with material properties such as quasi-incompressibility, nonlinearity and viscoelasticity, such as liver and kidney.

The above are only the preferred embodiments of the present invention. It should be noted that, for those skilled in the art, without departing from the principles of the present invention, several improvements and modifications can be made. These improvements and Retouching should also be regarded as the protection scope of the present invention.

Claims (1)

1. a nonlinear mass spring soft tissue deformation simulation method, is characterized in that, comprises the steps: Collect medical image data of soft tissue, segment the organ image, and obtain 3D surface mesh data by Delaunay triangular mesh division; On the basis of 3D mesh division, the mesh nodes are used as mass points, and the edges connected between nodes are spring components composed of structural springs and damping springs in parallel to form a mass-spring soft tissue model; Detect whether there is a collision between the surgical instrument and the mass spring soft tissue model. If there is a collision, calculate the force of each mass point, update the mass position and mass velocity, and calculate the deformation value of the mass point on the surface of the mass spring soft tissue model until the surgical instrument and the mass spring soft tissue model. until there is no collision between them; End the whole process after outputting the simulation results and visual rendering; Update the particle position and particle velocity according to the following explicit frame Verlet integration, and obtain the deformation value of the particle _xi on the surface of the particle spring soft tissue model;

f _i =fi ^ext _-fi ^s _-fi ^d _, f _i ^ext =m _i g+ ^fu ,

In the formula, m _i is the mass of the particle x _i , f _i is the resultant force on the particle x _i ; f _i ^ext is the external force on the particle x _i , including the particle gravity m _i g and the artificially applied external force f ^u ; f _i ^d is The damping resultant force on the particle _{xi, f i s} ^is the spring resultant force on the particle _xi ;

is the damping force between particles x _i and x _j ,

The damping coefficient of the spring, vi _represents the velocity of the particle _xi ;

is the spring force between the particles x _i and x _j ,

The spring coefficient of the spring, l ⁰ and l' represent the original length of the spring and the deformed spring length, respectively; The formula for the spring length l _(i,j) between two mass points x _i , x _j is as follows:

where a, b are two positive parameters, and the curvature κ _{(i, j)} is the average curvature of the particles x _i , x _j

to estimate the average; Use the Matllieu-Desbrun discretization method based on the mean curvature manifold to estimate the mean curvature of the particle x _i

in,

is the gradient operator of any particle x _i with respect to its coordinates (x, y, z), A is the total area of all triangular patches in the neighborhood of the particle x _i , α _j , β _j are the two included angles respectively, which are the particle points The diagonal of two adjacent triangular patches in the neighborhood of x _i ;

in,

is the edge vector from x _i to x _j ,

is the cross product, N _i is the estimated normal vector of the particle _xi ,

is the normal vector of the triangular patch around the particle _xi , and | _xi | is the total number of adjacent particles around _xi .

Images