CN109214082A - A kind of high resolution numerical simulation method of near field underwater blast wave load - Google Patents
A kind of high resolution numerical simulation method of near field underwater blast wave load Download PDFInfo
- Publication number
- CN109214082A CN109214082A CN201811018574.6A CN201811018574A CN109214082A CN 109214082 A CN109214082 A CN 109214082A CN 201811018574 A CN201811018574 A CN 201811018574A CN 109214082 A CN109214082 A CN 109214082A
- Authority
- CN
- China
- Prior art keywords
- grid
- zoning
- boundary
- formula
- explosive
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Pending
Links
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F30/00—Computer-aided design [CAD]
- G06F30/20—Design optimisation, verification or simulation
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- Computer Hardware Design (AREA)
- Evolutionary Computation (AREA)
- Geometry (AREA)
- General Engineering & Computer Science (AREA)
- General Physics & Mathematics (AREA)
- Management, Administration, Business Operations System, And Electronic Commerce (AREA)
Abstract
This invention proposes a kind of high resolution numerical simulation method of near field underwater blast wave load, belongs to underwater explosion numerical simulation technology field.The present invention carries out high-precision numerical discretization to the non-linear Eulerian equation of meter and strong compressibility using the golden method (RKDG) of interruption gal the Liao Dynasty and three rank TVD Runge-Kutta methods；The physical quantity of gas-liquid interface two sides is pre-processed using ghost fluid method (GFM), effectively cuts down non-physical concussion caused by being interrupted by physical quantity；Using the complex topology structure variation occurred on Level Set Method (LSM) treated substance interface, the accurate capture to moving interface is realized.Method proposed by the present invention can effectively simulate the generation and communication process of near field underwater blast wave load；Numerical simulation result and classical empirical equation, experimental result can quite well, demonstrate this method and handling strong compressible, strong discontinuity, there is some superiority in transient state strong nonlinearity problem.
Description
Technical field
This invention proposes a kind of high resolution numerical simulation method of near field underwater blast wave load, belongs to underwater
Explosion numerical simulation technology field.
Background technique
It is increasingly prominent with marine safety problem, countries in the world competitively develop high speed, big dose, precise guidance elder generation
Into Underwater Battery, great threat is constituted to large surface warships such as submarine and aircraft carriers.The shock wave that Underwater Battery explosion generates
Load, pressure peak may be up to GPa magnitude, globality will be caused to injure destruction to naval ship structure.Therefore, effectively simulation is underwater
Explosion wave load character is of great significance for naval vessels safety and Protection.In addition, Underwater Battery Performance Evaluation, ocean
Also there is great demand in the fields such as military installations construction to this technology.To sum up, method proposed by the present invention has important skill
Art background.
Traditional underwater explosion research mode mainly includes that theory analysis and experiment are tested.
Theory analysis is by establishing and solving the Closure equation including fluid motion equation, state of matter equation etc.
Group obtains the analytic solutions of related physical quantity.This method in the research of early stage have certain guidance meaning, but simultaneously there is also
More defect: (1) theoretical model is relatively simple, and research achievement can not solve actual engineering problem.(2) draw in derivation process
The artificial it is assumed that there are certain errors between actual condition of some simplicity is entered.(3) it is related to the mathematic(al) manipulation of many complexity,
Research difficulty is big, the period is long.
Experiment test is can the most directly to provide reliably for other correlative studys with accurate research method, result
Data basis, the empirical equation of many classics are also derived from the summary to abundant experimental results.The deficiency of experiment test is:
(1) duration of underwater explosion process, pressure peak reached as high as GPa magnitude in ms magnitude, and experiment measurement difficulty is very big.
(2) explosion test has greatly destructive, and the loss of measurement of correlation instrument is very big, and experimental cost is higher.Meanwhile testing participant
There is a certain security risk by member.(3) often there is uncertain factor in experiment, repeatability is poor.
With the continuous development of computer technology, numerical simulation has become the research essential method of underwater explosion problem
One of.The advantage of the technology is: (1) solve the problems, such as the high-efficient of Practical Project, it is adaptable.(2) compared with experiment, numerical value
That simulates is low in cost, and without security risk.(3) easy to operation.Currently, numerical simulation field of exploding under water, mainly
There are two kinds of research modes: first is that using existing business software, such as LS-DYNA, AUTODYN, ABAQUS, though this mode
Right simple possible, but its numerical result is often not accurate enough.Second is that the certain methods proposed using domestic and foreign scholars, such as finite element
Grid model, mesh free SPH method, ALE algorithm etc., achieve many positive achievements.However, there are still two at present for these methods
Item is insufficient: (1) the case where applicable far field is exploded mostly more lacks the research of near field underwater explosion；(2) many methods are adopted
With low order precision format, the accuracy of calculated result is affected.Precision guided weapon high speed development, naval vessels security protection demand under water
Under the overall background of surge, the above two o'clock deficiency is urgently made up.
For these deficiencies, the invention proposes the high resolution numerical simulation methods of near field underwater blast wave load.
Specifically, this method is a kind of based on the golden method (RKDG) of interruption gal the Liao Dynasty, while Level Set Method (LSM) and virtually is coupled
The unified algorithm of stream method (GFM).RKDG method can fluid motion equation carry out high-precision it is discrete；LSM method can have
Imitate the complex topology structure variation at treated substance interface；GFM can effectively near inhibiting substances interface Nonphysical Oscillation.Pass through
This method can carry out high-precision analog near field underwater blast wave load.
Summary of the invention
The invention proposes a kind of high resolution numerical simulation methods of near field underwater blast wave load.Its object is to
The relative absence of current near field underwater blast wave load value simulation field is made up, and improves the precision of calculated result.
Above-mentioned purpose that the invention is realized by the following technical scheme:
Step 1: determining the zoning of underwater explosion, establish rectangular coordinate system, be m × n grid by the region division
(wherein, m indicates that the direction x number of grid, n indicate the direction y number of grid).
Step 2: each physical quantity, correlation function in definition zoning under original state, and establish fluid governing equation
And state of matter equation, comprising:
(a) interface function of each substance, i.e. Level Set function:
Wherein, Γ (0) indicates that zero contour surface, i.e. material interface, d (x, y, Γ (0)) indicate that a bit (x, y) is arrived in computational domain
The distance of material interface, Ω1、Ω2The respectively region of material interface two sides.
(b) physical parameter of each substance: density, speed, the pressure of explosive, original state are lauched close under original state
Degree, speed, pressure.
(c) according to the physical characteristic of underwater explosion field, viscous adaptable nothing, irrotationality, meter and can fluid governing equation: be established
The fluid governing equation of compressibility:
Ut+ ▽ F (U)=S (U) (2)
Wherein,
(d) state equation (Tait equation) and relevant parameter of water:
Wherein, N, AW、BWIt is constant relevant to the physical property of water, ρW0Indicate the density of water in the initial state, pW、
ρWRespectively indicate the pressure and density of water.
(e) detonation products state equation (JWL equation) and relevant parameter:
Wherein, AE、BE、R1、R2, ω be constant relevant to explosive property, ρE0Indicate explosive in the initial state close
Degree, ρEIndicate the density of detonation product, PEIndicate the density of detonation product, E0Indicate the specific internal energy of detonation product.
Step 3: expanding zoning, boundary condition is set according to different operating conditions: since boundary lacks difference scheme institute
Required mesh point, therefore zoning is expanded to N layers outward.According to the format that this method uses, at least need to open up outward
Open up 3 layers of grid (i.e. N >=3).
For free boundary condition, have: UΓ=U-1,U-i=U-i+1, i=2N, wherein UΓIndicate zoning side
Random physical quantity on boundary's grid, U-iIndicate the random physical quantity on virtual grid.Specifically, first outside the boundary of zoning
Physical quantity on layer virtual grid is equal to the physical quantity on the boundary mesh of zoning；Second layer virtual net outside the boundary of zoning
Physical quantity on lattice is equal to the physical quantity on first layer virtual grid；Physics outside the boundary of zoning on third layer virtual grid
Amount is equal to the physical quantity on second layer virtual grid, and so on.
For wall boundary condition, have:Wherein, UΓIndicate zoning
Random physical quantity on boundary mesh in addition to speed, U-iIndicate the random physical quantity on virtual grid in addition to speed, VΓTable
Show the speed on the boundary mesh of zoning, V-iIndicate the speed on virtual grid.Specifically, first outside the boundary of zoning
Layer virtual grid on perpendicular in the speed and boundary mesh on boundary perpendicular to boundary velocity vector and be zero, other physical quantitys
It is equal；Outside the boundary of zoning on second layer virtual grid perpendicular to boundary speed and first layer virtual grid on perpendicular to side
The velocity vector on boundary and be zero, other physical quantitys are equal；Perpendicular to boundary on third layer virtual grid outside the boundary of zoning
On speed and second layer virtual grid perpendicular to boundary velocity vector and be zero, other physical quantitys are equal, and so on.
Step 4: choose CFL parameter and calculate time step: CFL parameter is the constant between 0 to 1；It is current to calculate
It is indicated with t=n Δ t total time used, n is to calculate step number, and time step Δ t can be obtained by following formula:
Δ x is the mesh spacing in the direction x in zoning；Δ y is the mesh spacing in the direction y in zoning；U is to calculate
Velocity component of the current grid in the direction x in region；Velocity component of the v for current grid in zoning in the direction y；C is meter
Calculate the velocity of sound in region at current grid.
Step 5: for the parameter assignment on zoning and virtual grid: as t=0, the calculating that will be provided in step 1
The physical parameter of each substance and interface function value in region, the parameter initial value as grid each in zoning；It, will as t > 0
The physical parameter of each substance and interface function value in the zoning being calculated in step 13, as net each in zoning
Initial value of the lattice parameter when next step calculates.For virtual grid, then according to the method established in step 3, pass through interpolation and extrapolation
Assignment is carried out to it.
Step 6: solving material interface normal direction: according to Level Set equation, to the normal of material interface each in zoning
Direction is solved, solution formula are as follows:
Step 7: solving Riemannian problem and to grid assignment near material interface: establishing explosive along material interface normal direction
Riemannian problem between water simultaneously solves the physical quantity between explosive and water at material interface.It is obtained with solution Riemannian problem
Physical quantity carries out assignment close to one layer of grid of material interface in explosive region.Meanwhile by extrapolate to explosive region other than
Grid carry out assignment, method particularly includes: at the flow field of treated substance interface explosive region side, by the material interface other side
Grid be set as Level Set method grid, the one of the pressure on these Level Set method grids and normal velocity and the side real units
It causes, and across the material interface interpolation and extrapolation of real fluid that tangential velocity and density then pass through explosive region acquires.
Step 8: by the discrete Eulerian equation of RKDG method, the numerical flux of positive pressure of shock wave is solved, specifically: to (2)
Formula both sides are obtained with multiplied by tentative function Φ (x, y), and in grid cell K upper integral:
After carrying out integration by parts to (7) formula, the line integral and Line Integral in Gauss quadrature formula solution formula are utilized.Due to net
Positive pressure of shock wave flux is discontinuous at lattice elementary boundary, therefore uses numerical fluxInstead of real fluxesUse approximate solutionInstead of accurately solving U, with the tentative function after numerical approximationInstead of tentative function Φ,
It can thus be concluded that:
Wherein, ω indicates weight coefficient, | e | indicate the length of integral domain grid cell boundary, | K | indicate integral domain
The area of grid cell.
For the numerical flux in (8) formulaIt is solved using (9) formula:
Wherein, α indicates Jacobi matrixThe maximum maximum value of characteristic value, U+(xel,yel, t) and it indicates
The U value of grid cell boundary, U are approached from positive direction-(xel,yel, t) and indicate the U that grid cell boundary is approached from negative direction
Value.For the pressure of explosive region detonation product, obtained by (4) formula in step 1.
For the approximate solution in (8) formulaIt is solved using (10) formula:
Wherein, for grid cell [xi-1/2,xi+1/2]×[yj-1/2,yj+1/2],Δxi=xi+1/2-xi-1/2, Δ yj=yj+1/2-yj-1/2。
The semi-discrete scheme that (9), (10) two formulas are brought into available (2) formula of (8) formula, is denoted as:
Ut=R U (11)
Wherein, the Discrete Operator of R U representation space derivative term.
Step 9: it is discrete by time term progress of the three rank TVD Runge-Kutta methods to (11) formula equation left side, it obtains
(2) the complete discrete scheme of formula, and then solve each physical quantity of detonation products in subsequent time zoning.Specifically
Formula are as follows:
Step 10: similar operations are carried out with step 7, solve Riemannian problem again and to grid assignment near material interface:
The Riemannian problem established between explosive and water along material interface normal direction simultaneously solves between explosive and water at material interface
Physical quantity.The physical quantity obtained with solution Riemannian problem carries out assignment close to one layer of grid of material interface in air section.
Meanwhile assignment is carried out to the grid other than air section by extrapolating, method particularly includes: in treated substance interfacial air region one
When the flow field of side, Level Set method grid is set by the grid of the material interface other side, the pressure on these Level Set method grids
It is consistent with the side real units with normal velocity, and tangential velocity and density then pass through the real fluid in explosive region across object
The extrapolation of matter interface interpolation acquires.
Step 11: on the basis of step 10, carrying out the operation similar with step 8, step 9.But the difference is that: step
In rapid 8, for the pressure of water area, obtained by (3) formula in step 1；What solution obtained in step 9 is that subsequent time calculates
Each physical quantity of water in region.
Step 12: the space derivation item of the Level Set equation of motion is carried out discrete:
For the lateral numerical flux of the interface function of explosive and water, can be solved by (14) formula:
Wherein,Indicate grid cell in zoningIn the interface function transverse direction of upper explosive and water
Numerical flux,Indicate the interface function numerical flux in lateral positive direction on the grid cell,Indicate the net
Interface function numerical flux in lateral negative direction on lattice unit；Indicate grid cell (i-r+s, j) in zoning
The flux of the lateral positive direction of the interface function of upper explosive and water；It indicates in zoning on grid cell (i-r+s, j)
The flux of the lateral negative direction of the interface function of explosive and water.Design parameter can be solved by (15), (16) formula:
Wherein,Indicate the interface function value of explosive and water on grid cell (i, j) in zoning.
For the Vertical Numerical flux of the interface function of explosive and water, method for solving is identical.
By step 12, the semi-discrete scheme of available (13) formula is denoted as:
Wherein,The Discrete Operator of representation space derivative term:
Step 13: by the time-derivative item of the three discrete Level Set equations of motion of rank TVD Runge-Kutta method,
The complete discrete scheme of (13) formula is obtained, and then solves the interface function value of each grid cell in subsequent time zoning.
Step 14: on the basis of step before, solve and export the pressure value of all grid cells in zoning: if
Grid cell is in detonation products region, then is acquired by (4) the formula calculating in step 1；If grid cell is in water
In region, then acquired by (3) formula in step 1.
Step 15: the current calculating total time t for calculating total time t=n Δ t and setting used of judgement*Relationship: if t
> t*, then terminate to calculate, export the pressure value of grid cell in zoning；If t > t*, then it is back to step 3, continues operation.
Beneficial effect
Underwater blast wave load high resolution numerical simulation method near field proposed by the present invention and traditional underwater explosion are rushed
It hits wave load method for numerical simulation to compare, have the advantage that
(1) present invention carries out Eulerian equation space derivation item using RKDG method discrete, using three rank TVD Runge-
Kutta method is discrete to the progress of Eulerian equation time-derivative item, has higher precision compared to conventional method；
(2) physical quantity discontinuous problem of the present invention using the method treated substance interface the DGFM two sides in GFM method group, height
Maintain to resolution ratio the discontinuous nature of material interface；DGFM method can effectively inhibit to solve compared to original GFM method
The non-physical reforming phenomena that interface generates when Eulerian equation, ensures the stable operation of algorithm；
(3) the Level Set method that the present invention uses can accurately capture moving interface, with traditional Lagrange circle
Face method for tracing is compared, and Level Set method more stable can effectively handle complex material interface and its topological structure change
Change, when the Nonlinear Large Deformations such as broken, tearing, fusion for especially occurring during solving the problems, such as underwater explosion, has obviously
Advantage.
(4) method proposed by the present invention can effectively simulate near field underwater blast wave load, accurately reflect underwater quick-fried
For fried shock wave in the propagation law near field, numerical result and experimental result meet ground preferably；With the tradition mainly for middle far field
Method for numerical simulation is compared, and the scope of application is bigger, more adaptable.
Detailed description of the invention
Fig. 1 is Riemannian problem schematic diagram
Fig. 2 is to construct local Riemannian problem schematic diagram
Fig. 3 is Level Set method schematic diagram
Fig. 4 is DGFM method schematic diagram
Fig. 5 is the initial geometric model figure of example
Fig. 6 is the pressure time-process curve graph in example at measuring point B
Fig. 7 is underwater blast wave simulation effect picture
Specific embodiment
The present invention is specifically described below with reference to example, the initial geometric model in example is presented in figure 5, is seen
Pressure time-process curve at measuring point B is presented in figure 6, and underwater blast wave simulation effect is presented in fig. 7.
Information needed for steps 1 and 2 is provided first: as shown in Fig. 5, using O point as coordinate origin, zoning size
A=1000mm, b=1800mm, transverse grid number are 500, and longitudinal grid number is 900, mesh spacing 2mm.The total time of setting
For 0.6ms, CFL number takes 0.2.Powder charge model is spherical charge, and centre coordinate is (0,0.9), radius r=11.15mm.Take measuring point
The coordinate of B is (0.5,0.9).Explosive area equation uses JWL state equation, and water area uses Tait state equation, specific object
The selection of reason parameter see the table below.
The JWL state equation parameter of 1 TNT explosive of table
| AE/Pa | BE/Pa | R1 | R2 | ω | ρE0/(kg/m3) | E0/(J/m3) |
| 3.712×1011 | 3.231×1109 | 4.15 | 0.95 | 0.3 | 1630 | 4.29×106 |
2 Tait state equation parameter of table
| N | Aw/Pa | Bw/Pa | ρW0/(kg/m3) |
| 7.15 | 1.0×105 | 3.31×108 | 1000 |
Steps 1 and 2 are completed as a result,.
Then, step 3 is carried out to step 15, obtains the numerical simulation result of underwater blast wave load.
Calculated result analysis:
Attached drawing 6 is the pressure time-process curve at observation point B, and solid line is the experiment value that pressure sensor measures, dotted line in figure
The pressure value calculated for numerical value.It can be seen that shock motion is to observation point B, and pressure increases suddenly at this time when t ≈ 0.30ms
Greatly to about 19.8MPa, the attenuation law of exponential form is presented in shock wave in water later.By correlation curve it is found that by we
The pressure trend and peak value that method is calculated match with experimental result, error is less than 5%, it was demonstrated that the numerical simulation side
Method is feasible effective.
Attached drawing 7 is that the underwater blast wave that this method is calculated simulates effect picture, is clearly labelled with impact in figure
Wave surface, air detonation product interface and cohesion rarefaction wave region.The biography as can be seen that wave surface and material interface explode under water
Good geometrical symmetry is maintained during broadcasting, interface svelteness meets objective physical rule, illustrates that this method can be quasi-
Really effectively capture underwater explosion material interface motion conditions.
The drawings and the specific embodiments two parts are combined above, and the present invention is elaborated.But it needs to state
Be: these pictures and text are merely illustrative, and any restrictions are not constituted to protection scope of the present invention, without departing from this hair
In the case where bright spirit and scope, those skilled in the art carry out a variety of improvement, equivalencing to the present invention and embodiments thereof
Or modification, these fall within the protection scope of the present invention.Protection scope of the present invention is subject to the appended claims.
Claims (4)
1. a kind of high resolution numerical simulation method of near field underwater explosion loading, it is characterised in that: its operating process are as follows:
Step 1: determining the zoning of underwater explosion, establish rectangular coordinate system, be m × n grid (its by the region division
In, m indicates that the direction x number of grid, n indicate the direction y number of grid)；
Step 2: each physical quantity, correlation function in definition zoning under original state, and establish fluid governing equation and object
Matter state equation, comprising:
(a) interface function of each substance, i.e. Level Set function:
Wherein, Γ (0) indicates that zero contour surface, i.e. material interface, d (x, y, Γ (0)) indicate that a bit (x, y) arrives substance in computational domain
The distance at interface, Ω1、Ω2The respectively region of material interface two sides；
(b) physical parameter of each substance: density, speed, the pressure of explosive under original state, density that original state is lauched, speed
Degree, pressure；
(c) according to the physical characteristic of underwater explosion field, viscous adaptable nothing, irrotationality, meter and compressible fluid governing equation: are established
The fluid governing equation of property:
Ut+ ▽ F (U)=S (U) (2)
Wherein,F (U)=[f (U), g (U)],
(d) state equation (Tait equation) and relevant parameter of water:
Wherein, N, AW、BWIt is constant relevant to the physical property of water, ρW0Indicate the density of water in the initial state, pW、ρWPoint
Not Biao Shi water pressure and density；
(e) detonation products state equation (JWL equation) and relevant parameter:
Wherein, AE、BE、R1、R2, ω be constant relevant to explosive property, ρE0Indicate the density of explosive in the initial state, ρETable
Show the density of detonation product, PEIndicate the density of detonation product, E0Indicate the specific internal energy of detonation product；
Step 3: expand zoning, according to different operating conditions be arranged boundary condition: by boundary lack difference scheme institute it is required
Mesh point, therefore zoning is expanded to N layers outward；According to the format that this method uses, at least need to expand 3 layers outward
Grid (i.e. N >=3)；
For free boundary condition, have: UΓ=U-1,U-i=U-i+1, i=2N, wherein UΓIndicate zoning boundary net
Random physical quantity on lattice, U-iIndicate the random physical quantity on virtual grid；Specifically, first layer is empty outside the boundary of zoning
Physical quantity on quasi- grid is equal to the physical quantity on the boundary mesh of zoning；Outside the boundary of zoning on second layer virtual grid
Physical quantity be equal to first layer virtual grid on physical quantity；Physical quantity etc. outside the boundary of zoning on third layer virtual grid
Physical quantity on second layer virtual grid, and so on；
For wall boundary condition, have:Wherein, UΓIndicate zoning boundary net
Random physical quantity on lattice in addition to speed, U-iIndicate the random physical quantity on virtual grid in addition to speed, VΓIt indicates to calculate
Speed on the grid of zone boundary, V-iIndicate the speed on virtual grid；Specifically, first layer is virtual outside the boundary of zoning
On grid perpendicular in the speed and boundary mesh on boundary perpendicular to boundary velocity vector and be zero, other physical quantitys are equal；
Outside the boundary of zoning on second layer virtual grid perpendicular to boundary speed and first layer virtual grid on perpendicular to boundary
Velocity vector and be zero, other physical quantitys are equal；Speed outside the boundary of zoning on third layer virtual grid perpendicular to boundary
With on second layer virtual grid perpendicular to boundary velocity vector and be zero, other physical quantitys are equal, and so on；
Step 4: choose CFL parameter and calculate time step: CFL parameter is the constant between 0 to 1；Used in current calculating
Indicated with t=n Δ t total time, n is to calculate step number, and time step Δ t can obtain by following formula:
Δ x is the mesh spacing in the direction x in zoning；Δ y is the mesh spacing in the direction y in zoning；U is zoning
Velocity component of the middle current grid in the direction x；Velocity component of the v for current grid in zoning in the direction y；C is to calculate area
Velocity of sound in domain at current grid；
Step 5: for the parameter assignment on zoning and virtual grid: as t=0, the zoning that will be provided in step 1
The physical parameter and interface function value of interior each substance, the parameter initial value as grid each in zoning；As t > 0, by step
The physical parameter of each substance and interface function value in the zoning being calculated in 13, as grid parameter each in zoning
Initial value when next step calculates；For virtual grid, then according in step 3 establish method, by interpolation and extrapolation to its into
Row assignment；
Step 6: solving material interface normal direction: according to Level Set equation, to the normal direction of material interface each in zoning
It is solved, solution formula are as follows:
Step 7: solving Riemannian problem and to grid assignment near material interface: establishing explosive and water along material interface normal direction
Between Riemannian problem and solve the physical quantity between explosive and water at material interface；The physics obtained with solution Riemannian problem
Amount carries out assignment close to one layer of grid of material interface in explosive region；Meanwhile by extrapolating to the net other than explosive region
Lattice carry out assignment, method particularly includes: at the flow field of treated substance interface explosive region side, by the net of the material interface other side
Lattice are set as Level Set method grid, and the pressure on these Level Set method grids is consistent with the side real units with normal velocity,
And across the material interface interpolation and extrapolation of real fluid that tangential velocity and density then pass through explosive region acquires；
Step 8: by the discrete Eulerian equation of RKDG method, the numerical flux of positive pressure of shock wave is solved, specifically: to (2) formula two
While being obtained with multiplied by tentative function Φ (x, y), and in grid cell K upper integral:
After carrying out integration by parts to (7) formula, the line integral and Line Integral in Gauss quadrature formula solution formula are utilized；Due to grid list
First boundary positive pressure of shock wave flux is discontinuous, therefore uses numerical fluxInstead of real fluxesUse approximate solutionInstead of accurately solving U, with the tentative function after numerical approximationInstead of tentative function Φ,
It can thus be concluded that:
Wherein, ω indicates weight coefficient, | e | indicate the length of integral domain grid cell boundary, | K | indicate integral domain grid
The area of unit；
For the numerical flux in (8) formulaIt is solved using (9) formula:
Wherein, α indicates Jacobi matrixThe maximum maximum value of characteristic value, U+(xel,yel, t) and it indicates from just
The U value of grid cell boundary, U are approached on direction-(xel,yel, t) and indicate the U value that grid cell boundary is approached from negative direction；It is right
In the pressure of explosive region detonation product, obtained by (4) formula in step 1；
For the approximate solution in (8) formulaIt is solved using (10) formula:
Wherein, for grid cell [xi-1/2,xi+1/2]×[yj-1/2,yj+1/2],Δ
xi=xi+1/2-xi-1/2, Δ yj=yj+1/2-yj-1/2；
The semi-discrete scheme that (9), (10) two formulas are brought into available (2) formula of (8) formula, is denoted as:
Ut=R U (11)
Wherein, the Discrete Operator of R U representation space derivative term；
Step 9: it is discrete by time term progress of the three rank TVD Runge-Kutta methods to (11) formula equation left side, obtain (2) formula
Complete discrete scheme, and then solve each physical quantity of detonation products in subsequent time zoning；Specific formula are as follows:
Step 10: carrying out similar operations with step 7, solve Riemannian problem again and to grid assignment near material interface: along object
Riemannian problem that matter interface normal direction is established between explosive and water simultaneously solves the physics between explosive and water at material interface
Amount；The physical quantity obtained with solution Riemannian problem carries out assignment close to one layer of grid of material interface in air section；Together
When, assignment is carried out to the grid other than air section by extrapolating, method particularly includes: in treated substance interfacial air region side
Flow field when, set Level Set method grid for the grid of the material interface other side, the pressure on these Level Set method grids with
Normal velocity is consistent with the side real units, and tangential velocity and density then pass through the real fluid in explosive region across substance
Interface interpolation extrapolation acquires；
Step 11: on the basis of step 10, carrying out the operation similar with step 8, step 9；But the difference is that: step 8
In, for the pressure of water area, obtained by (3) formula in step 1；What solution obtained in step 9 is subsequent time zoning
Each physical quantity of interior water；
Step 12: it is discrete to the space derivation item progress of the Level Set equation of motion, it can obtain:
For the lateral numerical flux of the interface function of explosive and water, can be solved by (14) formula:
Wherein,Indicate grid cell in zoningNumerical value in the interface function transverse direction of upper explosive and water
Flux,Indicate the interface function numerical flux in lateral positive direction on the grid cell,Indicate the grid list
Interface function numerical flux in lateral negative direction in member；It indicates in zoning on grid cell (i-r+s, j)
The flux of the lateral positive direction of the interface function of explosive and water；It indicates in zoning on grid cell (i-r+s, j)
The flux of the lateral negative direction of the interface function of explosive and water；Design parameter can be solved by (15), (16) formula:
Wherein,Indicate the interface function value of explosive and water on grid cell (i, j) in zoning；
For the Vertical Numerical flux of the interface function of explosive and water, method for solving is identical；
By step 12, the semi-discrete scheme of available (13) formula is denoted as:
Wherein,The Discrete Operator of representation space derivative term:
Step 13: by the time-derivative item of the three discrete Level Set equations of motion of rank TVD Runge-Kutta method, obtaining
(13) the complete discrete scheme of formula, and then solve the interface function value of each grid cell in subsequent time zoning；
Step 14: on the basis of step before, solving and export the pressure value of all grid cells in zoning: if grid
Unit is in detonation products region, then is acquired by (4) the formula calculating in step 1；If grid cell is in water area
It is interior, then it is acquired by (3) formula in step 1；
Step 15: the current calculating total time t for calculating total time t=n Δ t and setting used of judgement*Relationship: if t > t*,
Then terminate to calculate, exports the pressure value of grid cell in zoning；If t > t*, then it is back to step 3, continues operation.
By the operation of above-mentioned steps, the shock wave in underwater explosion zoning at each moment each grid cell can be obtained
Pressure value.
2. a kind of high resolution numerical simulation method of near field underwater blast wave load as described in claim 1, feature
It is: extrapolation assignment is carried out to grid cell near material interface using Riemann problem solution in its step 7, method particularly includes:
At the flow field of treated substance interface explosive region side, Level Set method net is set by the grid of the material interface other side
Lattice, the pressure and normal velocity on these Level Set method grids are consistent with the side real units, and tangential velocity and density
Then acquired by across the material interface interpolation and extrapolation of real fluid in explosive region.
3. a kind of high resolution numerical simulation method of near field underwater blast wave load as described in one of claims 1 or 2,
It is characterized by: passing through the joint of step 7, step 11 and step 12, have the two or more medium interaction problems of processing
Function.
4. a kind of high resolution numerical simulation method of near field underwater blast wave load as described in claims 1 to 3, special
Sign is: the discrete Eulerian equation of RKDG method is used in its step 8, it is only necessary to read grid cell self information and participate in calculating i.e.
The above precision of three ranks can be obtained.
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201811018574.6A CN109214082A (en) | 2018-09-03 | 2018-09-03 | A kind of high resolution numerical simulation method of near field underwater blast wave load |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201811018574.6A CN109214082A (en) | 2018-09-03 | 2018-09-03 | A kind of high resolution numerical simulation method of near field underwater blast wave load |
Publications (1)
| Publication Number | Publication Date |
|---|---|
| CN109214082A true CN109214082A (en) | 2019-01-15 |
Family
ID=64987167
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201811018574.6A Pending CN109214082A (en) | 2018-09-03 | 2018-09-03 | A kind of high resolution numerical simulation method of near field underwater blast wave load |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN109214082A (en) |
Cited By (15)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN109902376A (en) * | 2019-02-25 | 2019-06-18 | 北京理工大学 | A kind of fluid structurecoupling high resolution numerical simulation method based on Continuum Mechanics |
| CN110879919A (en) * | 2019-11-18 | 2020-03-13 | 中国人民解放军陆军防化学院 | Sectional type simulation method for poison diffusion under explosion action |
| CN111008492A (en) * | 2019-11-22 | 2020-04-14 | 电子科技大学 | High-order unit Euler equation numerical simulation method based on Jacobian-free matrix |
| CN111783365A (en) * | 2020-06-04 | 2020-10-16 | 三多(杭州)科技有限公司 | Virtual medium method, device and equipment applied to low-voltage interface processing |
| CN112395722A (en) * | 2019-08-15 | 2021-02-23 | 武汉理工大学 | Method for acquiring motion response of hull beam under action of underwater explosion and wave load |
| CN112861263A (en) * | 2021-02-22 | 2021-05-28 | 西北工业大学 | Calculation simulation method suitable for compressible two-phase flow |
| CN112861445A (en) * | 2020-12-29 | 2021-05-28 | 中国船舶重工集团公司第七研究院 | Simulation method of grid-free numerical model |
| CN113408168A (en) * | 2021-06-18 | 2021-09-17 | 北京理工大学 | High-precision numerical simulation method based on Riemann problem accurate solution |
| CN113515900A (en) * | 2021-04-07 | 2021-10-19 | 南京航空航天大学 | Numerical simulation method for simulating wake flow of continuous rotation detonation engine |
| CN113553739A (en) * | 2021-07-06 | 2021-10-26 | 西安近代化学研究所 | Method for calculating explosion output characteristics of mixed explosive |
| CN113591345A (en) * | 2021-07-08 | 2021-11-02 | 北京理工大学 | Explosion reaction flow high-precision prediction method based on generalized Riemann solution method device |
| CN114595603A (en) * | 2022-03-01 | 2022-06-07 | 哈尔滨工程大学 | Underwater structure impact resistance calculation method considering underwater explosion bubble jet flow slamming |
| CN114647990A (en) * | 2022-02-22 | 2022-06-21 | 哈尔滨工程大学 | Method for efficiently and highly accurately calculating pressure wavelet of air gun seismic source |
| CN116306364A (en) * | 2023-03-13 | 2023-06-23 | 武汉理工大学 | Wave-absorbing simulation method for explosion water mist in cabin |
| CN116992796A (en) * | 2023-09-27 | 2023-11-03 | 中国科学技术大学 | Self-adaptive low-dissipation SPH-HLLC Riemann solver coupling algorithm |
-
2018
- 2018-09-03 CN CN201811018574.6A patent/CN109214082A/en active Pending
Cited By (24)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN109902376A (en) * | 2019-02-25 | 2019-06-18 | 北京理工大学 | A kind of fluid structurecoupling high resolution numerical simulation method based on Continuum Mechanics |
| CN112395722B (en) * | 2019-08-15 | 2022-05-31 | 武汉理工大学 | Method for acquiring motion response of hull beam under action of underwater explosion and wave load |
| CN112395722A (en) * | 2019-08-15 | 2021-02-23 | 武汉理工大学 | Method for acquiring motion response of hull beam under action of underwater explosion and wave load |
| CN110879919A (en) * | 2019-11-18 | 2020-03-13 | 中国人民解放军陆军防化学院 | Sectional type simulation method for poison diffusion under explosion action |
| CN110879919B (en) * | 2019-11-18 | 2023-08-18 | 中国人民解放军陆军防化学院 | Sectional simulation method for poison diffusion under explosion effect |
| CN111008492A (en) * | 2019-11-22 | 2020-04-14 | 电子科技大学 | High-order unit Euler equation numerical simulation method based on Jacobian-free matrix |
| CN111008492B (en) * | 2019-11-22 | 2022-10-14 | 电子科技大学 | High-order unit Euler equation numerical simulation method based on Jacobian-free matrix |
| CN111783365A (en) * | 2020-06-04 | 2020-10-16 | 三多(杭州)科技有限公司 | Virtual medium method, device and equipment applied to low-voltage interface processing |
| CN111783365B (en) * | 2020-06-04 | 2023-09-22 | 三多(杭州)科技有限公司 | Virtual medium method, device and equipment applied to low-voltage interface processing |
| CN112861445A (en) * | 2020-12-29 | 2021-05-28 | 中国船舶重工集团公司第七研究院 | Simulation method of grid-free numerical model |
| CN112861263A (en) * | 2021-02-22 | 2021-05-28 | 西北工业大学 | Calculation simulation method suitable for compressible two-phase flow |
| CN112861263B (en) * | 2021-02-22 | 2024-02-13 | 西北工业大学 | Calculation simulation method suitable for compressible two-phase flow |
| CN113515900A (en) * | 2021-04-07 | 2021-10-19 | 南京航空航天大学 | Numerical simulation method for simulating wake flow of continuous rotation detonation engine |
| CN113408168A (en) * | 2021-06-18 | 2021-09-17 | 北京理工大学 | High-precision numerical simulation method based on Riemann problem accurate solution |
| CN113553739A (en) * | 2021-07-06 | 2021-10-26 | 西安近代化学研究所 | Method for calculating explosion output characteristics of mixed explosive |
| CN113591345B (en) * | 2021-07-08 | 2024-01-23 | 北京理工大学 | Explosion reaction flow high-precision prediction method based on generalized Riemann solver |
| CN113591345A (en) * | 2021-07-08 | 2021-11-02 | 北京理工大学 | Explosion reaction flow high-precision prediction method based on generalized Riemann solution method device |
| CN114647990A (en) * | 2022-02-22 | 2022-06-21 | 哈尔滨工程大学 | Method for efficiently and highly accurately calculating pressure wavelet of air gun seismic source |
| CN114595603A (en) * | 2022-03-01 | 2022-06-07 | 哈尔滨工程大学 | Underwater structure impact resistance calculation method considering underwater explosion bubble jet flow slamming |
| CN114595603B (en) * | 2022-03-01 | 2023-06-13 | 哈尔滨工程大学 | Underwater structure impact resistance calculation method considering underwater explosion bubble jet slamming |
| CN116306364A (en) * | 2023-03-13 | 2023-06-23 | 武汉理工大学 | Wave-absorbing simulation method for explosion water mist in cabin |
| CN116306364B (en) * | 2023-03-13 | 2024-03-22 | 武汉理工大学 | Wave-absorbing simulation method for explosion water mist in cabin |
| CN116992796B (en) * | 2023-09-27 | 2023-12-22 | 中国科学技术大学 | Self-adaptive low-dissipation SPH-HLLC Riemann solver coupling method |
| CN116992796A (en) * | 2023-09-27 | 2023-11-03 | 中国科学技术大学 | Self-adaptive low-dissipation SPH-HLLC Riemann solver coupling algorithm |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN109214082A (en) | A kind of high resolution numerical simulation method of near field underwater blast wave load | |
| Yan et al. | Experimental and numerical investigation of water impact on air-launched AUVs | |
| Cheng et al. | Towards the modeling of the ditching of a ground-effect wing ship within the framework of the SPH method | |
| Zhang et al. | SPH-BEM simulation of underwater explosion and bubble dynamics near rigid wall | |
| CN106124349A (en) | Determine the critical erosion rate of gas well, the method and device of flow | |
| Nguyen et al. | Multiphase flow simulation of water-entry and-exit of axisymmetric bodies | |
| Liang et al. | Numerical study of spherical blast-wave propagation and reflection | |
| Faucher et al. | High resolution adaptive framework for fast transient fluid-structure interaction with interfaces and structural failure–Application to failing tanks under impact | |
| Rebelo et al. | A Comparison between three air blast simulation techniques in LS-DYNA | |
| Ning et al. | Numerical analysis of the shaped charged jet with large cone angle | |
| Shenshen et al. | Two-phase flow characteristics of different charging structure in ignition and flame-spreading | |
| Ren et al. | Numerical studies of penetration problems by an improved particle method | |
| Mu et al. | Numerical simulation on the cavitation flow of high speed oblique water entry of revolution body | |
| Yang et al. | Characteristics and load reduction method of the deep-sea implosion of the ellipsoidal and egg-shaped pressure hulls | |
| Georgievskiy et al. | Shock focusing upon interaction of a shock with a cylindrical dust cloud | |
| Grządziela | Ship impact modeling of underwater explosion | |
| Hsiao et al. | Development of compressible-incompressible link to efficiently model bubble dynamics near floating body | |
| Zhu et al. | Modeling and analysis of numerical wave tank based on the ALE algorithm | |
| Pappalardo et al. | Comparison of CFD numerical approaches for the simulation of accidental gas release in energy applications | |
| Chen et al. | Transient Response of the Steel Plate to Underwater Explosion Bubble | |
| Lu et al. | Risk assessment method and protection goals of high concrete gravity dam subjected to far-field underwater nuclear explosion | |
| Jie et al. | Numerical simulation on underwater explosion in small-sized containers | |
| Charles et al. | A numerical study on cavity expansion in water: hydraulic ram under ballistic impacts | |
| Wei | Development of a new steam explosion model for the MC3D software | |
| Duan et al. | Design and Implementation of a Damage Assessment System for Large-Scale Surface Warships Based on Deep Learning |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| WD01 | Invention patent application deemed withdrawn after publication |
Application publication date: 20190115 |
|
| WD01 | Invention patent application deemed withdrawn after publication |