Modeling of non-reactive solute transport in fractured clayey till during variable flow rate and time

Fractures and biopores can act as preferential flow paths in clay aquitards and may rapidly transmit contaminants into underlying aquifers. Reliable numerical models for assessment of groundwater contamination from such aquitards are needed for planning, regulatory and remediation purposes. In this investigation, high resolution preferential water-saturated flow and bromide transport data were used to evaluate the suitability of equivalent porous medium (EPM), dual porosity (DP) and discrete fracture/matrix diffusion (DFMD) numerical modeling approaches for assessment of flow and non-reactive solute transport in clayey till. The experimental data were obtained from four large undisturbed soil columns (taken from 1.5 to 3.5 m depth) in which biopores and channels along fractures controlled 96-99% of water-saturated flow. Simulating the transport data with the EPM effective porosity model (FRACTRAN in EPM mode) was not successful because calibrated effective porosity for the same column had to be varied up to 1 order of magnitude in order to simulate solute breakthrough for the applied flow rates between 11 and 49 mm/day. Attempts to simulate the same data with the DP models CXTFIT and MODFLOW/MT3D were also unsuccessful because fitted values for dispersion, mobile zone porosity, and mass transfer coefficient between mobile and immobile zones varied several orders of magnitude for the different flow rates, and because dispersion values were furthermore not physically realistic. Only the DFMD modeling approach (FRACTRAN in DFMD mode) was capable to simulate the observed changes in solute transport behavior during alternating flow rate without changing values of calibrated fracture spacing and fracture aperture to represent the macropores.     

Analysis of the long-term behavior of solute transport with nonlinear equilibrium sorption using breakthrough curves and temporal moments

We analyzed the long-term behavior of breakthrough curves (BTCs) and temporal moments of a solute subjected to Freundlich equilibrium sorption (s = kc(n)). For one-dimensional transport in a homogeneous porous medium, we derived a power-law relation between travel time, tau, and solute displacement, chi, with the exponent being equal to the Freundlich n exponent. The mean solute velocity, derived from the first time moment, was found to change as tau(n-1). For n values larger than 0.66, the second time moment could be related to c chi(2/n), where c is a constant. An approach based on the use of a critical concentration was developed to estimate the presence of the asymptotic regime in the tail of the BTC. This approach was tested successfully using numerical case studies. One-dimensional numerical simulations with varying values of k, n and initial mass were run to verify the closed form analytical expressions for the large time behavior of temporal moments and the tailing part of breakthrough curves. Good agreement between the slope of the tailing part of log-log transformed BTCs and the predicted slope using asymptotic theory was found. Asymptotic theory in general underestimated the magnitude of the concentration in the tail. The quality of the estimated concentrations in the tail improved for small values of the dispersivity. Experimental BTCs of uranin and benazolin were analyzed in combination with sorption/desorption batch experiments using asymptotic theory. A good agreement between the value of n parameter derived from desorption experiment with benazolin and the value of the n parameter derived from the tail of the BTC was found.     

Choice of dispersion coefficients in reactive transport calculations on smoothed fields

Local dispersion dominates the mixing of compounds that are introduced separately into the subsurface and do not partition into any other than the aqueous phase. Thus, reactions between these compounds are controlled by dispersive mixing if they are not limited by kinetics. I quantify longitudinal dispersive mixing by the longitudinal effective dispersion coefficient of a conservative tracer introduced by a point-like injection [Water Resour. Res. 36 (12) (2000) 3591-3604]. In the upscaling of mixing-controlled reactive transport, I apply the mean velocity and the effective dispersion coefficient to the macroscopic transport calculations, whereas the reactive parameters on the macro-scale are identical to those on the local scale. The applicability of the approach is demonstrated for the transport of compounds undergoing a second-order irreversible bimolecular reaction. Ten realizations of a two-dimensional heterogeneous log-conductivity field are considered. Using the effective dispersion parameters, the overall mass balance is met in the ensemble average, whereas solute spreading is underestimated. I assess the lack of spreading by the difference between the expected macrodispersion and effective dispersion coefficients. I extend the approach to simulations on log-conductivity fields obtained by kriging of regularly spaced conductivity measurements. These fields contain the large-scale features of the true fields but do not resolve the small-scale variability. For the calculations on the kriged fields, the corresponding conditional covariance is substituted into the analytical expressions of effective dispersion, yielding a correction effective dispersion coefficient. The comparison between simulations on the fully resolved fields and on the kriged fields indicates that the approach is valid for wide plumes meeting the ergodicity condition. The high variability of mixing on small scales unresolved by kriging, however, leads to severe uncertainty when mixing-controlled reactions are predicted for narrow plumes.     

Variable-density groundwater flow and solute transport in heterogeneous porous media: approaches, resolutions and future challenges

In certain hydrogeological situations, fluid density variations occur because of changes in the solute or colloidal concentration, temperature, and pressure of the groundwater. These include seawater intrusion, high-level radioactive waste disposal, groundwater contamination, and geothermal energy production. When the density of the invading fluid is greater than that of the ambient one, density-driven free convection can lead to transport of heat and solutes over larger spatial scales and significantly shorter time scales than compared with diffusion alone. Beginning with the work of Lord Rayleigh in 1916, thermal and solute instabilities in homogeneous media have been studied in detail for almost a century. Recently, these theoretical and experimental studies have been applied in the study of groundwater phenomena, where the assumptions of homogeneity and isotropy rarely, if ever, apply. The critical role that heterogeneity plays in the onset as well as the growth and/or decay of convective motion is discussed by way of a review of pertinent literature and numerical simulations performed using a variable-density flow and solute transport numerical code. Different styles of heterogeneity are considered and range from continuously "trending" heterogeneity (sinusoidal and stochastic permeability distributions) to discretely fractured geologic media. Results indicate that both the onset of instabilities and their subsequent growth and decay are intimately related to the structure and variance of the permeability field. While disordered heterogeneity tends to dissipate convection through dispersive mixing, an ordered heterogeneity (e.g., sets of vertical fractures) allows instabilities to propagate at modest combinations of fracture aperture and separation distances. Despite a clearer understanding of the processes that control the onset and propagation of instabilities, resultant plume patterns and their migration rates and pathways do not appear amenable to prediction at present. The classical Rayleigh number used to predict the occurrence of instabilities fails, in most cases, when heterogeneous conditions prevail. The incorporation of key characteristics of the heterogeneous permeability field into relevant stability criteria and numerical models remains a challenge for future research.     

Operator-splitting errors in coupled reactive transport codes for transient variably saturated flow and contaminant transport in layered soil profiles

One possible way of integrating subsurface flow and transport processes with (bio)geochemical reactions is to couple by means of an operator-splitting approach two completely separate codes, one for variably-saturated flow and solute transport and one for equilibrium and kinetic biogeochemical reactions. This paper evaluates the accuracy of the operator-splitting approach for multicomponent systems for typical soil environmental problems involving transient atmospheric boundary conditions (precipitation, evapotranspiration) and layered soil profiles. The recently developed HP1 code was used to solve the coupled transport and chemical equations. For steady-state flow conditions, the accuracy was found to be mainly a function of the adopted spatial discretization and to a lesser extent of the temporal discretization. For transient flow situations, the accuracy depended in a complex manner on grid discretization, time stepping and the main flow conditions (infiltration versus evaporation). Whereas a finer grid size reduced the numerical errors during steady-state flow or the main infiltration periods, the errors sometimes slightly increased (generally less than 50%) when a finer grid size was used during periods with a high evapotranspiration demand (leading to high pressure head gradients near the soil surface). This indicates that operator-splitting errors are most significant during periods with high evaporative boundary conditions. The operator-splitting errors could be decreased by constraining the time step using the performance index (the product of the grid Peclet and Courant numbers) during infiltration, or the maximum time step during evapotranspiration. Several test problems were used to provide guidance for optimal spatial and temporal discretization.     

Virus transport in physically and geochemically heterogeneous subsurface porous media

A two-dimensional model for virus transport in physically and geochemically heterogeneous subsurface porous media is presented. The model involves solution of the advection-dispersion equation, which additionally considers virus inactivation in the solution, as well as virus removal at the solid matrix surface due to attachment (deposition), release, and inactivation. Two surface inactivation models for the fate of attached inactive viruses and their subsequent role on virus attachment and release were considered. Geochemical heterogeneity, portrayed as patches of positively charged metal oxyhydroxide coatings on collector grain surfaces, and physical heterogeneity, portrayed as spatial variability of hydraulic conductivity, were incorporated in the model. Both layered and randomly (log-normally) distributed physical and geochemical heterogeneities were considered. The upstream weighted multiple cell balance method was employed to numerically solve the governing equations of groundwater flow and virus transport. Model predictions show that the presence of subsurface layered geochemical and physical heterogeneity results in preferential flow paths and thus significantly affect virus mobility. Random distributions of physical and geochemical heterogeneity have also notable influence on the virus transport behavior. While the solution inactivation rate was found to significantly influence the virus transport behavior, surface inactivation under realistic field conditions has probably a negligible influence on the overall virus transport. It was further demonstrated that large virus release rates result in extended periods of virus breakthrough over significant distances downstream from the injection sites. This behavior suggests that simpler models that account for virus adsorption through a retardation factor may yield a misleading assessment of virus transport in "hydrogeologically sensitive" subsurface environments.     

Evaluating non-equilibrium solute transport in small soil columns

Displacement studies on leaching of bromide and two pesticides (atrazine and isoproturon) were conducted under unsaturated steady state flow conditions in 24 small undisturbed soil columns (5.7 cm in diameter and 10 cm long) each collected from two sites differing in soil structure and organic carbon content in North Germany. There were large and irregular variabilities in the characteristics of both soils, as well as in the shapes of breakthrough curves (BTCs) of different columns, including some with early breakthrough and increased tailing, qualitatively indicating the presence of preferential flow. It was estimated that one preferential flow column (PFC) at site A, and four at site B, contributed, respectively to 11% and 58% of the accumulated leached fraction and to more than 80% of the maximum observed standard deviation (SD) in the field-scale concentration and mass flux of pesticides at two sites. The bromide BTCs of two sites were analyzed with the equilibrium convection-dispersion equation (CDE) and a non-equilibrium two-region/mobile-immobile model. Transport parameters of these models for individual BTCs were determined using a curve fitting program, CXTFIT, and by the time moment method. For the CDE based equilibrium model, the mean values of retardation factor, R, considered separately for all columns, PFCs or non-preferential flow columns (NPFCs) were comparable for the two methods; significant differences were observed in the values of dispersion coefficients of two sites using the two estimation methods. It was inferred from the estimated parameters of non-equilibrium model that 5-12% of water at site A, and 12% at site B, was immobile during displacement in NPFCs. The corresponding values for PFCs of two sites were much larger, ranging from 25% to 51% by CXTFIT and from 24% to 72% by the moment method, suggesting the role of certain mechanisms other than immobile water in higher degrees of non-equilibrium in these columns. Peclet numbers in PFCs of both sites were consistently smaller than five, indicating the inadequacy of the non-equilibrium model to incorporate the effect of all forms of non-equilibrium in PFCs. Overall, the BTCs of individual NPFCs, PFCs and of field average concentration at the two sites were better reproduced with parameters obtained from CXTFIT than by the moment method. The moment method failed to capture the peak concentrations in PFCs, but tended to describe the desorption and tail branches of BTCs better than the curve fitting approach.     

Analytical solution of two-dimensional solute transport in an aquifer–aquitard system

This study deals with two-dimensional solute transport in an aquifer–aquitard system by maintaining rigorous mass conservation at the aquifer–aquitard interface. Advection, longitudinal dispersion, and transverse vertical dispersion are considered in the aquifer. Vertical advection and diffusion are considered in the aquitards. The first-type and the third-type boundary conditions are considered in the aquifer. This study differs from the commonly used averaged approximation (AA) method that treats the mass flux between the aquifer and aquitard as an averaged volumetric source/sink term in the governing equation of transport in the aquifer. Analytical solutions of concentrations in the aquitards and aquifer and mass transported between the aquifer and upper or lower aquitard are obtained in the Laplace domain, and are subsequently inverted numerically to yield results in the real time domain (the Zhan method). The breakthrough curves (BTCs) and distribution profiles in the aquifer obtained in this study are drastically different from those obtained using the AA method. Comparison of the numerical simulation using the model MT3DMS and the Zhan method indicates that the numerical result differs from that of the Zhan method for an asymmetric case when aquitard advections are at the same direction. The AA method overestimates the mass transported into the upper aquitard when an upward advection exists in the upper aquitard. The mass transported between the aquifer and the aquitard is sensitive to the aquitard Peclet number, but less sensitive to the aquitard diffusion coefficient.     

Modeling radium and radon transport through soil and vegetation

A one-dimensional flow and transport model was developed to describe the movement of two fluid phases, gas and water, within a porous medium and the transport of 226Ra and 222Rn within and between these two phases. Included in this model is the vegetative uptake of water and aqueous 226Ra and 222Rn that can be extracted from the soil via the transpiration stream. The mathematical model is formulated through a set of phase balance equations and a set of species balance equations. Mass exchange, sink terms and the dependence of physical properties upon phase composition couple the two sets of equations. Numerical solution of each set, with iteration between the sets, is carried out leading to a set-iterative compositional model. The Petrov-Galerkin finite element approach is used to allow for upstream weighting if required for a given simulation. Mass lumping improves solution convergence and stability behavior. The resulting numerical model was applied to four problems and was found to produce accurate, mass conservative solutions when compared to published experimental and numerical results and theoretical column experiments. Preliminary results suggest that the model can be used as an investigative tool to determine the feasibility of phytoremediating radium and radon-contaminated soil.     

Simple screening models of NAPL dissolution in the subsurface

Three simple screening models of nonaqueous phase liquid (NAPL) dissolution in the subsurface are proposed based on the NAPL mass conservation and the assumption of proportionality between the residual NAPL source zone concentration and the remaining residual NAPL mass. The purpose of the proposed models is to predict the solute concentration in the zone of the residual NAPL as a result of dissolution. The predicted source zone concentration decrease is used to simulate and account for the decrease of dissolution rate with time. The proposed simple NAPL dissolution models enable the pseudo-equilibrium formulation to be used and therefore the numerical simulations for field application problems can be simplified compared to the non-equilibrium counterpart. With proper choice of empirical parameters, the proposed simple screening models can work as well as more complex dissolution rate correlation models, such as that of Imhoff et al. [Water Resour. Res. 30 (1994) 307-320]. It is found that the proposed models are very good for quantifying non-equilibrium dissolution, which is characterized by tailing of breakthrough curves. The models are especially useful for situations of small residual NAPL saturation, which are typical for many field applications.     

Comparison of batch, stirred flow chamber, and column experiments to study adsorption, desorption and transport of carbofuran within two acidic soils

Different methods (batch, column and stirred flow chamber experiments) used for adsorption and desorption of carbofuran studies were compared. All tested methods showed that the carbofuran adsorption was higher in the soil with the higher organic matter content, whereas the opposite behaviour was observed for the percentage of carbofuran desorbed. However, different methods have revealed some discrepancies in carbofuran adsorption/desorption kinetics. Although batch method showed interesting data on equilibrium experiments, such as a low heterogeneity for the carbofuran adsorption sites independent of soil organic matter content, it had some disadvantages for carbofuran adsorption/desorption kinetic studies. The disadvantages were related with the excessive limitations of this method on kinetics, i.e., no difference could be detected between different soils. However, with column and stirred flow chamber methods the carbofuran adsorption/desorption kinetics of different soils could be compared. Moreover, the absolute values of carbofuran adsorption/desorption and its rate were higher in the stirred flow chamber than in the batch and column experiments. Using stirred flow chamber experiments the carbofuran desorption was significantly faster than its adsorption, whereas carbofuran using column experiments they were similar. These discrepancies should be considered when the results obtained only with one method is discussed.     

Numerical methods for improving sensitivity analysis and parameter estimation of virus transport simulated using sorptive-reactive processes

Using one- and two-dimensional homogeneous simulations, this paper addresses challenges associated with sensitivity analysis and parameter estimation for virus transport simulated using sorptive-reactive processes. Head, flow, and conservative- and virus-transport observations are considered. The paper examines the use of (1) observed-value weighting, (2) breakthrough-curve temporal moment observations, and (3) the significance of changes in the transport time-step size. The results suggest that (1) sensitivities using observed-value weighting are more susceptible to numerical solution variability, (2) temporal moments of the breakthrough curve are a more robust measure of sensitivity than individual conservative-transport observations, and (3) the transport-simulation time step size is more important than the inactivation rate in solution and about as important as at least two other parameters, reflecting the ease with which results can be influenced by numerical issues. The approach presented allows more accurate evaluation of the information provided by observations for estimation of parameters and generally improves the potential for reasonable parameter-estimation results.     

A novel two-dimensional model for colloid transport in physically and geochemically heterogeneous porous media

A two-dimensional model for colloid transport in geochemically and physically heterogeneous porous media is presented. The model considers patchwise geochemical heterogeneity, which is suitable to describe the chemical variability of many surficial aquifers with ferric oxyhydroxide-coated porous matrix, as well as spatial variability of hydraulic conductivity, which results in heterogeneous flow field. The model is comprised of a transient fluid flow equation, a transient colloid transport equation, and an equation for the dynamics of colloid deposition and release. Numerical simulations were carried out with the model to investigate the colloid transport behavior in layered and randomly heterogeneous porous media. Results demonstrate that physical and geochemical heterogeneities markedly affect the colloid transport behavior. Layered physical or geochemical heterogeneity can result in distinct preferential flow paths of colloidal particles. Furthermore, the combined effect of layered physical and geochemical heterogeneity may result in enhanced or reduced preferential flow of colloids. Random distribution of physical heterogeneity (hydraulic conductivity) results in a random flow field and an irregularly distributed colloid concentration profile in the porous medium. Contrary to random physical heterogeneity, the effect of random patchwise geochemical heterogeneity on colloid transport behavior is not significant. It is mostly the mean value of geochemical heterogeneity rather than its distribution that governs the colloid transport behavior.     

A multi-directional tracer test in the fractured Chalk aquifer of E. Yorkshire, UK

A multi-borehole radial tracer test has been conducted in the confined Chalk aquifer of E. Yorkshire, UK. Three different tracer dyes were injected into three injection boreholes and a central borehole, 25 m from the injection boreholes, was pumped at 330 m(3)/d for 8 days. The breakthrough curves show that initial breakthrough and peak times were fairly similar for all dyes but that recoveries varied markedly from 9 to 57%. The breakthrough curves show a steep rise to a peak and long tail, typical of dual porosity aquifers. The breakthrough curves were simulated using a 1D dual porosity model. Model input parameters were constrained to acceptable ranges determined from estimations of matrix porosity and diffusion coefficient, fracture spacing, initial breakthrough times and bulk transmissivity of the aquifer. The model gave equivalent hydraulic apertures for fractures in the range 363-384 microm, dispersivities of 1 to 5 m and matrix block sizes of 6 to 9 cm. Modelling suggests that matrix block size is the primary controlling parameter for solute transport in the aquifer, particularly for recovery. The observed breakthrough curves suggest results from single injection-borehole tracer tests in the Chalk may give initial breakthrough and peak times reasonably representative of the aquifer but that recovery is highly variable and sensitive to injection and abstraction borehole location. Consideration of aquifer heterogeneity suggests that high recoveries may be indicative of a high flow pathway adjacent, but not necessarily connected, to the injection and abstraction boreholes whereas low recoveries may indicate more distributed flow through many fractures of similar aperture.     

Virus retention and transport through Al-oxide coated sand columns: effects of ionic strength and composition

Knowledge of the factors that influence the fate and transport of viruses in porous media is very important for accurately determining groundwater vulnerability and for developing protective regulations. In this study, six saturated sand column experiments were performed to examine the effects of a positively charged Al-oxide, which was coated on sand particles, on the retention and transport of viruses (phiX174 and MS-2) in background solutions of different ionic strength and composition. We found that the Al-oxide coating on sand significantly removed viruses during their transport in a phosphate buffered saline (PBS) solution. Mass balance calculations showed that 34% of the input MS-2 was inactivated/irreversibly sorbed on the surface of Al-oxide coated sand whereas 100% of phiX174 was recovered. Results from this study also indicated that higher ionic strength facilitated the transport of both phiX174 and MS-2 through the Al-oxide coated sand. This was attributed to the effect of ion shielding, which at higher ionic strength decreased the electrostatic attraction between the viral particles and the sand surface and consequently decreased virus sorption. Strong effect of the ionic strength indicates that an outer-sphere complexation mechanism was responsible for the virus sorption on the Al-oxide coated sand. Ion composition of the background solutions was also found to be a significant factor in influencing virus retention and transport. Virus transport was enhanced in the presence of phosphate (HPO(4)(2-)) as compared to bicarbonate (HCO(3)(-)), and the effect of HPO(4)(2-) was more significant on MS-2 than on phiX174. The presence of bivalent cations (Ca(2+) and Mg(2+)) increased virus transport because the cations partially screened the negative charges on the viruses therefore decreased the electrostatic attraction between the positively charged sand surface and the negatively charged viruses. Mass recovery data indicated that bivalent cations gave rise to a certain degree of inactivation/irreversibly sorption of phiX174 on the surface of Al-oxide coated sand. On the contrary, the bivalent cations appeared to have protected MS-2 from inactivation/irreversibly sorption. This study provides some insights into the mechanisms responsible for virus retention and transport in porous media.     

Numerical modeling of coupled variably saturated fluid flow and reactive transport with fast and slow chemical reactions

The couplings among chemical reaction rates, advective and diffusive transport in fractured media or soils, and changes in hydraulic properties due to precipitation and dissolution within fractures and in rock matrix are important for both nuclear waste disposal and remediation of contaminated sites. This paper describes the development and application of LEHGC2.0, a mechanistically based numerical model for simulation of coupled fluid flow and reactive chemical transport, including both fast and slow reactions in variably saturated media. Theoretical bases and numerical implementations are summarized, and two example problems are demonstrated. The first example deals with the effect of precipitation/dissolution on fluid flow and matrix diffusion in a two-dimensional fractured media. Because of the precipitation and decreased diffusion of solute from the fracture into the matrix, retardation in the fractured medium is not as large as the case wherein interactions between chemical reactions and transport are not considered. The second example focuses on a complicated but realistic advective-dispersive-reactive transport problem. This example exemplifies the need for innovative numerical algorithms to solve problems involving stiff geochemical reactions.     

Modelling radionuclide transport for time varying flow in a channel network

Water flowrates and flow directions may change over time in the subsurface for a number of reasons. In fractured rocks flow takes place in channels within fractures. Solutes are carried by the advective flow. In addition, solutes may diffuse in and out of stagnant waters in the rock matrix and other stagnant water regions. Sorbing species may sorb on fracture surfaces and on the micropore surfaces in the rock matrix. We present a method by which solute particles can be traced in flowing water undergoing changes in flowrate and direction in a complex channel network where the solutes can also interact with the rock by diffusion in the rock matrix. The novelty of this paper is handling of diffusion in the rock matrix under transient flow conditions. The diffusive processes are stochastic and it is not possible to follow a particle deterministically. The method therefore utilises the properties of a probability distribution function for a tracer moving in a fracture where matrix diffusion is active. The method is incorporated in a three dimensional channel network model. Particle tracking is used to trace out a multitude of flowpaths, each of which consists of a large number of channels within fractures. Along each channel the aperture and velocity as well as the matrix sorption properties can vary. An efficient method is presented whereby a particle can be followed along the variable property flowpath. For stationary flow conditions and a network of channels with advective flow and matrix diffusion, a simple analytical solution for the residence time distribution along each pathway can be used. Only two parameter groups need to be integrated along each path. For transient flow conditions, a time stepping procedure that incorporates a stochastic Monte-Carlo like method to follow the particles along the paths when flow conditions change is used. The method is fast and an example is used for illustrative purposes. It is exemplified by a case where land rises due to glacial rebound. It is shown that the effects of changing flowrates and directions can be considerable and that the diffusive migration in the matrix can have a dominating effect on the results.     

A comparison of numerical solutions of partial differential equations with probabilistic and possibilistic parameters for the quantification of uncertainty in subsurface solute transport

Traditionally, uncertainty in parameters are represented as probabilistic distributions and incorporated into groundwater flow and contaminant transport models. With the advent of newer uncertainty theories, it is now understood that stochastic methods cannot properly represent non random uncertainties. In the groundwater flow and contaminant transport equations, uncertainty in some parameters may be random, whereas those of others may be non random. The objective of this paper is to develop a fuzzy-stochastic partial differential equation (FSPDE) model to simulate conditions where both random and non random uncertainties are involved in groundwater flow and solute transport. Three potential solution techniques namely, (a) transforming a probability distribution to a possibility distribution (Method I) then a FSPDE becomes a fuzzy partial differential equation (FPDE), (b) transforming a possibility distribution to a probability distribution (Method II) and then a FSPDE becomes a stochastic partial differential equation (SPDE), and (c) the combination of Monte Carlo methods and FPDE solution techniques (Method III) are proposed and compared. The effects of these three methods on the predictive results are investigated by using two case studies. The results show that the predictions obtained from Method II is a specific case of that got from Method I. When an exact probabilistic result is needed, Method II is suggested. As the loss or gain of information during a probability–possibility (or vice versa) transformation cannot be quantified, their influences on the predictive results is not known. Thus, Method III should probably be preferred for risk assessments.