A numerical Eulerian method of moment for solute transport in a nonstationary dual-porosity medium

An Eulerian perturbation approach was applied to develop a method of moment for solute transport in a nonstationary, fractured medium. The conceptualized fractured medium is described through a dual-porosity model. Stochastic governing equations for mean concentration and concentration covariance were analytically derived to the first-order accuracy of log-conductivity variance and solved with a numerical method––a finite difference method. The developed method is called a numerical Eulerian method of moment (NEMM). This method was compared with the stationary transport theory [Water Resour. Res. 36(7) (2000) 1665] for predicting mean concentration and its spatial moments. The comparison indicated that the two methods matched very well in predicting first and second spatial moments. NEMM solutions were also compared with Monte Carlo simulations for solute transport in stationary fractured media. The results of the two methods were consistent for calculating small log conductivity variance. The theory was then used to study effects of various parameters and nonstationarity of the medium on flow and transport processes. Results indicated that medium nonstationarity would significantly influence the solute transport process. The nonstationary transport theory relaxes many assumptions adopted in stationary theories and paves the way for applying the NEMM to many environmental projects, especially in analyzing uncertainty of solute transport.     

Stochastic study of solute transport in a nonstationary medium

A Lagrangian stochastic approach is applied to develop a method of moment for solute transport in a physically and chemically nonstationary medium. Stochastic governing equations for mean solute flux and solute covariance are analytically obtained in the first-order accuracy of log conductivity and/or chemical sorption variances and solved numerically using the finite-difference method. The developed method, the numerical method of moments (NMM), is used to predict radionuclide solute transport processes in the saturated zone below the Yucca Mountain project area. The mean, variance, and upper bound of the radionuclide mass flux through a control plane 5 km downstream of the footprint of the repository are calculated. According to their chemical sorption capacities, the various radionuclear chemicals are grouped as nonreactive, weakly sorbing, and strongly sorbing chemicals. The NMM method is used to study their transport processes and influence factors. To verify the method of moments, a Monte Carlo simulation is conducted for nonreactive chemical transport. Results indicate the results from the two methods are consistent, but the NMM method is computationally more efficient than the Monte Carlo method. This study adds to the ongoing debate in the literature on the effect of heterogeneity on solute transport prediction, especially on prediction uncertainty, by showing that the standard derivation of solute flux is larger than the mean solute flux even when the hydraulic conductivity within each geological layer is mild. This study provides a method that may become an efficient calculation tool for many environmental projects.     

A faster simulation method for the stochastic response of hysteretic structures subject to earthquakes

A semi-analytical forward-difference Monte Carlo simulation procedure is proposed for the determination of the lower order statistical moments and the joint probability density function of the stochastic response of hysteretic non-linear multi-degree-of-freedom structural systems subject to nonstationary gaussian white noise excitation, as an alternative to conventional direct simulation methods. The method generalizes the so-called Ermak-Allen algorithm developed for simulation applications in molecular dynamics to structural hysteretic systems. The proposed simulation procedure rely on an assumption of local gaussianity during each time step. This assumption is tantamount to various linearizations of the equations of motion. The procedure then applies an analytical convolution of the excitation process, hereby reducing the generation of stochastic processes and numerical integration to the generation of random vectors only. Such a treatment offers higher rates of convergence, faster speed and higher accuracy. The procedure has been compared to the direct Monte Carlo simulation procedure, which uses a fourth-order Runge-Kutta scheme with the white noise process approximated by a broad band Ruiz-Penzien broken line process. The considered system was a multi-dimenensional hysteretic shear frame, where the constitutive equation of the hysteretic shear forces are described by a bilinear hysteretic model. The comparisons show that significant savings in computer time and accuracy can be achieved.     

Sensitivity and moment analyses of head in variably saturated regimes

A numerical approach for approximating statistical moments of hydraulic heads of variably saturated flows in multi-dimensional porous media is developed. The approximation relies on a first-order Taylor series expansion of a finite element flow model and an adjoint state numerical method for variably saturated flows to evaluate sensitivities. This approach can be employed to analyze uncertainties associated with predictions of head of steady-state or transient flows in variably saturated porous media, with any type of boundary and initial conditions. Limitations of stochastic analytical methods such as spectral/perturbation approaches and the time-consuming Monte Carlo simulation technique are thus alleviated. An example is given to demonstrate the utility of the approach and to investigate the temporal evolution of head variances in a variably saturated flow regime. Results show that the fluctuation of the water table can have significant impacts on the propagation of the head variance.     

Nonstationary seismic responses of structure with nonlinear stiffness subject to modulated Kanai–Tajimi excitation

A simple structure under earthquake excitation is modeled as a single‐degree‐of‐freedom system with nonlinear stiffness subject to modulated Kanai–Tajimi excitation. The nonstationary responses including the nonstationary probability densities of the system responses and the statistical moments are obtained in semi‐analytical form. By applying the stochastic averaging method based on the generalized harmonic functions, the averaged Fokker–Planck–Kolmogorov(FPK) equation governing the nonstationary probability density of the amplitude is derived. Then, the solution of the FPK equation is approximately expressed by a series expansion in terms of a set of properly selected basis functions with time‐dependent coefficients. According to the Galerkin method, the time‐dependent coefficients are solved from a set of linear first‐order differential equations. Thus, the nonstationary probability densities of the amplitude and the state responses as well as the statistic moments of the amplitude are obtained. Finally, two types of the modulating functions, i.e. constant function and exponential function, are considered to give some semi‐analytical formulae. The proposed procedures are checked against the Monte Carlo simulation. The effects of the structure natural frequency and the intensity of the excitation as well as the ground stiffness on the system responses are discussed. It should be pointed out that the proposed method is good for broadband excitation and light damping. Copyright © 2011 John Wiley & Sons, Ltd.     

GPU加速的水体辐射传输Monte Carlo模拟模型研究

水体辐射传输方程是复杂的微积分方程,只能利用数值方法求解,如Monte Carlo光线追踪法、不变嵌入法、离散坐标法等,其中,Monte Carlo方法是目前解决水体水下光场三维问题的唯一有效方法.根据辐射传输理论,开发了水下光场的Monte Carlo模拟模型,主要包含大气、水-气界面、层化水体和水底边界4个模块.实现了模拟任意太阳角度、不同水体固有光学属性和任意深度条件下,考虑大气、粗糙水面和水底边界的水下光场,能够获取辐亮度、辐照度等辐射量的空间分布.该模型暂不考虑Raman散射、偏振、内部光源的影响.实现了GPU加速水下光场Monte Carlo模拟,并用Mobley等提出的海洋光学标准问题中的问题1~6进行验证.在两种计算环境下,通过对不同边界条件下的CPU、GPU运行时间及加速比的对比,发现GPU计算可以达到几百至上千倍的加速比.     

A comparative study of numerical approaches to risk assessment of contaminant transport

In risk analysis, a complete characterization of the concentration distribution is necessary to determine the probability of exceeding a threshold value. The most popular method for predicting concentration distribution is Monte Carlo simulation, which samples the cumulative distribution function with a large number of repeated operations. In this paper, we first review three most commonly used Monte Carlo (MC) techniques: the standard Monte Carlo, Latin Hypercube sampling, and Quasi Monte Carlo. The performance of these three MC approaches is investigated. We then apply stochastic collocation method (SCM) to risk assessment. Unlike the MC simulations, the SCM does not require a large number of simulations of flow and solute equations. In particular, the sparse grid collocation method and probabilistic collocation method are employed to represent the concentration in terms of polynomials and unknown coefficients. The sparse grid collocation method takes advantage of Lagrange interpolation polynomials while the probabilistic collocation method relies on polynomials chaos expansions. In both methods, the stochastic equations are reduced to a system of decoupled equations, which can be solved with existing solvers and whose results are used to obtain the expansion coefficients. Then the cumulative distribution function is obtained by sampling the approximate polynomials. Our synthetic examples show that among the MC methods, the Quasi Monte Carlo gives the smallest variance for the predicted threshold probability due to its superior convergence property and that the stochastic collocation method is an accurate and efficient alternative to MC simulations.     

An efficient probabilistic finite element method for stochastic groundwater flow

We present an efficient numerical method for solving stochastic porous media flow problems. Single-phase flow with a random conductivity field is considered in a standard first-order perturbation expansion framework. The numerical scheme, based on finite element techniques, is computationally more efficient than traditional approaches because one can work with a much coarser finite element mesh. This is achieved by avoiding the common finite element representation of the conductivity field. Computations with the random conductivity field only arise in integrals of the log conductivity covariance function. The method is demonstrated in several two- and three-dimensional flow situations and compared to analytical solutions and Monte Carlo simulations. Provided that the integrals involving the covariance of the log conductivity are computed by higher-order Gaussian quadrature rules, excellent results can be obtained with characteristic element sizes equal to about five correlation lengths of the log conductivity field. Investigations of the validity of the proposed first-order method are performed by comparing nonlinear Monte Carlo results with linear solutions. In box-shaped domains the log conductivity standard deviation σ_Y may be as large as 1.5, while the head variance is considerably influenced by nonlinear effects as σ_Y approaches unity in more general domains.     

Stochastic analysis of estuarine hydraulics

A new methodology is presented for the solution of the stochastic hydraulic equations characterizing steady, one-dimensional estuarine flow. The methodology is predicated on quasi-linearization, perturbation methods, and the finite difference approximation of the stochastic differential operators. Assuming Manning's roughness coefficient is the principal source of uncertainty in the model, stochastic equations are presented for the water depths and flow rates in the estuarine system. Moment equations are developed for the mean and variance of the water depths. The moment equations are compared with the results of Monte Carlo simulation experiments. The results confirm that for any spatial location in the estuary that (1) as the uncertainty in the channel roughness increases, the uncertainty in mean depth increases, and (2) the predicted mean depth will decrease with increasing uncertainty in Manning'sn. The quasi-analytical approach requires significantly less computer time than Monte Carlo simulations and provides explicit     

Scale effect and calibration of contaminant transport models

The influence on solute transport of the small-scale spatial variation of aquifer hydraulic conductivity (K) was analyzed by comparing results from fine-grid (2 m by 2 m) simulations of a synthetic heterogeneous aquifer to those from coarse-grid (8 m by 4 m) simulations of an equivalent homogeneous aquifer. Realizations of the K field of the heterogeneous aquifer were generated, using the Monte Carlo approach, from a lognormal distribution with mean log K of 2 (K in m/d) and three levels of log K variance of 0.1, 0.5, and 1.0. Numerical simulation results show that the average standard deviation of point concentrations increased from 1.21 to 5.78 when the value of log K variance was increased from 0.1 to 1.0. The average discrepancy between modeled concentrations (obtained from a coarse-grid deterministic numerical simulation) and the actual mean point concentrations (obtained from fine-grid Monte Carlo numerical simulations) increased from 0.91 to 4.23 with the increase in log K variance. The results from this study illustrate the uncertainty in predictions from contaminant transport models due to their inability to simulate the effects of heterogeneities at scales smaller than the model grid.     

Optimal allocation of computational resources in hydrogeological models under uncertainty

Flow and transport models in heterogeneous geological formations are usually large-scale with excessive computational complexity and uncertain characteristics. Uncertainty quantification for predicting subsurface flow and transport often entails utilizing a numerical Monte Carlo framework, which repeatedly simulates the model according to a random field parameter representing hydrogeological characteristics of the aquifer. The physical resolution (e.g. spatial grid resolution) for the simulation is customarily chosen based on recommendations in the literature, independent of the number of Monte Carlo realizations. This practice may lead to either excessive computational burden or inaccurate solutions. We develop an optimization-based methodology that considers the trade-off between the following conflicting objectives: time associated with computational costs, statistical convergence of the model prediction and physical errors corresponding to numerical grid resolution. Computational resources are allocated by considering the overall error based on a joint statistical–numerical analysis and optimizing the error model subject to a given computational constraint. The derived expression for the overall error explicitly takes into account the joint dependence between the discretization error of the physical space and the statistical error associated with Monte Carlo realizations. The performance of the framework is tested against computationally extensive simulations of flow and transport in spatially heterogeneous aquifers. Results show that modelers can achieve optimum physical and statistical resolutions while keeping a minimum error for a given computational time. The physical and statistical resolutions obtained through our analysis yield lower computational costs when compared to the results obtained with prevalent recommendations in the literature. Lastly, we highlight the significance of the geometrical characteristics of the contaminant source zone on the optimum physical and statistical resolutions.     

A numerical method of moments for solute transport in physically and chemically nonstationary formations: linear equilibrium sorption with random K<Subscript>d</Subscript>

A Lagrangian perturbation method is applied to develop a method of moments for solute flux through a three-dimensional nonstationary flow field. The flow nonstationarity stems from medium nonstationarity and internal and external boundaries of the study domain. The solute flux is described as a space-time process where time refers to the solute flux breakthrough through a control plane (CP) at some distance downstream of the solute source and space refers to the transverse displacement distribution at the CP. The analytically derived moment equations for solute transport in a nonstationarity flow field are too complicated to solve analytically, a numerical finite difference method is implemented to obtain the solutions. This approach combines the stochastic model with the flexibility of the numerical method to boundary and initial conditions. The developed method is applied to study the effects of heterogeneity and nonstationarity of the hydraulic conductivity and chemical sorption coefficient on solute transport. The study results indicate all these factors will significantly influence the mean and variance of solute flux.     

Development of fragility curves using high‐dimensional model representation

Fragility curves represent the conditional probability that a structure's response may exceed the performance limit for a given ground motion intensity. Conventional methods for computing building fragilities are either based on statistical extrapolation of detailed analyses on one or two specific buildings or make use of Monte Carlo simulation with these models. However, the Monte Carlo technique usually requires a relatively large number of simulations to obtain a sufficiently reliable estimate of the fragilities, and it is computationally expensive and time consuming to simulate the required thousands of time history analyses. In this paper, high‐dimensional model representation based response surface method together with the Monte Carlo simulation is used to develop the fragility curve, which is then compared with that obtained by using Latin hypercube sampling. It is used to replace the algorithmic performance‐function with an explicit functional relationship, fitting a functional approximation, thereby reducing the number of expensive numerical analyses. After the functional approximation has been made, Monte Carlo simulation is used to obtain the fragility curve of the system. Copyright © 2012 John Wiley & Sons, Ltd.     

Uncertainty estimation of pathlines in ground water models

A method is proposed to estimate the uncertainty of the location of pathlines in two-dimensional, steady-state confined or unconfined flow in aquifers due to the uncertainty of the spatially variable unconditional hydraulic conductivity or transmissivity field. The method is based on concepts of the semianalytical first-order theory given in Stauffer et al. (2002, 2004), which allows estimates of the lateral second moment (variance) of the location of a moving particle. However, this method is reformulated in order to account for nonuniform recharge and nonuniform aquifer thickness. One prominent application is the uncertainty estimation of the catchment of a pumping well by considering the boundary pathlines starting at a stagnation point. In this method, the advective transport of particles is considered, based on the velocity field. In the case of a well catchment, backtracking is applied by using the reversed velocity field. Spatial variability of hydraulic conductivity or transmissivity is considered by taking into account an isotropic exponential covariance function of log-transformed values with parameters describing the variance and correlation length. The method allows postprocessing of results from ground water models with respect to uncertainty estimation. The code PPPath, which was developed for this purpose, provides a postprocessing of pathline computations under PMWIN, which is based on MODFLOW. In order to test the methodology, it was applied to results from Monte Carlo simulations for catchments of pumping wells. The results correspond well. Practical applications illustrate the use of the method in aquifers.     

Numerical method of moments for solute transport in a nonstationary flow field conditioned on hydraulic conductivity and head measurements

The unconditional stochastic studies on groundwater flow and solute transport in a nonstationary conductivity field show that the standard deviations of the hydraulic head and solute flux are very large in comparison with their mean values (Zhang et al. in Water Resour Res 36:2107–2120, 2000; Wu et al. in J Hydrol 275:208–228, 2003; Hu et al. in Adv Water Resour 26:513–531, 2003). In this study, we develop a numerical method of moments conditioning on measurements of hydraulic conductivity and head to reduce the variances of the head and the solute flux. A Lagrangian perturbation method is applied to develop the framework for solute transport in a nonstationary flow field. Since analytically derived moments equations are too complicated to solve analytically, a numerical finite difference method is implemented to obtain the solutions. Instead of using an unconditional conductivity field as an input to calculate groundwater velocity, we combine a geostatistical method and a method of moment for flow to conditionally simulate the distributions of head and velocity based on the measurements of hydraulic conductivity and head at some points. The developed theory is applied in several case studies to investigate the influences of the measurements of hydraulic conductivity and/or the hydraulic head on the variances of the predictive head and the solute flux in nonstationary flow fields. The study results show that the conditional calculation will significantly reduce the head variance. Since the hydraulic head measurement points are treated as the interior boundary (Dirichlet boundary) conditions, conditioning on both the hydraulic conductivity and the head measurements is much better than conditioning only on conductivity measurements for reduction of head variance. However, for solute flux, variance reduction by the conditional study is not so significant.     

Stochastic analysis of velocity spatial variability in bounded rectangular heterogeneous aquifers

Transport of inert solutes in two-dimensional bounded heterogeneous porous media is investigated in a stochastic framework. After adopting a first-order approximation of the flow equations, analytical expressions are derived for the velocity covariances. Effects of the boundary conditions and aquifer size upon the statistical moments are analyzed. While the size of the domain is shown to have small influence on the covariances in most cases, the solutions are considerably modified by the boundaries. The results are compared with analytical solutions on infinite domains, and several discrepancies are demonstrated. For example, while the velocity variances on infinite domains are homogeneous, the present results are strongly non-stationary. Finally, the problem is solved numerically by the Monte Carlo simulation method. The results, including the behavior near the boundaries, are shown to be in close agreement with analytical solutions.     

基于Remez迭代算法求解差分系数的双相介质自适应有限差分正演模拟

利用传统有限差分方法对基于Biot理论的双相介质波动方程进行数值求解时,由于慢纵波的存在,数值频散效应较为明显,影响模拟精度.相对于声学近似方程及普通弹性波方程,Biot双相介质波动方程在同等数值求解算法和精度要求条件下,其地震波场正演模拟需要更多的计算时间.本文针对Biot一阶速度-应力方程组发展了一种变阶数优化有限差分数值模拟方法,旨在同时提高其正演模拟的精度和效率.首先结合交错网格差分格式推导Biot方程的数值频散关系式.然后基于Remez迭代算法求取一阶空间偏导数的优化差分系数,并用于Biot方程的交错网格有限差分数值模拟.在此基础上把三类波的平均频散误差参数限制在给定的频散误差阈值和频率范围内,此时优化有限差分算子的长度就能自适应非均匀双相介质模型中的不同速度区间.数值频散曲线分析表明:基于Remez迭代算法的优化有限差分方法相较传统泰勒级数展开方法在大波数范围对频散误差的压制效果更明显;可变阶数的优化有限差分方法能取得与固定阶数优化有限差分方法相近的模拟精度.在均匀介质和河道模型的数值模拟实验中将本文变阶数优化有限差分算法与传统泰勒展开算法、最小二乘优化算法进行比较,进一步证明其在复杂地下介质中的有效性和适用性.     

多点激励下结构随机地震反应分析的反应谱方法

基于随机振动理论,提出了多点激励作用下线性系统随机地震反应分析的均值反应谱方法,给出了结构峰值反应的均值、标准差以及反应平均频率的反应谱组合公式。这可以将反应谱方法推广应用到多点激励结构的抗震可靠度分析中。鉴于组合公式中谱参数和相关系数需要由烦琐的数值积分得到,本文进一步针对它们给出合理的简化计算式,从而使得建议的反应谱方法的计算效率大大增加。最后,以一个双塔斜拉桥为例,对本文方法进行了验证。基于建议方法的计算结果与Monte Carlo模拟结果吻合较好。与经典的多点激励反应谱方法(MSRS法)比较,本文方法具有其无法比拟的计算效率。     

起伏地表下基于改进BISQ模型双相介质中曲线网格有限差分法波场模拟

本文以基于改进BISQ模型的二维双相各向同性介质一阶速度-应力方程为基础,推导出了曲线坐标系下对应的方程,然后采用低频散、低耗散的同位网格MacCormack有限差分法来离散方程,并采用紧致的单边MacCormack差分格式结合牵引力镜像法来施加自由地表边界条件,实现了地震波场数值模拟.曲线网格有限差分法采用贴体网格来描述自由表面,地表的网格线紧贴地形,避免了台阶近似造成的数值散射.数值模拟结果表明,在双相介质起伏自由地表和分界面处,各类波型复杂的反射透射规律可以清晰展现,曲线网格有限差分法可以精确地解决地震波在含起伏地表的双相各向同性介质中的传播问题.