Anisotropic diffusion in a finite cylinder, with geochemical applications

Atomic diffusion in minerals may not be well represented by solutions to the diffusion equation for a sphere with a single-valued diffusivity, either because they have platy or elongated habits or because the energetics of diffusion is sensitive to crystallographic direction. In many cases, a cylinder having characteristic radial and axial diffusivities is arguably a better model, but rigorous solutions to the anisotropic diffusion equation for a finite cylinder have not been available. Here we develop general analytical solutions that capture both the internal distribution of diffusant as a function of time, C(r, z, t), and the fraction, F, of diffusant lost during a specified thermal history. These solutions are shown to conform with existing analytical expressions for limiting cases of diffusion in a slab or infinite cylinder. We present, in addition, a simple numerical (finite difference) approach that not only reproduces the results of our analytical expressions but also enables us to move beyond some of the limitations of the equations to simulate complex natural scenarios involving non-zero and time-dependent boundary conditions, arbitrary initial distribution of diffusant within the cylinder and simultaneous diffusion and radiogenic ingrowth. The complementary nature of the two approaches is emphasized and several illustrative applications to ‘real-world’ problems are described, including noble-gas thermochronometry and halogen-hydroxyl interdiffusion in apatite.     

A consistent algorithm for quasigeostrophic energy and energy flux analysis

The energy flux vector is familiar to geophysical fluid dynamics, but only in the abstract. Equations for computing the energy balance and the energy flux vector (S) for a quasigeostrophic numerical ocean model are provided in both analytic and finite difference form. The multilayer model includes wind, lateral and bottom friction, and bathymetry. The finite difference approximation to the analytic energy balance is consistent with the finite difference approximation to the quasigeostrophic potential vorticity equation. An application of the algorithms to a simple steady ocean gyre illustrates their use. This choice of the flux vector (from an infinity of possibilities) has several advantages: it is computationally convenient; it is capable of field measurement; it has nice mathematical properties in the limit of small amplitude motion; and it facilitates intuition about the dynamics.     

Using a shifted Grünwald-Letnikov scheme for the Caputo derivative to study anomalous solute transport in porous medium

We propose an extension of the shifted Grünwald-Letnikov method to solve fractional partial differential equations in the Caputo sense with arbitrary fractional order derivative α and with an advective term. The method uses the relation between Caputo and Riemann-Liouville definitions, the shifted Grünwald-Letnikov, and the traditional backward and forward finite difference method. The stability of the method is investigated for the implicit and explicit scheme with homogeneous boundary conditions, and a stability criterion is found for the advective-dispersive equation. An application of the method is used to solve contaminant diffusion and advective-dispersive problems. The numerical solution for the fractional diffusion and fractional advection-dispersion is compared with their respective analytical solutions for different time and space grid refinements. The diffusion simulation exhibited a good fit between the analytical and numerical solutions, with the explicit scheme going from stable to unstable as the time and space refinement changes. The fractional advection-dispersion application produced small deviations from the analytical solution. These deviations, however, are analogous to the numerical dispersions encountered in conventional finite difference solutions of the advection-dispersion equation. The new method is also compared with the traditional L2 method. Notably, an example that involves asymmetrical fractional conditions, a fractional diffusivity that depends on time, and a source term show how the methods compare. Overall, this study assesses the quality and easiness of use of the numerical method.     

Solution of the unsaturated soil moisture equation using repeated transforms

An alternative method of solution for the linearized ‘theta‐based’ form of the Richards equation of unsaturated flow is developed in two spatial dimensions. The Laplace and Fourier transformations are employed to reduce the Richards equation to an ordinary differential equation in terms of a transformed moisture content and the transform variables, s and ξ. Separate analytic solutions to the transformed equation are developed for initial states which are either in equilibrium or dis‐equilibrium. The solutions are assembled into a finite layer formulation satisfying continuity of soil suction, thereby facilitating the analysis of horizontally stratified soil profiles. Solution techniques are outlined for various boundary conditions including prescribed constant moisture content, prescribed constant flux and flux as a function of moisture change. Example solutions are compared with linearized finite element solutions. The agreement is found to be good. An adaptation of the method for treating the quasilinearized Richards equation with variable diffusivity is also described. Comparisons of quasilinear solutions with some earlier semi‐analytical, finite element and finite difference results are also favourable. Copyright © 2001 John Wiley & Sons, Ltd.     

Transient drainage volume characterization and flow simulation in reservoir models using the fast marching method

The Fast Marching Method (FMM) has been applied to characterize the transient drainage volume and to simulate flow as a function of time in porous media using the concept of the “diffusive time of flight” (DTOF). The DTOF (τ) provides a spatial coordinate which reduces the three dimensional pressure diffusivity equation to an equivalent one dimensional formulation. It is obtained from the solution to the Eikonal equation via the FMM. Previous applications of this approach have solved the flow equations numerically or by using an analytic asymptotic approximation. Both solution approaches rely upon three characteristics. (1) Accurate solution for the DTOF irrespective of the degree of heterogeneity within a reservoir model. (2) The approximation of the three dimensional pressure solution in terms of the τ coordinate. (3) Accurate representation and discretization of the drainage volume, on which the asymptotic and numerical flow simulations are based. The second and third of these characteristics are specific to reservoir engineering applications, and provide the focus of this study. Analysis of the drainage volume shows that the near well region requires special treatment, leading to a composite discretization for the drainage volume. This discretization has a direct impact upon the calculation of the well test pressure derivative when the asymptotic approximation is used for pressure transient interpretation. For flow simulation, the discretization directly impacts the calculation of the well index and the intercell transmissibility computed in the τ coordinate, and places additional constraints on the discretization of the drainage volume. These new results are validated by comparison with a commercial finite difference flow simulator, and are shown to be more accurate than earlier computational approaches.

     

Testing the predictive power of seagrass depth limit models

The power of equations predicting seagrass depth limit (Z_c) from light extinction (K _z) was tested on data on seagrass depth limits collected from the literature. The test data set comprised 424 reports of seagrass colonization depth and water transparency, including data for 10 seagrass species. This data set confirmed the strong negative relationship betweenZ _c andK _z. The regression equation in Duarte (1991) overestimated the realized seagrass colonization depths at colonization depths < 5 m, while there was no prediction bias above this threshold. These results indicated that seagrass colonizing turbid waters (K _z 0.27 m^-1) have higher apparent light requirements than those growing in clearer waters. The relationship between seagrass colonization depth and light attenuation shifts at a threshold of light attenuation of 0.27 m^-1, requiring separate equations to predictZ _c for seagrass growing in more turbid waters and clearer waters, and to set targets for seagrass restoration and conservation efforts.     

金红石Zr和锆石Ti含量地质温度计

作为近年来新提出的两种单矿物微量元素温度计(金红石Zr含量温度计和锆石Ti含量温度计),由于其简单实用,一经提出便引起了广泛注意,许多研究者尝试将温度计应用于各种不同类型的岩石中。到目前为止每种温度计都存在几个不同的计算公式、这些公式的适用范围和适用的地质情况目前已有统一认识,但是对于所测定温度的地质意义还存在争议。在对变质岩中金红石Zr含量温度计的应用研究中,一部分研究者发现这个温度计所得到的温度与造岩矿物阳离子配分温度计相吻合,因此可以指示峰期变质温度。然而,在对大别-苏鲁造山带超高压变质岩的研究中发现,金红石Zr含量温度计得到的温度比峰期变质温度明显偏低。通过对比国内外的研究分析,认识到不仅压力、活度、元素扩散、流体作用的参与导致的退变反应可能致使微量元素温度计所记录的温度偏低,而且矿物的不同生长世代或生长介质的不同都可能致使微量元素温度偏低。因此,在应用地质温度计时,要结合样品的岩相学、矿物包裹体和微量元素、U-Pb体系定年等方面予以综合考虑,并对矿物的形成环境和形成世代加以限定,从而为合理解释矿物中微量元素的分配及其记录的温度信息提供有效制约。     

涨落差法原理浅析及其在洪水预报中的应用

从河道洪水运动波方程的有限差分解入手,推导出了涨落差法的计算公式,进而证明该计算公式为河道洪水扩散波方程二阶精度的有限差分解。该洪水预报方法参数少,计算方法简单,有预见期。实际应用表明,涨落差法的特点是实用性强,而且具有较高的洪水预报精度,但同时也存在误差累积的缺点。     

Generalized poroelastic wellbore problem

This paper presents a novel analytical solution to the transient, z‐dependent, and asymmetric problem of an infinite wellbore drilled into a fluid‐saturated porous medium. The formulations are based on Biot's linear theory of poroelasticity, in which the dependency of poroelastic field variables to spatial coordinates as well as time domain is considered in the most general form. This gives flexibility to the solution in cases that cannot be analyzed using the conventional plane strain or symmetric models. One such case is when calculating the stress variations around an inclined wellbore where the far‐field stresses are acting over a finite vertical section. The results of our solution to this case with a three‐dimensional state of far‐field stress are used to analyze the stability of inclined wellbores passing through abnormally stressed formations. The presented solution is capable of finding expressions for fundamental solutions with stress or flow boundary conditions at the wellbore. These solutions are here adopted to analyze the pressure disturbances generated by multiprobe formation tester, a standard wireline device that is designed for downhole fluid sampling as well as estimating the directional permeabilities of subsurface earth formations. A comparison with the conventional solution for the relevant pressure diffusion equation indicates that the poroelastic effect is relatively significant in relation to the transient response of the pore pressure. Further, it is shown that the finite dimensions of sink probe would, to a great extent, contribute to the formation's pore pressure variations at its immediate proximity. Copyright © 2013 John Wiley & Sons, Ltd.     

Response of saturated porous media to cyclic thermal loading

The response of a semi‐infinite saturated porous medium subjected to a harmonic thermal loading on its free face is studied herein. The pressure diffusion equation that governs the fluctuation of the interstitial pressure is established. It allows us to obtain prevalent parameters, i.e. the thermal and fluid mass diffusivities and the coefficient of relative bulk variation. Closed‐form solutions of the maximum fluid pressure P_max and its location x_cr are derived. It is shown that the location x_cr of P_max is localized and depends on the diffusivity ratio and the frequency of the thermal loading while the magnitude of P_max depends on the diffusivity ratio and the thermal amplitude. Master curves for x_cr and P_max versus diffusivity ratio are built. It follows that three regimes can be distinguished: namely, thin spalling, thick spalling or in‐depth cracking and no cracking. Copyright © 2003 John Wiley & Sons, Ltd.     

Mg–Fe interdiffusion rates in wadsleyite and the diffusivity jump at the 410-km discontinuity

Mg–Fe interdiffusion rates have been measured in wadsleyite aggregates at 16.0–17.0 GPa and 1230–1530 °C by the diffusion couple method. Oxygen fugacity was controlled using the NNO buffer, and water contents of wadsleyite were measured by infrared spectroscopy. Measured asymmetric diffusion profiles, analyzed using the Boltzmann–Matano equation, indicate that the diffusion rate increases with increasing iron concentration and decreasing grain size. In the case of wadsleyite containing 50–90 weight ppm H₂O, the Mg–Fe interdiffusion coefficients at compositions of Mg/(Mg + Fe)=0.95 in the coarse-grained region (about 60 m) and 0.90 in the fine-grained region (about 6 m) were determined to be a D_Xmg = 0.95 (m² s^–1)=1.24 × 10^–9 exp[–172 (kJ mol^–1)/RT] and D_Xmg = 0.90 (m² s^–1)=1.77 × 10^–9 exp[–143 (kJ mol^–1)/RT], respectively. Grain-boundary diffusion rates were estimated to be about 4 orders of magnitude faster than the volume diffusion rate. Grain-boundary diffusion dominates when the grain size is less than a few tens of microns. Results for the nominally dry diffusion couple in the present study are roughly consistent with previous studies, taking into account differences in pressure and grain size, although water contents of samples were not clear in previous studies. We observed that the diffusivity is enhanced by about 1 order of magnitude in wadsleyite containing 300–2100 wt. ppm H₂O at 1230 °C, which is almost identical to the enhancement associated with a 300 °C increase in temperature. It is still not conclusive that a jump in diffusivity exists between olivine and wadsleyite because water contents of olivine in previous diffusion studies and effects of water on the olivine diffusivity are uncertain.     

Examination of the Chayes-Kruskal procedure for testing correlations between proportions

The Chayes-Kruskal procedure for testing correlations between proportions uses a linear approximation to the actual closure transformation to provide a null value,p _ij , against which an observed closed correlation coefficient,r _ij , can be tested. It has been suggested that a significant difference betweenp _ij andr _ij would indicate a nonzero covariance relationship between theith andjth open variables. In this paper, the linear approximation to the closure transformation is described in terms of a matrix equation. Examination of the solution set of this equation shows that estimation of, or even the identification of, significant nonzero open correlations is essentially impossible even if the number of variables and the sample size are large. The method of solving the matrix equation is described in the appendix.     

Mg tracer diffusion in aluminosilicate garnets at 750-850° C, 1 atm. and 1300° C, 8.5 GPa

 Tracer diffusion coefficients of Mg in natural aluminosilicate garnets of composition Alm₃₈Pyr₅₀Gr₁₀Sp₂ and Alm₇₃Pyr₂₁Gr₅Sp₁ have been measured at 1 bar, 750-850° C and at 8.5 GPa, 1300° C by chemically depositing a salt layer enriched in ²⁶Mg on the specially prepared surface of a garnet single crystal. Diffusion anneals at 1 atmosphere (101325 Pa) were carried out at a controlled f _{O 2} of ∼10^−17.5 bars maintained by a flowing gas mix of CO-CO₂. Annealing conditions were carefully chosen to avoid decomposition of garnet by redox reactions. High pressure anneals were carried out in a multianvil apparatus. Induced diffusion profiles (0.1–0.6 μm) were measured by an ion-microprobe with SIMS attachment. Diffusion coefficients at 1 atmosphere are in excellent agreement with extrapolation of data from high P-T experiments (Loomis et al. 1985; Chakraborty and Ganguly 1992) and also with the low temperature (750–900° C) dataset of Cygan and Lasaga (1985) if the diffusion coefficients are assumed to be proportional to f _{O 2} ^1/6. Such an f _{O 2} dependence, however, makes this dataset inconsistent with the recent dataset of Schwandt et al. (1995) on garnets of composition (Alm₁₅Pyr₇₂Gr₁₃Sp₀) unless a strong compositional dependence of Mg tracer diffusivity for Mg-rich garnets is invoked. The present experimental results show that such a compositional dependence is weak to non-existent for garnets with >38 mole percent almandine component. It is emphasized that the temperature dependence of diffusion coefficients at constant oxygen fugacities (activation energy ≈54 kcal/mol) are different from that along an oxygen fugacity buffer (activation energy ≈64.5 kcal/mol), as already pointed out by Chakraborty and Ganguly (1991). This distinction is of importance for modelling natural processes. The measurements at low temperatures either eliminate the need for, or greatly reduce the uncertainty of, extrapolation of laboratory data for modelling metamorphic processes. The high pressure results combined with those from Chakraborty and Ganguly (1992) and Loomis et al. (1985) indicate that pressure dependence of Mg tracer diffusivity in garnets is much stronger than that in forsterite (Chakraborty et al. 1994). This difference in pressure dependence of diffusivity may be caused by the difference in compressibility of the coordination polyhedra of Mg between olivines and garnets. Activation volumes of Mg tracer diffusion as high as 8 cm³/mol may be estimated using the present data in combination with earlier results. These data suggest that at a temperature of 1300° C, Mg tracer diffusion rates in garnets will decrease by an order of magnitude for every 100 km depth. The pressure effect will be stronger at lower temperatures. For calculations involving diffusion coefficients of garnets at high pressures (e.g. mantle xenoliths, eclogites) the pressure dependence of diffusivity must be taken into account. Received: 21 December 1994 / Accepted: 22 September 1995     

Interpretation of heterogeneity effects in synchrotron X-ray fluorescence microprobe data

Heterogeneity effects often limit the accuracy of synchrotron X-ray fluorescence microprobe elemental analysis data to ± 30%. The difference in matrix mass absorption at Kα and Kβ fluorescence energies of a particular element can be exploited to yield information on the average depth-position of the element or account for heterogeneity effects. Using this technique, the heterogeneous distribution of Cu in a simple layered sample could be resolved to a 2 × 2 × 10 (x, y, z, where z is the depth coordinate) micrometer scale; a depth-resolution limit was determined for the first transition metal series and several other elements in calcite and iron oxide matrices. For complex heterogeneous systems, determination of average element depth may be computationally limited but the influence of heterogeneity on fluorescence data may still be assessed. We used this method to compare solid-state diffusion with sample heterogeneity across the Ni-serpentine/calcite boundary of a rock from Panoche Creek, California. We previously reported that Ni fluorescence data may indicate solid state diffusion; in fact, sample heterogeneity in the depth dimension can also explain the Ni fluorescence data. Depth heterogeneity in samples can lead to misinterpretation of synchrotron X-ray microprobe results unless care is taken to account for the influence of heterogeneity on fluorescence data.     

Sources of radiation in the early Universe: The equation of radiative transfer and optical distances

We have derived the radiative-transfer equation for a point source with a specified intensity and spectrum, originating in the early Universe between the epochs of annihilation and recombination, at redshifts z_s = 10⁸?10⁴. The direct radiation of the source is separated from the diffuse radiation it produces. Optical distances from the source for Thomson scattering and bremsstrahlung absorption at the maximum of the thermal background radiation are calculated as a function of the redshift z. The distances grow sharply with decreasing z, approaching asymptotic values, the absorption distance increasing more slowly and reaching their limiting values at lower z. For the adopted z values, the optical parameters of the Universe can be described in a flat model with dusty material and radiation, and radiative transfer can be treated in a grey approximation.     

Solution of Kepler’s equation for rectilinear motion

The properties of the solution of the kinematic equation (Kepler’s equation) presented as a series in powers of a function of time are considered for the case of unperturbed, rectilinear elliptical and rectilinear hyperbolic motion. Kepler’s equation for unperturbed, rectilinear elliptical motion has the form E − sinE = z ³/6, where E is the eccentric anomaly determining the position in the orbit, z ³/6 is the mean anomaly, which is proportional to the time measured from an encounter, and its solution can be represented as a series in powers of z. It is established that the coefficients of the series are positive. The asymptotic for the coefficients in the region of convergence of the series is found, which covers the entire orbit. The series continues to converge over the entire boundary of the circle of convergence. The kinematic equation for unperturbed, rectilinear hyperbolic motion has the form sinhH − H = ζ ³/6, where, as before, ζ ³ is proportional to time. The substitution E = iH, z = iζ reduces one equation to the other. The series for the solution in the hyperbolic case differs from the series for the elliptical solution only in its alternating-sign coefficients. However, the region of convergence covers only part of the orbit in the hyperbolic case.     

Modelling Two-Phase Incompressible Flow in Porous Media Using Mixed Hybrid and Discontinuous Finite Elements

In this paper, we present a numerical model for simulating two-phase (oil–water and air–water) incompressible and immiscible flow in porous media. The mathematical model which is based on a fractional flow formulation is formed of two nonlinear partial differential equations: a mean pressure equation and a water saturation equation. These two equations can be solved in a sequential manner. Two numerical methods are used to discretize the equations of the two-phase flow model: mixed hybrid finite elements are used to treat the pressure equation, h-based Richards' equation and the diffusion term in the saturation equation, the advection term in the saturation equation is treated with the discontinuous finite elements. We propose a better way to calculate the nonlinear coefficients contained in our equations on each element of the discretized domain. In heterogeneous porous media, the saturation becomes discontinuous at the interface between two porous media. We show in this paper how to use the capillary pressure–saturation relationship in order to handle the saturation jump in the mixed hybrid finite element method. The two-phase flow simulator is verified against analytical solutions for some flow problems treated by other authors.     

The Individual Thicknesses of Stellar Disks of Galaxies Viewed at Arbitrary Angles

A method for determining the thickness of the stellar disk of a galaxy from a photometric image of the galaxy in the plane of the sky is proposed and justified. The method can be applied to determine the thickness of plane-parallel exponential disks with an arbitrary, radius-independent, luminosity distribution perpendicular to the plane of the disk J(r, z) = exp(-r/h)f(z). A special feature and advantage of the method is that it enables determination of the thicknesses of disks viewed at arbitrary angles to the plane of the sky (but not strictly edge-on or face-on). The key idea of the method is finding the true inclinations of galaxies viewed at arbitrary angles not from their isophotes, but instead from the azimuthal distribution of the exponential parameter h. The difference between the inclination determined in the traditional way using the isophotes and the true inclination enables estimation of the thickness of the disk. The effectiveness of the method for determining the inclinations of plane-parallel disks is confirmed using a sample of model isothermal galactic disks: I(r, z) = I₀ exp(-r/h)sech²(z/z₀). The inclinations of the planes of the model galaxies to the line of sight and the relative thicknesses in the model sample vary arbitrarily, making it possible to determine the limits of applicability of the method: z₀/h < 0.7 and 10° <i < 75°.A sample of 44 piecewise-exponential disks of galaxies of the southern sky clusters is used to illustrate the application of the technique to observational data. Comparing the distribution of inclinations calculated using the new method and the traditional isophote method shows that the new method yields a more uniform distribution of inclinations to the plane of the sky for the sample galaxies. The derived average disk thicknesses and the disk-thickness distributions are consistent with statistical estimates and observational data from the literature for samples of galaxies viewed edge-on

     

High‐order ADE scheme for solving the fluid diffusion equation in non‐uniform grids and its application in coupled hydro‐mechanical simulation

To improve the stability and efficiency of explicit technique, one proposed method is to use an unconditionally stable alternating direction explicit (ADE) scheme. However, the standard ADE scheme is only moderately accurate and restricted to uniform grids. This paper derives a novel high‐order ADE scheme capable of solving the fluid diffusion equation in non‐uniform grids. The new scheme is derived by performing a fourth‐order finite difference approximation to the spatial derivatives of the diffusion equation in non‐uniform grid. The implicit Crank‐Nicolson technique is then applied to the resulting approximation, and the subsequent equation is split into two alternating direction sweeps, giving rise to a new high‐order ADE scheme. Because the new scheme can be potentially applied in coupled hydro‐mechanical (H‐M) simulation, the pore pressure solutions from the new scheme are then sequentially coupled with an existing geomechanical simulator in the computer program Fast Lagrangian Analysis of Continua. This coupling procedure is called the sequentially explicit coupling technique based on the fourth‐order ADE scheme (SEA‐4). Verifications of well‐known consolidation problems showed that the new ADE scheme and SEA‐4 can reduce computer runtime by 46% to 75% to that of Fast Lagrangian Analysis of Continua's basic scheme. At the same time, the techniques still maintained average percentage error of 1.6% to 3.5% for pore pressure and 0.2% to 1.5% for displacement solutions and were still accurate under typical grid non‐uniformities. This result suggests that the new high‐order ADE scheme can provide an efficient explicit technique for solving the flow equation of a coupled H‐M problem, which will be beneficial for large‐scale and long‐term H‐M problems in geoengineering.     

Structure and properties of loaded silica contacts during pressure solution: impedance spectroscopy measurements under hydrothermal conditions

In order to investigate directly the structure and properties of grain boundaries in silicate materials undergoing pressure solution, in situ measurements of these properties are required. We report electrical impedance spectroscopy measurements, performed, under hydrothermal conditions, on individual glass–glass and glass-quartz contacts undergoing pressure solution. Resulting estimates of the average grain boundary diffusivity product ( Z = Dd_\textav C^* Z = D\delta_{\text{av}} C^{*} ) for silica transport and of the average grain boundary fluid film thickness ( d_\textav \delta_{\text{av}} ) fall in the ranges 6.3 ± 1.4 × 10⁻¹⁸ m³ s⁻¹ and 350 ± 210 nm, respectively. However, the average values for Z and d_\textav \delta_{\text{av}} obtained were likely influenced by cracking and irregular dissolution of the dissolving contact surfaces, rather than representing uniformly wetted grain boundary properties. Post-mortem SEM observations indicate that the contact surfaces were internally rough. Taken together, our data support the notion that during pressure solution of quartz, grain boundary diffusion is rapid, and interface processes (dissolution and precipitation) are more likely to be rate-limiting than diffusion.