On the upstream mobility scheme for two-phase flow in porous media

When neglecting capillarity, two-phase incompressible flow in porous media is modelled as a scalar nonlinear hyperbolic conservation law. A change in the rock type results in a change of the flux function. Discretizing in one-dimensional with a finite volume method, we investigate two numerical fluxes, an extension of the Godunov flux and the upstream mobility flux, the latter being widely used in hydrogeology and petroleum engineering. Then, in the case of a changing rock type, one can give examples when the upstream mobility flux does not give the right answer.Comment: A preprint to be published in Computational Geoscience

Multilevel Monte Carlo Finite Volume Methods for Random Conservation Laws with Discontinuous Flux

We consider conservation laws with discontinuous flux where the initial datum, the flux function, and the discontinuous spatial dependency coefficient are subject to randomness. We establish a notion of random adapted entropy solutions to these equations and prove well-posedness provided that the spatial dependency coefficient is piecewise constant with finitely many discontinuities. In particular, the setting under consideration allows the flux to change across finitely many points in space whose positions are uncertain. We propose a single- and multilevel Monte Carlo method based on a finite volume approximation for each sample. Our analysis includes convergence rate estimates of the resulting Monte Carlo and multilevel Monte Carlo finite volume methods as well as error versus work rates showing that the multilevel variant outperforms the single-level method in terms of efficiency. We present numerical experiments motivated by two-phase reservoir simulations for reservoirs with varying geological properties.Comment: 25 pages, 9 figures, 4 tables, major revision

A Parallel Implementation of the Newton's Method in Solving Steady State Navier-Stokes Equations for Hypersonic Viscous Flows. alpha-GMRES: A New Parallelisable Iterative Solver for Large Sparse Non-Symmetric Linear Systems

The motivation for this thesis is to develop a parallelizable fully implicit numerical Navier-Stokes solver for hypersonic viscous flows. The existence of strong shock waves, thin shear layers and strong flow interactions in hypersonic viscous flows requires the use of a high order high resolution scheme for the discretisation of the Navier-Stokes equations in order to achieve an accurate numerical simulation. However, high order high resolution schemes usually involve a more complicated formulation and thus longer computation time as compared to the simpler central differencing scheme. Therefore, the acceleration of the convergence of high order high resolution schemes becomes an increasingly important issue. For steady state solutions of the Navier-Stokes equations a time dependent approach is usually followed using the unsteady governing equations, which can be discretised in time by an explicit or an implicit method. Using an implicit method, unconditional stability can be achieved and as the time step approaches infinity the method approaches the Newton's method, which is equivalent to directly applying the Newton's method for solving the N-dimensional non-linear algebraic system arising from the spatial discretisation of the steady governing equations in the global flowfield. The quadratic convergence may be achieved by using the Newton's method. However one main drawback of the Newton's method is that it is memory intensive, since the Jacobian matrix of the non-linear algebraic system generally needs to be stored. Therefore it is necessary to use a parallel computing environment in order to tackle substantial problems. In the thesis the hypersonic laminar flow over a sharp cone at high angle of attack provides test cases. The flow is adequately modelled by the steady state locally conical Navier-Stokes (LCNS) equations. A structured grid is used since otherwise there are difficulties in generating the unstructured Jacobian matrix. A conservative cell centred finite volume formulation is used for the spatial discretisation. The schemes used for evaluating the fluxes on the cell boundaries are Osher's flux difference splitting scheme, which has continuous first partial derivatives, together with the third order MUSCL (Monotone Upwind Schemes for Conservation Law) scheme for the convective fluxes and the second order central difference scheme for the diffusive fluxes. In developing the Newton's method a simplified approximate procedure has been proposed for the generation of the numerically approximate Jacobian matrix that speeds up the computation and reduces the extent of cells in which the discretised physical state variables need to be used in generating the matrix element. For solving the large sparse non- symmetric linear system in each Newton's iterative step the ?-GMRES linear solver has been developed, which is a robust and efficient scheme in sequential computation. Since the linear solver is designed for generality it is hoped to apply the method for solving similar large sparse non-symmetric linear systems that may occur in other research areas. Writing code for this linear solver is also found to be easy. The parallel computation assigns the computational task of the global domain to multiple processors. It is based on a new decomposition method for the Nth order Jacobian matrix, in which each processor stores the non-zero elements in a certain number of columns of the matrix. The data is stored without overlap and it provides the main storage of the present algorithm. Corresponding to the matrix decomposition method any N-dimensional vector decomposition can be carried out. From the parallel computation point of view, the new procedure for the generation of the numerically approximate Jacobian matrix decreases the memory required in each processor. The alpha-GMRES linear solver is also parallelizable without any sequential bottle-neck, and has a high parallel efficiency. This linear solver plays a key role in the parallelization of an implicit numerical algorithm. The overall numerical algorithm has been implemented in both sequential and parallel computers using both the sequential algorithm version and its parallel counterpart respectively. Since the parallel numerical algorithm is on the global domain and does not change any solution procedure compared with its sequential counterpart, the convergence and the accuracy are maintained compared with the implementation on a single sequential computer. The computers used are IBM RISC system/6000 320H workstation and a Meiko Computer Surface, composed of T800 transputers

Depth-averaged and 3D Finite Volume numerical models for viscous fluids, with application to the simulation of lava flows

This Ph.D. project was initially born from the motivation to contribute to the depth-averaged and 3D modeling of lava flows. Still, we can frame the work done in a broader and more generalist vision. We developed two models that may be used for generic viscous fluids, and we applied efficient numerical schemes for both cases, as explained in the following. The new solvers simulate free-surface viscous fluids whose temperature changes are due to radiative, convective, and conductive heat exchanges. A temperature-dependent viscoplastic model is used for the final application to lava flows. Both the models behind the solvers were derived from mass, momentum, and energy conservation laws. Still, one was obtained by following the depth-averaged model approach and the other by the 3D model approach. The numerical schemes adopted in both our models belong to the family of finite volume methods, based on the integral form of the conservation laws. This choice of methods family is fundamental because it allows the creation and propagation of discontinuities in the solutions and enforces the conservation properties of the equations. We propose a depth-averaged model for a viscous fluid in an incompressible and laminar regime with an additional transport equation for a scalar quantity varying horizontally and a variable density that depends on such transported quantity. Viscosity and non-constant vertical profiles for the velocity and the transported quantity are assumed, overtaking the classic shallow-water formulation. The classic formulation bases on several assumptions, such as the fact that the vertical pressure distribution is hydrostatic, that the vertical component of the velocity can be neglected, and that the horizontal velocity field can be considered constant with depth because the classic formulation accounts for non-viscous fluids. When the vertical shear is essential, the last assumption is too restrictive, so it must relax, producing a modified momentum equation in which a coefficient, known as the Boussinesq factor, appears in the advective term. The spatial discretization method we employed is a modified version of the central-upwind scheme introduced by Kurganov and Petrova in 2007 for the classical shallow water equations. This method is based on a semi-discretization of the computational domain, is stable, and, being a high-order method, has a low numerical diffusion. For the temporal discretization, we used an implicit-explicit Runge-Kutta technique discussed by Russo in 2005 that permits an implicit treatment of the stiff terms. The whole scheme is proved to preserve the positivity of flow thickness and the stationary steady-states. Several numerical experiments validate the proposed method, show the incidence on the numerical solutions of shape coefficients introduced in the model and present the effects of the viscosity-related parameters on the final emplacement of a lava flow. Our 3D model describes the dynamics of two incompressible, viscous, and immiscible fluids, possibly belonging to different phases. Being interested in the final application of lava flows, we also have an equation for energy that models the thermal exchanges between the fluid and the environment. We implemented this model in OpenFOAM, which employs a segregated strategy and the Finite Volume Methods to solve the equations. The Volume of Fluid (VoF) technique introduced by Hirt and Nichols in 1981 is used to deal with the multiphase dynamics (based on the Interphase Capturing strategy), and hence a new transport equation for the volume fraction of one phase is added. The challenging effort of maintaining an accurate description of the interphase between fluids is solved by using the Multidimensional Universal Limiter for Explicit Solution (MULES) method (described by Marquez Damian in 2013) that implements the Flux-Corrected Transport (FCT) technique introduced by Boris and Book in 1973, proposing a mix of high and low order schemes. The choice of the framework to use for any new numerical code is crucial. Our contribution consists of creating a new solver called interThermalRadConvFoam in the OpenFOAM framework by modifying the already existing solver interFoam (described by Deshpande et al. in 2012). Finally, we compared the results of our simulations with some benchmarks to evaluate the performances of our model

High Resolution Schemes for Conservation Laws With Source Terms.

This memoir is devoted to the study of the numerical treatment of source terms in hyperbolic conservation laws and systems. In particular, we study two types of situations that are particularly delicate from the point of view of their numerical approximation: The case of balance laws, with the shallow water system as the main example, and the case of hyperbolic equations with stiff source terms. In this work, we concentrate on the theoretical foundations of highresolution total variation diminishing (TVD) schemes for homogeneous scalar conservation laws, firmly established. We analyze the properties of a second order, flux-limited version of the Lax-Wendroff scheme which avoids oscillations around discontinuities, while preserving steady states. When applied to homogeneous conservation laws, TVD schemes prevent an increase in the total variation of the numerical solution, hence guaranteeing the absence of numerically generated oscillations. They are successfully implemented in the form of flux-limiters or slope limiters for scalar conservation laws and systems. Our technique is based on a flux limiting procedure applied only to those terms related to the physical flow derivative/Jacobian. We also extend the technique developed by Chiavassa and Donat to hyperbolic conservation laws with source terms and apply the multilevel technique to the shallow water system. With respect to the numerical treatment of stiff source terms, we take the simple model problem considered by LeVeque and Yee. We study the properties of the numerical solution obtained with different numerical techniques. We are able to identify the delay factor, which is responsible for the anomalous speed of propagation of the numerical solution on coarse grids. The delay is due to the introduction of non equilibrium values through numerical dissipation, and can only be controlled by adequately reducing the spatial resolution of the simulation. Explicit schemes suffer from the same numerical pathology, even after reducing the time step so that the stability requirements imposed by the fastest scales are satisfied. We study the behavior of Implicit-Explicit (IMEX) numerical techniques, as a tool to obtain high resolution simulations that incorporate the stiff source term in an implicit, systematic, manner

High-Order Entropy Stable Finite Difference Schemes for Nonlinear Conservation Laws: Finite Domains

Developing stable and robust high-order finite difference schemes requires mathematical formalism and appropriate methods of analysis. In this work, nonlinear entropy stability is used to derive provably stable high-order finite difference methods with formal boundary closures for conservation laws. Particular emphasis is placed on the entropy stability of the compressible Navier-Stokes equations. A newly derived entropy stable weighted essentially non-oscillatory finite difference method is used to simulate problems with shocks and a conservative, entropy stable, narrow-stencil finite difference approach is used to approximate viscous terms