Particle-based simulation technique for medical applications

Abstract

This paper proposes a particle-based nonlinear elastic object simulation technique for virtual surgery. Particle-based techniques are used to model and simulate nonlinear elastic objects, such as the skin and internal organs. This enables the simulation to consider various factors, such as location, direction, and depth, when making incisions in the organs. However, the issue with this method is that it can only simulate precisely cut tissue during incision. Objects with elasticity, such as tissue, require the generation of complex debris during incision. This paper proposes a particle-based elastic object simulation technique to model the debris from torn tissue when making an incision in the organs. It can predict where the tissue will tear based on the maximum shear stress (MSS) theory and Tresca’s yield criterion when the force applied to the tissue exceeds the maximum stress. Newly generated particles at the predicted location are remeshed with nearby particles. We verified the superiority of our proposed method over traditional particle-based methods by accurately representing more complex debris and comparing the results of incisions in the same area of the body. This allows various incision types, such as stab wounds and lacerations, to be simulated.

KEYWORDS:

• Incision of the organs
• debris
• particle simulation

1. Introduction

The demand for healthcare services is increasing due to population aging and intensified competition, as well as psychological disorders. As a result, virtual reality (VR) technology is emerging as a solution for training medical professionals and providing alternative psychological treatments. VR-based applications were initially developed in the 1990s in the United States to treat post-traumatic stress disorder in combat veterans and are now widely used in the treatment of various phobias due to their affordability and convenience (Rizzo & Shilling, 2017). The range of applications has expanded from those for trauma therapy to those for surgery and rehabilitation. Virtual surgery is widely utilised in the medical field as an auxiliary tool for doctors (Shan et al., 2022), a learning tool for medical students, and a training or simulation tool for medical professionals. In order to provide virtual surgery that is as realistic as possible, technologies from multiple disciplines are utilised to accurately represent the environment and features of the surgery and handle the significant amount of data used (He & Tang, 2021). One important element in virtual surgery is the 3D human body model (Kim et al., 2021). To create a realistic surgical training experience, it is necessary to provide tactile feedback, resistance, viscosity, and other physical properties of the virtual body. Healthcare professionals, such as doctors and therapists, can use virtual environments and human body models to obtain useful information before performing direct surgeries on patients. Additionally, the realistic physical representation of soft tissue can be used in conjunction with haptic tactile feedback tools. Based on this, it can be utilised in medical applications such as surgical procedure planning, computer/robotic-assisted surgery, and palpation-based disease diagnosis.

FEM is a methodology for numerically solving differential equations (Baccouch, 2021). FEM can predict the state and changes in the force of elastic objects by combining configuration, equilibrium, and appropriate equations. However, complex nonlinear equations can be time-consuming to calculate. Although FEM produces physically accurate deformation results, it is unsuitable because it is too costly for applications requiring real-time calculations, such as VR-based applications.

Particle-based dynamics (PBD) simulates movement by considering the relative distances and angles between particles that make up an object to obtain the necessary information (Bender et al., 2014). Compared to FEM, this technology can more accurately represent the forces that particles experience, making it suitable for simulating elastic objects. However, if particles exceed their limits, the binding between particles may break, resulting in the mesh disappearing. To solve this issue, one of the closest edge particles is selected, and a division plane perpendicular to the edge direction is placed through that particle to divide it. All triangles on the division plane are assigned to the original particle, and all triangles below are assigned to the duplicate particle. This method allows the mesh to remain stable even under extreme conditions but does not accurately represent the debris.

When tissue damage occurs, various types of cells are activated at the site of the injury, leading to an inflammatory response. This inflammatory response is triggered by inflammatory cells that migrate to the site of injury through the bloodstream and blood vessels. These inflammatory cells interact with debris, such as tissue fragments, and initiate the process of new tissue formation and wound healing. The shape of tissue debris varies depending on the size and location of the wound, as well as the extent of tissue damage. Typically, debris can be small, round, or cylindrical and can be concentrated around the boundary of the wound. Moreover, debris in complex shapes can form near structures that receive more pressure at the site of the wound. These debris play an essential role in the process of generating new cells and tissues during wound healing.

We propose a method that can generate the shape of the debris resulting from surgical incisions. Specifically, this paper leverages stress theory, which can estimate the force applied to an object and its resulting effects. Conventional strength theory was used to predict the point where an object yield based on the calculations of the multidirectional forces and the shear stress applied to its surface. It was hypothesised that by calculating the angle and intensity of the force acting on the particles, stress theory could be utilised. To verify this, we utilised Tresca's yield criterion, which calculates the yield stress based on the maximum and minimum stress applied to the object rather than calculating the range of stress to which the object is subjected. The particle-based method determines the criterion for the force applied based on the distance between particles that comprise the object; thus, the maximum stress, minimum stress, and the resulting yield stress can also be determined based on the distance between particles that applied force. Tresca’s yield criterion was applied to the particle vectors of the maximum and minimum distances under stress to specify where new particles should be generated, and a random value ranging from −10 to 10% was assigned to create complex debris. We carried out tests using various percentages ranging from 3 to 50% to achieve results that resembled the debris generated during actual incision procedures. As a result, a 10% random range was deemed appropriate and applied for the experiment. The particles near the generated ones were remeshed. Through this, the shape of the debris can be realistically generated. This can be helpful in understanding the properties and potential impacts of tissues. Virtual surgery applications must realistically express the deformation of the tissue by tools. Also, since the deformation by interaction with the user is necessary, the cost must be reduced to make real-time calculation possible. The proposed algorithm makes it possible to provide a realistic real-time incision representation of the skin while satisfying these conditions. Real-time skin transformation representation allows users of virtual surgery to provide realistic feedback on using virtual surgical tools, further increasing the effect of learning and preparing for virtual surgery. Therefore, this algorithm can be used in virtual medical applications that require realistic tissue deformation representation.

This paper introduces related work on algorithms for simulating the elasticity of human tissues in the case of virtual surgery. It then describes a particle-based elastic object simulation algorithm based on Tresca’s yield criterion and an algorithm for generating particles and debris. Through comparative experiments using the proposed method and traditional PBD method, the effectiveness of the proposed method is demonstrated, and conclusions and future research directions are presented.

2. Related work

This section covers the research process of deformation techniques, including virtual surgery. Reznick et al. demonstrated that there was no difference in effectiveness between training with a bench station using actual surgery and human models and using a structured surgical evaluation method called the Objective Structured Assessment of Technical Skills (OSATS) (1997). Therefore, training systems to teach basic surgical techniques are being developed by simulating human models using VR technology. It is important to achieve real-time interactivity to enable the rapid repetition of perception and value judgment. To ensure that surgery is carried out in an accurate and systematic manner, it is necessary to be able to accurately track changes in elastic models, such as the human body. To create a realistic deformable elastic model, the mass-spring method, the linear elasticity-based FEM, and the PBD method are used. A mass-spring system is one way to implement a continuous physical model via points with mass and a virtual spring connecting them (Gibson & Mirtich, 1997). This method can be used to represent continuous models with some accuracy and speed unless scientifically accurate calculations are required. Mass-spring systems have been used to represent models with one-dimensional or two-dimensional structures, such as hair (Selle et al., 2008), cloth (Choi & Ko, 2002), and rigid bodies (Nealen et al., 2006), with a small degree of elasticity. However, the mass-spring system creates artificial anisotropy through mesh selection, making it difficult to describe the characteristics of soft tissue deformation, and it is also difficult to associate spring stability with material properties, such as Young’s modulus.

Linear elasticity-based FEM is a method of simulating a force-deformed object based on the assumption of small displacements (Bro-Nielsen & Cotin, 1996). Haptic feedback can be processed in real time by applying precomputation. However, there are two issues to consider when using this approach. The first occurs when a model that applies a precomputed response is used. Interactions that change the model’s topology, such as destruction, require changing the precomputed data structure, resulting in a problem of decreasing speed. The second is due to the substantial limitations of the precomputed response transformation. Linear elasticity-based FEM assumes that the deformation of an object is linear, but the deformation of an actual object becomes nonlinear as it increases. Thus, if the size of the deformation increases, or if nonlinear transformations such as rotation are applied to the deformation, this approach produces inaccurate results. Instead of updating the entire data structure as one way to solve the problem of using a model with a precomputed response, a method for changing the pre-inverted stiffness matrix that encompasses the geometric and physical specifications of an object was proposed (Lee et al., 2006). Instead of recalculating the entire matrix to improve efficiency, only some matrices that affect the deformation of the object must be recalculated in this method. Felippa proposed a solution called the co-rotational method to improve the results of deformation based on the second issue of linear elasticity-based methods based on the assumption of small displacements (2000). This method separates the dispositions of an object assumed to be a small strain into a rigid body and strained parts, and it calculates internal force using only the strained parts to calculate nonlinear variations (such as rotation) more accurately.

Another proposed method for representing a deformable model in real time is to use FEM to calculate variations based on an explicit integration scheme (Taylor et al., 2008). The advantage of this method is that it uses only mass matrices to calculate deformation. Mass matrices can be simplified into diagonal matrices through a process called mass lumping (Zienkiewicz et al., 2005). Then, each degree of freedom can be solved independently by decomposing the equations of motion. As a result, intuitive parallelisation becomes possible, and the speed in solving the equations increases (Comas et al., 2008). The explicit integration method is suitable for real-time nasal movements, such as brain deformation, because it can artificially increase mass to process different materials (Joldes et al., 2009).

Mozafary and Payvandy proposed an MSM system to simulate the draping of knitted fabrics (2018). Nedel and Thalmann proposed an angular spring-based MSM to simulate muscle deformation in real time (1998). Ren et al. built an intelligent simulation platform for intravenous surgery based on a mass-spring model (2017). The physical properties of the MSM are closely related to the attenuation coefficient and stiffness coefficient of the spring, so adjustment of the coefficients is very important. Bianchi et al. determined spring coefficients and mesh topologies based on linear elastic objects (2004). However, their scheme is not suitable for nonlinear soft tissue models because it can only be used when Poisson’s ratio is constant. San-Vicente et al. estimated the spring coefficient based on a cubic model tensile test (2012). Some researchers have introduced nonlinear behaviour in MSMs to improve accuracy. Omar et al. introduced conical springs into an MSM to provide nonlinear behaviour (2015). Li et al. proposed a new bending spring as a surface mass-spring model to correct the angle to improve the shape restoration performance and the accuracy of the proposed model (2018).

In order to improve the disadvantage that the calculation cost of the deformation algorithm is high, several studies using data-driven methods have also been conducted. Pezzementi et al. proposed training a 2D mass-spring system to behave similarly to a nonlinear finite element model to simulate the deformation due to the interaction between tool and tissue (2008). As a result, the 2D mass-spring system produced similar results to the FE model, which can be physically expressed more accurately. As deep learning models gradually develop, several methods have been proposed to maintain accuracy while reducing the high computational cost of FEM through deep learning networks. Mendizabal et al. proposed a deep learning-based technique based on the U-Net architecture to express the deformation of hyperelastic materials (2020). Lampen et al. developed a deep learning network based on the PointNet++ architecture that takes point cloud data and explicit boundary types as inputs instead of using FEM to predict how the face passively deforms with respect to bone movement in orthognathic surgery (2022).

A PBD method was proposed by Müller et al. to implement deformable objects in real time (2007). Most of the popular simulation methods used in computer graphics are force based, first calculating acceleration by accumulating internal and external forces based on Newton’s second law of motion. The speed of the object is obtained via the time integration of acceleration, and the position is finally calculated from this speed. PBD is a method of controlling the position of the constraint in addition to the process of obtaining speed in the force-based method. By directly controlling the position, it is possible to easily express an elastic object that is difficult to solve using the force-based method. Müller used a multi-grid method (2008) that assumes the original mesh to be the finest mesh before generating a coarser mesh (modified according to constraints) and interpolating the results to speed up the Gaussian–Seidel method used in error correction via the constraints in PBD. Bouaziz et al. proposed a method of using constructive dynamics to prevent the unintentional stiffening of targets as the number of iterations applying constraints increases (2014). To solve the problem where the simulation target becomes stiff (regardless of the stiffness parameter as the number of iterations from applying the constraints increases) and to calculate the constraint-force estimates that can be used in various applications, Macklin et al. proposed a method called XPBD, which modifies the multiplier used when changing positions via constraints (2016). Fratarcangeli et al. introduced a method to speed up Gauss-Seidel method to solve a sparse system of linear constraints by partitioning the set of equations governing the system and simultaneously solving the equations that are independent of each other (2016). Wang et al. proposed an accelerated algebraic multigrid framework for nonlinear cloth simulation using adaptive smoothing and demonstrated that it could be used with the Newton-Raphson method (2018).

Several methods have been studied to express the human body's soft tissue and elastic object using FEM or particle-based dynamics. It is necessary to consider the time required to calculate the displacement and the physical accuracy. Reducing computation time while maintaining accuracy has been proposed by using specific factors such as biomechanical characteristics of soft tissue or by combining a method with other simulation techniques. Kubiak et al. used PBD in a surgical simulation (2007). To express complex hairstyles, Rungjiratananon et al. used shape-matching approaches that implemented PBD (2010). Umetani et al. used a PBD model derived from Cosserat theory to simulate the complex bending and twisting of elastic rods (2015). Macklin and Müller used PBD to apply density constraints to particles to express fluids (2013). Segato et al. proposed a method of simulating the deformation of the brain white matter due to catheter insertion through PBD and using it for preoperative path planning and intraoperative guidance (2021).

In addition, a method of using modification of standard PBD distance constraints to express contraction and extension of fibres constituting muscle (Romeo et al., 2018), or a method of using tetrahedral mesh instead of triangular mesh in traditional PBD and adding volume constraints to express soft tissue deformation (Liu et al., 2020) have been proposed. Jayashudha and Kabadi proposed a method based on the Delaunay triangulated cube to implement soft tissue deformation and removal (2020). In this method, when the cube forming the soft tissue and the virtual object to be removed collide, the part where the epidermis is cut is calculated using the intersection points specified in advance. Xie et al. used an accelerated FEM method using a Kalman filter to represent soft tissue in surgical simulation (2020). Zhang et al. (2022) proposed a method combining FEM and step-variable fourth-order Runge-Kutta numerical calculation for the physical modelling of the lung. This method obtained higher single-step accuracy through higher-order equations and proposed using two descent methods to converge to accurate results quickly. Shi et al. expressed the soft tissue as a volume structure composed of two parts, a surface part and an internal part, implemented the cutting through FEM, and optimised the cutting process through a Bézier curve to smooth the cutting edge without increasing the computational load (2019). Tan et al. implemented organ incision by combining FEM and PBD (2021). The mesh representing the organ uses a hybrid mesh consisting of a coarse volume mesh and a fine surface mesh. After calculating the displacement of the coarse volume mesh with FEM, this value was put into one of the constraints of the PBD that calculates the displacement of the fine surface mesh to affect the displacement of the surface mesh and enable more accurate expression. Stomakhin et al. (2013) proposed the material point method (MPM) to observe the flow of particles within a grid, using equations based on the Cauchy stress tensor for an accurate representation of a force-applied object. Wolper et al. (2019) supplemented the traditional MPM with von Mises yield criterion and fracture theory to simulate the realistic process of an elastic object changing when it can no longer withstand the applied force. In August 2016, Dr. J. Bederson of Mount Sinai Hospital in New York successfully performed an arteriovenous malformation surgery using the surgery AR platform CaptiView developed jointly by Leica and Brainlab (Schwam et al., 2021). The surgeon did not need to turn his head because 2D/3D images of the patient’s brain taken before surgery were overlaid on the screen showing the microscopic image of the patient’s brain. Other medical information related to the surgery could also be called up on the screen, and a feature that automatically adjusted the focus by tracking where the surgeon was looking during the surgery was also included. Professor Makito Sasaki of the International University of Health and Welfare in Japan performed VR-assisted laparoscopic surgery and prostate cancer surgery (Yoshida et al., 2019).

3. A particle-based method for the simulation of the destruction of elastic objects

This section explains the method for modelling complex debris when an elastic object is destroyed. The point of debris formation is determined using Tresca's yield criterion, which considers the maximum and minimum forces exerted on the object. The concept of constraints used in traditional PBD is adopted to satisfy this condition, where forces are treated as constraints between particle distances or angles.

3.1. Elastic object simulation based on constraints using PBD

Most materials have both elasticity and plasticity, and which property is more dominant depends on the type of material, the magnitude of the external force, and the shape of the deformation. For example, a rubber band stretches can stretch easily and returns to its original shape. However, when affected by a strong external force, it breaks is cut immediately or is permanently deformed and does not return to its original shape.

In this paper, we introduce particle-based dynamic techniques to simulate the contraction and expansion of an elastic object subjected to external forces. Each particle was constrained to maintain the tension between itself and the other particles entering within a certain range connected to itself with a specific bond. (1) $C_{stretch}(p_1,p_2)=|p_1\text{--}{p}_2|\text{--}d$(1) In Equation (1), p₁ and p₂ are the positions of the two particles, and d is the distance that the two particles must maintain. Constraints were used to adjust the distance between particles that affected each other. Instead of indirect repositioning by force., PBD gave each particle a mass, position, and velocity and then simulated the particles by direct repositioning via constraints In PBD, each particle with mass, position and velocity is directly repositioned by constraints. (2) $\begin{array}{rl}{C}_{bend}(p_1,p_2,p_3,p_4)& =arccos((p_2\text{--}{p}_1)\times (p_3\text{--}{p}_1)/|(p_2\text{--}{p}_1)\times (p_3\text{--}{p}_1)|\cdot (p_2\text{--}{p}_1)\times (p_4\text{--}{p}_1)\\ & /|(p_2\text{--}{p}_1)\times (p_4\text{--}{p}_1)|)\text{--}{\varphi }_0\end{array}$(2) The bending constraint keeps the angle between consecutive triangles at a certain angle used the angles between successive triangles formed by the particles to keep the object from spreading beyond a certain angle. In Equation (2), (p₁, p₃, p₂) and (p₁, p₂, p₄) are triangles consisting of three particles each, and the two triangles are connected. The scalar ϕ₀ represents the dihedral angle between the two triangles. Angles that are not subject to excessive forces return to the dihedral angle ϕ₀ due to the reactive forces If the two triangles are subjected to a force that does not affect the elasticity, the angle between the two triangles returns to the dihedral angle ϕ₀ because of the reactive forces.

To accurately simulate elastic objects, stretch constraints that maintained a constant distance and bending constraints that maintained a constant bending angle were used, as shown in Figures 1 and 2. This also easily expressed the elasticity of an object and was utilised to allow the energised object to maintain its original shape The constraints also expressed an object's elasticity and were utilised to allow the object deformed under a force to maintain its original shape.

Figure 1. The constraint on the distance between particles.

Figure 2. The constraint on the angle between particles.

3.2. Representation of debris using Tresca’s yield criterion

In traditional particle-based dynamic simulations, the maintenance of interparticle bonds is judged by distances with characteristic thresholds. When the distance between the particles exceeds the threshold due to external force, the particles affected by external force are removed. As a result, an elastic object, such as tissue, is torn. However, real elastic objects leave debris when they tear above the threshold in real elastic objects torn above the threshold, debris remains. In traditional particle-based dynamic simulations, when the particle above the threshold is removed, all meshes formed by the particle do not accurately represent the object's surface after tearing. are not accurately represented. The representation of Ddebris not only gives a sense of reality during a simulation but also affects the results of elastic calculations.

This paper uses the yield criterion, which is a reference point where an object is deformed by an external force of a certain size or more. This paper uses the yield criterion, a reference point where an external force of a specific size or more permanently deforms an object. It was applied to predict the point at which debris from the destroyed elastic body object was generated. Typical hypotheses for yield criteria include maximum distortion energy theory and maximum shear stress theory (suitable for soft materials), maximum principal stress theory, and the Mohr–Coulomb theory (suitable for the destruction of brittle materials).

The maximum strain energy theory predicts that yielding occurs when the strain energy of an object is equal to or greater than the strain energy of the yield. g_e is the maximum stress, σ is the main stress, and τ is the shear stress in the maximum distortion energy yield criterion equation shown below: The maximum distortion energy yield criterion equation is shown in Equation (3). In Equation (3), g_e is the maximum stress, σ is the main stress, and τ is the shear stress. (3) $g_{{}_e}^{DE}=1/\surd 2[{({\sigma }_x\text{--}{\sigma }_y)}^2\text{--}{({\sigma }_y\text{--}{\sigma }_z)}^2\text{--}{({\sigma }_z\text{--}{\sigma }_x)}^2+6({{\tau }_{xy}}^2+{{\tau }_{yz}}^2+{{\tau }_{zx}}^2)](1/2)\ge 1$(3) The maximum shear stress theory argues that yielding begins when the shear stress an object receives is equal to or greater than the inherent yield criterion of the object. If σ₁ is the maximum main stress, σ₃ is the minimum main stress, and S_y is the yield strength. The maximum shear stress yield criterion is shown in Equation (4): (4) $g_e^{MSS}=({\sigma }_1\text{--}{\sigma }_3)/S_y\ge 1$(4) Our method uses Tresca’s yield criterion. Yielding occurs when the stress exerted on an object is equal to the maximum stress of the tensile stress test. When the bonds between the particles break under external force, it stores the vector applied to the particle that received the force. Then, the vector from which the particle is to be created based on Tresca’s yield criterion is inferred. Finally, particles are created using the inferred vectors, and nearby particles are remeshed to express the debris.

Figure 3 shows Tresca's yield criterion. It is related to principal stress, which means that stress occurs asmeaning stress occurs against an external force inside an object. It is expressed as half the difference between the maximum and minimum principal stress and can express the yield criterion of various objects (Cazacu et al., 2014). Tresca’s yield criterion in Equation (5) is expressed approximately: (5) ${\tau }_f=({\sigma }_1\text{--}{\sigma }_3)/2$(5) The yield criterion is half the difference between the maximum principal stress (σ₁) and the minimum principal stress (σ₃). When the pressure is rapidly increased by an external force When an external force rapidly increases the pressure, the elastic body object generates a reaction force to return to its original state. Loss occurs when the applied external force is greater than the inherent yield criterion of the object exceeds the object's inherent yield criterion. Tresca's yield criterion is useful in particle-based simulations because of its requirements. Traditional MSS theory requires calculations of calculating all forces applied to a unit area based on angles. However, particle-based simulations can be calculated based on physical quantities applied to each particle because the unit of an object is made up of particles. Tresca's yield criterion requires the largest and smallest principal stresses, so it which can be approximated by the particle with the strongest force and the particle with the smallest force.

Figure 3. An expression of stress and distance under Tresca’s yield criterion.

This paper assumed a threshold value for maintaining the connection between each particle as the maximum stress. It and defined the moment when it failed to withstand it, causing a defect in the yield criterion of the elastic object. It was assumed that the state with the minimum principal stress was a stabilised state, maintaining incompressibility. (6) $(|k{\rho }_{max}\text{--}k{\rho }_{min}|)/2={\tau }_f$(6) Equation (6) is the yield criterion equation for elasticity, in which ρ_max is the change in the distance of the particles under the maximum principal stress, and ρ_min is the change in the distance of the particles under the minimum principal stress. The elasticity coefficient k is the value that represents the elasticity of a specific body part. The particles are simulated based on PBD (Müller et al., 2007). Each particle has mass and velocity, and moves under the influence of gravity and external forces. In order to represent nonlinear motion, iterative techniques including the Newton-Rhapson method are used. When an external force is applied, the distance and angle between particles change. In order to maintain the values specified by constraints, the positions of the particles that are subject to the force are rearranged. Therefore, the distances and angles between particles are the results of the force, and this which can be used as a criterion for yield.

This paper proposes an algorithm to explain the changing of shapes via elasticity when objects are subjected to external forces. First, the distance between the particle P_i receiving the force and other particles connected through the edge E_j is calculated and compared with the initial value. If the distance between particles is greater than the maximum shear strength, then brittleness deformation is applied according to Tresca’s yield criterion to create particles half the length of the existing threshold and to reconstruct the mesh based on the results. When the distance between particles exceeds the maximum shear force, brittleness deformation is applied according to Tresca's yield criterion. New particles are generated at half the length of the threshold, and the mesh is reconstructed according to the result. The proposed algorithm is as follows: (Figure 4).

Figure 4. Pseudo code for the proposed algorithm.

In the simulated model, when the bond between particles is broken by external forces and the triangle primitive information between the connected particles is deleted, a new particle is created at the position predicted by Tresca's yield criterion, with a 10% random value applied to the position. Particles with broken bonds are bonded to the newly created particles. Finally, the distance constraint value for the newly created primitive is reconfigured to express the debris of the elastic bodyobject based on Tresca’s yield criterion.

4. Results

This experiment was conducted on a computer with an i5-8500 CPU and a GTX 3090 TI GPU. A 3D leg model with approximately 150,000 particles was used. Figure 5 shows the results using the proposed method proposed method results, in which particles forming an elastic body object increased to the threshold or more. There are various cases in which an object changes shape under force by forces the object receives. When pulling particles that constitute the flesh, both the traditional PBD method and the proposed method show identical results up to the maximum stress. The applied force is reflected in the particle as the same stress. When the applied force is removed from the particle, it moves in the opposite vector without exceeding the threshold, the particle moves in the opposite direction of the removed force due to the stress, and the particle bounces multiple times due to constraints and damping.

Figure 5. Simulation of the removal of particles that make up the skin on a human leg. In the traditional PBD method, the removed particles and the area of the mesh formed by the particles disappear without leaving any debris (a) and (c). When removing particles in the same region, the proposed method creates a new particle at the predicted position based on Tresca's yield criterion and remeshes it with nearby particles (b) and (d). (a) The result of removing the tissue on the leg in a rectangular shape using the traditional PBD method, zoomed out (left) and zoomed in (right). (b) The result of removing the tissue on the leg in a rectangular shape using the proposed method, zoomed out (left) and zoomed in (right). (c) The result of removing the tissue on the leg vertically using the traditional PBD method, zoomed out (left) and zoomed in (right). (d) The result of removing the tissue on the leg vertically using the proposed method, zoomed out (left) and zoomed in (right).

In contrast to this, if the flesh is torn in one direction, traditional PBD copies the forced particles and simulates the tearing of the tissue tearing as soon as the forced particles are removed. This method restores the minimum unit triangle that makes up the tissue. As a result, it is suitable for indicating the destruction of an object without elastic deformation because the affected tissue is cleanly removed. However, when an elastic object is destroyed, a stretched cut surface is formed due to the yielding criterion of the object, and rough debris is generated. The stress applied in the opposite vector is proportional to the position of the newly created particle. In the cases of (a) and (c), when about 200 particles were torn in a rectangular shape or vertically from a leg made of particles, the results expressed by the traditional PBD are shown. Conversely, (b) and (d) are the results of applying the proposed method. Under the same conditions same as (a) and (c), rough debris stretched as much as the red area due to the yield criterion.

4.1. Real-time performance

We verified real-time performance by comparing the processing time of the proposed idea versus that of the traditional PBD after causing loss by applying force to several particles simultaneously.

Figure 6 shows graphs that confirm whether the simulation of the proposed idea was method can be performed in real time compared to the traditional PBD method. We established 200 and 4,000 particles as benchmarks for the range of incisions. The use case of 200 particles may be suitable for representing small-scale treatments or incisions in a localised area. On the other hand, for a large-scale surgery, such as an incision on the front-facing skin area of the upper leg bone, the case of 4,000 particles would be needed can be used as a reference. There was no frame change according to the number of particles in Frame 0, in which where the object was disposed. In the two-seconds period during which loss and remaking proceeded due to the force applied to the particles, about 20% of the frame drops occurred under PBD. On the other hand, the proposed method had a frame drop of about 30%. The slightly higher computational cost is due to Tresca's yield criterion's addition to calculating where to generate a new particle between the two destroyed particles. With the proposed idea, about 30% of the frame drops occurred due to particle generation and remeshing by loss. However, but overall, a result close to 60 fps was reached, enabling real-time simulation even with more costs.

Figure 6. Comparison of the FPS of the traditional PBD and the proposed method when destroying 200 particles (a). Comparison of the FPS of the traditional PBD and the proposed method when destroying 4,000 particles (b).

4.2. Usability test

Based on eight non-experts who acquired background knowledge and two residents who had actual surgical experience, we conducted a surveyed with ten individuals labelled as respondent 1 through 10 in Table 1. Each respondent performed a incision simulation on a prepared leg model using both the traditional PBD method and the proposed method. The survey was designed based on three criteria to evaluate the effectiveness of the proposed method. First, whether the incision simulation can be performed without the inconvenience of speed. Second, whether the virtual surgical process was helpful for actual surgery. Third, whether the quality had increased compared to the traditional PBD method. The survey consisted of five questions to investigate the realism of the simulation system, and the respondents were asked to respond on a scale of 1 (strongly disagree) to 5 (strongly agree). The results are shown in Table 1.

Table 1. Survey on the realistic effect of simulation.

	Question contents
	①	②	③	④	⑤	⑥	⑦	⑧	⑨	⑩	Mean
Question number	Score (1: Strongly Disagree, 2: Disagree, 3: Neutral, 4: Agree, 5: Strongly Agree)
1	Is it possible to perform incision simulation without any inconvenience caused by the speed (frame rate)?
	5	4	5	5	4	5	5	4	4	4	4.5
2	Is there similarity between virtual incision of the organ surgery and actual incision of the organ process?
	5	4	4	4	3	3	4	4	3	4	3.8
3	Does the model with visible debris of incised organ provide a sense of realism compared to the model without them?
	4	5	5	4	5	4	5	5	4	4	4.5
4	Do you think each process of virtual surgery simulation will help trainees learn basic surgical skills?
	3	4	4	4	4	3	5	4	3	4	3.8

The first question, which measured whether simulation can be performed in real-time without inconvenience, received a score of 4.5. The second question, which measured the similarity between the simulations of cutting into virtual and actual organs, received a score of 3.8. Furthermore, the third question, which measured the realism of the process of cutting into the organs using the algorithm provided in this paper compared to traditional PBD, received a score of 4.5. The fourth question, which measured whether the virtual simulation was helpful for training, received a score of 3.8. According to the survey results, all items received scores above 3.5, confirming that the simulation system could be implemented for training in virtual surgery. Additionally, it was found that the process of cutting into organs using the proposed method had a higher degree of realism than the process of cutting into organs using traditional methods.

5. Conclusion

The proposed method in this paper allows for a more realistic representation of debris resulting from incisions made in virtual surgery. We propose a method that approximates Tresca's yield criterion to generate complex debris at the moment of destruction when an elastic object is subjected to strong forces. The traditional PBD method sets the position of a newly created particle as one of the positions of two particles when destruction between particles occurs. So, this method does not express the change of debris of an elastic object that changes due to stress when the distance between particles exceeds a certain threshold and is cut. As a result, when an object is torn, the cut area becomes relatively simple. On the other hand, the proposed method predicts the location of newly created particles through Tresca's yield criterion. The stress applied to the two particles determines the position of the new particle between the two particles, resulting in a rougher cut area when the object is torn. These results show the proposed method can express materials with complex cross-sections when cut. This method enables It allows the creation of 3D human body models with more realistic physical characteristics, which, in turn, allows for realistic diagnosis and surgical training. However, since the elasticity coefficient k of the human body varies depending on the body part, additional improvements to the equation are needed to apply Tresca's yield criterion to the human body as a wholewe have been using a specific value for the elasticity coefficient k of the human body. However, since the body's elasticity varies depending on the body part, additional improvements are needed to apply Tresca's yield criterion to the entire human body. It is necessary to apply different values of the elasticity coefficient k based on tensile experiments and depending on the specific body part. Furthermore, a more scientific approximation equation is required, as a 10% random range is a brute force value A 10% random range was also applied when creating new particles to generate complex debris. However, this is an arbitrary value, so a more scientifically derived approximation equation is needed. Additionally, future studies will proceed the potential of training artificial intelligence using actual patient surgical data to overcome hardware or algorithm limitations and to create more advanced simulations.

Disclosure statement

No potential conflict of interest was reported by the author(s).

Additional information

Funding

This work was supported by the National Research Foundation of Korea (NRF) grant funded by the Korean government (Nos. NRF-2022R1A2B5B01001553 and NRF-2022R1A4A1033549).