Real-time blood and water fusion simulation method based on improved Marching cubes algorithm accelerated by GPU

Abstract

The present invention discloses a real-time blood and water fusion simulation method based on an improved Marching cubes algorithm accelerated by GPU. To address the problems of low efficiency, poor real-time performance, and insufficient fidelity in reconstructing and rendering a large number of different types of particle surfaces using the Marching cubes algorithm, the present invention proposes first establishing an anisotropic kernel function to calculate implicit surfaces by considering the anisotropic characteristics of different liquid particles. Then, a data structure is established within the GPU thread for both water and blood particles to accelerate the search. Finally, the GPU is used to accelerate the construction of a triangular mesh using Marching cubes, and rendering is performed using OpenGL. This method can quickly and accurately reconstruct a triangular mesh from different mixed particle data, significantly improving the efficiency and accuracy of surface reconstruction using the traditional Marching cubes algorithm, and can better meet the real-time and realistic simulation requirements of blood and water in virtual surgery.

Description

Technical Field

[0001] The present invention belongs to the technical field of real-time fluid simulation, and in particular relates to a real-time blood and water fusion simulation method based on an improved Marchingcubes algorithm accelerated by GPU. Background Art

[0002] In real life, water is a Newtonian fluid, while blood, composed primarily of blood cells and plasma, is a non-Newtonian fluid with a non-constant viscosity coefficient. This means that the shear strain rate and shear stress do not conform to a linear relationship. This makes the fusion of blood and water much more complex than the physical properties of ordinary fluid mixing. Ignoring the non-Newtonian properties of blood and simplifying the treatment of blood as a Newtonian fluid significantly reduces the physical realism of the blood-water mixture. In virtual surgery, in addition to the requirement for realism, the simulation of blood-water fusion must also meet real-time requirements. Traditional physics-based fluid simulation methods, due to their computational complexity, result in low computational efficiency and struggle to meet real-time requirements.

[0003] Researchers at home and abroad have proposed numerous methods for simulating the mixing of two liquids. However, these methods often suffer from computational overhead, high storage costs, and coarse models for simulating the unique properties of non-Newtonian fluids such as blood and water, making it difficult to achieve both real-time and realistic simulations. Furthermore, the detailed changes in physical properties such as density and viscosity when mixing two liquids with different physical properties also determine the fidelity of the final rendering. Therefore, scholars at home and abroad have been diligently exploring new methods for simulating the mixing of different fluids.

[0004] Early fluid simulation methods mostly used two-dimensional heightfields to avoid extracting and rendering the fluid surface. Recently, this heightfield approach has been applied to simulate ocean waves and, in conjunction with particle systems, to simulate splashing waves. Fedkiw et al. successfully simulated smoke using grid-based surface tracking technology. However, rendering the flow of a true fluid surface requires not only methods based on fluid dynamics but also surface reconstruction and rendering based on different fluid dynamics models.

[0005] There are two main approaches to rendering liquid surfaces based on different fluid dynamics simulations – surface tracking methods based on the Euler method and surface rendering based on Lagrangian particles.

[0006] The surface tracking method based on the Euler method defines the implicit surface of the fluid by updating a scalar field. Generally speaking, this scalar field is discretized through a regular grid, and the scalar field is used to determine whether the fixed point of the Euler grid is inside or outside the fluid. The most classic method of fluid simulation based on the Lagrangian method is to use a large number of particles to simulate the motion state of the fluid, because the particles have already simulated the various physical properties of the fluid and will flow according to the physical properties of the fluid, which makes the motion of the fluid simple and does not cause stability problems. However, creating an ideal closed liquid surface by Lagrangian particles is a challenging task.

[0007] As research on mixed fluid rendering and simulation continues to deepen both domestically and internationally, the various rendering methods that support it have also continued to develop and advance. The initial challenges of rendering simulation have gradually shifted from stability and accuracy to the extraction and rendering of mixed fluid surfaces. Müller et al. proposed a Lagrangian particle-based surface rendering method that separates the fluid's interior from its exterior by defining a scalar field. The surface is then extracted and reconstructed using a Marching cubes algorithm.

[0008] Although there are currently multiple methods for fluid rendering simulation, the more classic methods include the ray tracing algorithm for volume rendering, the screen domain snowball algorithm, and the Marching Cubes algorithm for surface rendering, blood and the liquid formed by the mixture of blood and water have low transparency, do not require layering and light and shadow, and focus more on surface smoothness. In addition, as a module of the virtual surgery system, it has high real-time requirements. Therefore, the Marching Cubes algorithm based on surface rendering is more suitable for the present invention. Summary of the Invention

[0009] To address the existing challenges of balancing realism and real-time performance, this paper proposes a real-time blood-water fusion simulation method that improves computational efficiency and optimizes rendering realism. Experiments have shown that this method achieves excellent results in both visual and performance metrics.

[0010] The technical solution adopted in the present invention is:

[0011] A real-time blood-water fusion simulation method based on a GPU-accelerated Marching cubes algorithm, characterized by comprising the following steps:

[0012] S1. The weighted principal components and anisotropic characteristics of two different particles are analyzed, and an improved method of interpolating some kernels is used to calculate the implicit surface of the blood-water mixture;

[0013] S2. For a large number of different types of particles, a hash grid and spatial AABB search results are constructed separately for each particle in the GPU for search;

[0014] S3, using GPU to accelerate Marching cubes to construct triangular meshes of mixed liquid particles;

[0015] S4. Use OpenGL to render the final calculated triangle patch.

[0016] Specifically, step S1 is as follows: To address the problem that current implicit surface calculation methods may cause the extracted surface to be less smooth, resulting in less realistic rendering effects, the present invention uses an improved method for interpolating particle cores to calculate the implicit surface of a blood-water mixture, which can effectively extract a smoother surface. The specific method for calculating the implicit surface of a blood-water mixture based on this method is as follows:

[0017] The kernel function of a sphere is defined by the spatial position and the radius of the two particle spheres, as shown in the following formula (1):

[0018]

[0019] Where x i are the coordinates of the particle, h is the adjustable core radius, and x is the coordinate of the point in space. The influence range of each particle core is limited by h, which means that a particle only affects a fixed area in the space around itself.

[0020] After calculating the cumulative value on the grid, the Marching cubes algorithm is used to retrieve the triangular mesh from the implicit surface representation of the spatial grid. However, this method is affected by the particle sampling and the spherical shape of the kernel function, resulting in a concave and convex surface. Therefore, the principal component analysis method is used to extract the main distribution of a set of samples. By setting the sample to the position of a set of points, the coordinates of some main dimensions can be obtained using this method. The weighted covariance matrix of each particle kernel is calculated by adding the weighted covariance matrices of each adjacent particle. The weighted covariance matrix is ​​calculated using the distance from the particle to the adjacent particle. By performing SVD decomposition on the covariance matrix, the three axes of the elliptical kernel function can be calculated as follows (2):

[0021] [U,S,V]=svd(CovarianceM) (2)

[0022] In formula (2), Covariance M is the covariance matrix of each particle, U is a matrix whose columns are the eigenvectors of the covariance matrix, S is a diagonal matrix whose diagonal elements are the eigenvalues, and V is a matrix whose rows are the eigenvectors of the covariance matrix.

[0023] S is then modified so that the integral of W, which fills the adjacent particle kernels, is constant. Therefore, S needs to be normalized. Using the previously normalized S, U, and V matrices, a matrix G is calculated that transforms the original spherical kernel function into an elliptical kernel function. This matrix G can then be used to obtain the anisotropic kernel function.

[0024] W(r,G)=det(G)P(||Gr||) (3)

[0025] Using this kernel function to calculate the implicit surface, a smooth solution can be obtained.

[0026] Step S2 specifically involves the following: Due to the diverse and large number of particles, traditional search methods waste significant time. To address this issue, the present invention utilizes spatial hash grids and AABB data structures, and uses GPUs to accelerate the search, saving significant computational time. The most commonly used data structures for accelerated search are KD trees and spatial hash grids. For homogeneous particles, spatial hash grids with uniform node distribution offer advantages. For anisotropic particles like water and blood, each particle core has a different shape. However, normalization allows the core to be confined to a fixed-size circle. Using spatial hash grids, the particles can still be evenly distributed within the hash table's keyword table. For each anisotropic particle, searching for the exact neighboring particles requires calculating all spatial grids that intersect the ellipse. Direct calculation would be prohibitively expensive, so AABB boxes can be used to approximate the ellipse shape, reducing computational overhead and improving computational speed. Alternatively, the AABB boxes can be directly extracted from the G matrix. Transformations along the ellipse's three axes are performed, and the maximum and minimum x-coordinate values ​​are recorded. Repeating this process for all six directions yields the two AABB vertices. Since query operations on hash grids are read-only operations, if a hash grid has been built, it can be used for a single proximity grid search on a GPU thread.

[0027] Step S3 specifically includes the following steps: To address the low mesh traversal efficiency of the traditional Marching Cubes algorithm, the present invention proposes first adding weights of two particles in space to the mesh to obtain a discrete representation of an implicit surface space function, which is an isosurface. The Marching Cubes algorithm is then used to enumerate the possible weight configurations for each mesh, laying the foundation for building a triangular mesh. Because the Marching Cubes calculation for each mesh is independent, a histogram pyramid (HP) data structure can be used to allocate storage and assign indexes to triangle blocks. To count the number of triangles in each mesh, the number of triangles is first populated into the bottom storage space of each mesh. Then, 1 / 4 of the threads are run in parallel, and in each thread, the total number of four available triangles is assigned to the thread. This process continues until the final number of triangles is obtained. After calculating the total number of triangles, GPU memory should be allocated based on the number of triangles. In each mesh, the Marching Cubes algorithm is applied to obtain triangle blocks based on the mesh values. When storing the triangle information in a single buffer, a method for obtaining indexes from the histogram pyramid is required, which is also an iterative process. This method can effectively speed up the speed and accuracy of the Marching cubes algorithm for grid traversal.

[0028] The step S4 specifically includes: after completing the above three steps, rendering the calculated triangle facets using OpenGL.

[0029] The present invention uses relevant technologies of fluid dynamics simulation of the fusion of blood and water and surface reconstruction based on the improved Marching cubes algorithm accelerated by GPU. By improving the Marching cubes algorithm to establish implicit surfaces and accelerating the GPU, the surface reconstruction and rendering of the mixed liquid of blood and water can be achieved. This not only improves the speed of the Marching cubes algorithm in surface reconstruction of the mixed particles of blood and water and shortens the time, but also greatly improves the realism and accuracy of the rendering of the mixed liquid. BRIEF DESCRIPTION OF THE DRAWINGS

[0030] Figure 1 This is a flowchart of the GPU-accelerated Marching cubes algorithm of the present invention;

[0031] Figure 2 This is the principal component analysis diagram of blood and water particles according to the present invention;

[0032] Figure 3 This is an AABB box display diagram of the anisotropic core query of blood and water particles according to the present invention;

[0033] Figure 4 This is a diagram showing the effect of the present invention on the blood and water fusion model particles;

[0034] Figure 5 This is a diagram showing the effect of blood and water fusion according to the present invention. DETAILED DESCRIPTION

[0035] The present invention will be further described below with reference to the accompanying drawings.

[0036] like Figure 1 As shown in the figure, a real-time blood and water fusion simulation method based on the GPU-accelerated Marching cubes algorithm is proposed. First, an anisotropic kernel function for particles is established and used to calculate the implicit surface. Second, a data structure is created for both water and blood particles within the GPU thread to accelerate the search. Third, the GPU is used to accelerate the construction of the triangular mesh using Marching cubes, and OpenGL is used for rendering. The detailed process includes the following steps:

[0037] S1. The weighted principal components and anisotropic characteristics of two different particles are analyzed, and an improved method of interpolating some kernels is used to calculate the implicit surface of the blood-water mixture;

[0038] S2. For a large number of different types of particles, a hash grid and spatial AABB search results are constructed separately for each particle in the GPU for search;

[0039] S3, using GPU to accelerate Marching cubes to construct triangular meshes of mixed liquid particles;

[0040] S4. Use OpenGL to render the final calculated triangle patch.

[0041] Step S1 specifically includes: To address the problem that current implicit surface calculation methods may result in the extracted surface being not smooth enough, the present invention uses an improved method of interpolating kernels to calculate the implicit surface, which can extract a smoother surface. The present invention uses this method to calculate the implicit surface of a mixture of blood and water as follows:

[0042] The kernel function of a sphere is defined by the spatial position and the radius of the two particle spheres, as shown in the following formula (1):

[0043]

[0044] Where x i are the coordinates of the particle, h is the adjustable core radius, and x is the coordinate of the point in space. The influence range of each particle core is limited by h, which means that a particle only affects a fixed area in the space around itself.

[0045] After the cumulative value on the grid is calculated, the Marching cubes algorithm is used to retrieve the triangular mesh from the implicit surface representation of the spatial grid. However, this method is affected by the particle sampling and the spherical shape of the kernel function, which causes a concave and convex surface. Therefore, the principal component analysis method is used to extract the main distribution of a set of samples. By setting the samples to the positions of a set of points, the coordinates of some main dimensions can be obtained using this method. The shape is as follows Figure 2 .

[0046] The weighted covariance matrix of each particle kernel is calculated by adding the weighted covariance matrices of each adjacent particle. The weighted covariance matrix is ​​calculated using the distance from the particle to the adjacent particles. By performing SVD decomposition on the covariance matrix, the three axes of the elliptical kernel function can be calculated as follows (2):

[0047] [U,S,V]=svd(CovarianceM) (2)

[0048] In formula (2), Covariance M is the covariance matrix of each particle, U is a matrix whose columns are the eigenvectors of the covariance matrix, S is a diagonal matrix whose diagonal elements are the eigenvalues, and V is a matrix whose rows are the eigenvectors of the covariance matrix.

[0049] S is then modified so that the integral of W, which fills the adjacent particle kernels, is constant. Therefore, S needs to be normalized. Using the previously normalized S, U, and V matrices, a matrix G is calculated that transforms the original spherical kernel function into an elliptical kernel function. This matrix G can then be used to obtain the anisotropic kernel function.

[0050] W(r,G)=det(G)P(||Gr||) (3)

[0051] Using this kernel function to calculate the implicit surface, a smooth solution can be obtained.

[0052] Specifically, the step S2 is as follows: due to the different types and large number of particles, traditional search methods will waste a lot of time. To address this problem, the present invention utilizes some data structures and uses GPUs to accelerate the search, which can save a lot of computing time. For isomorphic particles, a spatial hash grid with uniform node distribution is more advantageous. For the two anisotropic particles of water and blood, the shape of each particle nucleus of the two particles is different, but the normalization operation allows the range of the nucleus to be limited to a circle of a fixed size, and then the particles can be evenly distributed in the keyword table of the hash table using the spatial hash grids separately. For each anisotropic particle, to query the precise adjacent particles, it is necessary to calculate all spatial grids that intersect with the ellipse. If calculated directly, the cost will be very high, so the AABB box can be used to approximate the shape of the ellipse to reduce the computing pressure and thereby increase the computing speed. In addition, the AABB box can be directly extracted from the structure of the G matrix; the transformation is performed along the three axes of the ellipse, and then the maximum and minimum values ​​of the x coordinate are recorded; this process is repeated for all six directions to obtain the two vertices of the AABB. As Figure 3 It shows how to use AABB for querying. Since query operations on hash grids are read-only operations, if the hash grid has been built, it can be used for a single proximity grid search on the GPU thread.

[0053] The step S3 is specifically as follows: in view of the low efficiency of the traditional Marching cubes algorithm in traversing the grid, by adding the weights of the two particles in space to the grid, a discrete representation of the spatial function of the implicit surface can be obtained, which is an isosurface. The Marching cubes algorithm enumerates the possible weight configurations of each grid, which is also the basis for constructing a triangular grid. The calculation of Marching cubes for each grid is independent, so a histogram pyramid (HP) data structure can be used to allocate storage and assign indexes to triangle blocks. Figure 4 In order to count the number of triangles in each mesh, the number of triangles is first filled into the bottom storage space of each mesh. Then 1 / 4 of the threads are run in parallel, and in each thread, the total number of 4 ready-made triangles is assigned to the thread. This process continues until the final number of triangles is obtained. After calculating the total number of triangles, GPU memory should be allocated according to the number of triangles. In each mesh, the Marching cubes algorithm is applied to obtain triangle blocks according to the mesh values. When storing triangle information in a single buffer, a method of obtaining indices from the histogram pyramid is required, which is also an iterative process.

[0054] The step S4 specifically includes: after completing the above three steps, rendering the calculated triangle facets using OpenGL.

[0055] like Figure 5 As shown in FIG. 1 , it is a rendering of the present invention. It can be seen that a more realistic effect of blood and water mixing will appear after rendering using the improved Marching cubs algorithm.

[0056] The above description merely represents the preferred embodiments of the present invention, and while the description is relatively specific and detailed, it should not be construed as limiting the scope of the present invention. It should be noted that a person skilled in the art would be able to make various modifications, improvements, and substitutions without departing from the spirit of the present invention, all of which fall within the scope of protection of the present invention. Therefore, the scope of protection of the present invention shall be determined by the appended claims.

Claims

1. A real-time blood-water fusion simulation method based on an improved Marching cubes algorithm accelerated by GPU, characterized by: The method comprises the following steps: S1. Analyze the weighted principal components and anisotropic characteristics of two different particles, and use an improved method of particle core interpolation to calculate the implicit surface of the blood-water mixture; S2. For a large number of different types of particles, a hash grid and spatial AABB search results are constructed separately for each particle in the GPU for search; S3, using GPU to accelerate Marching cubes to construct triangular meshes of mixed liquid particles; S4, using OpenGL to render the final calculated triangle patch; The step S1 is specifically as follows: S1-1. Define the kernel function of a sphere by the spatial position and the radius of two particle spheres; The kernel function is as follows: （1）； Where, are the coordinates of the particle, is the adjustable kernel radius, are the coordinates of a point in space; S1-2. After calculating the cumulative value on the grid, the Marching cubes algorithm is used to retrieve the triangular mesh from the implicit surface representation of the spatial grid, and the principal component analysis method is used to eliminate the concave and convex surfaces caused by the influence of particle sampling and the spherical shape of the kernel function; The weighted covariance matrix of each particle kernel is calculated by adding the weighted covariance matrices of each adjacent particle. The weighted covariance matrix is ​​calculated using the distance from the particle to the adjacent particles. By performing SVD decomposition on the covariance matrix, the three axes of the elliptical kernel function can be calculated as follows (2): （2）； In formula (2), Covariance M is the covariance matrix of each particle, U is a matrix whose columns are the eigenvectors of the covariance matrix, S is a diagonal matrix whose diagonal elements are eigenvalues, and V is a matrix whose rows are the eigenvectors of the covariance matrix; Then, S is modified so that the integral of W filled with adjacent particle nuclei is a constant, and S is standardized. For the previously standardized modified S, U, and V matrices, a matrix G is calculated that can transform the original spherical kernel function into an elliptical kernel function, and this matrix G can be used to obtain anisotropic kernel functions. (3)； S1-3. Use the kernel function of formula (3) to calculate the implicit surface and obtain a smooth solution; The step S3 is specifically as follows: S3-1. By adding the weights of the two particles in space to the grid, a discrete representation of the implicit surface space function is obtained as an isosurface. S3-2, enumerate the possible weight configurations of each grid through the Marching cubes algorithm, laying the foundation for building a triangular mesh; S3-3. Use the histogram pyramid data structure to allocate storage and assign indexes to triangle blocks to build a triangle mesh.

2. The real-time blood-water fusion simulation method based on the improved Marching cubes algorithm accelerated by GPU according to claim 1, characterized in that: Specifically, step S2 includes the following steps: for the two anisotropic particles of water and blood, each particle core has a different shape, and a normalization operation is performed to limit the return of the core to a circle of a fixed size. Then, a spatial hash grid is used for each of the two particles to evenly distribute the two particles into the keyword table of each hash table. For each anisotropic particle, an AABB box is used to reduce computational pressure, increase computing speed, and thus query accurate adjacent particles.

Images