Unstructured Tetrahedral Grids Research Articles

A three-dimensional numerical well-test model for pressure transient analysis in fractured horizontal wells with secondary fractures

During oil and gas reservoir development, multi-stage horizontal wells (MFHWs) and hydraulic fracturing techniques can effectively increase estimated ultimate recovery. However, there still lacks an understanding of the three-dimensional (3D) pressure transient behaviors of multi-stage fractured horizontal wells with secondary fractures. To narrow this gap, a three-dimensional numerical well-test model based on a discrete fracture model and unstructured tetrahedral grids is developed to study the pressure transient behaviors of MFHWs with secondary fractures. The pressure transient solutions of MFHWs with secondary fractures have been demonstrated by model verifications. The results show that the proposed model can accurately capture the complex transient flow around fractures, including early radial flow that is not easily captured by two-dimensional numerical well test models. The proposed model classifies the flow regimes of a MFHW as: wellbore storage and skin effects, early radial flow, bilinear flow, linear flow, elliptical flow, pseudo radial flow, and pseudo-boundary dominated flow. It is found that the fracture geometry has a relatively large effect on the shape of the pressure derivative curve in this work. The hydraulic fracture half-length has the greatest impact on the pressure transient behaviors of the MFHW, followed by fracture height and secondary fracture half-length, as found in this study. Additionally, fracture parameters are evaluated, and actual well testing data are interpreted, taking into account the fracture height. This work is meaningful to understand the three-dimensional pressure transient behaviors of MFHWs with secondary fractures.

Simulation of two-phase flow in 3D fractured reservoirs using a projection-based Embedded Discrete Fracture Model on Unstructured tetrahedral grids (pEDFM-U)

The modeling of fluid flow in heterogeneous, anisotropic and fractured porous media is relevant in many applications, including hydrocarbon and groundwater extraction, dispersion of contaminants, hydrogen or carbon dioxide (CO2) storage. Thus, accurate and scalable simulation of fluid flow through these formations continues to be a great challenge. The presence of fractures that are explicitly modeled, ranging from flow barriers to highly conductive channels, significantly increases the complexity of the numerical simulation. The Embedded Discrete Fracture Model (EDFM) produces results that can be as accurate as those obtained using equidimensional Discrete Fracture Models (DFM) at a much lower computational cost for highly conductive fractures, but it is not adequate for impermeable barriers. In contrast, the pEDFM (projection-based Embedded Discrete Fracture Model) produces accurate results for both, channels and barriers by using additional matrix-fracture connectivities. In its current versions, it is restricted to cartesian and corner-point grid geometries and to the classical Two Point Flux Approximation (TPFA) scheme. In this work, for the first time, the pEDFM is extended to handle unstructured tetrahedral meshes (pEDFM-U). In addition, interface fluxes are approximated by using the Multipoint Flux Approximation method with a Diamond stencil (MPFA-D). This allows for simulation of highly heterogeneous and anisotropic geo-models. The developed method is robust and can deal with full permeability tensors on arbitrary tetrahedral meshes. The advective terms are discretized considering an implicit First Order Upwind (FOU) method. Our method is implemented within the DARSim (Delft Advanced Reservoir Simulation) open-source simulator framework. Through several test cases, the proposed method showed to be robust and capable to accurately capture the effects of both high and low permeability fractures for general tetrahedral meshes, under arbitrary heterogeneous and anisotropic permeability tensors for the rock matrix.

An Efficient Parallel Iterative Solver for Controlled-source Electromagnetic 3-D Adaptive Forward Modeling in General Anisotropic Media

Controlled-source electromagnetic (CSEM) surveys are moving increasingly toward realistic and complicated scenarios where the survey region may contain undulated topography, complex geometries, and electrical anisotropy media. In this paper, the adaptive finite element numerical method is employed to discrete the total electric field equation, which can provide precise electromagnetic (EM) responses even with a coarse initial mesh. The unstructured tetrahedral grid is employed to effectively address arbitrary irregular geometry and mountain terrain. Then, the flexible generalized minimum residual solver (FGMRES), auxiliary-space Maxwell pre-conditioner, and grid division technique were used to solve the large-scale linear system of equations, which can stably solve ill-conditioned problems using fewer computing resources. Finally, we conducted a numerical experiment via our newly proposed forward modeling scheme on a synthetic model with multi-contrast electrical anisotropy. It validated that accurate EM fields could be obtained against the semi-analytic solutions, and this iterative solver has good robustness for various anisotropic media. As a result, we have developed a state-of-the-art 3D CSEM anisotropic forward modeling engine, which can quickly and accurately deal with large-scale and complex geo-models.

Three-Dimensional Separate and Joint Inversions of Multi-Component Frequency-Domain Airborne Electromagnetic Data: Synthetic Model Studies.

Airborne electromagnetic (AEM) surveys using airborne mobile platforms enable rapid and efficient exploration of areas where groundwork is difficult. They have been widely used in fields such as shallow resource exploration and environmental engineering. Three-dimensional AEM inversion is the main technique used in fine structural interpretation. However, most current methods focus on separate component data inversions, which limit the kinds of structures that can be recovered in the inversion results. To address this issue, a method for the robust 3D joint inversion of multicomponent frequency-domain AEM data was developed in this study. First, a finite element method based on unstructured tetrahedral grids was used to solve the forward problem of frequency-domain AEM data for both isotropic and anisotropic media. During inversion, a limited-memory quasi-Newton (L-BFGS) method was used to reduce the memory requirements and enable the joint inversion of large-scale multicomponent AEM data. The effectiveness of our algorithm was demonstrated using synthetic models for both isotropic and anisotropic cases, with 5% Gaussian noise added to the modeling data to simulate the measured data for separate and joint inversions. The results of the synthetic models show that joint inversion has advantages over separate inversion in that it enables the recovery of finer underground structures and provides a novel approach for the fine interpretation of frequency-domain AEM data.

Influence of the spatial grid type on the result of calculating the neutron fields in the nuclear power plant shielding

The paper considers the influence of the spatial grid type on the result of solving the equation of neutron transport in the nuclear power plant (NPP) shielding. Neutron fields were calculated in a realistic model of a liquid metal cooled fast neutron tank reactor with an integral equipment layout. Structured cubic and unstructured hexahedral grids (pmsnsys and FRIGATE codes) and unstructured tetrahedral and prismatic grids (RADUGA T code) are used. Limiting values of the group fluxes averaged over the material zones for refined grids have been obtained. It has been shown that the calculation results depend on the type of approximation for the curvilinear inner boundaries between the material zones rather than on the grid cell type (cube, hexahedron, tetrahedron, prism). Using “toothed” approximations for curvilinear boundaries leads to an increase in the area of the boundaries, as well as to the neutron flux refraction condition arising on them. These effects lead to an upward bias in the transport equation solution, and for all energy groups. Conclusion. When solving an equation of neutron transport in the NPP shielding by a grid technique, it is necessary to use grids other than leading to “toothed” approximations of the inner boundaries. Tetrahedral or prismatic grids, or grids of arbitrary hexahedrons can be recommended, as well as composite grids in which cubic cells are located inside the material zone, and hexahedron cells are located near the zone boundary.

3-D anisotropic modelling of geomagnetic depth sounding based on unstructured edge-based finite-element method

SUMMARY Geomagnetic depth sounding (GDS) is a geophysical electromagnetic (EM) method that studies the deep structure and composition of the Earth by using long-period EM signals from geomagnetic observatories and satellites. In this paper, a 3-D anisotropic GDS modelling algorithm is developed. The curl–curl equation is discretized using the edge-based finite-element method on unstructured tetrahedral grids. In order to solve the computationally demanding problem of EM modelling on a global scale, the complex linear system is first separated into the equivalent real linear systems and then the real system is iteratively solved by the flexible generalized minimum residual method with a block diagonal pre-conditioner. This will greatly reduce the condition number of the linear system and thus speed up the solution process. We verify the accuracy of the proposed algorithm by comparing our results with the existing methods. After that, we design a subduction zone model to simulate the EM responses under isotropic and anisotropic environments, respectively. The numerical results show the high efficiency of the proposed algorithm and the response differences between isotropic and anisotropic models. This research can provide theoretical and technical support for the high-accuracy and efficient inversion of GDS data for the geo-dynamic study.

Analyzing the effect of dipping anisotropy on 3D marine controlled-source electromagnetic response based on a goal-oriented adaptive finite element method

This paper presents a goal-oriented adaptive edge-based finite element (FE) algorithm for modelling marine controlled-source electromagnetic (CSEM) responses in three-dimensional (3-D) dipping anisotropic conductive media. We adopt the secondary field formulation for quasi-static Maxwell's equations to avoid the source singularity. The algorithm is implemented on unstructured tetrahedral grids, which allow for the simulation of complex model geometries, including tilted electrical interfaces and bathymetry. We validate the adaptive finite element code against a one-dimensional (1-D) layered model with vertical anisotropy. Further, we compare the 3-D FE numerical solutions with those obtained by a 2-D FE code for a two-dimensional dipping anisotropic model. The modelling results from the new 3-D FE algorithm are in good agreement with both the analytical solutions derived from the 1-D code and the numerical results from the anisotropic 2-D FE algorithm. We then use the algorithm to simulate the marine CSEM responses over different anisotropic 3-D models, including an anisotropic reservoir, anisotropic surrounding sediments, and an anisotropic overburden, and performed a detailed analysis of the effect of VTI, HTI, and TTI anisotropy on 3-D marine CSEM response. The modelling results illustrate that marine CSEM responses are considerably influenced by various anisotropic parameters, but to different extents. The research on different types of anisotropy within different media can help us avoid the misinterpretation of marine CSEM data when processing and interpreting of the measured geophysical electromagnetic data.

Resistance Evaluation for the Submerged Glider System using CFD Modelling

The underwater gliders are type of autonomous underwater vehicles, which are in their way to be used instead of traditional propellers or thrusters. They are used as submerged gliders or as a part connected to floating hull which is propelled by it. The underwater system depends on different number of hydrofoils (underwater wings) that allow it to glide forward while descending through the water. This paper concentrates on three critical design details: profile of the wings, multiple flapping wings and their angle of attack (AoA), and a method is focused on the mesh generation to predict calm water resistance for the submerged glider system while considering the profile and the maximum rotating angle of the wings and rudder. Flow around submerged wing system model has been calculated at different Froude numbers ranging from 0.1 to 0.6. The grid generator is established for meshing the computational domain with structured hexahedral and unstructured tetrahedral grid. Simulations are carried out using commercial CFD code ANSYS Fluent 19 to calculate calm water resistance of the submerged glider system with different number of hydrofoils. The results conclude that the cambered plate was chosen to design the rudder at 19o AOA and the hydrofoil NACA0012 can be practically applied in the design of the wings at 15o AOA from the resistance point of view. In addition, this investigation introduces a new application for CFD calculations in estimating the resistance of the submerged glider system.

Numerical simulation and experimental research of fiber motion in vortex spinning

The application of high speed air flow in the textile field represents the frontier development direction in this field. However, the coupling mechanism between air flow and fiber is extremely complicated which greatly limits the development and application of new textile technologies. For this reason, this research study is intended to focus on the common cutting-edge basic scientific problem of the air-flow-fiber coupling mechanism, the vortex spinning technology is taken as a breakthrough point for research. Three-dimensional numerical models of air flow in vortex spinning nozzle and a free-end fiber were established. The hybrid grids including the structured hexahedral grids and the unstructured tetrahedral grids were used to divide the air-flow computational region, the realizable κ- ε model was used to solve the turbulent characteristics of air flow, and the wall function method was used to solve the air flow in the near wall region. The displacement and deformation of a single free-end fiber under the action of the air-flow force was solved by the second-order nonlinear double asymptotic method and interpolation method. The dynamic twist process of a free-end fiber was simulated in three-dimensional space by taking into account the fiber characteristics, such as straightening, bending, and torsion. In addition, a high magnification microscope was used to observe the free-end fiber motion during the spinning process. The results show that the numerical simulation is consistent with the experimental observation, proving that dynamic numerical simulation can solve difficult problems in the unobserved dynamic twisting process, and this deepens the understanding of the spinning mechanism in vortex spinning.

3D joint inversion of airborne ZTEM and ground MT data using the finite element method with unstructured tetrahedral grids

As an airborne electromagnetic method induced by natural sources, the Z-axis tipper electromagnetic (ZTEM) system can primarily recover near-surface shallow structures, due to band-limited frequencies (usually 30–720 Hz) of the airborne survey and high sample rate acquisition along the terrain. In contrast, traditional ground magnetotellurics (MT) allows better recovery of deep structures as the data acquired are typical of large site intervals (usually higher than 1 km) and lower frequencies (usually lower than 400 Hz). High-resolution MT surveys allow for shallow small and deep large anomalies to be adequately interpreted but need large site intervals and broadband frequency range, which are seldom used as they are quite costly and laborious. ZTEM data are tippers that relate local vertical to orthogonal horizontal fields, measured at a reference station on the ground. As the 1D structures produce zero vertical magnetic fields, ZTEM is not sensitive to background resistivity. Thus, in general, ZTEM can only reveal relative resistivities and not real resistivities. A combination of the ZTEM and MT methods can be an effective technique, alleviating the shortcomings of the individual methods. At present, complex underground structures and topography introduce difficulties for data inversion and interpretation, as conventional ZTEM and MT forward modeling are generally used on structured grids with limited accuracy. To effectively recover complex underground structures with topography, we developed a 3D framework for joint MT and ZTEM inversion with unstructured tetrahedral grids. The finite element method is used for the forward problem because of its flexibility with unstructured tetrahedral meshes. The limited-memory quasi-Newton algorithm (L-BFGS) for optimization is used to solve the joint inverse problem, which saves memory and computational time by avoiding the explicit calculation of the Hessian matrix. To validate our joint inversion algorithm, we run numerical experiments on two synthetic models. The first synthetic model uses two conductive anomalous bodies of different sizes and depths. At the same time, a simple quadrangular is used for comparing the inversions with and without topography. In contrast, the second synthetic model represents a realistic topography with two different conductivities and the same depth. Both single-domain and joint inversions of the ZTEM and MT data are carried out for the two synthetic models to demonstrate the complementary advantages of joint inversion, while the second model is also used to test the adaptability of the joint inversion to complex topography. The results demonstrate the effectiveness of the finite element method with unstructured tetrahedral grids and the L-BFGS method for joint MT and ZTEM inversion. In addition, the inversion results prove that joint MT and ZTEM inversion can recover deep structures from the MT data and small near-surface structures from the ZTEM data by alleviating the weaknesses of the individual methods.

Hybrid mesh for magnetotelluric forward modeling based on the finite element method

Unstructured tetrahedral grids have been applied in magnetotelluric (MT) forward modeling using the finite element (FE) method because of their adaptability to complex anomalies. However, high-quality results require an extreme refinement of the near-surface area, which leads to excessive meshes and an increased degree of freedom (DoF) of the governing equation of the finite element system. To reduce the computational cost, we have developed a hybrid mesh based on triangular prisms and tetrahedrons. The required elements in the near-surface area are reduced because the quality of the triangular prism is not limited by the element aspect ratio. The deep area is discretized by tetrahedral elements to ensure the flexibility of the unstructured grids. The superiority of this hybrid mesh has been tested on a layered model, the DTM1 model and terrain relief models. The results show that the modeling efficiency has been improved, especially for high-frequency data. The accuracy of the modeling using the hybrid mesh is significantly higher than that of the tetrahedral mesh with a similar DoF. Usage of the hybrid mesh can be easily adapted to complex geoelectric models with strong terrain fluctuations, which requires less computational cost than using conventional unstructured elements.

3D minimum-structure inversion of controlled-source EM data using unstructured grids

In this study, a minimum-structure inversion procedure for controlled-source electromagnetic data that parameterizes the Earth model in terms of an unstructured tetrahedral grid is presented. To minimize the objective function defining the inversion problem Gauss–Newton iterations are used. Iterative solvers are employed for the solution of the system of equations at each Gauss–Newton step. These solvers do not need the matrix to be explicitly formed to solve the linear system, and hence, these kinds of solvers ensure memory efficiency of the developed algorithm. For the forward-modelling problem, a potentials formulation, in which the electric field is expressed in terms of vector and scalar potentials and the relevant partial differential equations are the Helmholtz and conservation of charge equations, is discretized by the finite-element method. This decomposition results in a well conditioned matrix equation compared to the E-field equation and gives a more natural description of electromagnetic phenomena. The linear system of equations for the forward-modelling problem is solved by a direct solver. The advantage of this type of solver is that the factorization that is generated, which takes up the bulk of the solution time for this system, can be reused when solving for multiple right-hand sides. This is important for the efficient calculation of the matrix–vector products in the above-mentioned Gauss–Newton systems. Two synthetic examples are given to illustrate the capabilities and efficiency of the presented procedure. In the first example, inversion of a dataset generated over a geometric shaped conductive body for a flat Earth’s surface is presented. In the second example, a realistic scenario is given. This example shows that topography and a realistic ore body can be naturally and effectively handled by unstructured grids.

A Study on Magnetotelluric Characteristics of Magmatic Geothermal Systems

The Magnetotelluric (MT) method is a widely used and effective method of exploring geothermal resources because it can reveal geological information at a great depth and is cost effective. In order to further improve the reliability and rationality of MT data interpretation, MT responses for a typical hydrothermal system and a hot dry rock (HDR) and partial melting system are investigated by a finite-element (FE) forward modeling approach based on unstructured tetrahedral grids that can handle with complex-shaped geothermal systems. These two geothermal models, designed by the 3ds Max software, are comprised of a clay cap, a reservoir, and a heat source, and are discretized into tetrahedral elements by TetGen software. The results show that the apparent resistivities at the broadband of frequencies are mainly affected by the shallow low-resistivity clay cap due to its strong shielding effects, and the induction arrows effectively reflect the boundary of the clay cap. The conductive heat-conducting path and heat source (1 S/m) considered in the models of the HDR and partial melting system could cause significant changes in the apparent resistivities and induction arrows at low frequencies. It suggests that in addition to the apparent resistivities, the induction arrows should be taken into consideration in MT data processing and inversion for better lateral resolution when exploring geothermal resources.

Joint Modeling of Wave Phenomena by Applying the Grid-Characteristic Method and the Discontinuous Galerkin Method

The aim of this work is to develop a hybrid computational method that combines the grid-characteristic method on regular structured grids with the discontinuous Galerkin method on unstructured tetrahedral grids. The proposed method makes it possible to describe integration domains with complex-shaped boundaries and contact boundaries and to calculate seismic fields taking into account the topography of the Earth’s surface, while saving computational resources. This modification of the method in the three-dimensional case has been proposed in this paper for the first time. Examples of using the developed method for calculating elastic wave phenomena arising during seismic prospecting are given. By way of testing, a comparison is made with results produced by the grid-characteristic method on curvilinear structured grids. The proposed hybrid method can be used not only for seismic prospecting, but also for calculating wave phenomena in other objects of complex shape.

Production Performance Simulation of the Fractured Horizontal Well Considering Reservoir and Wellbore Coupled Flow in Shale Gas Reservoirs

Unmodified original shale gas reservoirs have been demonstrated to have extremely low permeability, and the long horizontal drilling and hydraulic fracturing techniques are widely used to obtain economic gas production. However, these techniques can result in a nonlinear flow of fluids during transport through the multi-scale pore system and the complex fracture network around the horizontal well, as well as a two-phase gas-water flow during shale production. To solve these issues, in this paper, first, the complex fracture network around the horizontal well is characterized using the discrete main hydraulic fractures combined with composite regions with different properties. Then, the dual porosity media and discrete-fracture model are implemented to describe the nonlinear flow behavior in the multi-scale pore system, and a wellbore pressure drop model is established to characterize the transient two-phase gas–liquid flow in the different sections of the fractured horizontal well. Finally, the coupled reservoir wellbore flow model is fully implicitly solved using the control volume finite element method and unstructured tetrahedral grids. The results of the simulation model of a multi-stage fractured horizontal well in the shale of the Longmaxi Formation will help field engineers deeply understand gas and water production performance and the dynamics of wellbore and reservoir pressure, as well as optimize development schemes.

A parallel adaptive finite-element approach for 3-D realistic controlled-source electromagnetic problems using hierarchical tetrahedral grids

SUMMARY To effectively and efficiently interpret or invert controlled-source electromagnetic (CSEM) data which are recorded in areas with the kind of complex geological environments and arbitrary topography that are typical, 3-D CSEM forward modelling software that can quickly solve large-scale problems, provide accurate electromagnetic responses for complex geo-electrical models and can be easily incorporated into inversion algorithms are required. We have developed a parallel goal-oriented adaptive mesh refinement finite-element approach for frequency-domain 3-D CSEM forward modelling with hierarchical tetrahedral grids that can offer accurate electromagnetic responses for large-scale complex models and that can efficiently serve for inversion. The approach uses the goal-oriented adaptive vector finite element method to solve the total electric field vector equation. The geo-electrical model is discretized by unstructured tetrahedral grids which can deal with complex underground geological models with arbitrary surface topography. Different from previous adaptive finite element software working on unstructured tetrahedral grids, we have utilized a novel mesh refinement technique named the longest edge bisection method to generate hierarchically refined grids. As the refined grids are nested into the coarse grids, the refinement technique can precisely map the electrical parameters of inversion grids onto the forward modelling grids so that the extra numerical errors generated by the inconsistency of electrical parameters between inversion grids and forward modelling grids are eliminated. In addition, we use the parallel domain-decomposition technique to further accelerate the computations, and the flexible generalized minimum residual solver (FGMRES) with an auxiliary Maxwell solver pre-conditioner to solve the final large-scale system of linear equations. In the end, we validate the performance of the proposed scheme using two synthetic models and one realistic model. We demonstrate that accurate electromagnetic fields can be obtained by comparison with the analytic solutions and that the code is highly scalable for large-scale problems with millions or even hundreds of millions of unknowns. For the synthetic 3-D model and the realistic model with complex geometry, our solutions match well with the results calculated by an existing 3-D CSEM forward modelling code. Both synthetic and realistic examples demonstrate that our newly developed code is an effective, efficient forward modelling engine for interpreting CSEM field data acquired in areas of complex geology and topography.

Massively parallel FDTD full-band Monte Carlo simulations of electromagnetic THz pulses in p-doped silicon at cryogenic temperatures

In semiconductor samples with a low density of free charge carriers the current response to an intense electromagnetic pulse is weak and can be treated as a perturbation in the finite-difference-time-domain calculation of the electromagnetic field. Up to first order the electric field that drives the particles does not depend on their density and the particles can be simulated independently. This opens the door for massively parallel Monte Carlo simulations that can be processed efficiently on a computer cluster and a large number of particles can be simulated resulting in a very low noise level. By this full-band Monte Carlo finite-difference-time-domain approach the generation of higher harmonics in p-doped silicon at cryogenic temperatures is investigated. It is found that the higher harmonics are mostly due to the warped valence bands, which are discretized on adaptive unstructured tetrahedral grids in the k-space to avoid numerical artifacts.

Three‐Dimensional Magnetotelluric Inversion for Arbitrarily Anisotropic Earth Using Unstructured Tetrahedral Discretization

AbstractWe propose a new three‐dimensional anisotropic inversion scheme for magnetotelluric (MT) data. In this method, the earth is discretized into unstructured tetrahedral grids that can fit complex structures well, such as the earth topography and coastline. We use a 3 × 3 tensor to describe the anisotropic conductivity in the governing equation for MT inversions that can be simplified to three principal conductivities and three rotation angles. The edge‐based finite‐element (FE) method is adopted to do the forward and adjoint modeling to compute full impedance and sensitivity. To make the inversion flexible in the restoration of anisotropic and isotropic structures tendentiously, we add an additional regularization in the objective function, with its weight controlled by a distribution function. Before solving the inverse problem using L‐BFGS optimization, we add upper‐bound and lower‐bound constraints to the parameters to suppress the solution space. Four synthetic experiments are designed to test the resolvability and resolution for the six anisotropic parameters in the inversion. The results show that our inversion method can well recover the resistivity in the principal x and y directions and the strike angle γ. Finally, we apply our inversion algorithm to the USArray MT data with the aim to reveal anisotropic structure under Oregon State.

A two-grid simulation framework for fast monitoring of fault stability and ground deformation in multiphase geomechanics

We present a simulation framework for fast monitoring of fault stability and ground deformation in multiphase geomechanics and demonstrate its efficacy for a joint CO2 sequestration and enhanced oil recovery case study. Fault stability is estimated by tracking the time evolution of the change in Coulomb Failure Function at critical points on the fault. The staggered solution algorithm for the coupled problem is augmented with a feature that allows for the flow and geomechanics sub-problems to be solved on two different unstructured tetrahedral grids. We demonstrate the accuracy of the two-grid method on the classical Mandel's problem. For the field scale problem, the geomechanics grid extends to the ground surface while the flow grid is truncated at a depth above which no well activity occurs and the flow field can be assumed to remain unperturbed. This framework reduces the computational burden associated with initialization of pressure and stress fields in the overburden, allows for a study of the critical interaction between overburden and faults, allows for fast renditions of flow-induced ground deformation, and allows choosing flow and geomechanics grid resolutions independently to capture disparate length scales of the underlying physics.

Conservative Finite Element Modeling of EEG and MEG on Unstructured Grids.

For interpretation of electroencephalography (EEG) and magnetoencephalography (MEG) data, multiple solutions of the respective forward problems are needed. In this paper, we assess performance of the mixed-hybrid finite element method (MHFEM) applied to EEG and MEG modeling. The method provides an approximate potential and induced currents and results in a system with a positive semi-definite matrix. The system thus can be solved with a variety of standard methods (e.g. the preconditioned conjugate gradient method). The induced currents satisfy discrete charge conservation law making the method conservative. We studied its performance on unstructured tetrahedral grids for a layered spherical head model as well as a realistic head model. We also compared its accuracy versus the conventional nodal finite element method ( P1 FEM). To avoid modeling singular sources, we completed our computations with a subtraction approach; the derived expression for the MEG response different from earlier published and involves integration of finite quantities only. We conclude that although the MHFEM is more computationally demanding than the P1 FEM, its use is justified for EEG and MEG modeling on low-resolution head models where P1 FEM loses accuracy.