Soft tissue modeling method in virtual operation based on biomechanical characteristics

Abstract

A soft tissue modeling method in a virtual operation based on biomechanical characteristics comprises the following steps: 1) setting a mass point spring system, and establishing a kinetic equation for each mass point; 2) moving least square method approximation of a displacement field; 3) simulation updating: calculating a new displacement field by using MLS approximation, and ensuring physical constraints; 4) solving a kinetic equation, and solving the equation by using an implicit Euler method to obtain system response changing along with time; and (5) iteration and adjustment are carried out, an implicit Euler method is used for solving the equation, and system response changing along with time is obtained. According to the method, the biomechanical characteristics of the soft tissue are effectively simulated by adding multiple functional springs, then the MLS method is integrated into the mass point spring model, the capacity of the model for processing complex deformation and large displacement events is enhanced, meanwhile, the overall numerical stability and accuracy are improved, and the method is suitable for large-scale popularization and application. And the established soft tissue model is more suitable for the requirements of a virtual operation system.

Description

Technical Field

[0001] The invention belongs to the technical field of virtual reality medical treatment, and in particular relates to a soft tissue modeling method in a virtual surgery simulation system. Background Art

[0002] At present, soft tissue models can be divided into geometric models and physical models according to whether the physical properties of organ soft tissue are considered.

[0003] The geometric modeling method only needs to consider the basic geometric features of soft tissues, and has a low computational cost. It is the most commonly used method in the early days of soft tissue modeling, such as surface modeling and three-dimensional modeling based on medical image data. Based on the surface modeling method, researchers further invented the free-form deformation method (FFD). Sederberg (Sederberg TW, Parry SR. Free-form deformation of solid geometric models [C]. Proceedings of the 13th annual conference on Computer graphics and interactive techniques. 1986: 151-160.) optimized the free-form deformation method based on the ternary Bernstein polynomial and proposed a method for free deformation on solid geometric models. This method uses the concepts of control points and control bodies. By adjusting the shape and position of the control body, the deformation operation of the solid model can be achieved.

[0004] Compared with geometric modeling, physical modeling is a model established based on the biomechanical properties of soft tissues and the physical properties of real soft tissues. Therefore, the physical model can accurately reflect the deformation process of soft tissues and can more realistically reproduce the morphological characteristics of soft tissues during surgery. Currently, the commonly used physical modeling methods mainly include Mass-Spring Model (MSM), Finite Element Method (FEM), Mesh-Free Methods (MFM), etc.

[0005] At present, there are still some problems with these soft tissue modeling methods, which can be mainly divided into two categories. The first category is the real-time problem of soft tissue modeling and simulation. When using traditional high-precision finite element models, the calculation cost is high and the real-time performance of the system is relatively low. For example, Dick (Dick C, Georgii J, Westermann RA real-time multigrid finitehexahedra method for elasticity simulation using CUDA[J]. Simulation Modelling Practice and Theory, 2011, 19(2):801-816.) et al. proposed a multigrid method for real-time simulation of elastic deformable objects on Nvidia's CUDA GPU architecture. This method is based on finite element discretization of deformable objects using hexahedrons. In order to accurately simulate large deformations, it considers the co-rotating strain formula.

[0006] The second category is the authenticity of the soft tissue model. When the mass spring model with a lower computational cost is used, it cannot accurately reflect the real human organ soft tissue due to its unique biomechanical properties. For example, Mollemans et al. (Mollemans W, Schutyser F, Van Cleynenbreugel J, et al. Tetrahedral mass spring model for fast soft tissue deformation [C]. International Symposium on Surgery Simulation and Soft Tissue Modeling. Berlin, Heidelberg: Springer Berlin Heidelberg, 2003: 145-154) et al. proposed a new extended mass spring model based on tetrahedral geometry using tetrahedral structure and local mesh deformation method.

[0007] Regarding the meshless methods mentioned above, there are also people such as De (De S, Kim J, Lim YJ, et al. The point collocation-based method of finite spheres (PCMFS) for real time surgery simulation [J]. Computers & structures, 2005, 83 (17-18): 1515-1525.) who developed an innovative meshless method called the Point Collocation-based Method of Finite Spheres (PCMFS), which is specially designed for real-time simulation of soft tissue in the context of minimally invasive surgery.

[0008] None of the above methods can guarantee both authenticity and real-time performance of simulation. Therefore, it is necessary to study a suitable method to improve the authenticity and real-time performance of soft tissue modeling and simulation. Summary of the invention

[0009] The purpose of the present invention is to provide a soft tissue modeling method in virtual surgery based on biomechanical characteristics, so as to improve the stability and authenticity of organ soft tissue.

[0010] The present invention is achieved through the following technical solutions.

[0011] The soft tissue modeling method in virtual surgery based on biomechanical characteristics of the present invention comprises the following steps:

[0012] Step 1, mass-spring system setting. At the beginning of the simulation, a network consisting of mass points and springs is first set up to represent the physical structure of the simulation. The mass points represent discrete nodes, while the springs represent the connections and interaction forces between these nodes. According to Newton's second law, a dynamic equation is established for each mass point. The displacement, velocity, and acceleration of each mass point are determined by its mass and the force generated by the spring connected to it.

[0013] As attached Figure 1 As shown, in view of the biomechanical properties of soft tissue, the present invention adds three functional springs to simulate the nonlinearity, viscoelasticity and incompressibility of soft tissue respectively.

[0014] 1. Nonlinear spring

[0015] For the nonlinear characteristics of soft tissue, that is, the variation trend between strain and stress, there are two stages: the nonlinear stage with smaller external force and the linear stage with larger external force. When the displacement of a particle under the action of external force is large, the elastic force on the particle should be significantly enhanced to pull the particle back to its original position, and along with the linear increase of strain, the stress also shows a linear increase trend. i and a j , the calculation formula of the nonlinear spring force is as follows:

[0016]

[0017] Among them, k ij represents point a i and a j The elastic coefficient of the nonlinear spring between, k1 and k2 are the nonlinear strength coefficient and linear strength coefficient of the spring respectively, and represents the spring displacement vector between the mass points, represents the initial length of the spring, Δl a It is the standard for dividing the linear force and nonlinear force of the spring into two stages. From the above formula, it can be seen that when the displacement of the spring is Δl ij When Δl is small, the change in spring length and force is a third-order polynomial relationship. ij Gradually increase to greater than Δl a When , it will be transformed into a linear relationship, which is consistent with the nonlinear biomechanical characteristics of soft tissue.

[0018] 2. Virtual Spring

[0019] For the simulation of incompressibility, that is, the deformation of soft tissue will be kept within a certain range, it is not only necessary to ensure that the surface nodes of the model are not subject to too much change, but also to consider the external deformation caused by the interaction between the internal particles of the model. When an external force is applied to the soft tissue to cause it to deform, the particles will pull the virtual spring, thereby generating a virtual spring force in the opposite direction of the displacement to avoid excessive changes in the tetrahedral unit. Therefore, the model of the present invention can effectively achieve constant volume and shape recovery after the external force disappears. The calculation formula of the virtual spring force is as follows:

[0020] F V =k3Δs i γ i

[0021] Where k3 represents the elastic coefficient of the virtual spring, Δs i represents the deformation of the virtual spring, γ i Represents the unit direction vector in the dynamic length direction.

[0022] When external force is applied to soft tissue to cause it to deform, the mass point will pull the virtual spring, thereby generating a virtual spring force in the opposite direction of the displacement to avoid excessive changes in the tetrahedron unit. This can avoid excessive deformation of the tetrahedron and simulate the incompressibility of soft tissue. The calculation formula of the virtual volume force can be further obtained:

[0023]

[0024] Among them, Z i is the particle a i The set of connected surrounding mass points, k3 is the virtual spring coefficient, V j and V j 0 These are the initial volume and the changed volume of the tetrahedral unit used in the present invention.

[0025] 3. Damping spring

[0026] The damping spring is mainly set up to meet the viscoelastic requirements of soft tissues, and its performance characteristics are creep and relaxation. i and a j The damping force formula of the damping spring between is as follows:

[0027]

[0028] In the above formula, d0 and d1 are two different constants. For particle a i The initial position of Represents point a i and a j The relative speed between them.

[0029] Step 2: Moving Least Squares (MLS) approximation of the displacement field. The MLS method is used to calculate the displacement of each mass point and its surrounding points, generating a continuous displacement field at the node position. When the mass point is displaced under the action of external forces, MLS helps determine its displacement trajectory to ensure the continuity and accuracy of the overall model. By defining a weight function, MLS can minimize the difference between the actual displacement of the mass point and the predicted displacement field.

[0030] Step 3, simulation update. As the simulation time progresses, the position, velocity, and acceleration of each particle are continuously updated. The MLS approximation is used to calculate the new displacement field and ensure physical constraints.

[0031] Step 4, solve the dynamic equation. Use the implicit Euler method to solve the equation and obtain the system response that changes with time. According to the implicit solution method in the Euler method, the differential equation is discretized by approximating the difference quotient between the antecedent and the latter term as the derivative. Assuming that the implicit integration method is used, the final expression formula for the displacement of the particle force is:

[0032] v(t+Δt)=v(t)+Δtm -1 F(s(t)

[0033] s(t+Δt)=s(t)+Δtv(t+Δt)

[0034] Among them, s(t+Δt) and s(t) are the positions of the particle at t+Δt and t respectively, v(t+Δt) and v(t) are the movement speed of the particle at t+Δt and t+Δt respectively, and Δt is the iteration step per unit time. The formula can be used to calculate the new coordinates of each particle in space in each unit time step, thereby simulating the phenomenon of soft tissue deformation.

[0035] Step 5, iteration and adjustment. After each particle position update, a weighted average calculation is performed on each particle according to the weight value set for each particle to adjust the weight function and other parameters of the MLS to provide the best displacement field approximation under the current configuration.

[0036] Compared with the prior art, the beneficial effects of the present invention are: first, by adding a variety of functional springs, the biomechanical properties of soft tissue are effectively simulated; and then, by integrating the MLS method into the particle-spring model, the model's ability to handle complex deformation and large displacement events is enhanced, while improving the overall numerical stability and accuracy, making the established soft tissue model more suitable for the requirements of the virtual surgery system. BRIEF DESCRIPTION OF THE DRAWINGS

[0037] Figure 1 It is a schematic structural diagram of the improved mass spring of the present invention.

[0038] Figure 2 Schematic diagram of liver soft tissue established in the present invention.

[0039] Figure 3 This is the extrusion test of the model of the present invention, wherein (a) is the initial extrusion, (b) is a small degree of extrusion, and (c) is a large degree of extrusion.

[0040] Figure 4 It is a flow chart of the present invention. DETAILED DESCRIPTION

[0041] The present invention will be further described through the following embodiments in conjunction with the accompanying drawings so that those skilled in the art can implement the invention with reference to the description.

[0042] A smoke simulation method for virtual surgery according to an embodiment of the present invention comprises the following steps:

[0043] Step 1, mass-spring system setting. At the beginning of the simulation, a network consisting of mass points and springs is first set up to represent the physical structure of the simulation. The mass points represent discrete nodes, while the springs represent the connections and interaction forces between these nodes. According to Newton's second law, a dynamic equation is established for each mass point. The displacement, velocity, and acceleration of each mass point are determined by its mass and the force generated by the spring connected to it.

[0044] As attached Figure 1 As shown, in view of the biomechanical properties of soft tissue, the present invention adds three functional springs to simulate the nonlinearity, viscoelasticity and incompressibility of soft tissue respectively.

[0045] 1. Nonlinear spring

[0046] For the nonlinear characteristics of soft tissue, that is, the variation trend between strain and stress should have two stages: the nonlinear stage with smaller external force and the linear stage with larger external force. When the displacement of a particle under the action of external force is large, the elastic force on the particle should be significantly enhanced to pull the particle back to its original position, and along with the linear increase of strain, the stress also shows a linear increase trend. i and a j , the calculation formula of the nonlinear spring force is as follows:

[0047]

[0048]

[0049] Among them, k ij represents point a i and a j The elastic coefficient of the nonlinear spring between , k1 and k2 are the nonlinear strength coefficient and linear strength coefficient of the spring respectively, k1 and k2 are set to 0.2 and 0.85 respectively. represents the spring displacement vector between the mass points, represents the initial length of the spring, Δl a It is the standard for dividing the linear force and nonlinear force of the spring into two stages. From the above formula, it can be seen that when the displacement of the spring is Δl ij When Δl is small, the change in spring length and force is a third-order polynomial relationship. ij Gradually increase to greater than Δl a When Δl a Set to 2mm.

[0050] 2. Virtual Spring

[0051] For the simulation of incompressibility, that is, the deformation of soft tissue will be kept within a certain range, it is not only necessary to ensure that the surface nodes of the model are not subject to too much change, but also to consider the external deformation caused by the interaction between the internal particles of the model. When an external force is applied to the soft tissue to cause it to deform, the particles will pull the virtual spring, thereby generating a virtual spring force in the opposite direction of the displacement to avoid excessive changes in the tetrahedral unit. Therefore, the model of the present invention can effectively achieve constant volume and shape recovery after the external force disappears. The calculation formula of the virtual spring force is as follows:

[0052] F V =k3Δs i γ i

[0053] Where k3 represents the elastic coefficient of the virtual spring, which is set to 0.04 here; Δs i represents the deformation of the virtual spring, and γi represents the unit direction vector in the dynamic length direction.

[0054] When external force is applied to soft tissue to cause it to deform, the mass point will pull the virtual spring, thereby generating a virtual spring force in the opposite direction of the displacement to avoid excessive changes in the tetrahedron unit. This can avoid excessive deformation of the tetrahedron and simulate the incompressibility of soft tissue. The calculation formula of the virtual volume force can be further obtained:

[0055]

[0056] Among them, Z i is the particle a i The set of connected surrounding mass points, k3 is the virtual spring coefficient, V j and V j 0 These are the initial volume and the changed volume of the tetrahedral unit used in the present invention.

[0057] 3. Damping spring

[0058] The damping spring is mainly set up to meet the viscoelastic requirements of soft tissues, and its performance characteristics are creep and relaxation. i and a j The damping force formula of the damping spring between is as follows:

[0059]

[0060] In the formula, d0 and d1 are two different constants. For particle a i The initial position of Represents point a i and a j The relative speed between them.o and d1 are set to 0.05 and 0.04, respectively.

[0061] Step 2: MLS approximation of displacement field. The MLS method is used to calculate the displacement of each mass point and its surrounding points, generating a continuous displacement field at the node position. When the mass point is displaced under the action of external forces, MLS helps determine its displacement trajectory to ensure the continuity and accuracy of the overall model. By defining a weight function, MLS can minimize the difference between the actual displacement of the mass point and the predicted displacement field. The general weight function w I (x) has two commonly used forms, namely, Gaussian weight function and cubic spline weight function. The weight function selected in this embodiment is the cubic spline weight function, and its formula is:

[0062]

[0063] in, Represents the distance vector, d = xx I is the distance between the point to be fitted and the node, d max Indicates the influence range of the node, that is, the area of ​​action of each weight function. In the calculation process of the MLS method, in order to As the best approximation of the desired function y(x), different fitting coefficients b are selected i (x), and The weighted square sum K of the errors at all nodes is minimized to obtain:

[0064]

[0065] The basis functions Most of the time, monomials are chosen, such as in one-dimensional problems:

[0066] Linear basis: P T (x) = [1, x], m = 2

[0067] Quadratic basis: P T (x) = [1, x, x 2 ],m=3

[0068] Cubic basis: P T (x) = [1, x, x 2 , x 3 ],m=4

[0069] Step 3, simulation update. As the simulation time progresses, the position, velocity, and acceleration of each particle are continuously updated. The MLS approximation is used to calculate the new displacement field and ensure physical constraints.

[0070] Step 4, solve the dynamic equation. Use the implicit Euler method to solve the equation and obtain the system response that changes with time. According to the implicit solution method in the Euler method, the differential equation is discretized by approximating the difference quotient between the antecedent and the latter term as the derivative. Assuming that the implicit integration method is used, the final expression formula for the displacement of the particle force is:

[0071] v(t+Δt)=v(t)+Δtm -1 F(s(t)

[0072] s(t+Δt)=s(t)+Δtv(t+Δt)

[0073] Among them, s(t+Δt) and s(t) are the positions of the particle at t+Δt and t respectively. v(t+Δt) and v(t) are the movement speeds of the particle at t+Δt and t+Δt respectively. Δt is the iteration step per unit time. The formula can be used to calculate the new coordinates of each particle in space per unit time step, thereby simulating the phenomenon of soft tissue deformation.

[0074] Step 5, iteration and adjustment. After each particle position update, the weight function or other parameters of the MLS may need to be adjusted to provide the best displacement field approximation for the current configuration.

[0075] Table 1 shows the error comparison between the present invention (M-MSM) and the traditional MSM and FEM models.

[0076] Table 1

[0077]

[0078]

[0079] As attached Figure 2 The structure diagram of the liver soft tissue is shown in Figure 1, which is composed of 4384 triangular basic units, with a length, height and width of 41 cm, 34 cm and 34 cm respectively. The surgical instrument is composed of 3303 triangular basic units, with a length, height and width of 25 cm, 2 cm and 2 cm respectively.

[0080] Attached Figure 3 is the effect of the extrusion test of the established model, where Figure 3 (a) To start extrusion, Figure 3 (b) small degree of compression, Figure 3 (c) A greater degree of compression. It can be seen that the present invention can effectively simulate the details of real soft tissues with a high degree of realism. Figure 4 It is a flow chart of the present invention.

Claims

1. A soft tissue modeling method in virtual surgery based on biomechanical characteristics, characterized by The following steps are involved: Step 1, mass-spring system setting: At the beginning of the simulation, first set up a network consisting of mass points and springs to represent the physical structure of the simulation; the mass points represent discrete nodes, and the springs represent the connections and interaction forces between these nodes; according to Newton's second law, establish a dynamic equation for each mass point; the displacement, velocity and acceleration of each mass point are determined by its mass and the force generated by the spring connected to it; Three functional springs are used to simulate the nonlinearity, viscoelasticity and incompressibility of soft tissues: (1) Nonlinear spring For two adjacent mass points a in the model i and a j , the calculation formula of the nonlinear spring force is as follows: Among them, k ij represents point a i and a j The elastic coefficient of the nonlinear spring between, k1 and k2 are the nonlinear strength coefficient and linear strength coefficient of the spring respectively, and represents the spring displacement vector between the mass points, represents the initial length of the spring, Δl a It is the standard for dividing the linear force and nonlinear force of the spring into two stages; when the displacement of the spring is Δl ij When Δl is small, the change in spring length and force is a third-order polynomial relationship. ij Gradually increase to greater than Δl a When , it will be transformed into a linear relationship, which is consistent with the nonlinear biomechanical characteristics of soft tissue; (2) Virtual spring The calculation formula of virtual spring force is as follows: F V =k3Δs i ·c i Where k3 represents the elastic coefficient of the virtual spring, Δs i represents the deformation of the virtual spring, γ i Represents the unit direction vector in the dynamic length direction; When external force is applied to soft tissue to cause it to deform, the mass point will pull the virtual spring, thereby generating a virtual spring force in the opposite direction of the displacement to avoid excessive changes in the tetrahedron unit, thereby avoiding excessive deformation of the tetrahedron, simulating the incompressibility of soft tissue, and further deriving the calculation formula for the virtual volume force: Among them, Z i is the particle a i The set of connected surrounding mass points, k3 is the virtual spring coefficient, V j and V j 0 is the initial volume and the changed volume of the tetrahedral unit used; (3) Damping spring For mass point a i and a j The damping force formula of the damping spring between is as follows: In the formula, d0 and d1 are two different constants. For point a i The initial position of Represents point a i and a j The relative speed between Step 2, MLS approximation of displacement field: Use the MLS method to calculate the displacement of each mass point and its surrounding points, and generate a continuous displacement field at the node position; when the mass point is displaced under the action of external force, MLS helps determine its displacement trajectory to ensure the continuity and accuracy of the overall model; by defining a weight function, MLS can minimize the difference between the actual displacement of the mass point and the predicted displacement field; Step 3, simulation update: As the simulation time moves forward, the position, velocity and acceleration of each particle are continuously updated. The MLS approximation is used to calculate the new displacement field and ensure physical constraints. Step 4, solve the dynamic equation: According to the implicit solution method in the Euler method, the differential equation is discretized; assuming that the implicit integration method is used, the final expression formula for the force displacement of the particle is: v(t+Δt)=v(t)+Δtm -1 F(s(t) s(t+Δt)=s(t)+Δtv(t+Δt) Among them, s(t+Δt) and s(t) are the positions of the particle at t+Δt and t respectively, v(t+Δt) and v(t) are the movement speeds of the particle at t+Δt and t+Δt respectively, Δt is the iteration step per unit time, and the new coordinates of each particle in space per unit time step are calculated by the formula to simulate the phenomenon of soft tissue deformation; Step 5, iteration and adjustment: After each particle position update, a weighted average calculation is performed on each particle according to the weight value set for each particle to adjust the weight function and other parameters of the MLS to provide the best displacement field approximation under the current configuration.

Images