Analysis of dynamic spherical cavity expansion in undrained modified Cam Clay soil

In order to capture the influence of the cavity expansion velocity, this paper presents a semianalytical solution for dynamic spherical cavity expansion in modified Cam Clay (MCC) soil. The key problem is solving the six coupled partial differential equations (PDEs) of cavity expansion, in which the dynamic term is considered in the stress equilibrium equation. The similarity transformation technique is used to transform the PDEs into ordinary differential equations (ODEs). Subsequently, the numerical method using the function “ODE45” in MATLAB is selected to solve the ODEs, which allows the stress and excess pore pressure around the expanding spherical cavity wall to be obtained. The proposed semianalytical solution for dynamic spherical cavity expansion was validated by comparting the degenerate solution with the published quasistatic solution for the MCC model. Parametric study was then conducted to capture the influence of the cavity wall velocity on the cavity expansion response. The proposed solution has potential application to geotechnical problems such as dynamic pile driving, the dynamic cone penetration test, and so forth.     

Multiscale finite-volume formulation for multiphase flow in porous media: black oil formulation of compressible,three-phase flow with gravity

Most practical reservoir simulation studies are performed using the so-called black oil model, in which the phase behavior is represented using solubilities and formation volume factors. We extend the multiscale finite-volume (MSFV) method to deal with nonlinear immiscible three-phase compressible flow in the presence of gravity and capillary forces (i.e., black oil model). Consistent with the MSFV framework, flow and transport are treated separately and differently using a sequential implicit algorithm. A multiscale operator splitting strategy is used to solve the overall mass balance (i.e., the pressure equation). The black-oil pressure equation, which is nonlinear and parabolic, is decomposed into three parts. The first is a homo geneous elliptic equation, for which the original MSFV method is used to compute the dual basis functions and the coarse-scale transmissibilities. The second equation accounts for gravity and capillary effects; the third equation accounts for mass accumulation and sources/ sinks (wells). With the basis functions of the elliptic part, the coarse-scale operator can be assembled. The gravity/capillary pressure part is made up of an elliptic part and a correction term, which is computed using solutions of gravity-driven local problems. A particular solution represents accumulation and wells. The reconstructed fine-scale pressure is used to compute the fine-scale phase fluxes, which are then used to solve the nonlinear saturation equations. For this purpose, a Schwarz iterative scheme is used on the primal coarse grid. The framework is demonstrated using challenging black-oil examples of nonlinear compressible multiphase flow in strongly heterogeneous formations.     

Natural convection of compressible and incompressible gases in undeformable porous media under cold climate conditions

A numerical model for convective heat and mass transport of compressible or incompressible gas flows with soil-water phase change is presented. In general, the gaseous phase is considered as compressible and the model accounts for adiabatic processes of compression heating and expansion cooling. The inherently compressible gaseous phase may nevertheless be considered as incompressible by adopting the Oberbeck–Boussinesq approximations. The numerical method used to solve the equations that describe natural convection is based on a Galerkin finite element formulation with adaptive mesh refinement and dynamic time step control. As most existing numerical studies have focused on the behavior of incompressible fluids, model substantiation examines the influence of fluid compressibility on two-widely used benchmarks of steady-state convective heat and mass transport. The relative importance of the effect of pressure-compressibility cooling is shown to increase as the thermal gradient approaches the magnitude of the adiabatic gradient. From these results, it may be concluded that pore-air compressibility cannot be neglected in medium to large-sized enclosures at small temperature differentials. After demonstrating its ability to solve fairly complex transient problems, the model is used to further our understanding of the thermal behavior of the toe drain at the LA2-BSU dam in the province of Quebec, Canada.     

A method of alternating characteristics with application to advection-dominated environmental systems

We present a numerical method for solving a class of systems of partial differential equations (PDEs) that arises in modeling environmental processes undergoing advection and biogeochemical reactions. The salient feature of these PDEs is that all partial derivatives appear in linear expressions. As a result, the system can be viewed as a set of ordinary differential equations (ODEs), albeit each one along a different characteristic. The method then consists of alternating between equations and integrating each one step-wise along its own characteristic, thus creating a customized grid on which solutions are computed. Since the solutions of such PDEs are generally smoother along their characteristics, the method offers the potential of using larger time steps while maintaining accuracy and reducing numerical dispersion. The advantages in efficiency and accuracy of the proposed method are demonstrated in two illustrative examples that simulate depth-resolved reactive transport and soil carbon cycling.     

A global method for coupling transport with chemistry in heterogeneous porous media

Modeling reactive transport in porous media, using a local chemical equilibrium assumption, leads to a system of advection–diffusion PDEs coupled with algebraic equations. When solving this coupled system, the algebraic equations have to be solved at each grid point for each chemical species and at each time step. This leads to a coupled non-linear system. In this paper, a global solution approach that enables to keep the software codes for transport and chemistry distinct is proposed. The method applies the Newton–Krylov framework to the formulation for reactive transport used in operator splitting. The method is formulated in terms of total mobile and total fixed concentrations and uses the chemical solver as a black box, as it only requires that one be able to solve chemical equilibrium problems (and compute derivatives) without having to know the solution method. An additional advantage of the Newton–Krylov method is that the Jacobian is only needed as an operator in a Jacobian matrix times vector product. The proposed method is tested on the MoMaS reactive transport benchmark.     

A parallel global-implicit 2-D solver for reactive transport problems in porous media based on a reduction scheme and its application to the MoMaS benchmark problem

In this article, an approach for the efficient numerical solution of multi-species reactive transport problems in porous media is described. The objective of this approach is to reformulate the given system of partial and ordinary differential equations (PDEs, ODEs) and algebraic equations (AEs), describing local equilibrium, in such a way that the couplings and nonlinearities are concentrated in a rather small number of equations, leading to the decoupling of some linear partial differential equations from the nonlinear system. Thus, the system is handled in the spirit of a global implicit approach (one step method) avoiding operator splitting techniques, solved by Newton’s method as the basic algorithmic ingredient. The reduction of the problem size helps to limit the large computational costs of numerical simulations of such problems. If the model contains equilibrium precipitation-dissolution reactions of minerals, then these are considered as complementarity conditions and rewritten as semismooth equations, and the whole nonlinear system is solved by the semismooth Newton method.     

Similarity solution for cavity expansion in thermoplastic soil

This paper presents an analytical solution for cavity expansion in thermoplastic soil considering non‐isothermal conditions. The constitutive relationship of thermoplasticity is described by Laloui's advanced and unified constitutive model for environmental geomechanical thermal effect (ACMEG‐T), which is based on multi‐mechanism plasticity and bounding surface theory. The problem is formulated by incorporating ACMEG‐T into the theoretical framework of cavity expansion, yielding a series of partial differential equations (PDEs). Subsequently, the PDEs are transformed into a system of first‐order ordinary differential equations (ODEs) using a similarity solution technique. Solutions to the response parameters of cavity expansion (stress, excess pore pressure, and displacement) can then be obtained by solving the ODEs numerically using mathematical software. The results suggest that soil temperature has a significant influence on the pressure‐expansion relationships and distributions of stress and excess pore pressure around the cavity wall. The proposed solution quantifies the influence of temperature on cavity expansion for the first time and provides a theoretical framework for predicting thermoplastic soil behavior around the cavity wall. The solution found in this paper can be used as a theoretical tool that can potentially be employed in geotechnical engineering problems, such as thermal cone penetration tests, and nuclear waste disposal problems.     

A discontinuous Galerkin method for two-phase flow in a porous medium enforcing H(div) velocityand continuous capillary pressure

We consider the slightly compressible two-phase flow problem in a porous medium with capillary pressure. The problem is solved using the implicit pressure, explicit saturation (IMPES) method, and the convergence is accelerated with iterative coupling of the equations. We use discontinuous Galerkin to discretize both the pressure and saturation equations. We apply two improvements, which are projecting the flux to the mass conservative H(div)-space and penalizing the jump in capillary pressure in the saturation equation. We also discuss the need and use of slope limiters and the choice of primary variables in discretization. The methods are verified with two- and three-dimensional numerical examples. The results show that the modifications stabilize the method and improve the solution.     

Godunov格式下高精度二维水流-输运耦合模型

针对复杂流体运动中物质输运方程的数值求解面临地形复杂、数值阻尼过大以及数值振荡等难题,建立了Godunov格式下求解二维水流-输运方程的高精度耦合数学模型,提出了集成输运对流项的HLLC (Harten-Lax-van Leer-Contact)型近似黎曼算子,可同时计算水流通量及输运通量,不仅有效模拟了复杂地形上水流运动,而且解决了输运方程中对流项产生的数值阻尼过大和不稳定振荡等难题。采用水深-水位加权重构技术和Minmod限制器,提高了模型处理复杂混合流态的能力,同时结合Hancock预测-校正方法,使模型具有时空二阶精度。算例结果表明,模型精度高、稳定性好,能有效抑制数值阻尼,适合模拟实际复杂流体运动中物质的输运过程,具有较好的推广应用价值。     

Numerical Analysis of Space–Time Finite Element Formulations for Miscible Displacements

A finite element formulation is proposed to approximate a nonlinear system of partial differential equations, composed by an elliptic subsystem for the pressure–velocity and a transport equation (convection–diffusion) for the concentration, which models the incompressible miscible displacement of one fluid by another in a rigid porous media. The pressure is approximated by the classical Galerkin method and the velocity is calculated by a post-processing technique. Then, the concentration is obtained by a Galerkin/least-squares space–time (GLS/ST) finite element method. A numerical analysis is developed for the concentration approximation. Then, stability, convergence and numerical results are presented confirming the a priori error estimates.     

Finite element modelling of transport of organic chemicals through soils

Aquifer contamination by organic chemicals in subsurface flow through soils due to leaking underground storage tanks filled with organic fluids is an important groundwater pollution problem. The problem involves transport of a chemical pollutant through soils via flow of three immiscible fluid phases: namely air, water and an organic fluid. In this paper, assuming the air phase is under constant atmospheric pressure, the flow field is described by two coupled equations for the water and the organic fluid flow taking interphase mass transfer into account. The transport equations for the contaminant in all the three phases are derived and assuming partition equilibrium coefficients, a single convective – dispersive mass transport equation is obtained. A finite element formulation corresponding to the coupled differential equations governing flow and mass transport in the three fluid phase porous medium system with constant air phase pressure is presented. Relevant constitutive relationships for fluid conductivities and saturations as function of fluid pressures lead to non-linear material coefficients in the formulation. A general time-integration scheme and iteration by a modified Picard method to handle the non-linear properties are used to solve the resulting finite element equations. Laboratory tests were conducted on a soil column initially saturated with water and displaced by p-cymene (a benzene-derivative hydrocarbon) under constant pressure. The same experimental procedure is simulated by the finite element programme to observe the numerical model behaviour and compare the results with those obtained in the tests. The numerical data agreed well with the observed outflow data, and thus validating the formulation. A hypothetical field case involving leakage of organic fluid in a buried underground storage tank and the subsequent transport of an organic compound (benzene) is analysed and the nature of the plume spread is discussed.     

Godunov Mixed Methods on Triangular Grids for Advection–Dispersion Equations

A time-splitting approach for advection–dispersion equations is considered. The dispersive and advective fluxes are split into two separate partial differential equations (PDEs), one containing the dispersive term and the other one the advective term. On triangular elements a triangle-based high resolution Finite Volume (FV) scheme for advection is combined with a Mixed Hybrid Finite Element (MHFE) technique to solve dispersion. This approach introduces an error proportional to the time step and the overall scheme is only first order accurate if special care is not taken in the definition of the numerical flux approximation for advection. By incorporating the diffusive effects into the definition of this numerical flux, near second order accuracy (up to a logh factor) can be proved theoretically and validated by numerical experiments in both one- and two-dimensional cases.     

Numerical Analysis of Space–Time Finite Element Formulations for Miscible Displacements

A finite element formulation is proposed to approximate a nonlinear system of partial differential equations, composed by an elliptic subsystem for the pressure–velocity and a transport equation (convection–diffusion) for the concentration, which models the incompressible miscible displacement of one fluid by another in a rigid porous media. The pressure is approximated by the classical Galerkin method and the velocity is calculated by a post-processing technique. Then, the concentration is obtained by a Galerkin/least-squares space–time (GLS/ST) finite element method. A numerical analysis is developed for the concentration approximation. Then, stability, convergence and numerical results are presented confirming the a priori error estimates.     

非饱和渗流模拟中非均匀空间网格的改进方法

Richards方程在非饱和渗流模拟及其他相关领域应用广泛。在数值求解过程中,可以采用有限差分方法进行数值离散并迭代求解,为了获得较可靠的数值解,常规的均匀网格空间步长往往是较小的。在一些不利数值条件下,如入渗于干燥土壤,迭代计算费时甚至精度也不能得到很好改善。因此,文章提出Chebyshev空间网格改进方法,结合有限差分方法对Richards方程进行数值离散以获得线性方程组,并通过经典的Picard迭代方法进行迭代求解线性方程组以得到Richards方程的数值解。通过均质土和分层土2个不利情况下的非饱和渗流算例,又结合模型解析解和软件Hydrus-1D,对比研究了改进网格方法与均匀网格方法获得数值解的精度。结果表明,提出的Chebyshev网格方法相较于传统的均匀网格,可以在较少的节点数下获得较高的数值精度,又具有较小的计算开销,有较好的应用前景。     

Hybrid Gaussian-cubic radial basis functions for scattered data interpolation

Scattered data interpolation schemes using kriging and radial basis functions (RBFs) have the advantage of being meshless and dimensional independent; however, for the datasets having insufficient observations, RBFs have the advantage over geostatistical methods as the latter requires variogram study and statistical expertise. Moreover, RBFs can be used for scattered data interpolation with very good convergence, which makes them desirable for shape function interpolation in meshless methods for numerical solution of partial differential equations. For interpolation of large datasets, however, RBFs in their usual form, lead to solving an ill-conditioned system of equations, for which, a small error in the data can cause a significantly large error in the interpolated solution. In order to reduce this limitation, we propose a hybrid kernel by using the conventional Gaussian and a shape parameter independent cubic kernel. Global particle swarm optimization method has been used to analyze the optimal values of the shape parameter as well as the weight coefficients controlling the Gaussian and the cubic part in the hybridization. Through a series of numerical tests, we demonstrate that such hybridization stabilizes the interpolation scheme by yielding a far superior implementation compared to those obtained by using only the Gaussian or cubic kernels. The proposed kernel maintains the accuracy and stability at small shape parameter as well as relatively large degrees of freedom, which exhibit its potential for scattered data interpolation and intrigues its application in global as well as local meshless methods for numerical solution of PDEs.     

A semi-analytical solution for the transient response of one-dimensional unsaturated single-layer poroviscoelastic media

This paper develops a semi-analytical solution for the transient response of an unsaturated single-layer poroviscoelastic medium with two immiscible fluids by using the Laplace transformation and the state-space method. Using the elastic–viscoelastic correspondence principle, we first introduce the Kelvin–Voigt model into Zienkiewicz’s unsaturated poroelastic model. The vibrational response for unsaturated porous material can be obtained by combining these two models and assuming that the wetting and non-wetting fluids are compressible, the solid skeleton and solid particles are viscoelastic, and the inertial and mechanical couplings are taken into account. The Laplace transformation and state-space method are used to solve the basic equations with the associated initial and boundary conditions, and the analytical solution in the Laplace domain is developed. To evaluate the responses in the time domain, Durbin’s numerical inverse Laplace transform method is used to obtain the semi-analytical solution. There are three compressional waves in porous media with two immiscible fluids. Moreover, to observe the three compressional waves clearly, we assume the two immiscible fluids are water and oil. Finally, several examples are provided to show the validity of the semi-analytical solution and to assess the influences of the viscosity coefficients and dynamic permeability coefficients on the behavior of the three compressional waves.     

基于不规则三角形网格和有限体积法的物理性流域水文模型

传统的栅格离散方式不能很好反映流域水文过程的边界特征,且难以实现流域水文过程的多尺度模拟.采用有限体积法构建了基于不规则三角形网格的物理性水文模型,将物理性描述的偏微分方程组在控制体积内积分得到空间半离散的常微分方程组,保证数值求解中的水量平衡,并可与概念性描述部分水文过程(如截留、填洼等)的常微分方程组更好地耦合;建立了数值求解方案,采用Triangle对计算区域进行离散,并在沁河上游流域进行了验证,结果表明模型具有较高的模拟精度和良好的应用前景.     

Efficient solution algorithms for multiaxial probabilistic elasto‐plastic constitutive simulations of soils

A Fokker‐Planck‐Kolmogorov (FPK) equation approach has recently been developed to probabilistically solve any elastic‐plastic constitutive equation with uncertain material parameters by transforming the nonlinear, stochastic constitutive rate equation into a linear, deterministic partial differential equation (PDE) and thereby simplifying the numerical solution process. For an uniaxial problem, conventional numerical techniques, such as the finite difference or finite element methods, may be used to solve the resulting univariate FPK PDE. However, for a multiaxial problem, an efficient algorithm is necessary for tractability of the numerical solution of the multivariate FPK PDE. In this paper, computationally efficient algorithms, based on a Fourier spectral approach, are presented for solving FPK PDEs in (stress) space and (pseudo) time, having space‐independent but time‐dependent coefficients and both space‐ and time‐dependent coefficients, that commonly arise in probabilistic elasto‐plasticity. The algorithms are illustrated by probabilistically simulating 2 common laboratory constitutive experiments in geotechnical engineering, namely, the unconfined compression test and the unconsolidated undrained triaxial compression test.     

Streamline method for resolving sharp fronts for complex two-phase flow in porous media

In this paper, we present a fast streamline-based numerical method for the two-phase flow equations in high-rate flooding scenarios for incompressible fluids in heterogeneous and anisotropic porous media. A fractional flow formulation is adopted and a discontinuous Galerkin method (DG) is employed to solve the pressure equation. Capillary effects can be neglected in high-rate flooding scenarios. This allows us to present an improved streamline approach in combination with the one-dimensional front tracking method to solve the transport equation. To handle the high computational costs of the DG approximation, domain decomposition is applied combined with an algebraic multigrid preconditioner to solve the linear system. Special care at the interior interfaces is required and the streamline tracer has to include a dynamic communication strategy. The method is validated in various two- and three-dimensional tests, where comparisons of the solutions in terms of approximation of flow front propagation with standard fully implicit finite-volume methods are provided.