CN102789650A - Sea surface track parallel simulation method based on particle system - Google Patents

Abstract

The invention relates to a sea surface track parallel simulation method based on a particle system. A sea surface track is generated by the method combining a real grid and a virtual particle, the parallel operation of the grid and a particle model is realized on the basis of CUDA (Compute Unified Device Architecture), and sea surface and track rendering in GPU (Programmable Graphics Processing Unit) is realized by GLSL (OpenGL Shading Language). According to the method disclosed by the invention, all the core computations and rendering are carried out in the GPU, the model computation and the rendering constitute a whole process in the GPU, the communication between the internal memory and the video memory is reduced, and the parallel processing capacity of the GPU is fully exerted. The experiment of the inventor proves that the simulation efficiency can be obviously improved under the condition of ensuring a good simulation effect.

Description

A kind of sea flight path parallelization emulation mode based on particIe system

Technical field

The invention belongs to Yu Haiyang visualization technique field, be specifically related to a kind of GPU of utilization parallel processing, based on the sea track Simulation method of particIe system.

Background technology

Ocean visual (Ocean Visualization) is the utilization computer graphics techniques, in visual form, and the simulation real-time change of ocean surface on time and space and the visual effect that is demonstrated.

The sea track Simulation is that the simulation naval vessel navigates by water the sea spray that evokes across the sea and rolls and diffusion phenomena.From simulated effect, the whole simulation content comprises three parts: the one, and the common zone, sea that flight path does not have influence on; The 2nd, the rolling and diffusion phenomena of flight path; The 3rd, flight path is to the disturbing phenomenon on normal sea.

For the common sea zone emulation that flight path does not have influence on, key problem is the height field that makes up the sea random fluctuation.The highest emulation mode of present existing fidelity is based on the method for physical model; Level and smooth particle flux body dynamics (SPH) method that adopts like people such as Guo Qiulei based on the GPU framework (Guo Qiulei, Tang Yi it, Liu Shiqiu; The clear .2011 of Li Gui; The full GPU implementation framework of a SPH fluid real time modelling. computer utility and software, 28 (11), 69-72.); The method of finding the solution the Navier-Stokes equation (the Aziz Belmiloudi that people such as Aziz Belmiloudi adopt; Francoise Brossier.1997; Regularity Results for a Navier-Stokes Type Problem Related to Oceangraphy.Acta Applicandae Mathematicae; 48,299-316.) etc.Though this emulation mode fidelity based on physical model is the highest, its internal physical rule computing is too complicated, though through optimizing and quickening, still can't be applicable to the real-time simulation field at present.The higher method of another kind of validity is based on the method for wave spectral model.As Zhu mustached using spectral methods based Philips (Zhu Cai .2010, dynamic marine environment simulation studies a number of key technologies, Huazhong University), Liang Xiufeng, who used a method based Jonswap spectrum (Liang Xiufeng, Yang Jianmin, Li Jun , XIAO Long-fei, LI Xin .2011, marine engineering applications for numerical wave tank. Chinese Science: Physics, Mechanics & Astronomy ,41,1123-122) and so on.These methods depend on mathematical model, but less with respect to operand based on the method for physics, after taking some optimized Measures, relatively are applicable to sea emulation in real time on a large scale.

For rolling and diffusion phenomena emulation of flight path, key problem is the motion of the realistic simulation flight path spray.By the method that the Kelvin ship wave is theoretical and particIe system combines that is based on of extensive employing, roll and diffusion phenomena at present with the tail spray that produces in the spray particle simulation ship navigation.Michael E.Goss has just adopted this method (Michael E.Goss.1997; A real time particle system for display of ship wakes.Computer Graphics and Applications; IEEE, 48:299-316.).Use particIe system Simulated Water wave divergence comparatively realistically, but along with the increasing of number of particles, the real-time of system receives very big influence.In order to promote simulation efficiency, people such as Zhao Xin have proposed a kind of LOD-GPU and have mixed the method for quickening (justice is kept in war for Zhao Xin, Li Fengxia; The quick .2008 of Chen Hong is based on the accelerated method of particIe system realization ship track emulation. Maritime Affairs University Of Dalian's journal, and 34 (1), 54-57.); People such as Zhang Qin have proposed some and have reduced method (Zhang Qin, Zhang Jian, the Min Jianping .2003 of computation complexity; Improve the method research of particIe system real-time. computer engineering, 29 (18), 46-48.).

For the disturbing phenomenon emulation of flight path to normal sea, traditional realization method is the stack of simple visual effect, like the three-dimensional naval vessel flight path method for visualizing (Ma Tian of people such as Ma Tian proposition; The yellow foundation; You ram the king, Zhang Qunfei .2011, the naval vessel flight path three-dimensional visualization new method of application ripples disperse characteristic. XI AN JIAOTONG UNIVERSITY Subject Index; 45 (8), 49-53.).This mode makes the sea grid almost not have interactivity with the naval vessel flight path in itself, so simulated effect is true inadequately.Sun Hui proposed the mutual ripples real time modelling method of a kind of grid and particle (Sun Hui. based on the water Behavior modeling .2010 of particle and ripple; East China Normal University .); But by the restriction of its interactive algorithm, increase when increasing with number of particles in the grid fineness, system real time is had a strong impact on.

In sum, existing sea track Simulation method mainly focuses on two problems, and the one, the sense of reality of enhancing realistic model, the 2nd, the arithmetic speed of lifting realistic model.The existing methods defective is, can not simulate the disturbing phenomenon of flight path to normal sea comparatively realistically on the one hand, can not reach the good coordination of authenticity and real-time on the other hand.

Parallel processing based on GPU is the effective ways that promote the model arithmetic speed.The full name of GPU is Programmable Graphics Processing Unit, GPU promptly able to programme.Have more ALU (Arithmetic Logic Unit, ALU) in the hardware configuration of GPU and be used for data processing, therefore be fit to very much intensive data are carried out parallel processing.Grid computing and particle computing in the track Simulation of sea all have comparatively intensive operand, and the internal arithmetic rule is identical, are fit to very much utilize GPU to carry out concurrent operation and handle.

CUDA is the parallel processing framework of a cover of NVIDIA company release based on GPU.CUDA as the data parallel computing equipment, is doing significant improvement with GPU aspect the GPU general-purpose computations: the one, adopted the Unified Treatment framework, and can effectively utilize the computational resource of summit renderer and pixel rendering device more; The 2nd, introduced shared storage in the sheet, support communicating by letter between random writing and parallel thread.In the CUDA programming model, CPU is called main frame (Host), and GPU is called equipment (Device).The CUDA parallel computation function that operates on the GPU is called Kernel (kernel function), and Kernel is that unit carries out with Block, but Kernel can be on a plurality of Block executed in parallel, and a plurality of Thread executed in parallel are arranged in each Block.

Summary of the invention

In order to reach the balance coordination of the sense of reality and real-time, the present invention proposes a kind of sea flight path parallelization emulation mode based on particIe system.Adopt grid to generate the sea flight path for the method that empty, actual situation combine,, play up in conjunction with sea and flight path among the GLSL realization GPU based on the parallelization computing that CUDA realizes grid and particle model for real, particle.This method is carried out with playing up all core calculations in GPU; Model Calculation is constituted an overall flow with playing up in GPU; Reduced the communication between calculator memory and the video memory; Give full play to the GPU parallel processing capability, under the situation that guarantees good simulated effect, significantly promote simulation efficiency.

Technical scheme of the present invention is following:

A kind of sea flight path parallelization emulation mode based on particIe system may further comprise the steps (Fig. 1):

(1) in CPU, carries out some initialized work; Specifically comprise: initialization CUDA equipment, application video memory, establishment FFT Plan, establishment VBO, initialization particIe system;

(2) in GPU, carry out the parallelization computing of sea grid; Specifically comprise: find the solution the Phillips wave frequency spectrum of a certain particular moment the sea level height field, carry out CUFFT find the solution Phillips wave frequency spectrum the wave diffusion field, result of calculation is outputed to the Global Memory in the video memory;

(3) in GPU, carry out the parallelization computing of flight path particIe system; Specifically comprise: calculate and upgrade particle property, with particle position be mapped to the sea grid system, calculate particle to the disturbance of grid around the position go forward side by side disturbance stack computing in the line range, operation result is outputed to Global Memory.

(4) in GPU Shader, use GLSL that environment, illumination etc. are carried out in the sea and play up, form dynamic sea flight path effect.

In each above-mentioned step; The parallelization computing of sea grid, the parallelization computing of flight path particIe system, particle are to disturbance stack computing in the scope of position; This three is that the inventor comes out through the research and design of hardships, below will describe above three's method for designing in detail.

The parallelization fortune of sea grid is what on the basis of Phillips spectral model, to be carried out.The Phillips frequency spectrum is the wave frequency domain amplitude under wind action that goes out through Research on Statistics and Analysis, possesses wave model features such as randomness, periodicity and multifrequency be forthright.On the Phillips frequency spectrum, add random perturbation and time driving conversion, can simulate the sea fluctuation under the wind action well.Frequency spectrum function after the conversion comprise T constantly the sea level height field find the solution with whole grid scope in the wave diffusion field find the solution two parts.

Find the solution the sea level height field constantly for T, each sea net point all is independently, therefore can utilize the two dimensional form of sea grid to carry out parallel processing.The mesh-density of supposing simulating area is meshW*meshH, splits according to two dimension, and the thread dimension in the Block is made as BlockDim (x, y, 1), and the Block number after then splitting is:

Wherein, meshH is the grid number that simulating area length is divided, and meshW is the grid number that the simulating area width is divided, and BlockDim.x and BlockDim.y are respectively the component of thread dimension BlockDim (x, y, 1) on two dimensions in the Block.

Thread number in each Block is:

munThreads＝BlockDim.x*BlckDim.y

Just corresponding net point computing unit of each thread like this.In the Kernel function, (Fig. 2) can be located through the index BlockIdx of Block and the dimension BlockDim of Block in the position of this thread.

Find the solution for wave diffusion field in the gamut grid, its essence is frequency-domain function is transformed into the spatial domain.Therefore can adopt fast two-dimensional fourier transformation (FFT) to find the solution.Adopt NVIDIA in the CUDA framework, to provide based on the GPU Algorithm FFT transformation library CUFFT that divides and rules at this and find the solution, basic procedure is:

1) allocation space, create a CUFFT handle:

cufftPlan2d(&fftPlan,meshW,meshH,CUFFT_C2R)

2) carry out the FFT conversion according to the Plan that creates:

cufftExecC2R(fftPlan,(cufftComplex*)d_ht,g_hptr)

3) calculate the back that finishes and destroy handle:

cufftDestroy(fftPlan)

The present invention is a diffusion process of simulating the flight path ripples with particIe system.The parallelization computing of flight path particIe system also is in the Kernel of CUDA function, to find the solution.The attribute update algorithm that one big characteristics of particIe system are each particles is identical, and mutual between the particle takes place between adjacent element.Therefore the present invention is to be that parallel computation is carried out in the unit with each particle.

Because particle number is unrestricted, the maximum Block number of system then is a fixed value, and therefore, in order to improve the hardware utilization factor, the present invention designs the dimension of Block and Thread.Suppose that particle number is N in the initialized particIe system, the maximum Block number of system is nBlockSize, and so, after the two dimension fractionation, the Thread Count among each Block is:

numThreads＝min(nBlockSize,N)

Total Block number is:

In Kernel, it mainly is that the locus is upgraded that the attribute of ripples particle upgrades.The present invention calculates the initial position of particle according to the peak line equation of locus of Kelvin ship wave, and concrete computing method are following:

${\mathrm{<figure-callout\; id="0"\; label="X"\; filenames="HDA00001907810200021.png"\; state="\{\{state\}\}">X</figure-callout>}}_0=\mathrm{\&Delta;X}+V*T*\cos \mathrm{\&theta;}*(1-\frac{1}{2}{\mathrm{<figure-callout\; id="2"\; label="cos"\; filenames="HDA00001907810200021.png"\; state="\{\{state\}\}">cos</figure-callout>}}^2\mathrm{\&theta;})+\mathrm{Random}\left(\right)$

${\mathrm{<figure-callout\; id="0"\; label="Y"\; filenames="HDA00001907810200021.png"\; state="\{\{state\}\}">Y</figure-callout>}}_0=\mathrm{\&Delta;Y}+\frac{1}{2}*V*T*{\mathrm{<figure-callout\; id="2"\; label="cos"\; filenames="HDA00001907810200021.png"\; state="\{\{state\}\}">cos</figure-callout>}}^2\mathrm{\&theta;}*\sin \mathrm{\&theta;}+\mathrm{Random}\left(\right)$

Z ₀＝ΔZ+Random()

Wherein, (Δ X, Δ Y, Δ Z) is the initial offset position of particle, and V is the ship headway, and T is a hours underway, and θ is the line put on particle source and the peak line and the angle of particle source direction of motion (being the shipping agency direction).The motion conditions of particle after emission is by initial velocity, acceleration, and parameter determining such as duration, these parameters be at x, y, the value that three different directions of z are had nothing in common with each other.

In addition, also need be in Kernel through the life cycle disappearance and the generation that come simulation particle of control particle.The present invention has added a vital values and a dead value among particle property, can simulate the dynamic change (Fig. 3) of whole particIe system easily through control particle life and death.

As shown in Figure 2, when each particle is dead, the survival sign of particle is made as death, vital values makes zero, and gives a death value at random for it, begins dead timing, brings in constant renewal in dead the value.When dead timing finished, particle activated again, and the time-to-live makes zero; For it gives a vital values at random; And all properties of initialization particle, getting into a new life circulation, the duration of particle survival is also restarted timing; Particle position and vital values also begin to upgrade, and be dead until particle.

Along with the continuous variation of particle existence and dead state, the continuous conversion of survival particle number in the whole particIe system, and also the process of this variation is at random, so in follow-up drafting, just formed random effect more at random.Simultaneously, adopt this kind mode to upgrade, all particle numbers of whole particIe system are fixed, and the video memory space that occupies is fixed, and update algorithm is also fixed, and need not carry out alternately with the external world.

Among the present invention, particIe system does not need explicit drafting, but as a kind of virtual media of Simulated Water wave energy diffusion process, finally is reflected in the disturbance to sea grid height field.The energy of particle can be by the Z value reflection of its locus, and each particle can both make its present located grid produce the height that departs from of Δ H.Along with the renewal of particle position, the grid and the energy intensity of its influence constantly change, and be dead until particle.When all particles in the particIe system during, just formed the effect of flight path ripples diffusions on the whole all with this kind mode disturbance sea grid.

In order to make particle more level and smooth to the disturbance of grid, the grid of not only disturbance current location of each particle also will produce the disturbance that dies down gradually to the grid in the certain limit on every side.Suppose that the range of disturbance radius is L, then disturbance energy is obtained maximal value at the net point nearest apart from particle, and along with decay level and smooth with the increase of particle distance (Fig. 4).

Theoretical according to the linear wave of Gerstner; Wave is that the trigonometric function swash property a plurality of various amplitudes, cycle, angular frequency and first phase is formed by stacking; Therefore; The present invention uses cosine function to simulate the process of this level and smooth decay, and cosine function can obtain well to seamlessly transit (Fig. 5) in the edge of attenuation range L in the tangential direction level of wave crest and trough.

If dis is (p ₀, be the absolute value of the interior point of grid t t) to the horizontal range of some P0, λ is the wavelength of cosine function, then puts departing from highly of t acquisition P0 point particle energy disturbance generation to be:

${\mathrm{\&Delta;H}}_t=\mathrm{\&Delta;H}*(1+\cos \left(\frac{\mathrm{dis}({\mathrm{<figure-callout\; id="0"\; label="p"\; filenames="HDA00001907810200021.png"\; state="\{\{state\}\}">p</figure-callout>}}_0,t)}{L}\mathrm{\&lambda;\&pi;}\right)),\mathrm{dis}({\mathrm{<figure-callout\; id="0"\; label="p"\; filenames="HDA00001907810200021.png"\; state="\{\{state\}\}">p</figure-callout>}}_0,t)\&le;L$

ΔH _t＝0,dis(p ₀,t)＞L

Wherein, Δ H is the maximum deviation height, Δ H _tDepart from highly for what some t place's disturbance in the grid produced.

For make particle to the impact effect of grid on every side on the whole more at random, all λ is carried out random valued within the specific limits one time after new life cycle of the every entering of particle, make the whole disturbance effect of appearance more at random, more true to nature like this.

After having calculated the disturbance of particle to grid, need be with the sea grid height field that highly is added to that departs from of disturbance generation.Because the CUDA computing kernel of grid and particIe system adopts different thread dividing modes, the grid computing is that each net point is mapped to a kernel thread, and the particIe system computing then is that each particle is mapped to a kernel thread.Therefore, the stack computing of the two must be by Global memory, same video memory data of operation in different kernels.

Among the present invention, stack computing kernel adopts the kernel thread dividing mode of particIe system, and the disturbance that in kernel, each particle is produced is mapped to mesh coordinate.If the horizontal coordinate of particle be position (x, y), mesh-density is meshW*meshH, the length of net region is Width, width is Height, then the mesh coordinate that is mapped to of this particle is:

$\mathrm{GridX}=[\frac{\mathrm{position}.x*\mathrm{meshW}}{\mathrm{Width}}+\frac{\mathrm{meshW}}{2}]$

$\mathrm{GridY}=[\frac{\mathrm{position}.y*\mathrm{meshH}}{\mathrm{Height}}+\frac{\mathrm{meshH}}{2}]$

According to coordinate (GridX, GridY) can be in the Global memory of storage sea grid height field network for location lattice point rapidly, the disturbance of stacked particles highly gets final product on the altitude information of this net point.When all particles all execute this operation to all net points in the range of disturbance after, just formed the disturbance effect of flight path to the sea grid.

The playing up the stage of flight path parallelization emulation on the sea, the present invention uses GLSL (OpenGL Shading Language) that sea and flight path are unified to play up.Summit buffer object VBO (Vertex Buffer Object) through with OpenGL is mapped to the interoperability that the CUDA address space is accomplished CUDA kernel and OpenGL.Concrete steps are:

1) configuration CUDA equipment.

2) the OpenGL buffer object is registered to CUDA.

3) in Kernel, the video memory address space is mapped to the OpenGL buffer object.

4) can realize read-write through the video memory address space of mapping to the OpenGL buffer object.

In Vertex Shader, obtain after the vertex data of sea level height field, according to illumination and observation point position calculation surface normal and normalization, painted in Fragment Shader afterwards.

Advantage of the present invention and good effect:

At first, the present invention has designed the parallelization method for solving based on the Phillips ocean spectral model of CUDA framework.Secondly, the present invention proposes the parallelization ripples particIe system construction method a kind of CUDA framework under, can on GPU, utilize limited video memory space to carry out particle property efficiently and upgrade, not with CPU and internal memory carry out unnecessary alternately.In addition; The invention allows for ripples particle under a kind of CUDA framework to the Disturbance Model of sea grid; Designed the Disturbance Model of ripples particle to the sea grid; With the grid is that reality, particle are method empty, that actual situation combines, and it is the disturbance that decays gradually in the center that the sea grid is produced with the particle, and whole disturbance field can smoothly be docked with basic sea level height field.

The present invention is with all calculating and play up all and in GPU, carry out, and Model Calculation is constituted an overall flow with playing up in GPU, reduces the communication between internal memory and the video memory, gives full play to the GPU performance, significantly promotes simulation efficiency.Through above method, reach the good coordination of the effect of visualization sense of reality and the emulation real-time of sea track Simulation.

The inventor uses the method, in the described hardware environment of table 1, has carried out sea track Simulation experiment, the design sketch (Fig. 6) of the sea track Simulation operation that obtains.The inventor carries out same algorithm on CPU, the contrast situation that obtains operational efficiency is as shown in table 2.

Table 1:CUDA track Simulation experimental configuration

Table 2:CPU and GPU carry out the efficient contrast

Visible by table 1 and table 2, even if under the situation that hardware performance is not given prominence to, the parallel processing method among the present invention also considerably beyond traditional CPU serial processing method, has good real time performance on performance.Particularly under the thinner situation of sea grid granularity, because the calculated amount of net point real-time operation is bigger, the advantage of method for parallel processing of the present invention shows more outstandingly.

Description of drawings

Fig. 1 is based on the sea flight path parallelization emulation mode basic step synoptic diagram of particIe system among the present invention;

Fig. 2 is that the kernel task of grid computing in sea among the present invention splits synoptic diagram;

Fig. 3 is existence of GPU particIe system and the dead constitutional diagram that designs among the present invention;

Fig. 4 is the particIe system that designs among the present invention disturbing influence synoptic diagram (top view) to grid;

Fig. 5 is the particIe system that designs among the present invention disturbing influence synoptic diagram (side view) to grid;

Fig. 6 is the resulting experiment effect figure of method of the present invention.

Embodiment

The overall coverage range in zone, emulation sea is Length*Width; Length is an emulation sea zone length; Width is an emulation sea peak width, and mesh-density is meshH*meshW, and meshH is the grid number that simulating area length is divided; MeshW is the grid number that the simulating area width is divided, and the particle number in the flight path particIe system is ParticleNumber.Maximum Thread (CUDA kernel thread) number in system's video card hardware in each Block is BlockSize.

The first step: in CPU, carry out some initialized work.Comprise: initialization CUDA equipment, application video memory, establishment FFT Plan (the execution parameter configuration of Fast Fourier Transform (FFT)), establishment VBO, initialization particIe system.

Second step: the sea grid is carried out the parallelization computing.

A) find the solution the sea level height field of the Phillips wave frequency spectrum of a certain particular moment.

The wave wavelength coverage that the sea grid need cover is:

$(\frac{\mathrm{Length}}{\mathrm{meshH}},\frac{\mathrm{Width}}{\mathrm{meshW}})~(\mathrm{Length},\mathrm{Width})$

Expression wave wavelength coverage on the length and width both direction (upper limit is exactly the length and the width in zone, whole sea, is maximal value) in simulating area in the bracket

The sea level height field Superposition Formula that is obtained by Phillips wave spectral change is:

$\mathrm{\&Delta;y}(P,t)=\underset{K}{\mathrm{\&Sigma;}}A\left(K\right)*e^{\mathrm{j\&Delta;Hs}\left(K\right)}*e^{\mathrm{j\&Delta;Ht}\left(K\right)}$

$=\underset{m=-\frac{\mathrm{<figure-callout\; id="2"\; label="meshH"\; filenames="HDA00001907810200021.png"\; state="\{\{state\}\}">meshH</figure-callout>}}{2}}{\overset{\frac{\mathrm{<figure-callout\; id="2"\; label="meshH"\; filenames="HDA00001907810200021.png"\; state="\{\{state\}\}">meshH</figure-callout>}}{2}}{\mathrm{\&Sigma;}}}\underset{n=-\frac{\mathrm{<figure-callout\; id="2"\; label="meshW"\; filenames="HDA00001907810200021.png"\; state="\{\{state\}\}">meshW</figure-callout>}}{2}}{\overset{\frac{\mathrm{meshW}}{2}}{\mathrm{\&Sigma;}}}\left[A\left(K\right)e^{\mathrm{j\&Delta;Hs}\left(K\right)}\right]*e^{\mathrm{j\&Delta;Ht}\left(K\right)}---(3-10)$

$=\underset{m=-\frac{\mathrm{<figure-callout\; id="2"\; label="meshH"\; filenames="HDA00001907810200021.png"\; state="\{\{state\}\}">meshH</figure-callout>}}{2}}{\overset{\frac{\mathrm{<figure-callout\; id="2"\; label="meshH"\; filenames="HDA00001907810200021.png"\; state="\{\{state\}\}">meshH</figure-callout>}}{2}}{\mathrm{\&Sigma;}}}\underset{n=-\frac{\mathrm{<figure-callout\; id="2"\; label="meshW"\; filenames="HDA00001907810200021.png"\; state="\{\{state\}\}">meshW</figure-callout>}}{2}}{\overset{\frac{\mathrm{meshW}}{2}}{\mathrm{\&Sigma;}}}[A\left(K\right)e^{j*\mathrm{\&omega;t}\left(K\right)}+A(-K)e^{j*\mathrm{\&omega;t}(-K)}]*e^{j(K\&CenterDot;P)}$

Wherein:

1) K is the propagation vector of wave, and its direction cosine equal the direction of propagation of wave, with the point on the frequency domain sea to the vector of initial point as the wave propagation direction, so just can contain the wave on all directions.P is the horizontal coordinate position of the net point of current calculating.

$K=(\frac{2\mathrm{\&pi;n}}{\mathrm{Length}},\frac{2\mathrm{\&pi;m}}{\mathrm{Width}}),$ $P=(\frac{n*\mathrm{Length}}{\mathrm{meshW}},\frac{m*\mathrm{Width}}{\mathrm{meshH}})$

2) A is self-defining whole amplitude amount of zoom, is used for dynamically controlling the whole height of wave.

3) m, n represent wave current wavelength on the length and width both direction in simulating area respectively.E is the mathematics natural constant.

4) Ae ^{J Δ Hs}And Ae ^{J Δ Ht}It is fourier basis functions.

5) time-angle frequencies omega t (K) is tried to achieve by Space Angle frequency (being the mould of K).ω t (K) has different computing method on different seas, and establishing g is acceleration of gravity, and ω s is the Space Angle frequency, and tan () is a tan, then:

The degree of depth is in the water of d: $\mathrm{\&omega;\; t}\left(K\right)=\sqrt{g*\mathrm{\&omega;\; s}*\mathrm{Tan}(\mathrm{\&omega;\; s}*d)},$

In the deep-sea,, therefore remove tan (ω s*d), become because the vibration frequency of wave does not receive effect of depth: $\mathrm{\&omega;\; t}\left(K\right)=\sqrt{g*\mathrm{\&omega;\; s}}.$

Shallow water surface:

(L is the parameter that is used for controlling the tension degree of the water surface).

Can increase the surperficial High-frequency water ripple of low frequency wave through such conversion, form the more water surface of horn of plenty of effect.

B) carry out the wave diffusion field that CUFFT finds the solution Phillips wave frequency spectrum.

1) allocation space, create a CUFFT (CUDA Fast Fourier Transform (FFT)) handle:

CufftPlan2d (&fftPlan; MeshW, meshH, the CUFFT_C2R) (configurations such as execution parameter of the Fast Fourier Transform (FFT) that fftPlan representes to create; MeshH is the grid number that simulating area length is divided; MeshW is the grid number that the simulating area width is divided, and CUFFT_C2R is the enumerated variable of building in the CUDA, the execution pattern of expression Fast Fourier Transform (FFT))

2) carry out the FFT conversion according to the Plan that creates:

CufftExecC2R (fftPlan; (cufftComplex*) d_ht; G_hptr) (the configurations such as execution parameter of the Fast Fourier Transform (FFT) that fftPlan representes to create; CufftComplex representes that Fast Fourier Transform (FFT) is to carry out with complex data type, and d_ht is the data before the conversion, and g_hptr is the data after the conversion)

3) calculate the back that finishes and destroy handle:

cufftDestroy(fftPlan)

C) result of calculation is outputed to the Global Memory in the video memory.

The 3rd step: carry out the parallelization computing of flight path particIe system.

A) calculating and renewal particle property.Each particle all carries out position renewal according to following initial position:

${\mathrm{<figure-callout\; id="0"\; label="X"\; filenames="HDA00001907810200021.png"\; state="\{\{state\}\}">X</figure-callout>}}_0=\mathrm{\&Delta;X}+V*T*\cos \mathrm{\&theta;}*(1-\frac{1}{2}{\mathrm{<figure-callout\; id="2"\; label="cos"\; filenames="HDA00001907810200021.png"\; state="\{\{state\}\}">cos</figure-callout>}}^2\mathrm{\&theta;})+\mathrm{Random}\left(\right)$

${\mathrm{<figure-callout\; id="0"\; label="Y"\; filenames="HDA00001907810200021.png"\; state="\{\{state\}\}">Y</figure-callout>}}_0=\mathrm{\&Delta;Y}+\frac{1}{2}*V*T*{\mathrm{<figure-callout\; id="2"\; label="cos"\; filenames="HDA00001907810200021.png"\; state="\{\{state\}\}">cos</figure-callout>}}^2\mathrm{\&theta;}*\sin \mathrm{\&theta;}+\mathrm{Random}\left(\right)$

Z ₀＝ΔZ+Random()

Wherein, (Δ X, Δ Y, Δ Z) is the initial offset position of particle, and V is a ship navigation speed, and T is a hours underway, and θ is the line put on particle source and the peak line and the angle of shipping agency direction.Random () is the function that obtains random number, and x0, y0, z0 represent the current particle position of trying to achieve respectively.Position renewal is all simulated with accelerated motion on X, Y, three directions of Z.

B) particle position is mapped to the sea grid system.Concrete grammar is:

If the horizontal coordinate of particle be position (x, y), mesh-density is meshW*meshH, the length of net region is Width, width is Height, then the mesh coordinate that is mapped to of this particle is:

$\mathrm{GridX}=[\frac{\mathrm{position}.x*\mathrm{meshW}}{\mathrm{Width}}+\frac{\mathrm{meshW}}{2}]$

$\mathrm{GridY}=[\frac{\mathrm{position}.y*\mathrm{meshH}}{\mathrm{Height}}+\frac{\mathrm{meshH}}{2}]$

According to coordinate (GridX, GridY) can be in the Global memory of storage sea grid height field network for location lattice point rapidly.

C) calculate particle to the disturbance of grid around the position disturbance stack computing in the line range of going forward side by side.

If dis is (p ₀, be that some t is to the absolute value of the horizontal range of some P0 in the grid t), λ is the wavelength of cosine function, and L is the maximum perturbation scope, and Δ H is the maximum deviation height.Then putting departing from highly of t acquisition P0 point particle energy disturbance generation is:

${\mathrm{\&Delta;H}}_t=\mathrm{\&Delta;H}*(1+\cos \left(\frac{\mathrm{dis}({\mathrm{<figure-callout\; id="0"\; label="p"\; filenames="HDA00001907810200021.png"\; state="\{\{state\}\}">p</figure-callout>}}_0,t)}{L}\mathrm{\&lambda;\&pi;}\right)),\mathrm{dis}({\mathrm{<figure-callout\; id="0"\; label="p"\; filenames="HDA00001907810200021.png"\; state="\{\{state\}\}">p</figure-callout>}}_0,t)\&le;L$

ΔH _t＝0,dis(p ₀,t)＞L

For make particle to the impact effect of grid on every side on the whole more at random; All its wavelength X that influences forcing function is carried out random valued within the specific limits one time after new life cycle of the every entering of particle; Make the whole disturbance effect of appearance more at random, more true to nature like this.

D) operation result is outputed to Global Memory.

The 4th step: in GPU Shader, use GLSL that environment, illumination etc. are carried out in the sea and play up painted.

A) the Fresnel effect computing method played up of illumination are:

TFresnel＝r+(1-r)×pow(1.0-dot(viewDir,normal),5.0)

Wherein, r is the reflection coefficient of air to water.Dot be build in the GLSL (OpenGL Shading Language) ask the dot-product function, pow () is the exponentiation function of building in the GLSL.The final water surface color that obtains is:

FinalColor＝waterColor+reflectColor×TFresnel

+refractColor×(1-TFresnel)

Wherein, the meaning of each parametric representation is:

1) TFresnel representes that the Fresnel coefficient (is the design factor of Fresnel effect.The Fresnel effect is meant: between air and water surface two media, the angle of the observer and the water surface is more little, and reflecting effect is obvious more; Angle is big more, and refraction effect is obvious more);

2) viewDir is the direction of eye-observation;

3) normal is a normal vector;

4) FinalColor is the color of final water;

5) waterColor is the vertex color of initial water;

6) reflectColor is the water-reflected color;

7) refractColor is the overwater refraction color.

B) last, handle the color transition on the spray and sea.

The wave vertex color is carried out alpha mix, make the place that wave is high more, the color of water is white more.So just can render the effect of wave whiting on the whole.Computing method are:

FinalColor＝mix(FinalColor,whiteColor,alpha×hightRate)

Wherein, FinalColor is the final water surface color after play up the front, and whiteColor is the white colour of wave, and alpha is a transparency; HightRate is the height ratio of wave, and mix () is the mixed function of building in GLSL (the OpenGL Shading Language) language.

Claims (7)

1. sea flight path parallelization emulation mode based on particIe system is characterized in that may further comprise the steps:

(1) in CPU, carries out some initialized work: specifically comprise: initialization CUDA equipment, application video memory, establishment FFT Plan, establishment VBO, initialization particIe system;

(2) in GPU, carry out the parallelization computing of sea grid: specifically comprise: find the solution the Phillips wave frequency spectrum of a certain particular moment the sea level height field, carry out CUFFT find the solution Phillips wave frequency spectrum the wave diffusion field, result of calculation is outputed to the Global Memory in the video memory;

(3) in GPU, carry out the parallelization computing of flight path particIe system: specifically comprise: calculate and upgrade particle property, with particle position be mapped to the sea grid system, calculate particle to the disturbance of grid around the position go forward side by side disturbance stack computing in the line range, operation result is outputed to Global Memory;

(4) in GPU Shader the utilization GLSL to the sea carry out environment, illumination is played up, and forms dynamic sea flight path effect.

2. the sea flight path parallelization emulation mode based on particIe system as claimed in claim 1; It is characterized in that; Among the initialization particIe system, particle property comprises a vital values and a dead value in the said step (1), and initial value all is endowed the random integer value of a positive number.

3. the sea flight path parallelization emulation mode based on particIe system as claimed in claim 1 is characterized in that, the method for calculating in the said step (3) and renewal particle property is:

When each particle is dead, the survival sign of particle is made as death, and gives a death value at random for it, begin dead timing, bring in constant renewal in dead the value; When dead timing finished, particle activated again, for it gives a vital values at random, got into a new life circulation, and the duration of particle survival is also restarted timing.

4. the sea flight path parallelization emulation mode based on particIe system as claimed in claim 1 is characterized in that, the CUDA thread fractionation mode that the sea level height field of Phillips wave frequency spectrum is found the solution in the parallelization in the said step (2) is:

The mesh-density of supposing simulating area is meshW*meshH, splits according to two dimension, and the thread dimension in the Block of CUDA kernel is made as BlockDim (x, y, 1), and the Block number after the fractionation is:

Wherein, meshH is the grid number that simulating area length is divided, and meshW is the grid number that the simulating area width is divided, and BlockDim.x and BlockDim.y are respectively the component of thread dimension BlockDim (x, y, 1) on two dimensions in the Block;

Thread number in each Block is:

munThreads＝BlockDim.x*BlockDim.y。

5. the sea flight path parallelization emulation mode based on particIe system as claimed in claim 1 is characterized in that the particle position in the said step (3) is mapped to sea mesh coordinate point methods and is:

$\mathrm{GridX}=[\frac{\mathrm{position}.x*\mathrm{meshW}}{\mathrm{Width}}+\frac{\mathrm{meshW}}{2}]$

$\mathrm{GridY}=[\frac{\mathrm{position}.y*\mathrm{meshH}}{\mathrm{Height}}+\frac{\mathrm{meshH}}{2}]$

Wherein, the horizontal coordinate of particle be position (x, y), mesh-density is meshW*meshH, the length of net region is Width, width is Height, the net point coordinate that obtains of mapping for (GridX, GridY).

6. the sea flight path parallelization emulation mode based on particIe system as claimed in claim 1 is characterized in that, in the said step (3), calculates particle and the disturbance of grid around the position is gone forward side by side the method for disturbance stack computing is in the line range:

ParticIe system is not carried out explicit drafting, but, finally be reflected in the disturbance to sea grid height field as a kind of virtual media of Simulated Water wave energy diffusion process; The energy of particle is by the Z value reflection of its locus, and each particle can both make its present located grid produce the height that departs from of Δ H; Along with the renewal of particle position, the grid and the energy intensity of its influence constantly change, and be dead until particle;

Particle is to the departing from highly of disturbance generation of grid is on every side:

${\mathrm{\&Delta;H}}_t=\mathrm{\&Delta;H}*(1+\cos \left(\frac{\mathrm{dis}(p_0,t)}{L}\mathrm{\&lambda;\&pi;}\right)),\mathrm{dis}(p_0,t)\&le;L$

ΔH _t＝0,dis(p ₀,t)＞L

Wherein, dis (p ₀, be that some t arrives a some P in the grid t) ₀The absolute value of horizontal range, λ is the wavelength of cosine function, L is the maximum perturbation scope, Δ H is the maximum deviation height, Δ H _tDepart from highly for what some t place's disturbance in the grid produced; And, after new life cycle of the every entering of particle, all λ is carried out random valued within the specific limits one time.

7. the sea flight path parallelization emulation mode based on particIe system as claimed in claim 1 is characterized in that, in the said step (4), the color transition on the spray and sea is played up disposal route is:

FinalColor＝mix(FinalColor,whiteColor,alpha×hightRate)

Wherein, FinalColor is the final water surface color after play up the front, and whiteColor is the white colour of wave, and alpha is a transparency, and hightRate is the height ratio of wave, and mix () is the mixed function of building in the GLSL language.

Images