旋转轴对称弱各向异性介质模型的射线追踪方法

Ray tracing is a rapid and effective method for wave field calculation. Not only in the field of seismic-wave theory, but also in the field of seismic inversion and migration imaging,the seismic ray tracing method has become one of the most important methods. In anisotropic media, group velocity and phase velocity have different propagation directions. The seismic wave propagates along the direction of group velocity , it does not depend on the direction of phase velocity. Ray angle is a complex function with respect to phase angle, it is difficult to measure and calculate. But most rocks are weak anisotropic, so the expression of phase velocity can be simplified greatly. Based on the approximate expression of phase velocity this thesis for rotating axisymmetric weak anisotropic media deduces an expression of the partial derivative of phase velocity and an expression of group velocity with the method of linear approximation. This paper uses the fourth order Runge-Kutta method together with the two-dimensional interpolation and linear interpolation to obtain the parameters of the physical locations. At last the paths of seismic wave in rotating axisymmetric weak anisotropic media are computed. According to the analysis of the computational results, it indicates that the method developed in this paper has strong adaptability, high computational efficiency and high accuracy for rotating axisymmetric weak anisotropic media.

复杂介质中地震波场数值模拟研究——基于广义正交多项式迭积微分算子法

In this paper, we propose a new numerical modeling method – Convolutional Forsyte Polynomial Differentiator (CFPD), aimed at simulating seismic wave propagation in complex media with high efficiency and accuracy individually owned by short-scheme finite differentiator and general convolutional polynomial method. By adjusting the operator length and optimizing the operator coefficient, both global and local informations can be easily incorporated into the wavefield which is important to invert the undersurface geological structure. The key issue in this paper is to introduce the convolutional differentiator based on Forsyte generalized orthogonal polynomial in mathematics into the spatial differentiation of the first velocity-stress equation. To match the high accuracy of the spatial differentiator, this method in the time coordinate adopts staggered grid finite difference instead of conventional finite difference to model seismic wave propagation in heterogeneous media. To attenuate the reflection artifacts caused by artificial boundary, Perfectly Matched Layer (PML) absorbing boundary is also being considered in the method to deal with boundary problem due to its advantage of automatically handling large-angle emission. The PML formula for acoustic equation and first-order velocity-stress equation are also derived in this paper. There is little difference to implement the PML boundary condition in all kind of wave equations, but in Biot media, special attenuation factors should be taken. Numerical results demonstrate that the PML boundary condition is better than Cerjan absorbing boundary condition which makes it more suitable to hand the artificial boundary reflection. Based on the theories of anisotropy, Biot two-phase media and viscous-elasticity, this paper constructs the constitutive relationship for viscous-elastic and two-phase media, and further derives the first-order velocity-stress equation for 3D viscous-elastic and two-phase media. Numerical modeling using CFPD method is carried out in the above-mentioned media. The results modeled in the viscous-elastic media and the anisotropic pore elastic media can better explain wave phenomena of the true earth media, and can also prove that CFPD is a useful numerical tool to study the wave propagation in complex media.

离散裂缝诱导各向异性粘弹性地震波模拟与特征分析

The real media always attenuate and distort seismic waves as they propagate in the earth. This behavior can be modeled with a viscoelastic and anisotropic wave equation. The real media can be described as fractured media. In this thesis, we present a high-order staggered grid finite-difference scheme for 2-D viscoelastic wave propagation in a medium containing a large number of small finite length fractures. We use the effective medium approach to compute the anisotropic parameters in each grid cell. By comparing our synthetic seismogram by staggered-grid finite-difference with that by complex-ray parameter ray tracing method, we conclude that the high-order staggered-grid finite-difference technique can effectively used to simulate seismic propagation in viscoelastic-anisotropic media. Synthetic seismograms demonstrate that strong attenuation and significant frequency dispersion due to viscosity are important factors of reducing amplitude and delaying arrival time varying with incidence angle or offset. On the other hand, the amount of scattered energy not only provides an indicator of orientation of fracture sets, but can also provide information about the fracture spacing. Analysis of synthetic seismograms from dry- and fluid-filled fractures indicates that dry-filled fractures show more significant scattering on seismic wavefields than fluid-filled ones, and offset-variations in P-wave amplitude are observable. We also analyze seismic response of an anticlinal trap model that includes a gas-filled fractured reservoir with high attenuation, which attenuates and distorts the so-called bright spot.

反射地震转换波成像方法研究

Multi-waves and multi-component get more and more attentions from oil industry. On the basis of existent research results, My research focuses on some key steps of OBC 4C datum processing. OBC datum must be preprocessed quite well for getting a good image. We show a flow chart of preprocess including attenuation of noise on multi-component datum、elimination ghost by summing P and Z and rotation of horizontal components. This is a good foundation for the coming steps about OBC processing. How to get exact converted point location and to analyze velocity are key points in processing reflection seismic converted wave data. This paper includes computing converted point location, analyzing velocity and nonhyperbolic moveout about converted waves. Anisotropic affects deeply the location of converted wave and the nonhyperbolic moveout. Supposed VTI, we research anisotropic effect on converted wave location and the moveout. Since Vp/Vs is important, we research the compute method of Vp/Vs from post-stack data and pre-stack data. It is a part of the paper that inversing anisotropic parameter by traveltime. Pre-stack time migration of converted wave is an focus, using common-offset Kirchhoff migration, we research the velocity model updating in anisotropic media. I have achieved the following results: 1) using continued Fractions, we proposed a new converted point approximate equation, when the offset is long enough ,the thomsen’s 2 order equation can’t approximate to the exact location of converted point, our equation is a good approximate for the exact location. 2) our new methods about scanning nonhyperbolic velocity and Vp/Vs can get a high quality energy spectrum. And the new moveout can fit the middle and long offset events. Processing the field data get a good result. 3) a new moveout equation, which have the same form as Alkhalifah’s long offset P wave moveout equation, have the same degree preciseness as thomsen’s moveout equation by testing model data. 4) using c as a function of the ratio offset to depth, we can uniform the Li’s and thomsen’s moveout equation in a same equation, the model test tell us choice the reasonable function C can improve the exact degree of Li’s and thomsen’s equation. 5) using traveltime inversion ,we can get anisotropic parameter, which can help to flat the large offset event and propose a model of anisotropic parameter which will useful for converted wave pre-stack time migration in anisotropic media. 6）using our pre-stack time migration method and flow, we can update the velocity model and anisotropic parameter model then get good image. Key words: OBC, Common converted Point (CCP), Nonhyperbolic moveout equation, Normal moveout correction, Velocity analysis, Anisotropic parameters inversion, Kirchhoff anisotropic pre-stack time migration, migration velocity model updating

旋转轴对称介质中的qP波一阶射线追踪

Elastic anisotropy is a very common phenomenon in the Earth’s interior, especial for sedimentary rock as important gas and oil reservoirs. But in the processing and interpretation of seismic data, it is assumption that the media in the Earth’s interior is completely elastic and isotropic, and then the methods based on isotropy are used to deal with anisotropic seismic data, so it makes the seismic resolution lower and the error on images is caused. The research on seismic wave simulation technology can improve our understanding on the rules of seismic wave propagation in anisotropic media, and it can help us to resolve problems caused by anisotropy of media in the processing and interpretation of seismic data. So researching on weakly anisotropic media with rotated axis of symmetry, we study systematically the rules of seismic wave propagation in this kind of media, simulate the process with numerical calculation, and get the better research results. The first-order ray tracing (FORT) formulas of qP wave derived can adapt to every anisotropic media with arbitrary symmetry. The equations are considerably simpler than the exact ray tracing equations. The equations allow qP waves to be treated independently from qS waves, just as in isotropic media. They simplify considerably in media with higher symmetry anisotropy. In isotropic media, they reduce to the exact ray tracing equations. In contrast to other perturbation techniques used to trace rays in weakly anisotropic media, our approach does not require calculation of reference rays in a reference isotropic medium. The FORT-method rays are obtained directly. They are computationally more effective than standard ray tracing equations. Moreover the second-order travel time corrections formula derived can be used to reduce effectively the travel time error, and improve the accuracy of travel time calculation. The tensor transformation equations of weak-anisotropy parameters in media with rotated axis of symmetry derived from the Bond transformation equations resolve effectively the problems of coordinate transformation caused by the difference between global system of coordinate and local system of coordinate. The calculated weak-anisotropy parameters are completely suitable to the first-order ray tracing used in this paper, and their forms are simpler than those from the Bond transformation. In the numerical simulation on ray tracing, we use the travel time table calculation method that the locations of the grids in the ray beam are determined, then the travel times of the grids are obtained by the reversed distance interpolation. We get better calculation efficiency and accuracy by this method. Finally we verify the validity and adaptability of this method used in this paper with numerical simulations for the rotated TI model with anisotropy of about 8% and the rotated ORTHO model with anisotropy of about 20%. The results indicate that this method has better accuracy for both media with different types and different anisotropic strength. Keywords: weak-anisotropy, numerical simulation, ray tracing equation, travel time, inhomogeneity

叠前时间偏移：各向异性、多次波衰减及应用

At present, in order to image complex structures more accurately, the seismic migration methods has been developed from isotropic media to the anisotropic media. This dissertation develops a prestack time migration algorithm and application aspects for complex structures systematically. In transversely isotropic media with a vertical symmetry axis (VTI media), the dissertation starts from the theory that the prestack time migration is an approximation of the prestack depth migration, based on the one way wave equation and VTI time migration dispersion relation, by combining the stationary-phase theory gives a wave equation based VTI prestack time migration algorithm. Based on this algorithm, we can analytically obtain the travel time and amplitude expression in VTI media, as while conclude how the anisotropic parameter influence the time migration, and by analyzing the normal moveout of the far offset seismic data and lateral inhomogeneity of velocity, we can update the velocity model and estimate the anisotropic parameter model through the time migration. When anisotropic parameter is zero, this algorithm degenerates to the isotropic time migration algorithm naturally, so we can propose an isotopic processing procedure for imaging. This procedure may keep the main character of time migration such as high computational efficiency and velocity estimation through the migration, and, additionally, partially compensate the geometric divergence by adopting the deconvolution imaging condition of wave equation migration. Application of this algorithm to the complicated synthetic dataset and field data demonstrates the effectiveness of the approach. In the dissertation we also present an approach for estimating the velocity model and anisotropic parameter model. After analyzing the velocity and anisotropic parameter impaction on the time migration, and based on the normal moveout of the far offset seismic data and lateral inhomogeneity of velocity, through migration we can update the velocity model and estimate the anisotropic parameter model by combining the advantages of velocity analysis in isotropic media and anisotropic parameter estimation in VTI media. Testing on the synthetic and field data, demonstrates the method is effective and very steady. Massive synthetic dataset、2D sea dataset and 3D field datasets are used for VTI prestack time migration and compared to the stacked section after NMO and prestack isotropic time migration stacked section to demonstrate that VTI prestack time migration method in this paper can obtain better focusing and less positioning errors of complicated dip reflectors. When subsurface is more complex, primaries and multiples could not be separated in the Radon domain because they can no longer be described with simple functions (parabolic). We propose an attenuating multiple method in the image domain to resolve this problem. For a given velocity model，since time migration takes the complex structures wavefield propagation in to account, primaries and multiples have different offset-domain moveout discrepancies, then can be separated using techniques similar to the prior migration with Radon transform. Since every individual offset-domain common-reflection point gather incorporates complex 3D propagation effects, our method has the advantage of working with 3D data and complicated geology. Testing on synthetic and real data, we demonstrate the power of the method in discriminating between primaries and multiples after prestack time migration, and multiples can be attenuated in the image space considerably.

通过第一性原理计算研究内地核的结构与成分

The determination of the composition and structure of the Earth’s inner core has long been the major subject in the study of the Earth’s deep interior. It’s widely believed that the Earth’s core is formed by iron with a fraction of nickel. However, light elements must exist in the inner core because the earth core is less dense than pure iron-nickel alloy (~2-3% in the solid inner core and ~6-7% in the liquid outer core). The questions are what and how much light element is there in the iron-nickel alloy. Besides the composition, the crystal structure of the iron with or without light element is also not well known. According to the seismological observations, the sound waves propagate 3-4% faster along the spin axis than in the equatorial plane. That means the inner core is anisotropic. The densest structure of iron-nickel alloy should be h.c.p structure under the very high pressures. However, the h,c,p structure does not propagate waves anisotropic ally. Then what is the structure of the iron-nickel alloy or the iron-nickle-light element alloy. In this study, we tried to predict the composition and the structure of the inner core through ab initio calculation of the Gibbs free energy, which is a function of internal energy, density and entropy. We conclude that the h.c.p structure is more stable than the b.c.c structure under high pressure and 0 K, but with the increase of temperature, the free energy of the b.c.c structure is decreasing much faster than the h.c.p structure caused by the vibration of the atomics, so the b.c.c structure is more stable at high temperatures. With the addition of light elements (S or Si or both), the free energy of b.c.c. decreases even faster, about 3at% of Si not only explains why the inner core is about 2-3 % lighter than the iron-nickle alloy, but also reasons why the inner core is anisotropic, since the b.c.c. structure becomes more stable than the h.c.p structure at 5500-6000K and b.c.c. is anisotropic in propagating seismic waves. Therefore, we infer that the inner core of the earth is formed by b.c.c iron and a fraction of nickel plus ~3at.% Si, with a temperature higher than 5500K, which is consistent with the studies from other approaches.

纵波方位属性裂缝检测技术综合研究

Fractured oil and gas reservoir is an important type of oil and gas reservoir, which is taking a growing part of current oil and gas production in the whole world. Thus these technologies targeted at exploration of fractured oil and gas reservoirs are drawing vast attentions. It is difficult to accurately predict the fracture development orientation and intensity in oil and gas exploration. Focused on this problem, this paper systematically conducted series study of seismic data processing and P-wave attributes fracture detection based on the structure of ZX buried mountain, and obtained good results. This paper firstly stimulated the propagation of P-wave in weak anisotropic media caused by vertical aligned cracks, and analyzed the rule of P-wave attributes’ variation associated with observed azimuth, such as travel-time, amplitude and AVO gradient and so on, and quantitatively described the sensitive degree of these attributes to anisotropy of fracture medium. In order to further study the sensitive degree of these attributes to anisotropy of fractures, meanwhile, this paper stimulated P-wave propagation through different types and different intensity anisotropic medium respectively and summarized the rule of these attributes’ variation associated with observed azimuth in different anisotropic medium. The results of these studies provided reliable references for predicting orientation, extensity and size of actual complicated cracked medium by P-wave azimuth attributes responses. In the paper, amounts of seismic data processing methods are used to keep and recover all kinds of attributes applied for fracture detection, which guarantee the high accurate of these attributes, thus then improve the accurate of fracture detection. During seismic data processing, the paper adopted the three dimensional F-Kx-Ky field cone filter technique to attenuate ground roll waves and multiple waves, then enhances the S/N ratio of pre-stack seismic data; comprehensively applying geometrical spread compensation, surface consistent amplitude compensation, residual amplitude compensation to recover amplitude; common azimuth processing method effectively preserves the azimuthal characteristics of P-wave attributes; the technique of bend ray adaptive aperture pre-stack time migration insures to obtain the best image in each azimuth. Application of these processing methods guaranteed these attributes’ accuracy, and then improved the accuracy of fracture detection. After comparing and analyzing a variety of attributes, relative wave impedance (relative amplitude) attribute is selected to inverse the orientation of fracture medium; attenuation gradient and corresponding frequency of 85% energy are selected to inverse the intensity of fracture medium; then obtained the fracture distribution characteristics of lower Paleozoic and Precambrian in ZX ancient buried mountains. The results are good accord with the characteristics of faults system and well information in this area.

复杂裂隙岩体统计断裂本构关系与变形参数研究

Rock mass is widely recognized as a kind of geologic body which consists of rock blocks and discontinuities. The deformation and failure of rock mass is not only determined by rock block，but also by discontinuity which is virtually more important. Mutual cutting and combination of discontinuities controlled mechanical property of rock mass. The complex cutting of discontinuities determine the intense anisotropy on mechanical property of rock mass，especially under the effect of ground stress. Engineering practice has show that the brittle failure of hard rock always occurs when its working stress is far lower than the yield strength and compressive strength，the failure always directly related to the fracture propagation of discontinuities. Fracture propagation of discontinuities is the virtue of hard rock’s failure. We can research the rock mass discontinuous mechanical properties precisely by the methods of statistical analysis of discontinuities and Fracture Mechanics. According to Superposition Principle in Fracture Mechanics，A Problem or C Problem could be chosen to research. Problem A mainly calculates the crack-tip stress field and displacement field on internal discontinuities by numerical method. Problem C calculate the crack-tip stress field and displacement field under the assumption of that the mainly rock mass stress field has been known. So the Problem C avoid the complex mutual interference of stress fields of discontinuities，which is called crack system problem in Fracture Mechanics. To solve Problem C, field test on stress field in the rock mass is needed. The linear Superposition of discontinuities strain energies are Scientific and Rational. The difference of Fracture Mechanics between rock mass and other materials can mostly expression as：other materials Fracture Mechanics mostly face the problem A，and can’t avoid multi-crack puzzle, while the Rock mass Fracture Mechanics answer to the Problem C. Problem C can avoid multi-discontinuities mutual interference puzzle via the ground stress test. On the basis of Problem C, Fracture Mechanics could be used conveniently in rock mass. The rock mass statistics fracture constitutive relations, which introduced in this article, are based on the Problem C and the Discontinuity Strain Energy linear superposition. This constitutive relation has several merits: first, it is physical constitutive relation rather than empirical; second, it is very fit to describe the rock mass anisotropy properties; third, it elaborates the exogenous factors such as ground stress. The rock mass statistics fracture constitutive relation is the available approach to answer to the physical, anisotropic and ground stress impacted rock mass problems. This article stand on the foundation of predecessor’s statistics fractures constitutive relation, and improved the discontinuity distributive function. This article had derived the limitation of negative exponential distribution in the course of regression analysis, and advocated to using the two parameter negative exponential distribution for instead. In order to solve the problems of two-dimension stability on engineering key cross-sectional view in rock mass, this article derived the rock mass planar flexibility tensor, and established rock mass two-dimension penetrate statistics fracture constitutive relation on the basis of penetrate fracture mechanics. Based on the crack tip plasticity research production of penetrate fracture, for example the Irwin plasticity equifinality crack, this article established the way to deal with the discontinuity stress singularity and plastic yielding problem at discontinuity tip. The research on deformation parameters is always the high light region of rock mass mechanics field. After the dam foundation excavation of XiaoWan hydroelectric power station, dam foundation rock mass upgrowthed a great deal of unload cracks, rock mass mechanical property gotten intricacy and strong anisotropy. The dam foundation rock mass mostly upgrowthed three group discontinuities： the decantation discontinuity, the steep pitch discontinuity, and the schistosity plane. Most of the discontinuities have got partial unload looseness. In accordance with ground stress field data, the dam foundation stress field greatly non-uniform, which felled under the great impaction of tectonic stress field, self-weight stress field, excavation geometric boundary condition, and excavation, unload. The discontinuity complexity and stress field heterogeneity, created the rock mass mechanical property of dam foundation intricacy and levity. The research on the rock mass mechanics, if not take every respected influencing factor into consideration as best as we can, major errors likely to be created. This article calculated the rock mass elastic modulus that after Xiao Wan hydroelectric power station dam foundation gutter excavation finished. The calculation region covered possession monolith of Xiao Wan concrete double-curvature arch dam. Different monolith were adopted the penetrate fracture statistics constitutive relation or bury fracture statistics constitutive relation selectively. Statistics fracture constitutive relation is fit for the intensity anisotropy and heterogeneity rock mass of Xiao Wan hydroelectric power station dam foundation. This article had contrastive analysis the statistics fracture constitutive relation result with the inclined plane load test actual measurement elastic modulus and RMR method estimated elastic modulus, and find that the three methods elastic modulus have got greatly comparability. So, the statistics fracture constitutive relations are qualified for trust. Generally speaking，this article had finished following works based on predecessors job: “Argumentation the C Problems of superposition principle in Fracture Mechanics, establish two-dimension penetrate statistics fracture constitutive relation of rock mass, argue the negative exponential distribution limitation and improve it, improve of the three-dimension berry statistics fracture constitutive relation of rock mass, discontinuity-tip plastic zone isoeffect calculation, calculate the rock mass elastic modulus on two-dimension cross-sectional view”. The whole research clue of this article inherited from the “statistics rock mass mechanics” of Wu Faquan（1992）.

三维非均匀各向异性介质中的纵波走时反演研究

Many observations show that seismic anisotropy is very common in the crust and upper mantle of the Earth. Seismic anisotropy can provide some clue about the changing and transporting process inside the earth. in recent years, abundant earthquake travel time data are accumulated, computers become more powerful, and these make the inversion of earthquake travel time data practical. In this thesis we studied the theory of elastic wave in anisotropic media, some formule for travel time inversion were derived. We present an iterative procedure to determine 21 elastic parameters from qP wave travel times. No a priori assumptions about heterogeneity and anisotropy of the model are made. The procedure is suitable for the case when we know nothing about the symmetry of anisotropy of the media, as well as for the case of earthquake travel time inversion which may contain various symmetry of anisotropy. The procedure is tested with a synthetic multiple-source offset VSP experiment. The results proved that the formulae are correct, and the procedure is practical. The results and the related theory indicate that the anisotropic inversion needs more rays than isotropic case. For a 2-D weak anisotropic (WA) medium, we need at least 5 rays in different directions to retrieve the elastic parameters on one grid point, and for a 3-D WA medium we need at least 15 rays in different directions to retrieve the elastic parameters on one grid point. The results also indicate that the starting background velocity has no influence on the final results, at least for the model we specified. Our results also show that insufficient illumination coverage will slow down the convergence rate, and make the results more sensitive to noise. We apply the procedure to a set of field travel time data. The data is from an artificial seismic observation. This observation is for locating micro-seismic events around a tunnel, its purpose is to find out if the digging process and the stress condition around the tunnel can generate micro-cracks. The size of this area is around 100m. The anisotropy derived from qP travel times is the same as the anisotropy showed by apparent velocities, and is also consistent with the anisotropy derived from S-wave splitting phenomena.

各向同性与各向异性介质中地震波旅行时正演模拟与反演成像

By seismic tomography, interesting results have been achieved not only in the research of the geosphere with a large scale but also in the exploration of resources and projects with a small scale since 80'. Compared with traditional inversion methods, seismic tomography can offer more and detailed information about subsurface and has been being paid attention by more and more geophysicists. Since inversion based on forward modeling, we have studied and improved the methods to calculate seismic traveltimes and raypaths in isotropic and anisotropic media, and applied the improved forward methods to traveltime tomography. There are three main kinds of methods to calculate seismic traveltime field and its ray path distribution, which are ray-tracing theory, eikonal equation by the finite-difference and minimum traveltime tree algorithm. In ray tracing, five methods are introduced in the paper, including analytic ray tracing, ray shooting, ray bending, grid ray tracing and rectangle grid ray perturbation with three points. Finite-difference solution of eikonal equation is very efficient in calculation of seismic first-break, but is awkward in calculation of reflection traveltimes. We have put forward a idea to calculate traveltimes of reflected waves using a combining way of eikonal equation method and other one in order to improve its capability of dealing with reflection waves. The minimum traveltime tree algorithm has been studied with emphases. Three improved algorithms are put forward on the basis of basic algorithm of the minimum traveltime tree. The first improved algorithm is called raypath tracing backward minimum traveltime algorithm, in which not only wavelets from the current source but also wavelets from upper source points are all calculated. The algorithm can obviously improve the speed of calculating traveltimes and raypaths in layered or blocked homogeneous media and keep good accuracy. The second improved algorithm is raypath key point minimum traveltime algorithm in which traveltimes and raypaths are calculated with a view of key points of raypaths (key points of raypths mean the pivotal points which determine raypaths). The raypath key point method is developed on the basis of the first improved algorithm, and has better applicability. For example, it is very efficient even for inhomogeneous media. Another improved algorithm, double grid minimum traveltime tree algorithm, bases upon raypath key point scheme, in which a model is divided with two kinds of grids so that the unnecessary calculation can be left out. Violent undulation of curved interface often results in the phenomenon that there are no reflection points on some parts of interfaces where there should be. One efficacious scheme that curved interfaces are divided into segments, and these segments are treated respectively is presented to solve the problem. In addition, the approximation to interfaces with discrete grids leads to large errors in calculation of traveltimes and raypaths. Noting the point, we have thought a new method to remove the negative effect of mesh and to improve calculation accuracy by correcting the traveltimes with a little of additional calculation, and obtained better results.

泥岩裂缝储层地震识别的理论和方法

During the exploration of fractured reservoirs, worldwide difficult problems will be encountered: how to locate the fractured zones, how to quantitatively determine the azimuth, density, and distribution of the fractures, and how to compute the permeability and porosity of the fractures. In an endeavor to solve these problems, the fractured shale reservoir in SiKou area of ShengLi oil field was chosen as a study area. A study of seismic predictive theory and methods for solving problems encountered in fractured reservoir exploration are examined herein. Building on widely used current fractured reservoir exploration techniques, new seismic theories and methods focusing on wave propagation principles in anisotropic medium are proposed. Additionally, integrated new seismic data acquisition and processing methods are proposed. Based on research and application of RVA and WA methods from earlier research, a new method of acoustic impedance varying with azimuth (IPVA) creatively is put forth. Lastly combining drilling data, well log data, and geologic data, an integrated seismic predictive method for cracked reservoir bed was formed. A summary of the six parts of research work of this paper is outlined below. In part one, conventional geologic and geophysical prediction methods etc. for cracked reservoir exploration are examined, and the weaknesses of these approaches discussed. In part two, seismic wave propagation principles in cracked reservoirs are studied. The wave equation of seismic velocity and attenuation factor in three kinds of fracture mediums is induced, and the azimuth anisotropy of velocity and attenuation in fracture mediums is determined. In part three, building on the research and application of AVA and WA methods by a former researcher, a new method of acoustic impedance creatively varying with azimuth (IPVA) is introduced. A practical software package utilizing this technique is also introduced. In part four, Base on previously discussed theory, first a large full azimuth 3d seismic data (70km~2) was designed and acquired. Next, the volume was processed with conventional processing sequence. Then AVA, WA, and IPVA processing was applied, and finally the azimuth and density of the fractures were quantitatively determined by an integrated method. Predictions were supported by well data that indicate the approach is highly reliable. in part five, geological conditions contributing to cracked reservoir bed formation are analyzed in the LuoJia area resulting in the discovery that the main fractured zones are related to fault distribution in the basin, that also control the accumulation of the oil and gas, the generation mechanisms and types of fractured shale reservoirs are studied. Lastly, by using full 3D seismic attributes, azimuth and density of cracked reservoir zones are successfully quantitative predicted. Using an integrated approach that incorporates seismic, geologic and well log data, the best two fractured oil prospects in LouJia area are proposed. These results herein represent a break through in seismic technology, integrated seismic predictive theory, and production technology for fractured reservoirs. The approach fills a void that can be applied both inside China, and internationally. Importantly, this technique opens a new exploration play in the ShengLi oil field that while difficult has substantial potential. Properly applied, this approach could play an important role toward stabilizing the oil field' production. In addition, this technique could be extended fracture exploration in other oil fields producing substantial economic reward.

弱各向异性介质中传播的地震波及其速度重建

We begin our studies to make the best of information of seismic data and carry out the description of cracks parameters by extracting anisotropic information. The researching contents are: (1) velocity and polarization anomaly of seismic wave (qP and qSV wave) in weak anisotropic media; (2) reflection seismic synthetic record in anisotropic media; (3) multiple scattering induced by cracks; (4) anisotropic structure inversion and velocity reconstruction with VSP (Vertical Seismic Profile) data; (5) multi-parameters analysis of anisotropy in time-domain and depth-domain. Then we obtain results as follows: (1) We achieve approximate relation of qP and qSV wave's velocity and polarization property in weak anisotropic media. At the same time, we calculate anisotropic velocity factors and polarization anomaly of several typical sedimentary rocks. The results show there are different anisotropic velocity factors and polarization anomaly in different rocks. It is one of the primary theoretical foundation which is expected to identify lithology; (2) We calculate reflection seismic synthetic record with theoretical model; (3) We simulate scattering induced by cracks with Boundary Element Method. Numerical studies show that in the presence of cracks; spatial and scale-length distributions are important and cannot be ignored in modeling cracked solids; (4) From traveltimes information of VSP data, we study the velocity parameter inversion of seismic wave under isotropic and anisotropic models, and its result indicate that the inversion imaging under anisotropic model will not destroy the original features of isotropic model, but it will bring on some bigger error if we adopt the method of isotropic model for anisotropic model data. Further more, basing on the study we develop the CDP mapping technology of reflecting structure under isotropic and anisotropic models, and we process real data as a trial of the methods; (5) We study the problem of initial model reconstruction of anisotropic parameters structure represented by Anderson parameter in depth domain for surface data.

各向异性介质中的地震面波和S波分裂

Surface wave propagation in the anisotropic media and S-wave splitting in China mainland are focused in this M.S. dissertation. We firstly introduced Anderson parameters in the research of surface wave propagation in the anisotropic media were deduced, respectively. By applying the given initial model to the forward calculation of Love wave, we compared dispersion curves of Love wave in the anisotropic media with the one in the isotropic media. the results show that, although the two kind of results are similar with each other, the effect of anisotropy can not be neglected. Furthermore, the variation of anisotropy factors will result in the variation of dispersion curves, especially for high-mode one. The method of grid dispersion inversion was then described for further tectonic inversion. We also deduced inversion equation on the condition that the layered media is anisotropic, and calculated the phase-velocity partial derivatives with respect to the model parameters, P- and S-wave velocities, density, anisotropic parameters for Rayleigh wave and Love wave. Having analyzed the results of phase-velocity partial derivatives, we concluded that the derivatives within each period decreased with the depth increasing, the phase-velocity of surface wave is sensitive to the S-wave velocities and anisotropic factors and is not sensitive to the densities of layers. Dispersion data of Love wave from the events occurred during the period from 1991 to 1998 around the Qinghai and Tibet Plateau, which magnitudes are more than 5.5, have been used in the grid dispersion inversion. Those data have been preprocessed and analyzed in the F-T domain. Then the results of 1°*1° grid dispersion inversion, the pure path dispersion data, in the area of Qianghai and Tibet Plateau were obtained. As an example, dispersion data have been input for the tectonic inversion in the anisotropic media, and the results of anisotropic factors under the region of Qianghai and Tibet Plateau were initially discussed. As for the other part of this dissertation. We first introduced the phenomena of S-wave splitting and the methods for calculation the splitting parameters. Then, We applied Butterworth band-pass filter to S-wave data recorded at 8 stations in China mainland, and analyzed S-wave splitting at different frequency bands. The results show the delay time and the fast polarization directions of S-wave splitting depend upon the frequency bands. There is an absence of S-wave splitting at the station of Wulumuqi (WMQ) for the band of 0.1-0.2Hz. With the frequency band broaden, the delay time of S-wave splitting decreases at the stations of Beijing (BJI), Enshi (ENH), Kunming (KMI) and Mudanjiang (MDJ); the fast polarization direction at Enshi (ENH) changes from westward to eastward, and eastward to westward at Hailaer (HIA). The variations of delay time with bands at Lanzhou (LZH) and qiongzhong (QIZ) are similar, and there is a coherent trend of fast polarization directions at BJI, KMI and MDJ respectively. Initial interpretations to the results of frequency band-dependence of S-wave splitting were also presented.

三维油藏数值模拟正反演算法研究

The dynamic prediction of complex reservoir development is one of the important research contents of dynamic analysis of oil and gas development. With the increase development of time, the permeabilities and porosities of reservoirs and the permeability of block reservoir at its boundaries are dynamically changing. How to track the dynamic change of permeability and porosity and make certain the permeability of block reservoir at its boundary is an important practical problem. To study developing dynamic prediction of complex reservoir, the key problem of research of dynamic prediction of complex reservoir development is realizing inversion of permeability and porosity. To realize the inversion, first of all, the fast forward and inverse method of 3-dimension reservoir simulation must be studied. Although the inversion has been widely applied to exploration and logging, it has not been applied to3-dimension reservoir simulation. Therefore, the study of fast forward and inverse method of 3-dimension reservoir simulation is a cutting-edge problem, takes on important realistic signification and application value. In this dissertation, 2-dimension and 3-dimension fluid equations in porous media are discretized by finite difference, obtaining finite difference equations to meet the inner boundary conditions by Peaceman's equations, giving successive over relaxation iteration of 3-dimension fluid equations in porous media and the dimensional analysis. Several equation-solving methods are compared in common use, analyzing its convergence and convergence rate. The alternating direction implicit procedure of 2-dimension has been turned into successive over relaxation iteration of alternating direction implicit procedure of 3-dimension fluid equations in porous media, which possesses the virtues of fast computing speed, needing small memory of computer, good adaptability for heterogeneous media and fast convergence rate. The geological model of channel-sandy reservoir has been generated with the help of stochastic simulation technique, whose cross sections of channel-sandy reservoir are parabolic shapes. This method makes the hard data commendably meet, very suit for geological modeling of containing complex boundary surface reservoir. To verify reliability of the method, theoretical solution and numerical solution are compared by simplifying model of 3-dimension fluid equations in porous media, whose results show that the only difference of the two pressure curves is that the numerical solution is lower than theoretical at the wellbore in the same space. It proves that using finite difference to solve fluid equations in porous media is reliable. As numerical examples of 3-dimension heterogeneous reservoir of the single-well and multi-well, the pressure distributions have been computed respectively, which show the pressure distributions there are clearly difference as difference of the permeabilities is greater than one order of magnitude, otherwise there are no clearly difference. As application, the pressure distribution of the channel-sandy reservoir have been computed, which indicates that the space distribution of pressure strongly relies on the direction of permeability, and is sensitive for space distributions of permeability. In this dissertation, the Peaceman's equations have been modified into solving vertical well problem and horizontal well problem simultaneously. In porous media, a 3D layer reservoir in which contain vertical wells and horizontal wells has been calculated with iteration. For channel-sandy reservoir in which there are also vertical wells and horizontal wells, a 3D transient heterogeneous fluid equation has been discretized. As an example, the space distribution of pressure has been calculated with iteration. The results of examples are accord with the fact, which shows the modification of Peaceman's equation is correct. The problem has been solved in the space where there are vertical and horizontal wells. In the dissertation, the nonuniform grid permeability integration equation upscaling method, the nonuniform grid 2D flow rate upscaling method and the nonuniform grid 3D flow rate upscaling method have been studied respectively. In those methods, they enhance computing speed greatly, but the computing speed of 3D flow rate upscaling method is faster than that of 2D flow rate upscaling method, and the precision of 3D flow rate upscaling method is better than that of 2D flow rate upscaling method. The results also show that the solutions of upscaling method are very approximating to that of fine grid blocks. In this paper, 4 methods of fast adaptive nonuniform grid upscaling method of 3D fluid equations in porous media have been put forward, and applied to calculate 3D heterogeneous reservoir and channel-sandy reservoir, whose computing results show that the solutions of nonuniform adaptive upscaling method of 3D heterogeneous fluid equations in porous media are very approximating to that of fine grid blocks in the regions the permeability or porosity being abnormity and very approximating to that of coarsen grid blocks in the other region, however, the computing speed of adaptive upscaling method is 100 times faster than that of fine grid block method. The formula of sensitivity coefficients are derived from initial boundary value problems of fluid equations in porous media by Green's reciprocity principle. The sensitivity coefficients of wellbore pressure to permeability parameters are given by Peaceman's equation and calculated by means of numerical calculation method of 3D transient anisotropic fluid equation in porous media and verified by direct method. The computing results are in excellent agreement with those obtained by the direct method, which shows feasibility of the method. In the dissertation, the calculating examples are also given for 3D reservoir, channel-sandy reservoir and 3D multi-well reservoir, whose numerical results indicate: around the well hole, the value of the sensitivity coefficients of permeability is very large, the value of the sensitivity coefficients of porosity is very large too, but the sensitivity coefficients of porosity is much less than the sensitivity coefficients of permeability, so that the effect of the sensitivity coefficients of permeability for inversion of reservoir parameters is much greater than that of the sensitivity coefficients of porosity. Because computing the sensitivity coefficients needs to call twice the program of reservoir simulation in one iteration, realizing inversion of reservoir parameters must be sustained by the fast forward method. Using the sensitivity coefficients of permeability and porosity, conditioned on observed valley erosion thickness in wells (hard data), the inversion of the permeabilities and porosities in the homogeneous reservoir, homogeneous reservoir only along the certain direction and block reservoir are implemented by Gauss-Newton method or conjugate gradient method respectively. The results of our examples are very approximating to the real data of permeability and porosity, but the convergence rate of conjugate gradient method is much faster than that of Gauss-Newton method.