CN103425849B - The emulation mode of the lower fluid of a kind of dynamic barrier impact - Google Patents

Abstract

The present invention is the emulation technology that the mode based on physics is simulated the lower fluid behavior of dynamic barrier impact. By advective term, external force term and the pressure term in computational fluid dynamics equation (Navier-Stokes equation) successively on GPU, upgrade the physical quantity such as speed, pressure of fluid, complete the simulation calculation of fluid. In the time that pressure term is calculated, consider influencing each other between fluid and barrier, set up the coupling linear equation between fluid and barrier according to the lagrange equation of motion of barrier and the incompressibility of fluid, upgrade pressure term by solving this linear equation. Because the simulation calculation of fluid realizes in GPU, and the simulation of barrier is carried out in CPU, in order to realize the coupling of fluid and barrier, the present invention carries out voxelization process to dynamic barrier in simulation process, be mapped on the emulation graticule of fluid the calculating of the fluid pressure that combined with GPU by CPU.

Description

The emulation mode of the lower fluid of a kind of dynamic barrier impact

Technical field

The invention belongs to virtual reality technology field, relate in particular to the dynamic barrier that CPU combines with GPU based on physics modeHinder the emulation mode of the lower fluid of thing impact.

Background technology

In actual life, can be found everywhere about the phenomenon of fluid, such as the kitchen smoke rising up slowly, the rivers that flow ahead fast ceaselessly,Flame blazing etc. In fluid, exist various barriers toward contact meeting, have and also have dynamically static state, such as free-swimming fish in river, smog passes through forest etc. At natural phenomena computer simulation field, fluidEmulation is one of them important research direction, and it similarly is a kind of useful building puzzle, is a lot of natural phenomenas of emulationBasis. But due to its huge assessing the cost, it all only belongs to the field that high-quality off-line is played up all the time. In real timeIn application, the effect of 3D fluid is all often to depend on particIe system to realize, and has very much attraction although particIe system can produceThe result of power, but they and do not meet true outward appearance and the behavior of fluid, the emulation of fluid in real time remains a kind of challenge, outstandingIt is the analogue simulation that also exists various barriers in fluid.

Summary of the invention

For the real motion behavior that shows the lower fluid of dynamic barrier impact that can be real-time, the invention provides onePlant the emulation mode of the lower fluid of dynamic barrier impact. The present invention adopts the physical message of 3D texture storage fluid, as speed, pressurePower etc., are used 3D texture that simulation calculation can be realized completely on GPU, have improved greatly simulation efficiency. In advection mistakeIn the middle of carrying out two in journey, semi-Lagrange step adds two-way error compensation and aligning step, and advection precision is improved by single orderTo second order, reduce the numerical dissipation in advection process. Improve the method for solving of pressure term in fluid motion equation, consideredInfluencing each other of fluid and dynamic barrier, carries out voxelization process to dynamic barrier, is mapped in fluid grid, according to barrierHinder the lagrange equation of motion of thing and the incompressibility of fluid to set up the coupling linear equation between barrier and fluid,Combine with the GPU pressure of Fluid Computation of CPU, synchronously upgrades the motion state of fluid and barrier, has improved dynamic barrierThe authenticity of the lower fluid motion of impact. In GPU, when computational speed divergence and barometric gradient, sampled point is offset, utilizesThe account form of MAC grid is sampled to 3D texture, has reduced the error in computational process, has improved computational accuracy.

The behavior emulation mode of the lower fluid of a kind of dynamic barrier impact of the present invention, the method comprises the following steps:

Step 1, the physical quantity (as speed, pressure etc.) of initialization fluid, is stored in the middle of 3D texture, is convenient in GPUCalculate.

Step 2, uses Semi Lagrangian scheme, the midrange speed field of advection Fluid Computation; For reducing in advection processDissipate, obtain higher advection precision, further improve Semi Lagrangian scheme, by carrying out in advance two centre half glugsBright day advection step, forward advection and a reverse advection (that is, the time is contrary), then carry out two-way error compensationWith aligning step, obtain more high-precision advection value.

Step 3, all dynamic barriers of voxelization, obtain 3D barrier texture and 3D boundary position texture. To dynamicallyWhen barrier carries out voxelization, can further comprise again:

Step 3.1, utilizes the decrement operations that increases progressively of template buffer memory to carry out inside and outside voxelization to all dynamic barriers, markNote barrier unit is SOLID, and element of fluid is FLUID, is saved in 3D barrier texture;

Step 3.2, barrier model, (Geometry in geometric coloration are played up in each grid section one timeShader) the each triangle surface in barrier is carried out to the crossing operation of cutting into slices with grid, use bilinear interpolation to calculatePosition of intersecting point coordinate, and unit, corresponding intersection point place is labeled as to BOUNDARY, represent the boundary element of barrier. In pixelIn tinter, label information is rendered in barrier texture, corresponding position of intersecting point is rendered in boundary position texture.

Step 4, according to 3D barrier texture, dyscalculia thing boundary normal, is kept in 3D boundary normal texture. PinTo boundary element, judge the label information of its neighbours unit according to barrier texture, adjacent at boundary element and element of fluidThe local normal that calculates, is saved in corresponding boundary normal texture. Wherein, normal direction is pointed to boundary element from element of fluid.

Step 5, reads 3D barrier texture, records the index position of obstacles borders in texture.

Step 6, the index position according to obstacles borders in texture, extracts boundary position and boundary normal and is saved inIn corresponding list of locations and normal list.

Step 7, according to the incompressibility of the barrier equation of motion and fluid, sets up between fluid and barrierCoupling linear equation. The coupling linear equation of setting up according to the incompressibility of the equation of motion of barrier and fluid is as public affairsShown in formula (1):

（1）

Wherein,For the transition matrix between generalized coordinates system and cartesian coordinate system,(k borderUnit number);For the mass matrix of barrier;For the speed under generalized coordinates before barrier and fluid coupling.

Step 8, CPU is combined with GPU item, solves coupling linear equation, the pressure field of Fluid Computation. Solving coupling lineProperty when equation, the method that adopts GPU to combine with CPU, makes、, first according to the glug of barrierBright daily motion equation calculates A and B in CPU, then respectively A and B is write to each self-corresponding 3D according to obstacles borders indexIn texture, in GPU, complete the calculating of terminal pressure. In the time of Fluid Computation pressure, and though solve be coupling linear equation alsoBe Poisson's equation, all need the divergence of computational speed in advance.

Step 9, according to the gradient of pressure field, goes on the midrange speed field of projection fluid is long without divergence to it, obtains fluidFinal speed field. In this step, need to carry out projection calculating to middle velocity field according to the gradient of pressure field. Calculating speedWhen degree divergence and barometric gradient, adopt the account form of MAC grid (a kind of staggered-mesh, corresponding with central gridding) to 3D lineReason is sampled, and in three directions of 3D texture, is offset separately 0.5 unit and samples, and can reduce in computational processError, improve computational accuracy.

The invention has the beneficial effects as follows: the physical message of fluid is kept in the middle of 3D texture, emulation is obtained on GPUTo realize, improve simulation efficiency, simultaneously the mutual corresponding simulation calculation process of also having simplified of 3D texture and Eulerian mesh. At meterWhen calculation speed divergence and barometric gradient, improve the sample mode of 3D texture, in three directions of 3D texture, sampled point has been dividedBe not offset 0.5 unit, utilize the account form of MAC grid to sample to 3D texture, reduced the error in computational process,Improve computational accuracy. The voxelization process that improves dynamic barrier, adds a boundary voxel process, in geometric colorationIn the each triangle surface in barrier is carried out and the crossing operation of target slice, calculate intersection point by bilinear interpolation(intersection point is obstacles borders position), outputs to pixel coloring device by these points, is labeled as by pixel coloring deviceBOUNDARY outputs in the middle of 3D barrier texture, the position of these intersection points is saved in 3D border texture simultaneously. Due to barrierHinder the motion of thing to be simulated in CPU, and the motion of fluid is simulated in GPU, is coupled at fluid and barrierWhen calculating, CPU is combined with GPU, the boundary information that voxelization obtains according to barrier, counts the pressure of boundaryCalculate, upgrade the motion state of fluid and barrier by this Simultaneous Pressure, not only improved the efficiency of emulation, also improved obstacleThe authenticity of the lower fluid motor behavior of thing impact.

Brief description of the drawings

Fig. 1 shows the emulation mode flow chart of the lower fluid behavior of a kind of dynamic barrier impact of the present invention.

Fig. 2 shows semi-Lagrange advection method.

Fig. 3 shows the semi-Lagrange advection step after improvement.

Fig. 4 shows advection error analysis and rectification building-out.

Fig. 5 shows the storage mode of MAC grid in the present invention.

Fig. 6 shows the derivation of the linear equation that is coupled between fluid and barrier.

Specific embodiments

Below in conjunction with accompanying drawing and embodiment, the preferential embodiment of the present invention is further illustrated:

Flow chart shown in Fig. 1 has provided the detailed process of the whole enforcement of the present invention:

Step 1, the 3D texture of initialization fluid grid establishment and the corresponding equal size of grid, preserves initial shapeThe physical quantity of fluid under state (speed, pressure etc.). 3D texture can be regarded as a kind of array of three-dimensional, for GPU computational processMiddle save data reads the value in texture by sampling operation in GPU.

Step 2, the process that is transmitted fluid itself and the physical quantity in fluid by speed is called advection, based on putting down of GPUFlowing a calculating is to use Semi Lagrangian scheme, from each grid cell inverse time tracking racetrack, obtains its position in the past,And to the grid cell starting, copy the physical quantity (such as speed) of that position, as shown in figure (2), the inverse time is followed the trail of velocity fieldTo X position, value in the grid of the most close X position is carried out to bilinear interpolation, draw required speed, and be stored in 3D speedIn the middle of degree texture. For reducing the numerical dissipation in advection process, improve advection precision, advection process can further be improved, and passes throughSemi-Lagrange step in the middle of carrying out in advance two, a forward and reverse (, the time is contrary). Carrying out halfWhen Lagrangian method, if first carry out a forward carry out again one reverse, must return in theory before executionInitial value, but the calculating of Semi Lagrangian scheme can exist certain error, supposes that this error is e, is carrying out soAfter complete forward and reverse semi-Lagrange, will produce the error of 2e with respect to initial value,, can obtain error. During in the end a step semi-Lagrange advection calculates, two-way error compensation and aligning step, use in advanceCurrency deducts error and carries out semi-Lagrange advection again, and so operation can reduce the numerical dissipation causing because of error,Improve advection precision. The flow process that whole advection calculates is as shown in figure (3), and error analysis process is as shown in figure (4).

Step 3, carries out voxelization process to all dynamic barriers, is mapped on fluid grid, and concrete process is as follows:

Step 3.1, carries out inside and outside voxelization to dynamic barrier. This embodiment adopts rectangular projection the barrier of inputHinder object model to be played up once in each section of object 3D texture structure. In the time of rendering model, utilize template buffer memory (with sectionDimension identical) increase progressively decrement operations, the back side increases progressively, successively decrease in front. So, all voxels that are positioned at model inside will obtainTo the stencil value of a non-zero, outside all voxels are still zero. Finally, do once finally traversal, stencil value is copied toIn 3D texture, wherein barrier unit is labeled as SOLID, and element of fluid is labeled as FLUID.

Step 3.2, to the voxelization of dynamic barrier exercise boundary, the boundary position of dyscalculia thing. This embodimentIn, boundary voxel process is similar with inside and outside voxelization, and model is played up in each grid section one time, but does not this time use template slowDeposit, but in geometrical processor, judge the intersection point of each triangle surface with target slice, bilinear interpolation is calculated intersection point positionPut, this position is boundary position. At pixel coloring device, active cell is labeled as to BOUNDARY and outputs to 3D barrier textureIn, export active cell position in corresponding 3D border texture simultaneously.

Step 4, plays up fluid grid one time, in pixel coloring device, judges that according to 3D barrier texture current cell isNo is boundary element (marking whether as BOUNDARY), dyscalculia thing boundary normal. The computational methods of boundary normal are as follows:Detecting the neighbours unit of current cell before and after up and down according to 3D barrier texture is barrier unit (SOLID), obstacleThing boundary element (BOUNDARY) or element of fluid (FLUID). Such as, if its left neighbours (x negative direction) are element of fluid,The normal of active cell points to barrier for 1(normal direction by fluid on directions X so). In other direction, be same behaviourDo, finally normal information is outputed in corresponding 3D boundary normal texture.

Step 5, creates CPU addressable 3D texture, by the information reproduction in 3D barrier texture to the line just having createdIn reason. Access texture, records the index of position in the cell place texture that is labeled as BOUNDARY, is saved in boundary position ropeDraw in list.

Step 6, equally creates the addressable 3D texture of CPU with step 5, by boundary position texture and boundary normal textureInformation copy to respectively in respective texture. Access texture, the index of trying to achieve according to step 5 reads corresponding boundary positionBe saved in boundary position list and boundary normal list with boundary normal.

Step 7, builds in normal direction according to the incompressibility of the lagrange equation of motion of barrier and fluidVertical coupling linear equation. The lagrange equation of motion of barrier is as shown in Equation (2):

（2）

Wherein,Represent the mass matrix of barrier,Represent suffered coriolis force;WithRepresentThe generalized force that barrier is inside and outside.

The incompressibility of fluid, the divergence of fluid velocity field is zero, as shown in Equation (3):

（3）

The process of establishing of coupling linear equation is as shown in figure (6), and the coupling linear equation of foundation is as formula (4)

Shown in:

（4）

Wherein, S is selection matrix, and for selecting the boundary element of barrier, coupling linear equation is just single for borderUnit's normal direction.

Step 8, the pressure field of Fluid Computation, can be further divided into two kinds of situations.

Situation 1, is labeled as the boundary element of BOUNDARY, in its normal direction by the iterative linear equation that is coupledCalculating pressure. Order、, can hold very much according to the boundary position of barrier and boundary normal informationEasy A and the B of calculating, is then filled into the value of A and B in the middle of 3D texture, carries out calculation of pressure in pixel coloring device. IfThe normal direction of boundary element is the positive direction towards X, first divergence is carried out in formula (4) equal sign the right and calculates, divergence calculating sideFormula as shown in Equation (5), is then calculated the force value of active cell according to formula (6).

（5）

（6）

Situation 2, is labeled as in the tangential direction of BOUNDARY and is labeled as the element of fluid part of FLUID, by repeatedlyIn generation, solves pressure Poisson's equation calculating pressure, and Poisson's equation as shown in Equation (7).

（7）

The mode that solves of Poisson's equation as shown in formula (8), formula (9),

（8）

（9）

Wherein, plus-minus 1/2 is illustrated in corresponding coordinate or to bearing or to 1/2 unit of positive direction skew, 3D texture being carried outSampling, as shown in figure (5).

Step 9, the pressure field required according to previous step, calculates its Grad, adopts Helmholtz-Hodge to decompose(any vector field can be decomposed into two other vector field sums to theorem: one is without divergence vector field, and another is markThe gradient of amount field, is shown in shown in formula (10). ) the midrange speed field of fluid is projected to it get on without Divergence Field, obtain fluidFinal speed field, ensures that emulation front and back all meet the incompressibility of fluid.

（10）

Wherein,Represent the velocity field that divergence is zero;Represent the velocity field of divergence non-zero.

Claims (3)

1. an emulation mode for the lower fluid of dynamic barrier impact, is characterized in that, the method comprising the steps of:

Step 1, the physical quantity of initialization fluid, is stored in the middle of 3D texture;

Step 2, the midrange speed field of Semi Lagrangian scheme advection Fluid Computation; Described Semi Lagrangian scheme advection calculatesThe midrange speed field of fluid, is specially: by carrying out in advance semi-Lagrange advection steps in the middle of two, forward advection andA reverse advection, then carries out two-way error compensation and aligning step, obtains more high-precision advection value;

Step 3, all dynamic barriers of voxelization, obtain 3D barrier texture and 3D boundary position texture;

Step 4, according to 3D barrier texture, computation bound normal, is kept in 3D normal map;

Step 5, reads 3D barrier texture, records the index position of obstacles borders in texture;

Step 6, the index position according to obstacles borders in texture, extracts corresponding boundary position and boundary normal and preservesIn corresponding list of locations and normal list;

Step 7, according to the incompressibility of the barrier equation of motion and fluid, sets up the coupling between fluid and barrierLinear equation; Described coupling linear equation is as shown in Equation (1):

（1）

Wherein,For the transition matrix between generalized coordinates system and cartesian coordinate system,, k is that border is singleUnit's number;For the mass matrix of barrier;For the speed under generalized coordinates before barrier and fluid coupling;

Step 8, CPU combines with GPU, solves coupling linear equation, the pressure field of Fluid Computation;

Step 9, according to the gradient of pressure field, the midrange speed field of projection fluid without on Divergence Field, obtains the final of fluid to itVelocity field.

2. the emulation mode of the lower fluid of dynamic barrier impact according to claim 1, is characterized in that described stepThe all dynamic barriers of voxelization in 3, are specially:

Step 3.1, utilizes the decrement operations that increases progressively of template buffer memory to carry out inside and outside voxelization to all dynamic barriers, mark barrierHindering thing unit is SOLID, and element of fluid is FLUID, is saved in 3D barrier texture;

Step 3.2, barrier model is played up in the section of each grid one time, in geometric coloration to the each triangle in barrierShape dough sheet is carried out the crossing operation with grid section, uses bilinear interpolation to calculate position of intersecting point coordinate, and by corresponding intersection pointUnit, place is labeled as BOUNDARY, represents the boundary element of barrier; In pixel coloring device, label information is rendered into barrierHinder in thing texture, corresponding position of intersecting point is rendered in boundary position texture.

3. the emulation mode of the lower fluid of dynamic barrier impact according to claim 1, is characterized in that described step 9In, according to the gradient of pressure field, middle velocity field is carried out to projection calculating, in the time of computational speed divergence and barometric gradient, adopt MACThe account form of grid is sampled to 3D texture, in three directions of 3D texture, is offset separately 0.5 unit and adoptsSample.