Nonlinear Multigrid Methods for Numerical Solution of the Unsaturated Flow Equation in Two Space Dimensions

Picard and Newton iterations are widely used to solve numerically the nonlinear Richards’ equation (RE) governing water flow in unsaturated porous media. When solving RE in two space dimensions, direct methods applied to the linearized problem in the Newton/Picard iterations are inefficient. The numerical solving of RE in 2D with a nonlinear multigrid (MG) method that avoids Picard/Newton iterations is the focus of this work. The numerical approach is based on an implicit, second-order accurate time discretization combined with a second-order accurate finite difference spatial discretization. The test problems simulate infiltration of water in 2D unsaturated soils with hydraulic properties described by Broadbridge–White and van Genuchten–Mualem models. The numerical results show that nonlinear MG deserves to be taken into consideration for numerical solving of RE.     

A Meshless Solution Technique for the Solution of 3D Unsaturated Zone Problems,Based on Local Hermitian Interpolation with Radial Basis Functions

Recent developments in meshless numerical methods have led to algorithms that can be used to solve arbitrarily large problems without the support of a connected mesh, and without the computational cost and numerical ill-conditioning issues usually associated with such solution techniques. This work applies the Local Hermitian Interpolation (LHI) method, based on local interpolation with Radial Basis Functions (RBFs), to the solution of 3D unsaturated porous media problems. The proposed implementation is capable of handling real soil properties, provided either as an analytical function or as a series of pointwise measurements. The technique is implemented with implicit and explicit timestepping, and is validated against two transient Richards’ equation models, of which one has a known analytical solution. In addition, a real-world infiltration problem based on a saturated–unsaturated formulation is modelled, using a realistic variation of soil properties with water-pressure.     

Computing flow-induced stresses of injection molding based on the Phan–Thien–Tanner model

A numerical approach is introduced to solve the viscoelastic flow problem of filling and post-filling in injection molding. The governing equations are in terms of compressible, non-isothermal fluid, and the constitutive equation is based on the Phan–Thien–Tanner model. By introducing some hypotheses according to the characteristics of injection molding, a quasi-Poisson type equation about pressure is derived with part integration. Besides, an analytical form of flow-induced stress is also generalized by using the undermined coefficient method. The conventional Galerkin approach is employed to solve the derived pressure equation, and the ‘upwind’ difference scheme is used to discrete the energy equation. Coupling is achieved between velocity and stress by super relax iteration method. The flow in the test mold is investigated by comparing the numerical results and photoelastic photos for polystyrene, showing flow-induced stresses are closely related to melt temperatures. The filling of a two-cavity box is also studied to investigate the viscoelastic effects on real injection molding.     

A Numerical Model of Transient Unsaturated Flow Under Centrifugation Based on Mass Balance

Numerical modeling of unsaturated flow in porous media under centrifugation is studied. A precise and numerically efficient approximation is presented for the mathematical model, based on Richards’ nonlinear and degenerate equation expressed in terms of effective saturation using the Van Genuchten–Mualem ansatz. The main difference with other methods is the utilization of a nonlocal condition based on mass balance. The method is suitable for determination of soil parameters, including the saturated conductivity, via the solution of an inverse problem in an iterative way. First, the fully saturated sample is centrifugated with a free outflow boundary during some time interval. Next, the output boundary is sealed and the sample is centrifugated for a prescribed time interval, or up to the creation of an equilibrium. Finally, the centrifugation is continued with a free outflow boundary. This procedure can be repeated to increase the information to drive the inverse problem. The application of the present method requires only non-intrusive, cheap measurements: rotational momentum and/or gravitational center of the sample, and optionally, the amount of expelled water.     

Finite Element-Based Characterization of Pore-Scale Geometry and Its Impact on Fluid Flow

We present a finite element (FEM) simulation method for pore geometry fluid flow. Within the pore space, we solve the single-phase Reynold’s lubrication equation—a simplified form of the incompressible Navier–Stokes equation yielding the velocity field in a two-step solution approach. (1) Laplace’s equation is solved with homogeneous boundary conditions and a right-hand source term, (2) pore pressure is computed, and the velocity field obtained for no slip conditions at the grain boundaries. From the computed velocity field, we estimate the effective permeability of porous media samples characterized by section micrographs or micro-CT scans. This two-step process is much simpler than solving the full Navier–Stokes equation and, therefore, provides the opportunity to study pore geometries with hundreds of thousands of pores in a computationally more cost effective manner than solving the full Navier–Stokes’ equation. Given the realistic laminar flow field, dispersion in the medium can also be estimated. Our numerical model is verified with an analytical solution and validated on two 2D micro-CT scans from samples, the permeabilities, and porosities of which were pre-determined in laboratory experiments. Comparisons were also made with published experimental, approximate, and exact permeability data. With the future aim to simulate multiphase flow within the pore space, we also compute the radii and derive capillary pressure from the Young–Laplace’s equation. This permits the determination of model parameters for the classical Brooks–Corey and van-Genuchten models, so that relative permeabilities can be estimated.     

Derivation of a Darcy’s Law for a Porous Medium Composed of Two Solid Phases Saturated by a Single- Phase Fluid: A Homogenization Approach

The objective of this article is to derive a macroscopic Darcy’s law for a fluid-saturated moving porous medium whose matrix is composed of two solid phases which are not in direct contact with each other (weakly coupled solid phases). An example of this composite medium is the case of a solid matrix, unfrozen water, and an ice matrix within the pore space. The macroscopic equations for this type of saturated porous material are obtained using two-space homogenization techniques from microscopic periodic structures. The pore size is assumed to be small compared to the macroscopic scale under consideration. At the microscopic scale the two weakly coupled solids are described by the linear elastic equations, and the fluid by the linearized Navier–Stokes equations with appropriate boundary conditions at the solid–fluid interfaces. The derived Darcy’s law contains three permeability tensors whose properties are analyzed. Also, a formal relation with a previous macroscopic fluid flow equation obtained using a phenomenological approach is given. Moreover, a constructive proof of the existence of the three permeability tensors allows for their explicit computation employing finite elements or analogous numerical procedures.     

Finite volume element method for analysis of unsteady reaction–diffusion problems

A finite volume element method is developed for analyzing unsteady scalar reaction-diffusion problems in two dimensions. The method combines the concepts that are employed in the finite volume and the finite element method together. The finite volume method is used to discretize the unsteady reaction-diffusion equation, while the finite element method is applied to estimate the gradient quantities at cell faces. Robustness and efficiency of the combined method have been evaluated on uniform rectangular grids by using available numerical solutions of the two-dimensional reaction-diffusion problems. The numerical solutions demonstrate that the combined method is stable and can provide accurate solution without spurious oscillation along the high-gradient boundary layers.     

Numerical Solvers for Radiation and Conduction in High Temperature Gas Flows

In this paper, the authors introduce a robust numerical technique for radiation–conduction heat transfer in the high temperature fields of gas turbine combustors. The conduction and radiation effects are analyzed by a differential and an integral equation, respectively. Using discrete ordinates for the angular discretization of the integral equation for the radiation effects and a Galerkin discretization for the heat equation, the authors propose a fast multilevel algorithm to solve the fully discretized problem. The algorithm uses the same mesh hierarchy for both radiation and conduction effects, but with two different smoothing operators. Numerical results are shown for test problems in three space dimensions, and comparisons to other methods are also given.     

A method of determining the piezoelectric modulus of a nonuniformly polarized rod

A method is proposed for determination of the depolarization function of a rod of piezoelectric ceramics for a specified amplitude-frequency characteristic of the current. The problem is reduced to a nonlinear integral equation solved by means of a combination of Tikhonov’s linearization and regularization methods. The uniqueness of the solution is shown and a series of numerical experiments is carried out with the aim to determine the polarization law. Translated from Prikladnaya Mekhanika i Tekhnicheskaya Fizika, Vol. 40, No. 3, pp. 204–210, May–June, 1999.     

Comparison of Three FE-FV Numerical Schemes for Single- and Two-Phase Flow Simulation of Fractured Porous Media

We benchmark a family of hybrid finite element–node-centered finite volume discretization methods (FEFV) for single- and two-phase flow/transport through porous media with discrete fracture representations. Special emphasis is placed on a new method we call DFEFVM in which the mesh is split along fracture–matrix interfaces so that discontinuities in concentration or saturation can evolve rather than being suppressed by nodal averaging of these variables. The main objective is to illustrate differences among three discretization schemes suitable for discrete fracture modeling: (a) FEFVM with volumetric finite elements for both fractures and porous rock matrix, (b) FEFVM with lower dimensional finite elements for fractures and volumetric finite elements for the matrix, and (c) DFEFVM with a mesh that is split along material discontinuities. Fracture discontinuities strongly influence single- and multi-phase fluid flow. Continuum methods, when used to model transport across such interfaces, smear out concentration/saturation. We show that the new DFEFVM addresses this problem producing significantly more accurate results. Sealed and open single fractures as well as a realistic fracture geometry are used to conduct tracer and water-flooding numerical experiments. The benchmarking results also reveal the limitations/mesh refinement requirements of FE node-centered FV hybrid methods. We show that the DFEFVM method produces more accurate results even for much coarser meshes.     

The Analytical Solution of the Boussinesq Equation for Flow Induced by a Step Change of the Water Table Elevation Revisited

The classical problem of flow induced by a sudden change of the piezometic head in a semi-infinite aquifer is re-examined. A new analytical solution is derived, by combining an expression describing the water table elevation upstream, obtained by the Adomian’s decomposition approach, to an existing polynomial expression (Tolikas et al. in Water Resour Res 20:24–28, 1984), adequate for the downstream region; the parameters of both approximations are computed by matching the two solutions at the inflection point of the water table. Although several analytical solutions are available in the literature, we demonstrate that the expression we have developed in this issue is the most accurate for strong or moderate non-linear flows, where the degree of non-linearity is defined as the ratio of the piezometric head elevation at the origin to the initial water table elevation. For this type of flows the perturbation-series solution of Polubarinova-Kochina, characterized by previous studies as the best available analytical solution provides physically unacceptable results, while the analytical solution of Lockington (J Irrig Drain Eng 123:24–27, 1997), used to check the accuracy of numerical schemes, underestimates the penetration distance of the recharging front.     

Numerical stress-strain analysis of a nonthin elastic plate

A numerical solution to elastic-equilibrium problems for nonthin plates is proposed. The solution is obtained by using the curvilinear-mesh method in combination with Vekua’s method. The efficiency (rapid convergence and accuracy) of this approach is demonstrated by solving test problems for thick plates that can also be solved exactly or approximately by other methods. A numerical solution is obtained to the bending problem for orthotropic nonthin plates of constant and varying thickness __________ Translated from Prikladnaya Mekhanika, Vol. 42, No. 3, pp. 119–126, March 2006.     

An unsteady point vortex method for coupled fluid–solid problems

A method is proposed for the study of the two-dimensional coupled motion of a general sharp-edged solid body and a surrounding inviscid flow. The formation of vorticity at the body’s edges is accounted for by the shedding at each corner of point vortices whose intensity is adjusted at each time step to satisfy the regularity condition on the flow at the generating corner. The irreversible nature of vortex shedding is included in the model by requiring the vortices’ intensity to vary monotonically in time. A conservation of linear momentum argument is provided for the equation of motion of these point vortices (Brown–Michael equation). The forces and torques applied on the solid body are computed as explicit functions of the solid body velocity and the vortices’ position and intensity, thereby providing an explicit formulation of the vortex–solid coupled problem as a set of non-linear ordinary differential equations. The example of a falling card in a fluid initially at rest is then studied using this method. The stability of broadside-on fall is analysed and the shedding of vorticity from both plate edges is shown to destabilize this position, consistent with experimental studies and numerical simulations of this problem. The reduced-order representation of the fluid motion in terms of point vortices is used to understand the physical origin of this destabilization.      

A phenomenological and extended continuum approach for modelling non-equilibrium flows

This paper presents a new technique that combines Grad’s 13-moment equations (G13) with a phenomenological approach to rarefied gas flows. This combination and the proposed solution technique capture some important non-equilibrium phenomena that appear in the early continuum-transition flow regime. In contrast to the fully coupled 13-moment equation set, a significant advantage of the present solution technique is that it does not require extra boundary conditions explicitly; Grad’s equations for viscous stress and heat flux are used as constitutive relations for the conservation equations instead of being solved as equations of transport. The relative computational cost of this novel technique is low in comparison to other methods, such as fully coupled solutions involving many moments or discrete methods. In this study, the proposed numerical procedure is tested on a planar Couette flow case, and the results are compared to predictions obtained from the direct simulation Monte Carlo method. This test case highlights the presence of normal viscous stresses and tangential heat fluxes that arise from non-equilibrium phenomena, which cannot be captured by the Navier–Stokes–Fourier constitutive equations or phenomenological modifications.      

Application of a modified Lindstedt–Poincaré method in coupled TDOF systems with quadratic nonlinearity and a constant external excitation

An analytical approach is developed for the nonlinear oscillation of a conservative, two-degree-of-freedom (TDOF) mass-spring system with serial combined linear–nonlinear stiffness excited by a constant external force. The main idea of the proposed approach lies in two categories, the first one is the transformation of two nonlinear differential equations of a two-mass system using suitable intermediate variables into a single nonlinear differential equation. Another is the treatment a quadratic nonlinear oscillator (QNO) by the modified Lindstedt–Poincaré (L-P) method presented recently by the authors. The first-order and second-order analytical approximations for the modified L-P method are established for the QNOs with satisfactory results. After solving the nonlinear differential equation, the displacements of two-mass system can be obtained directly from the governing linear second-order differential equation. Unlike the common perturbation method, the modified L-P method is valid for weak as well as strong nonlinear oscillation systems. On the other hand, the new approach yields simple approximate analytical expressions valid for small as well as large amplitudes of oscillation. In short, this new approach yields extended scope of applicability, simplicity, flexibility in application, and avoidance of complicated numerical integration as compared to the previous approaches such as the perturbation and classical harmonic balance methods. Two examples of nonlinear TDOF mass-spring systems excited by a constant external force are selected and the approximate solutions are verified with the exact solutions derived from the Jacobi elliptic function and also the numerical fourth-order Runge–Kutta solutions.     

Numerical simulation of the interaction of gravity waves with elastic ice floes

The self-consistent motion of a fluid and elastically oscillating plates partially covering the fluid is simulated numerically in the linear approximation. The problem is reduced to the simultaneous solution of the Laplace equation for the fluid and the equation of elastic plate oscillations for the ice. The numerical and analytical solutions, the latter obtained from an integral equation containing the Green’s function, are compared. To solve the problem numerically, the boundary element method for the Laplace equation and the finite element method for the equation describing the elastic plate are proposed. The coefficients of transmission and reflection of surface gravity waves from the floating plates are calculated. It is shown that the solution may be quasi-periodic with characteristics determined by the initial values of the wave and ice-floe parameters. The ice floes may exert a filtering effect on the surface wave spectrum, essentially reducing its most reflectable components. Sankt-Peterburg. Translated from Izvestiya Rossiiskoi Akademii Nauk, Mekhanika Zhidkosti i Gaza, No. 3, pp. 123–131, May–June, 2000.     

A comparison of three adsorption equations and sensitivity study of parameter uncertainty effects on adsorption refrigeration thermal performance estimation

This paper presents isosteric-based adsorption equilibrium tests of three activated carbon samples with methanol as an adsorbate. Experimental data was fitted into Langmuir equation, Freundlich equation and Dubinin-Astakov (D–A) equation, respectively. The fitted adsorption equations were compared in terms of agreement with experimental data. Moreover, equation format’s impacts on calculation of the coefficient of performance (COP) and refrigeration capacity of an adsorption refrigeration system was analyzed. In addition, the sensitivity of each parameter in each adsorption equation format to the estimation of cycle’s COP and refrigeration capacity was investigated. It was found that the D–A equation is the best form for presenting the adsorptive property of a carbon-methanol working pair. The D–A equation is recommended for estimating thermal performance of an adsorption refrigeration system because simulation results obtained using the D–A equation are less sensitive to errors of experimentally determined D–A equation’s parameters.     

A control volume-based discretization of the Reynolds equation for the numerical solution of elastohydrodynamic lubrication problems

This paper deals with the discretization of the one-dimensional Reynolds equation coupled with the film shape equation, that is used for the numerical solution of elastohydrodynamically lubricated contacts. The derivation of the developed discretization formula is based on the control volume approach. To reduce the discretization error caused by the upwind expression of the Couette (velocity) term, non-symmetric control volumes are used for discretization of the Reynolds equation, while for the elasticity equation the standard approach is used. A numerical method for the solution of the pressure and the film thickness profiles of elastohydrodynamically lubricated isothermal line contacts is presented. Results are presented for chosen typical parameters of a highly loaded contact. To show the formula efficiency, the convergence speed of both the presented discretization formula and a chosen comparative discretization formula (A.A. Lubrecht, Ph.D. Thesis, University of Twente, The Netherlands, 1987 and C.H. Venner, Ph.D. Thesis, University of Twente, The Netherlands, 1991) are checked. The results show that the presented formula gives better approximations of film thickness values for a given number of equidistant grid nodes. Moreover, the presented approach is probably suitable for more sophisticated cases, such as transient situations and elliptical contacts. © 1998 John Wiley & Sons, Ltd.     

An unstructured finite‐volume method for coupled models of suspended sediment and bed load transport in shallow‐water flows

The aim of this work is to develop a well‐balanced finite‐volume method for the accurate numerical solution of the equations governing suspended sediment and bed load transport in two‐dimensional shallow‐water flows. The modelling system consists of three coupled model components: (i) the shallow‐water equations for the hydrodynamical model; (ii) a transport equation for the dispersion of suspended sediments; and (iii) an Exner equation for the morphodynamics. These coupled models form a hyperbolic system of conservation laws with source terms. The proposed finite‐volume method consists of a predictor stage for the discretization of gradient terms and a corrector stage for the treatment of source terms. The gradient fluxes are discretized using a modified Roe's scheme using the sign of the Jacobian matrix in the coupled system. A well‐balanced discretization is used for the treatment of source terms. In this paper, we also employ an adaptive procedure in the finite‐volume method by monitoring the concentration of suspended sediments in the computational domain during its transport process. The method uses unstructured meshes and incorporates upwinded numerical fluxes and slope limiters to provide sharp resolution of steep sediment concentrations and bed load gradients that may form in the approximate solutions. Details are given on the implementation of the method, and numerical results are presented for two idealized test cases, which demonstrate the accuracy and robustness of the method and its applicability in predicting dam‐break flows over erodible sediment beds. The method is also applied to a sediment transport problem in the Nador lagoon.Copyright © 2013 John Wiley & Sons, Ltd.     

Buoyancy-driven motion of a two-dimensional bubble or drop through a viscous liquid in the presence of a vertical electric field

The effect of an electric field on the buoyancy-driven motion of a two-dimensional gas bubble rising through a quiescent liquid is studied computationally. The dynamics of the bubble is simulated numerically by tracking the gas–liquid interface when an electrostatic field is generated in the vertical gap of the rectangular enclosure. The two phases of the system are assumed to be perfect dielectrics with constant but different permittivities, and in the absence of impressed charges, there is no free charge in the fluid bulk regions or at the interface. Electric stresses are supported at the bubble interface but absent in the bulk and one of the objectives of our computations is to quantify the effect of these Maxwell stresses on the overall bubble dynamics. The numerical algorithm to solve the free-boundary problem relies on the level-set technique coupled with a finite-volume discretization of the Navier–Stokes equations. The sharp interface is numerically approximated by a finite-thickness transition zone over which the material properties vary smoothly, and surface tension and electric field effects are accounted for by employing a continuous surface force approach. A multi-grid solver is applied to the Poisson equation describing the pressure field and the Laplace equation governing the electric field potential. Computational results are presented that address the combined effects of viscosity, surface tension, and electric fields on the dynamics of the bubble motion as a function of the Reynolds number, gravitational Bond number, electric Bond number, density ratio, and viscosity ratio. It is established through extensive computations that the presence of the electric field can have an important effect on the dynamics. We present results that show a substantial increase in the bubble’s rise velocity in the electrified system as compared with the corresponding non-electrified one. In addition, for the electrified system, the bubble shape deformations and oscillations are smaller, and there is a reduction in the propensity of the bubble to break up through increasingly larger oscillations.