We present a study on using a multigrid scheme for a preconditioner for acoustic frequency-domain finite-difference (FDFD) forward modeling. To achieve fourth-order accuracy, we employ the perfectly matched layer (PML) approach. The multi-grid preconditioner combined with the bi-conjugate gradient stabilized (BiCGStab) iterative solver has been successfully applied to solve acoustic forward problems by using the radiation boundary condition, the sponge layer absorbing boundary condition, and the first-order PML boundary condition. The condition number of the discrete Helmholtz equation is highly dependent on the boundary conditions. The high-order PML enables us to match the accuracy of the boundary condition with the accuracy of the finite-difference (FD) scheme in the interior domain, especially for frequency-domain approaches. The geometric multigrid is employed to construct a preconditioner for a fourth-order FDFD forward solver equipped with the PML boundary condition. For efficiency this preconditioner is constructed using a second-order FD scheme with a negligible attenuation function inside the PML domain. The preconditioner is used for accelerating the convergence of the FDFD forward solver for cases where the discretization grids are over-sampled (the number of discretization points per minimum wavelength is greater than 10). The number of multi-grid levels is also chosen adaptively depending on the number of discretization points. Our study shows that the multigrid preconditioner can speed-up the convergence the BiCGStab solver total computational time by a factor of three for cases with over-sampled discretization grids. We also observe that the BiCGStab solver using an accurate PML boundary condition does not have convergence problems as one may encounter with iterative solvers using the radiation boundary condition.
We propose two new boundary conditions to regulate coherent reflections from the model boundaries in numerical solutions of wave equations. Both boundary conditions have the common feature that the boundary condition is varied with respect to time. The first boundary condition expands or contracts the computational model during a modeling simulation. The effect is to cause a Doppler shift in the reflected wavefield that can be used to shift energy outside a frequency band of interest. Additionally, when the computational domain is expanding, the range of possible incidence angles on the boundary is restricted. This can be used to increase the effectiveness of many existing absorbing boundary conditions that are more effective for incidence angles close to normal. The second boundary condition is an extension of random boundaries. By carefully changing the realization of a random boundary over time, a more diffusive wavefield can be simulated. We show results with 2D numerical simulations of the scalar wave equation for both these boundary conditions. While the first boundary condition has application to modeling, both of these boundary conditions have potential applications within algorithms that rely upon modeling kernels, such as reverse-time migration and full-waveform inversion.
We demonstrate that implicit time integration methods for second-order-in-time wave equations can be derived from rational expansions to the cosine function of pseudo-differential operators. Using the scalar wave equation as an example, we give complete stability condition and grid dispersion analysis results for general implicit time integration methods. Furthermore, we propose an optimization method to develop unconditionally stable implicit time stepping schemes.
We have developed new heterogeneous 3D second and fourth order staggered-grid finite-difference schemes for modeling seismic wave propagation in the Laplace-Fourier domain. Recent interest in full waveform inversion in the Laplace-Fourier domain has been the motivation for the development for these types of wave-field simulators. Our approach is based on the principles of an integral equation approximation technique for the velocity-stress formulation in the Cartesian coordinate system. The fourth-order scheme is obtained by the combination of integral identities for the two elementary volumes — “small” and “large” around nodes where the wave variables are defined. The final matrix formulation for the fourth or second-order scheme is performed only for velocity components and the resulting linear system is solved by Krylov iterative methods. We have applied these simulators for the investigation of the wave fields of the SEG/EAGE model in the Laplace-Fourier domain along with other test models. Our eventual goal is to embed it in an inversion scheme for joint seismic-electromagnetic imaging.
Xu, Feng (SOUTHWEST PETROLEUM UNIVERSITY) | Guo, Xiao (the State Key Lab of Oil and Gas Reservoir Geology and Exploitation) | Wang, Wanbin (the State Key Lab of Oil and Gas Reservoir Geology and Exploitation) | Zhang, Nan (the State Key Lab of Oil and Gas Reservoir Geology and Exploitation) | Jia, Sa (the State Key Lab of Oil and Gas Reservoir Geology and Exploitation) | Wang, Xiaoqin
Low permeability reservoir occupies an important position in China's petroleum industry, as an important oil industry resource over a period of time in future, which helps to increase reserves and production. Surfactant flooding is one of the most effective ways to improve development effect in the low permeability reservoirs, which can reduce the injection pressure and increase injection rate, thereby enhancing oil recovery. Common mechanisms that how the surfactants enhance oil recovery are discussed in the paper. On the basis of assumptions, a seepage model of surfactant flooding was built to describe the complex process. Experiments was done to determine the surfactant's performance, such as the relationship between surfactant concentration and oil/water interfacial tension, the relationship between surfactant concentration and rock adsorbance, and the relationship between surfactant concentration and water viscosity. Take the Yanchang low permeability reservoir for example, where the study of surfactant flooding is conducted after the optimal water injection plan to evaluate the effect of different surfactant concentration on the development results in the low permeability reservoir. The simulation results show that surfactant flooding has the role of enhancing oil recovery, and that the optimal surfactant concentration is 2%, which can enhance oil recovery by the percentage of 0.22.
This paper deals with a time stepping method for the finite element simulation of two phase flow hydraulic and mechanical (H2M) coupled processes in porous media, which is a common phenomenon in geological applications such as CO2 storage facilities. The computation task arising from the numerical modeling of H2M coupled processes in such real geological application is intensive. Therefore, the high performance computing is of interest to the corresponding researchers. In the present study, we present a time stepping method with PI (proportional and integral feedback) automatic control to improve the computation efficiency of the modeling of H2M coupled processes. We apply the PI control to solve the nonlinear coupled partial differential equations with a first order finite difference scheme for time discretization and the Picard method for linearization. The efficiency of the present method is demonstrated by applying it to a CO storage benchmark.
To reduce anthropogenic greenhouse gas emissions into the atmosphere, the carbon dioxide capture and storage (CCS) concept is introduced by some researchers as an emerging transition technology [1,2]. The study of CCS is therefore under active consideration recently. According to various studies, deep saline aquifers provide the most substantial carbon dioxide storage capacity [3,4,5,6,7], and are often located near possible CO2 sources such as coal-fired power plants.
To ascertain migration and trapping of CO2 in the formations and assess the capacity and the safety (possible leakage) of the reservoir, the numerical simulation of injection and spreading of carbon dioxide in the underground is essential for understanding the physical and chemical processes at different length and time scales. In the numerical analysis of time dependent thermo-hydraulic processes in porous media, the time stepping is a crucial issue for numerical stability and computational efficiency. Practically, the fixed time step size does not often satisfy the stability and efficiency requirements in solving problems that exhibit complexity in geometry and nonlinearity in material properties. Therefore, adaptive time stepping methods including high-order integration have been developed and are widely applied . Among the available adaptive time stepping methods, the well-known techniques for prediction of the time step size h are e.g. Courant number approach based on Courant-Friedrichs- Lewy condition  for the finite difference method, primary variable based prediction (e.g. those presented in ref. [10,11] and local error control methods (cf.[8,12,13]). The local error control methods especially those based on theoretical control ideas are problem independent for any numerical methods for ODEs [8,12,13,14]. For nonlinear equations, the theory based automatic controls such as P (proportional feedback) or PI (proportional and integral feedback) permit stable and efficient time stepping [8,14] for numerical solver. The present work is subjected to apply the adaptive time stepping with automatic control for the finite element modeling of two phase flow hydraulic and mechanical coupled processes in CO2 storage facilities.
To this purpose, we present an approach of PI (proportional and integral feedback)  automatic time stepping for modeling the problems with different coupled physical processes. Within the context of the presented time stepping approach, each process uses the time step size predicted by the PI control itself to guarantee the stability of the simulation of each process under coupling.
Controlled blasting techniques are used to control overbreak and to aid in the stability of the remaining rock formation. The less competent the rock mass itself is, the more care has to be taken in avoiding damage. Presplitting is one of the most common methods which is used in many open pit mining and surface blast design. The purpose of presplitting is to form a fracture plane across which the radial cracks from the production blast cannot travel. Presplitting should be thought of as a protective measure to keep the final wall from being damaged by the production blasting. The purpose of this study is to investigate of effect of presplitting on the generation of a smooth wall in a rock domain under a surface blast process. The 2D distinct element code was used for simulation of presplitting in a rock slope. The blast load history as a function of time applied to the inner wall of each blasthole. Important parameters that were considered in the analysis were stress tensor and fracturing pattern. The blast loading magnitude and blasthole spacing were found to be very significant in the final results.
Drilling and blasting continues to be an important method of block production and block splitting. Drill and blast technique has a disadvantage that sometimes it produces cracks in uncontrolled manner and also produces micro cracks in the block as well as in remaining rock, if not carefully carried out. Recovery by this method is low as compared to other methods. Therefore, attempts have been made to develop controlled growth of crack in the desired direction. The control of fractures in undamaged brittle materials is of considerable interest in several practical applications including rock fragmentation and overbreak control in mining [1–3]. One way of achieving controlled crack growth along specific directions and inhibit growth along other directions is to generate stress concentrations along those preferred directions. Several researchers have suggested a number of methods for achieving fracture plane control by means of blasting. Fourney et al.  suggested a blasting method which utilizes a ligamented split-tube charge holder. Nakagawa et al.  examined the effectiveness of the guide hole technique by model experiments using acrylic resin plates and concrete blocks having a charge hole and circular guide holes. Katsuyama et al.  suggested a controlled blasting method using a sleeve with slits in a borehole. Mohanty [7,8] suggested a fracture plane control technique using satellite holes on either side of the central pressurized hole, and demonstrated its use through laboratory experiments and field trials in rock. Nakamura et al.  suggested a new blasting method for achieving crack control by utilizing a charge holder with two-wedge-shaped air cavities. Nakamura  performed model experiments to examine the effectiveness of the guide hole with notches. Cho et al.  performed experiments using a notched charge hole to visualize fracturing and gas flow due to detonation ofexplosives
Quasi-static mechanical testing of granular structures (representing rock) is simulated using a Generalized Interpolation Material Point method (GIMP). The numerical analysis is carried out by representing the rock as a dense packing of non-uniformly-sized particles that are cemented in the vicinity of points of contact. Mechanical behavior of the cement is imparted by cohesive zone elements. With progressively increasing application of far-field stresses, local stresses increase to the point where the cohesive zone elements lose load-bearing capacity. The sensitivity of this decohesion to the number of cohesive zone segments was investigated (i.e., the discretization of cemented zones into conjoined cohesive elements). A parametric study also varied the traction strengths for these cohesive features. Cementitious micro-properties were estimated to represent cohesive characteristics of typical porous rock materials and to result in representative peak-loading behavior for the aggregate sample. Uniaxial compression and Brazilian tests were numerically performed on granular assemblies with rock-like properties. INTRODUCTION
Experimental testing of rocks under various loading conditions has revealed numerous complex mechanisms occurring at macro- and microscopic levels. The results of these phenomena are manifested in an emerging mechanical response and ultimately failure of the rock [1-8]. This preliminary program has been carried out to develop a model that will simulate these micro- and/or macroscale events and replicate experimental results from laboratory testing of rock samples. In this investigation, a cohesive zone model in the Generalized Interpolation Material Point Method has been implemented to simulate representative laboratory loading protocols.
The behavior of granular materials has been extensively studied by various researchers. Cundall & Strack  performed numerical calculations for interparticle forces within assemblies of discs and spheres. Chang et.al  derived closed-form solution for random packing under low levels of deviatoric stress. Potyondy and Cundall  proposed the bonded-particle model for describing the mechanical behavior of packing of non-uniform sized circular or spherical particles. Particle-type models of rocks have been used to address various phenomenona in rock in geomechanics [12-13]. Numerical/analytical concepts of cohesion in granular situations are not new. Barenblatt  used the concept of equivalent traction on contacting surfaces to describe localized separation under progressive loading. A number of researchers have used the cohesive layer modeling approach to represent localized failures in concrete [15-18]. These models offer advantage of incorporating non-linear material properties with mesh independent post-localization behavior.
As is no surprise, the candidate samples are represented by an assemblage of grains immobilized by frictional resistance and by cementation at contact points (as well as any pore filling cement). In these demonstrations, the grains were treated as spherical or with realistic structure determined from appropriate tomography. Cementation between the grains was described by cohesive zone elements. These cohesive zone elements were characterized using published traction-displacement relationships with peak strength. Numerical simulations of uniaxial compressive testing and laboratory loading are shown to demonstrate deformation and strength of typical granular assemblies.
In this simulation exercise, porous rocks have been represented as cemented, granular media with following assumptions:
Rock masses consist of intact rock and discontinuities such as faults, joints and bedding planes. The presence of such discontinuities in rock masses dominates the response of jointed rock masses to static and dynamic loading. This paper focuses on the propagation and dynamic effects of blast waves in faulted rock masses. In order to investigate the effect of faults, a numerical simulation was conducted. The 2D distinct element code (UDEC) was used to model fault effect on rock failure and stress distribution through the rock mass due to blast wave propagation. The blast loading history was simplified and applied to the blasthole walls. Accordingly, the interaction of explosive energy transferred to the rock mass from the blasthole pressure was examined as a function of distance to the fault plane from the blasthole. A Mohr-Coulomb material model was used for host rock to allow for plastic failure calculations. The conducted numerical study describes the role of fault in blasting in a qualitative manner. On the other hand, a free face boundary was considered as a common blast operation which is conducted in surface mining.
The process of rock fragmentation by blasting is a complicated phenomenon which is controlled by many variables and parameters. Considering all this parameters in a single analysis is not possible at the present time, especially when some of them are not clearly understood yet and the effect of others is difficult to quantify. In most blasting practices, empirical or semi-empirical techniques are used for blast design and fragmentation analysis. These techniques are based on information obtained for certain range of rock types and blasting conditions and cannot be generalized for all blasting conditions. With regard to the limitations of empirical methods, numerical methods are viable tools to further understand and illustrate the fragmentation process. Application of numerical methods in blasting allows for consideration of complex boundary conditions, material non-linearity, dynamic material behavior, geometric non-linearity and complexity associated with blasting operations. The nature and degree of heterogeneity of the rock mass is very important in blast design. That is, discontinuities such as joints, bedding planes, faults, and soft seams can allow the explosive''s energy to be wastefully dissipated rather than perform the work intended. In some cases, the discontinuities can dominate the fracture pattern produced by the explosive, and the influence of the structural geology often overshadows that of the rock''s mechanical and physical properties. Best fragmentation is usually obtained where the face is parallel to the major discontinuity set [1 J. The last few decades have seen a variety of numerical studies on the blast-induced waves and their propagation in rock masses with much efforts being placed on the study of dynamic responses of continuous rock masses under blast loading [2,3 J. However, rock masses encountered in reality generally contain geological discontinuities (e.g. joints, faults and bedding p1anes).The properties of rock masses are determined by both the properties of the intact rock and the discontinuities [4,5].
Reliable methods must be used to perform a realistic assessment of mine roof support requirement and address the geotechnical risks associated with longwall mining. Dependable tools provide a safe working environment, increased production, efficient management of resources and reduce environmental impacts of mining. Although various methods, for example, analytical, experimental and empirical are being adopted in mining, in recent days numerical tools are becoming popular due to the advancement in computer hardware and numerical methods. Empirical rules based on past experiences do provide a general guide, however due to the heterogeneous nature of mine geology (i.e. none of the mine sites are identical), numerical simulations of mine site specific conditions would lend better insights into some underlying issues. The paper highlights the use of a continuum mechanics based tool in coal mining with a mine scale model. The continuum modeling can provide close to accurate stress fields and deformation. The paper describes the use of existing mine data to calibrate and validate the model parameters, which then are used to assess geotechnical issues related with installing a new high capacity longwall mine at the mine site. A variety of parameters, for example, chock convergences, caveability of top coal and overlying sandstones have been estimated.
The paper describes numerical simulations of longwall mining. Geotechnical specifications depend on the geological information of mines, so none of the mines can be considered to be identical and the experience from any previous mine is not completely portable. In some respect the past experience can be used as a preface step in planning the new mine. The detailed investigation of site specific issues, for example chock capacity requirement, strata caving behaviour, roof displacement and top coal caveability (in case of longwall top coal caving, LTCC) should be evaluated appropriately in order to gain confidence in mine design. The roof support system is critical in successful mining operations . The requirement and selection of roof support systems depend on site specific geotechnical conditions , for example nature and strength of overlaying strata, strength, orientation and height of mining, panel geometry and layouts, and the capital cost. These factors, which may vary in various mine sites, necessitate selecting an optimum roof support system. So the understanding of response of support capacity under local geotechnical conditions is vital for a site specific mining environment. Poor ground condition combined with poor choice of support may cause face instability [3,4]. Understanding the causes of face instability specifically determining whether they were simply the result of insufficient support capacity or whether the deployment of better mining techniques or more appropriate support could have avoided the events is important. Support design is based primarily on the replacement of the extracted coal with mechanized support capable of controlling the deformation in the immediate roof, the design being of sufficient capacity to allow effective mechanized extraction at the desired production rate. An optimum capacity of roof support systems provides safe as well as efficient mining operations .