# Mathematical formulation

In this section we will present the mathematical formulation to represent multiphase flows through porous media in one of the most general ways available (i.e., with the least number of assumptions). As mentioned before, we will be volume-averaging the Navier-Stokes equations to obtain the so-called VARANS equations. For the sake of brevity, not all the steps of the mathematical derivation will be shown here. The reader is referred to Higuera (2015) for the complete details.

## Definitions

Before presenting the governing equations, we need to understand the physical processes involved and define the mathematical operations to average the RANS equations.

When flow passes through a porous medium it loses momentum. Energy is dissipated via viscous effects (boundary layer around all the individual components that comprise the medium) and wakes (laminar or turbulent). The size, shape and arrangement of the solids inside the medium also play an important role. For example, a uniformly distributed and regular pattern of solids may present preferential flow directions, while a random distribution of the same elements can yield higher tortuosity, thus enhancing the dissipative effects. All these factors need to be considered when representing the effects that the porous media exert on the flow.

In the volume-averaging technique, a volume of fixed size (*V*) acts as an averaging filter for the porous media while swept all over the domain of interest. There are three requirements for the volume-averaging volume. First, *V* needs to be large enough to capture the internal structure of the porous media. At the same time, the volume needs to be small enough so that physical properties can be assumed constant within it, as defined in Gray (1975) and Whitaker (1996). Finally, *V* also needs to be small enough to capture the scale of the relevant flow features.

The volume averaging operator (( )) for a given field / inside *V* is called extended or superficial averaging, and is defined as follows:

in which *V _{v}* is the volume of voids, i.e., the volume outside the solid phases (capable of containing fluid), located within the averaging volume, as sketched in Figure 1.

The volume averaging operator defined in Equation 1 includes the variation of *V _{v}* in time and space. While

*V*is fixed, it can be placed at different locations in which the solid geometry enclosed may vary. This means that the extended volume-average of a constant

Figure 1: Sketch of the VARANS domains.

field may present space or time gradients due to the variations of *V _{v}* itself. To correct this behaviour, the intrinsic volume average ({ )^) is defined:

In order to link the extended and intrinsic averages we need to introduce the variable porosity *(ф),* which is defined as the quotient between the volume of voids and the control volume: *ф =* Therefore, the two average types are linked in this way: *(a) = ф (a)^.*

The following variable decomposition will be adopted:

in which / is the real value (obtained by measurement or via Direct Numerical Simulation). The over-bar denotes that / is a Reynolds-averaged variable. (/)/ is the intrinsic volume- averaged value and *f"* is the spatial fluctuation. The spatial fluctuation term represents the small-scale features of the flow, lost when applying volume-averaging. Moreover, *f"* has a spatial mean equal to zero by definition, similar to the turbulent fluctuation in the Reynolds decomposition, which has a zero time mean, namely *=* 0. However, for two arbitrary

fields, generally (*,f"g")^ ф* 0, thus, the spatial fluctuation components will originate terms that cannot be solved, and need to be modelled to close the equations.

For convenience the over-bars will be dropped from now on, understanding that volumeaveraging is applied on the Reynolds-averaged variables.

There are a number of techniques and theorems that are an important part of the derivations. The reader is referred to Whitaker (1967) and Slattery (1967) for the theorems for the local volume average of a gradient and the theorem for the local volume average of a time derivative, and to the appendices in Higuera (2015) for additional materials.

A single porosity *(ф)* only suffices to represent a system with a unique fluid phase and an obstacle. However, for the sake of completeness, we will propose a general framework for multiphase systems and different porous materials. As shown in Figure 2, the fluid phases (immiscible fluids) are denoted by ft. Solids that have a rigid matrix which cannot deform or move are denoted by 7*j,* whereas moving solids are denoted by ft.

A new variable is required to link all the materials together: volume fraction *(e).* This variable is also known as saturation if the material is a fluid. By definition, volume fraction is the unit volume of a given material within a control volume, therefore:

Figure 2: Sketch of the VARANS phases.

Given this definition, porosity can now be expressed as:

Furthermore, it is convenient to separate the materials that will not deform or move from those that will, and derive an independent porosity term for each of them. This will allow us, for example, to distinguish between the effects of the materials forming the armour layers of breakwaters and solids that may go into suspension, as sediment grains. We will call these static porosity and dynamic porosity, respectively, and define them as follows:

Finally, the total porosity can also be posed in terms of static and dynamic porosities, namely: