CN103425881A - Method for certainty numerical simulation of crack medium seismic wave response - Google Patents
Method for certainty numerical simulation of crack medium seismic wave response Download PDFInfo
- Publication number
- CN103425881A CN103425881A CN2013103408431A CN201310340843A CN103425881A CN 103425881 A CN103425881 A CN 103425881A CN 2013103408431 A CN2013103408431 A CN 2013103408431A CN 201310340843 A CN201310340843 A CN 201310340843A CN 103425881 A CN103425881 A CN 103425881A
- Authority
- CN
- China
- Prior art keywords
- crack
- node
- finite element
- fracture
- nodes
- 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
Images
Landscapes
- Geophysics And Detection Of Objects (AREA)
Abstract
The invention discloses a method for certainty numerical simulation of crack medium seismic wave response. The method for certainty numerical simulation of the crack medium seismic wave response comprises the steps that a crack in a geologic model is designed to be two faces in an intermittent contact relation and finite element unit division is conducted on the model with the crack serving as a boundary; according to the contact relation of the two faces of the crack, the simulation processing mode of finite element nodes where the crack passes is determined and wave field numerical value simulation is conducted on the model. According to the method for certainty numerical simulation of the crack medium seismic wave response, simulation calculation of crack medium low and medium frequency seismic waves is achieved.
Description
Technical field
The present invention relates to the forward simulation technology field containing seismic event in fracture medium, relate in particular to a kind of determinacy method for numerical simulation of fracture medium seismic event response.
Background technology
Along with the fast development of oil-gas exploration and development cause, slit formation hydrocarbon-bearing pool will become one of the efforts will be concentrated on prospecting object in this century.Because fractured model is extremely important for the hydrocarbon-bearing pool of prediction slit formation, therefore, the research of fractured model receives much concern.
Crack is the important reservoir space of reservoir of oil and gas and main seepage channel, and the growth in crack and distribution are the key factors of controlling hydrocarbon occurrence.The growth in crack, not only make the important rock of several classes, as carbonatite, magmatite, metamorphosed rock and mud stone etc. can become reservoir, and also plays an important role in low permeability sandstone reservoir.At present, the output of whole world slit formation hydrocarbon-bearing pool has accounted for the over half of petroleum gas total production.Therefore, how quantitative forecast Reservoir Fracture systematic parameter accurately and effectively, be an important research content in exploration of oil and gas field and exploitation.Yet the theoretical foundation of subterranean fracture system prediction is the communication theory of seismic event in fracture medium, must accurately hold the propagation law of seismic event in fracture medium, could develop Fracture System Forecasting Methodology accurately and effectively.
Research, containing the response problem of seismic event in fracture medium, mainly contains three kinds of means: i.e. rock physics experiment, EFFECTIVE MEDIUM THEORY and numerical simulation.Rock physics experiment can only be carried out ultrasonic measuring to the fritter sample in laboratory usually, is difficult to obtain the fluctuation response data of exploration frequency range, and sampling and test environment are difficult to represent actual conditions.EFFECTIVE MEDIUM THEORY is a kind of important means of Characteristics of Seismic Wave Propagation in the understanding fracture medium, it can set up contacting between seismic wave attribute and middle sight or micro crack parameter on a macro scale, is to utilize Seismic Data directly to obtain the basis of subterranean fracture system information.But, due to the difficulty on physics and mathematics manipulation, the crack EFFECTIVE MEDIUM THEORY of current popular, as Hudson(1980) model, Schoenberg(1995) model, Chapman(2003) model etc., they are all also perfect not to the utmost.These models for be all the low-density Fracture System basically, and at aspects such as fracture shape, interstitial interaction, fluid effects, many approximate simplification are arranged, be difficult to reliable treatments high density complex fracture system problem.In numerical simulation study, two analoglike methods are arranged.One is based on the method for numerical simulation of EFFECTIVE MEDIUM THEORY, and it has inherited defect and the deficiency of EFFECTIVE MEDIUM THEORY.Another kind of is fracture medium to be regarded as to the determinacy method for numerical simulation (Wu etc., 2005) of the inhomogeneous composite be comprised of matrix and crack stuff.In theory, this class determinacy method for numerical simulation can be simulated earthquake wave propagation in any complex fracture medium, has EFFECTIVE MEDIUM THEORY and the incomparable advantage of rock physics experimental study.But, because the scale of actual subterranean fracture gap is generally less, particularly relatively low frequency seismic exploration wave-wave is long for the aperture of crack gap, almost can regard zero width as, this just requires, in simulation, model is adopted to a large amount of mesh generations, thereby produce calculated amount and the memory requirements of magnanimity, cause under the current computer condition, also can't effectively carry out the determinacy numerical simulation of the fracture medium seismic event response of low-frequency earthquake waves field (exploration frequency range).
Summary of the invention
Because the problem of the determinacy numerical simulation of the fracture medium seismic event response that can't effectively carry out the low-frequency earthquake waves field that prior art exists, fundamental purpose of the present invention is to provide a kind of determinacy method for numerical simulation of fracture medium seismic event response, wherein:
According to the determinacy method for numerical simulation of the fracture medium seismic event of embodiment of the present invention response, comprise:
Crack in geologic model is set to have two faces of interrupted contact relation, take crack, as border, described model is carried out to the finite element element subdivision; Determine the simulation process mode of the finite element node passed through in crack according to the contact relation of two faces in crack, and described model is carried out to wave field numerical.
Wherein, described method also comprises: by high-order limited first method, a plurality of nodes are set on crack, described node comprises: low order node and/or high-order encryption node.
Wherein, the contact relation of two faces in described crack comprises: separation relation or contact relation; Described method also comprises: the type of determining node according to the contact relation of two faces in crack: if the fracture plane of Nodes be separation relation this node be discontinuous node, if the fracture plane of Nodes be contact relation this node be continuous nodes.
Wherein, the described contact relation according to two faces in crack is determined the step of the simulation process mode of the finite element node passed through in crack, comprise: for discontinuous node, according to two isolated nodes, carry out simulation process, wherein, described two isolated nodes are identical but two isolated nodes not acting on mutually of physical location; Carry out simulation process for continuous nodes according to continuous medium.
Wherein, described step of described model being carried out to wave field numerical comprises: calculate mass matrix and the stiffness matrix of finite element unit, and the outer force vector of computing node; Assemble a plurality of finite elements unit and form the overall equation of elastic wave finite element; According to the requirement of default focal shock parameter and receive mode and simulation duration, the overall equation of finite element is carried out to iterative computing node wave field.
Wherein, calculate mass matrix and the stiffness matrix of finite element unit by following formula:
Wherein, M means mass matrix, and K means stiffness matrix, and n means the element number of whole regional subdivision, and ρ means the density of medium, and N means the row vector of shape function, and B means gradient matrix; D means elastic matrix.
Wherein, by the outer force vector of following formula computing node:
Wherein, U means the nodal displacement column vector,
Mean the second derivative of U to the time; F means the suffered external force column vector of node.
Wherein, by following formula computing node wave field: U
T+ Δ t=Δ t
2M
-1(F-KU
t)+2U
t-U
T-Δ t
Wherein, Δ t is time step.
Wherein, the type in described crack comprises: orthogonal fracture, dipping fracture.
According to technical scheme of the present invention, by being modeled as, subterranean fracture there is interrupted contact relation and distance for infinitesimal two faces, usining crack is split into finite element grid to model as inner boundary, decide the assembling intergration model of cell node on fracture plane according to the actual contact relation of two faces in crack, realized the numerical simulation containing the multiple dimensioned complex fracture medium medium and low frequency seismic wave field of many groups.With determinacy method for numerical simulation of the prior art, compare, technical scheme of the present invention makes more than the number of meshes of fracture medium model subdivision reduced by two or three orders of magnitude, corresponding analog computation amount makes the analog computation of fracture medium medium and low frequency seismic event be achieved more than also having reduced by two or three orders of magnitude.
The accompanying drawing explanation
Accompanying drawing described herein is used to provide a further understanding of the present invention, forms the application's a part, and schematic description and description of the present invention the present invention does not form inappropriate limitation of the present invention for explaining.In the accompanying drawings:
Fig. 1 is the process flow diagram according to the determinacy method for numerical simulation of the fracture medium seismic event response of the embodiment of the present invention;
Fig. 2 is the schematic diagram according to the fractured model of the embodiment of the present invention;
Fig. 3 A and Fig. 3 B are according to the fractured model subdivision of the embodiment of the present invention and the schematic diagram of finite element Node configuration;
Fig. 4 is the schematic diagram according to the geologic model of the embodiment of the present invention;
Fig. 5 is the schematic diagram that ground wave script holder that application the present invention simulates acquisition records the x component;
Fig. 6 is the schematic diagram that ground wave script holder that application the present invention simulates acquisition records the z component.
Embodiment
For making the purpose, technical solutions and advantages of the present invention clearer, below in conjunction with drawings and the specific embodiments, the present invention is described in further detail.
According to embodiments of the invention, provide a kind of determinacy method for numerical simulation of fracture medium seismic event response.
Fig. 1 is that as shown in Figure 1, the method comprises according to the process flow diagram of the determinacy method for numerical simulation of the fracture medium seismic event response of the embodiment of the present invention:
Step S102, the crack in geologic model is set to have two faces of interrupted contact relation, take crack, as border, described model is carried out to the finite element element subdivision.
Step S104, determine the simulation process mode of the finite element node passed through in crack, and described model carried out to wave field numerical according to the contact relation of two faces in crack.
Wherein, the contact relation of two faces in crack comprises: separation relation or contact relation.If the fracture plane of Nodes be separation relation this node be discontinuous node, if the fracture plane of Nodes be contact relation this node be continuous nodes.Carry out simulation process for discontinuous node according to two isolated nodes, wherein, described two isolated nodes are identical but two isolated nodes not acting on mutually of physical location; Carry out simulation process for continuous nodes according to continuous medium.
By above-described embodiment, realized sign and the simulation of fracture in a large geologic model, and the simulation calculation of realizing seismic event in fracture medium.
Below in conjunction with Fig. 2 to Fig. 3, the embodiment of the present invention is described in further detail.
1, as shown in Figure 2, a series of crack is set on request in geologic model.The type in described crack comprises: orthogonal fracture and any dipping fracture, the present invention is not limited for the type in crack, for the sake of simplicity at Fig. 2 and hereinafter take orthogonal fracture and describe as example.
2, geologic model is carried out to element subdivision, it is to make crack become elementary boundary that subdivision requires, and the finite element grid node on it comprises contact point and burble point.
Particularly with reference to figure 3.As shown in Figure 3A, crack (thick black line in figure) regarded as and had interrupted contact relation and, apart from being infinitesimal two faces, using crack as inner boundary, model is split into to finite element grid.As shown in Figure 3 B, the finite element grid node that crack is passed through comprises: contact point (as A, B point), and burble point (as C, D point), burble point has independently two nodes (as the C point has C1 and C2, the D point has D1 and D2) at same position.
In one embodiment of the invention, can also more nodes be set on crack by high-order limited first method, mean the contact relation of both sides, crack meticulouslyr (as large black joint A, B, C, D point in figure are the low order node, the little black node between A, B, C, D is the high-order encryption node).
3, the finite element node is sorted out, divided continuity point (being contact point) and point of discontinuity (being burble point), and mark the crack elementary boundary and the node that pass through.
At each discontinuous node location, arrange on the fracture plane both sides and there are identical but two isolated nodes that do not act on mutually of physical location.
4, by following formula computing unit mass matrix and stiffness matrix, and calculate outer force vector, assemble all unit and form the overall equation of elastic wave finite element:
In formula (1), M, K are respectively mass matrix and stiffness matrix; U is the nodal displacement column vector;
For the second derivative of U to the time; F is the suffered external force column vector of node.The expression formula of M and K is respectively:
In formula (2), the unit number that n is whole regional subdivision, the integral operation after the summation symbol is carried out on each unit;
The density that ρ is medium;
The row vector that N is shape function, the Interpolation that it adopts with element displacement is relevant;
B is gradient matrix;
D is elastic matrix, and it is difference with the difference of medium elastic property.
5, according to the focal shock parameter and the receive mode that design, and according to the requirement of simulating duration, to the iterative that carries out of the differential equation (1).Solve the employing centered Finite Difference Methods.Its recursive form is
U
T+ Δ t=Δ t
2M
-1(F-KU
t)+2U
t-U
T-Δ tFormula (3)
In formula (3), Δ t is time step (being sampling time interval).
Below adopt the embodiment of the present invention to carry out the field wave numerical simulation.Fig. 4 is a geologic model with certain practical significance, and it has three-decker, and each layer of physical parameter is in Table 1.Middle layer is containing each 350 of vertical and horizontal fractures, and the fracture length scope is 5m-20m, and length and location all determines at random, explosive source be positioned in the middle of earth's surface below the 5m depths, dominant frequency is 50Hz.Fig. 5 records the x component for the ground wave script holder that simulation obtains, and Fig. 6 records the z component for the ground wave script holder that simulation obtains.Analog result shows, applies the result of calculation of the inventive method except having obtained the clear record of reflection wave accurately, has also obtained the stronger scattered wave field of subterranean fracture to seismic wave field, and this is that fracture medium causes one of key factor of seismic wave attenuation.
Table 1 model physical parameter
| ? | Velocity of longitudinal wave Vp(m/s) | Shear wave velocity Vs(m/s) | Density p (kg/m3) |
| Ground floor | 3500 | 2000 | 2150 |
| The second layer | 3800 | 2100 | 2200 |
| The 3rd layer | 4300 | 2350 | 2300 |
According to technical scheme of the present invention, by being modeled as, subterranean fracture there is interrupted contact relation and distance for infinitesimal two faces, usining crack is split into finite element grid to model as inner boundary, decide the assembling intergration model of cell node on fracture plane according to the actual contact relation of two faces in crack, realized the numerical simulation containing the multiple dimensioned complex fracture medium medium and low frequency seismic wave field of many groups.With determinacy method for numerical simulation of the prior art, compare, technical scheme of the present invention makes more than the number of meshes of fracture medium model subdivision reduced by two or three orders of magnitude, corresponding analog computation amount makes the analog computation of fracture medium medium and low frequency seismic event be achieved more than also having reduced by two or three orders of magnitude.
The foregoing is only embodiments of the invention, be not limited to the present invention, for a person skilled in the art, the present invention can have various modifications and variations.Within the spirit and principles in the present invention all, any modification of doing, be equal to replacement, improvement etc., within all should being included in claim scope of the present invention.
Claims (9)
1. the determinacy method for numerical simulation of a fracture medium seismic event response, is characterized in that, comprising:
Crack in geologic model is set to have two faces of interrupted contact relation, take crack, as border, described model is carried out to the finite element element subdivision;
Determine the simulation process mode of the finite element node passed through in crack according to the contact relation of two faces in crack, and described model is carried out to wave field numerical.
2. method according to claim 1, is characterized in that, also comprises:
By high-order limited first method, a plurality of nodes are set on crack, described node comprises: low order node and/or high-order encryption node.
3. method according to claim 1, is characterized in that, the contact relation of two faces in described crack comprises: separation relation or contact relation; Described method also comprises:
Determine the type of node according to the contact relation of two faces in crack: if the fracture plane of Nodes be separation relation this node be discontinuous node, if the fracture plane of Nodes be contact relation this node be continuous nodes.
4. method according to claim 3, is characterized in that, the described contact relation according to two faces in crack is determined the step of the simulation process mode of the finite element node passed through in crack, comprising:
Carry out simulation process for discontinuous node according to two isolated nodes, wherein, described two isolated nodes are identical but two isolated nodes not acting on mutually of physical location;
Carry out simulation process for continuous nodes according to continuous medium.
5. method according to claim 4, is characterized in that, described step of described model being carried out to wave field numerical comprises:
Calculate mass matrix and the stiffness matrix of finite element unit, and the outer force vector of computing node;
Assemble a plurality of finite elements unit and form the overall equation of elastic wave finite element;
According to the requirement of default focal shock parameter and receive mode and simulation duration, the overall equation of finite element is carried out to iterative computing node wave field.
6. method according to claim 5, is characterized in that, calculates mass matrix and the stiffness matrix of finite element unit by following formula:
Wherein, M means mass matrix, and K means stiffness matrix, and n means the element number of whole regional subdivision, and ρ means the density of medium, and N means the row vector of shape function, and B means gradient matrix; D means elastic matrix.
7. method according to claim 5, is characterized in that, by the outer force vector of following formula computing node:
Wherein, U means the nodal displacement column vector,
Mean the second derivative of U to the time; F means the suffered external force column vector of node.
8. method according to claim 5, is characterized in that, by following formula computing node wave field:
U
t+Δt=Δt
2M
-1(F-KU
t)+2U
t-U
t-Δt
Wherein, Δ t is time step.
9. method according to claim 1, is characterized in that, the type in described crack comprises: orthogonal fracture, dipping fracture.
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN2013103408431A CN103425881A (en) | 2013-08-07 | 2013-08-07 | Method for certainty numerical simulation of crack medium seismic wave response |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN2013103408431A CN103425881A (en) | 2013-08-07 | 2013-08-07 | Method for certainty numerical simulation of crack medium seismic wave response |
Publications (1)
| Publication Number | Publication Date |
|---|---|
| CN103425881A true CN103425881A (en) | 2013-12-04 |
Family
ID=49650610
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN2013103408431A Pending CN103425881A (en) | 2013-08-07 | 2013-08-07 | Method for certainty numerical simulation of crack medium seismic wave response |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN103425881A (en) |
Cited By (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN105160093A (en) * | 2015-08-28 | 2015-12-16 | 重庆交通大学 | Method for analyzing change of force applied to fractured rock slope along depth under heavy seismic load |
| CN107656308A (en) * | 2017-08-18 | 2018-02-02 | 中国科学院地质与地球物理研究所 | A kind of common scattering point pre-stack time migration imaging method based on time depth scanning |
| CN108460212A (en) * | 2018-03-02 | 2018-08-28 | 北京交通大学 | Underground structure earthquake motion three-dimensional numerical simulation method |
Citations (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| JP2011203834A (en) * | 2010-03-24 | 2011-10-13 | Toshiba Corp | Method and apparatus for electromagnetic field simulation, and manufacturing method of semiconductor device |
| CN102253415A (en) * | 2011-04-19 | 2011-11-23 | 中国石油大学(华东) | Method for establishing earthquake response mode based on fracture equivalent medium model |
-
2013
- 2013-08-07 CN CN2013103408431A patent/CN103425881A/en active Pending
Patent Citations (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| JP2011203834A (en) * | 2010-03-24 | 2011-10-13 | Toshiba Corp | Method and apparatus for electromagnetic field simulation, and manufacturing method of semiconductor device |
| CN102253415A (en) * | 2011-04-19 | 2011-11-23 | 中国石油大学(华东) | Method for establishing earthquake response mode based on fracture equivalent medium model |
Non-Patent Citations (3)
| Title |
|---|
| 张美根 等: "各向异性弹性波场的有限元数值模拟", 《地球物理学进展》, 30 September 2002 (2002-09-30) * |
| 张美根 等: "复杂裂缝介质中地震波场数值模拟研究", 《中国地球物理2010》, vol. 17, no. 3, 31 December 2010 (2010-12-31) * |
| 汪文帅 等: "基于多辛结构谱元法的保结构地震波场模拟", 《地球物理学报》, vol. 55, no. 10, 31 October 2012 (2012-10-31) * |
Cited By (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN105160093A (en) * | 2015-08-28 | 2015-12-16 | 重庆交通大学 | Method for analyzing change of force applied to fractured rock slope along depth under heavy seismic load |
| CN107656308A (en) * | 2017-08-18 | 2018-02-02 | 中国科学院地质与地球物理研究所 | A kind of common scattering point pre-stack time migration imaging method based on time depth scanning |
| CN107656308B (en) * | 2017-08-18 | 2019-05-31 | 中国科学院地质与地球物理研究所 | A kind of common scattering point pre-stack time migration imaging method based on time depth scanning |
| CN108460212A (en) * | 2018-03-02 | 2018-08-28 | 北京交通大学 | Underground structure earthquake motion three-dimensional numerical simulation method |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| US9665537B2 (en) | Method for generating a fractured reservoir mesh with a limited number of nodes in the matrix medium | |
| Lei et al. | A new approach to upscaling fracture network models while preserving geostatistical and geomechanical characteristics | |
| Wu et al. | Integrated characterization of the fracture network in fractured shale gas reservoirs—stochastic fracture modeling, simulation and assisted history matching | |
| Bogdanov et al. | Effective permeability of fractured porous media in steady state flow | |
| Odling et al. | Permeability scaling properties of fault damage zones in siliclastic rocks | |
| US9103194B2 (en) | Method for constructing a fracture network grid from a Voronoi diagram | |
| EP3018502A2 (en) | Modeling fluid-conducting fractures in reservoir simulation grids | |
| Jiang et al. | Seepage flow with free surface in fracture networks | |
| US20120116740A1 (en) | Method for characterizing the fracture network of a fractured reservoir and method for developing it | |
| US10534877B2 (en) | Adaptive multiscale multi-fidelity reservoir simulation | |
| Ouenes et al. | Integrated fractured reservoir modeling using both discrete and continuum approaches | |
| Yao et al. | Tuning fractures with dynamic data | |
| Fachri et al. | Sensitivity of fluid flow to deformation-band damage zone heterogeneities: A study using fault facies and truncated Gaussian simulation | |
| US11041976B2 (en) | Method and system for creating and using a subsurface model in hydrocarbon operations | |
| Blöcher et al. | Evaluation of three exploitation concepts for a deep geothermal system in the North German Basin | |
| Beskardes et al. | Modelling DC responses of 3-D complex fracture networks | |
| Salinas et al. | Dynamic mesh optimisation for geothermal reservoir modelling | |
| Li et al. | Development of 3‐D curved fracture swarms in shale rock driven by rapid fluid pressure buildup: Insights from numerical modeling | |
| CN103425881A (en) | Method for certainty numerical simulation of crack medium seismic wave response | |
| Hjeij et al. | Comparing advanced discretization methods for complex hydrocarbon reservoirs | |
| Huang et al. | Connectivity evaluation for three‐dimensional fracture network in support‐based model: A case study in the Ordos Basin, China | |
| Limaye et al. | A vector‐based method for bank‐material tracking in coupled models of meandering and landscape evolution | |
| Biver et al. | Modeling uncertainties in carbonate karstic reservoirs: presentation of a tool and its application to a real field case in Russia | |
| Swyer et al. | Permeability potential modeling of geothermal prospects combining regional crustal strain rates with geomechanical simulation of fault slip and volcanic center deformation: a case study for Washington State geothermal play fairways | |
| CN115964962A (en) | Three-dimensional seam network modeling digital-analog integration method based on geological engineering integration idea |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| C06 | Publication | ||
| PB01 | Publication | ||
| C10 | Entry into substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| RJ01 | Rejection of invention patent application after publication |
Application publication date: 20131204 |
|
| RJ01 | Rejection of invention patent application after publication |