Modelling aerated flows is a complex application of computational fluid dynamics (CFD) since the interfaces between air and water change rapidly. In this work, the simulation of aerated flows with the smoothed particle hydrodynamics (SPH) method is investigated with a focus towards the application in engineering practice. To prove the accuracy of the method, the processes of air entrainment and rising air bubbles are studied. Through monitoring the evolution of the bubble contours it is shown that the novel approach of adding artificial repulsion forces at the interface does not alter the dynamics but stabilizes the flow. Building on these fundamental processes we extend the discussion to practical applications with a special focus on forced aeration. Since the employment of a detailed SPH model to practical problems remains out of bounds due to the high computational demand, we propose a combined experimental and numerical study where experimental bubble characteristics are imposed on the numerical simulation. Based on the data of the conducted bubble column experiment, the computational demand is significantly decreased such that the oxygen consumption due to biokinetic processes can be modelled. The future perspective is to apply SPH to urban water systems, e.g., for simulating detailed processes in wastewater treatment and sewer hydraulics.

NOTATION

  • = bond number (−)

  • c

    = speed of sound ()

  • D

    = diameter of bubble column ()

  • = distance between orifices ()

  • = representative field variable at position of particle i (−)

  • ɛ

    = particle repulsion strengthness parameter (−)

  • g

    = gravitational acceleration

  • h

    = smoothing length

  • H

    = water surface level

  • mi

    = particle mass

  • Pi

    = particle pressure

  • = reference density

  • = particle density

  • = dimensionless kernel distance (−)

  • = aeration rate

  • ri

    = particle position vector in space

  • = bubble radius

  • = Reynolds number (−)

  • = particle sampling distance

  • = elapsed simulation time (s)

  • = reference temperature (°C)

  • vi

    = particle velocity

  • = particle volume

  • = kinematic viscosity

  • = smoothing kernel

Subscripts

  • = reference particle index

  • = neighbouring particle index

  • = water phase subscript

  • = air phase subscript

INTRODUCTION

If air is entrained in water, the flow dynamics are governed by multi-phase gas-liquid fluid. Since both phases are immiscible with different densities, the overall flow becomes compressible. Aerated flows are frequently encountered in hydraulics, most importantly in connection with weirs, chutes and spillways – thus contributing to energy dissipation. In environmental engineering, forced aeration of water by means of air bubble entrainment is essential to provide the oxygen needed for biokinetic processes (see, e.g., Henze et al. 1987). In the field of wastewater treatment particularly, aeration is decisive for both performance and cost of the treatment process (Bischof et al. 1996).

Chanson (2013) provides a broad introduction to the hydraulics of aerated flows, emphasizing its importance to the field of engineering. At a micro-scale the equations of fluid motion can be separately derived for both phases and then coupled with each other. It is crucial to supply an explicit model to account for the evolution of the interphases between air and water. However, due to the high computational effort, such direct numerical simulations can only be obtained in a small number of cases. As an alternative, averaged forms of the conservation equations can be used, resulting in two-phase Navier–Stokes equations. However, this simplified form of the governing equations gives rise to additional constraints and relations which again make the numerical solution challenging (Chanson 2013). In this paper, we introduce the smoothed particle hydrodynamics (SPH) method to study these problems. The advantages and drawbacks of the SPH method for modelling aerated flows of engineering applications are discussed and a novel hybrid solution is proposed.

SPH, which was introduced by Gingold & Monaghan (1977) and Lucy (1977), is a fully Lagrangian meshless computational fluid dynamics (CFD) method initially developed for astrophysical simulations. Since its foundation SPH has been successfully applied to fluid mechanics to simulate, e.g., free surface flows (Monaghan 1994; Gomez-Gesteira et al. 2010), transport phenomena (Tartakovsky et al. 2007) and multi-phase problems (Colagrossi & Landrini 2003; Hu & Adams 2007; Grenier et al. 2009). Owing to its Lagrangian nature, SPH has several advantages in comparison with Eulerian CFD methods and other numerical methods like the particle method of characteristics (see, e.g., Hwang et al. 2013): for example, advection is accounted for exactly and mass, momentum and energy are conserved (Monaghan 2005). With SPH the flow dynamics are modelled by dividing the fluid into a discrete set of particles which are evolved based on a weighted influence of their neighbourhood. By specifying corresponding densities and viscosities for particles of different phases, the method intrinsically supports multi-fluid simulations in the case of comparable phase densities. However, a drawback encountered by Colagrossi & Landrini (2003) for high density ratios (as given for air and water) is that the phase interface becomes unstable. In this work, we develop a novel multi-fluid SPH algorithm resembling standard weakly compressible SPH (Monaghan 2013). The key deviation is to extend the momentum equation by an artificial repulsion force which is applied between particles of different phases. While this procedure is known to cure the instability at the phase interface, it has not previously been applied to aerated flows.

Instead of focussing on a single detailed small-scale study, the purpose of the paper is to provide a practically oriented overview of the capabilities and limitations of the SPH method in modelling aerated flows. We begin by summarizing the present multi-fluid weakly compressible SPH (WCSPH) algorithm used. The accuracy of the method is then demonstrated by modelling the fundamental processes of air entrainment and rising air bubbles. The correct evolution of the bubble contours proves that the artificial repulsion force stabilizes the phase interface while otherwise not affecting the flow dynamics. Building on these fundamental examples the challenges and extensions required for advancing to practical applications are discussed. Referring to forced aeration (e.g., in wastewater treatment), the key issue of modelling continuous air inlets is addressed. We then raise the case study of a bubble column simulation as an example to identify the lack of computational effectiveness as a drawback of the method. Even though improvements like graphics processing unit (GPU) computing enhance the computing time, the detailed modelling of the aeration-induced mixing and the interfacial oxygen transfer remains out of bounds if conventional multi-phase SPH algorithms are employed. As a novel method of resolution we propose a combined experimental and numerical study to account for the homogeneous flow regime of a bubble column. Notwithstanding its limitations, this hybrid approach significantly reduces the computational cost such that the simulation of aerated flows in engineering applications becomes reasonable.

PRINCIPLES

The principle of the SPH method is to simulate the hydrodynamics by dividing continuous matter into a discrete set of points which are referred to as particles. Each particle carries physical parameters like density , pressure P, velocity , volume V and constant mass m. Even though these interpolation points are initially distributed on a rectangular grid, SPH is a meshless Lagrangian method where the observer follows the movement of individual particles instead of retaining a fixed coordinate system. The time evolution of the particles' physical parameters is governed by the continuum equations of fluid dynamics which are numerically solved by interpolating from neighbouring particles with a smoothing kernel . Thereby is the position vector in space and h, which controls the width of the kernel, is referred to as the smoothing length. Required kernel properties are convergence to the -function in the limit of vanishing smoothing length, symmetry and partition of unity if a Shepard filter correction is used.

If subindices denote the value of a field variable at the position vector of particle i such that , then with SPH this value is computed by weighting the influence of neighbouring particles by 
formula
1

This fundamental interpolation formula allows the governing equations to be expressed as a set of discrete differential equations whose derivatives of variables are solely computable by kernel operations. This allows for an efficient numerical solution of the equations in a highly parallel framework.

Governing equations

The governing equations can either be derived by integrating the fundamental interpolation formula by parts (see, e.g., Monaghan 2005) or from the Eckart Lagrangian by following a variational principle (Monaghan 2013). The latter derivation ensures that the SPH method reflects the symmetry properties of the Hamiltonian, e.g., linear and angular momentum are conserved. First, the differential form of the SPH continuity equation describes the rate of change of a reference particle's density by weighting the physical parameters of its neighbourhood 
formula
2
Even though the subindex j includes all particles except for the reference particle, only neighbouring particles contribute due to the finite influence radius of the smoothing kernel. The double subindex abbreviates the difference between vectors such that except for the kernel where .
Second, particle acceleration is governed by the momentum equation which includes a multi-fluid stabilization term proposed by (Monaghan 2013) 
formula
3
The first term approximates the pressure gradient and directly derives from the Lagrangian. The second term is attributed to the standard single fluid viscosity, although the averaging of density and viscosity is altered to support fluid phases with a large difference in viscosities. Conversely, the third term is not physically derivable, but empirically added to ensure stability of the method. Even though WCSPH is fundamentally multi-phase, for large density ratios such as for air and water the discontinuity in the density distribution causes an unstable phase interface which is initially repressed by adding an unphysically high surface tension (Colagrossi & Landrini 2003). Later corrections include the introduction of auxiliary functions, e.g., a particle number density (Hu & Adams 2007) or a particle volume distribution to avoid the discontinuity (Grenier et al. 2009). Whereas these algorithms require an extensive reformulation of the method, the present method is strongly related to standard SPH except for the term in the momentum equation which stabilizes the phase interface by applying the original artificial repulsion force by Monaghan (2013) between particles of different phases 
formula
4

Based on a sensitivity analysis we, however, adapt the repulsion strengthness parameter to since the higher value proposed by Monaghan (2013) restricts the evolution of the bubble contours. Minor deviations do not affect the stability of the method but, e.g., alter the evolution of the contours of a rising bubble. Note that this correction, which is controlled by the phase density difference and the parameter , slightly increases the particles' pressure in the neighbourhood band of the phase interface.

Following the update of velocities, in a Lagrangian method like SPH, particles are advected with the local fluid velocity such that 
formula
5
Third, the set of governing equations is closed by relating density to pressure by an equation of state of the form 
formula
6

To limit the density fluctuations to 1% as required to model an effectively incompressible fluid, the liquid's speed of sound is set to (Morris et al. 1997). For and if the values of other parameters are also replaced, an equivalent equation of state holds for the gaseous phase. Note that for air-water flow, the speed of sound ratio is which results in small time steps. The particles' physical parameters for water as liquid and air as gaseous phase are summarized in Table 1.

Table 1

Physical parameters of the fluids at a reference temperature of T = 20 °C

Phase ρ0 (kg m−3ν (mm2 s−1γ (−) 
Water X (simulation) 998.207 1.004 
Air Y (simulation) 1.204 15.11 1.4 
Phase ρ0 (kg m−3ν (mm2 s−1γ (−) 
Water X (simulation) 998.207 1.004 
Air Y (simulation) 1.204 15.11 1.4 
Using the Wendland kernel of order 4 with the ratio as the weighting function , the set of governing equations is integrated in time by a second-order velocity-Verlet scheme with a dynamically controlled step size (see, e.g., Meister & Rauch 2014) according to a Courant–Friedrichs–Lewy (CFL) condition. The kernel's expression as a function of the dimensionless variable is 
formula
7

Wall boundary conditions and viscous interaction

In contrast to Monaghan (2013), solid boundaries are modelled by a multi-phase formulation of generalized wall boundary conditions which impose a local force balance between fluid and solid wall particles as initially proposed by Adami et al. (2012). For time integration we recommend the second-order velocity–Verlet scheme since for the boundary model used the computationally expensive evaluation of particle accelerations can then be limited to once per time step. Later, in this work, we propose a combined experimental and numerical bubble column model, where bubbles are assigned a solid ellipsoidal shape and experimental rising velocities. In this study, the following shear viscous forces are substituted in the momentum equation 
formula
8
To account for the correct viscous interaction it is important to ensure that the ellipsoidal bubbles carry sufficient layers of air particles such that the influence radius of the kernel is exceeded.

Implementation

The results of this work are produced with the in-house developed SPH code titled SPHASE which abbreviates Smoothed Particle Hydrodynamics Activated Sludge Engine (SPHASE). The code is subject to open source licensing by the end of the research project. The name SPHASE is attributed to the aim of establishing SPH as a numerical CFD engine in the field of wastewater treatment by coupling it to the description of biokinetic processes. Still, the code is very general and allows application to a variety of problems, including aerated flows. The key feature of SPHASE is the real-time visualization of the computation, which is performed in a highly parallel framework. Open multi-processing (OpenMP) is the employed application programing interface that provides multithreading of C + + code on various operating systems including Linux, Mac OS X and Windows. Constructs for workload balancing, thread creation and synchronization are provided by OpenMP on all parallelizable code parts. In SPHASE these parts correspond to the majority of the code which is supported by the highly parallel structure of the SPH method. The particle neighbourhood search, which is the limiting structure for performance, is optimized by performing the radix sorting algorithm provided by the thrust library on a uniform grid which leads to a cost of (Merrill & Grimshaw 2011) with n denoting the number of particles. To fully exploit the advantage of the instant visualization, real-time user interaction is an important feature of this code. Physical properties of particles are accessible at all times and inflow rates are adjustable during the simulation. For visualization the cross-platform application QT is employed. The code further complies to the H5Part format standard for exporting simulation data. Finally, SPHASE can be independently used as a post-processing tool for SPH data available in this format. Throughout this work we reference simulation runtimes for an Intel(R) Xeon(R) CPU E5-2695 v2 processor @ 2.4 GHz × 24 with a cache size of 30 MB and 32 GB of DDR3-1600 RAM. The operating system is Arch Linux 64-bit.

Fundamental processes

Entrained air bubbles and pockets are an important feature of aerated flows and occur if the shear forces at the surface exceed the resisting surface tension. In fluid dynamics, the modelling of entrained air is a key issue since it significantly affects the flow due to its intense interaction with the liquid. In nature, entrained air enhances the interfacial oxygen transfer which is required for the natural self-purification and mixing of a water body. Even though aerated flows are subject to extensive studies, detailed CFD simulations are recommended to advance future understanding (Chanson 2013). However, due to the complex highly dynamic processes involved the modelling of entrained air is challenging with Eulerian CFD methods, whereas the present SPH algorithm is a powerful method to capture these multi-fluid interactions. This is demonstrated by two examples: first, the modelling of entrained air is illustrated by the breaking of a wave. We then study the rising bubble problem to prove that the present method allows for an accurate simulation of the bubble contours and the coalescence with the free surface. This test case is particularly suitable for demonstrating that the artificial repulsion force successfully stabilizes but otherwise does not affect the flow dynamics.

Breaking wave

The standard example to illustrate the power of SPH in simulating highly dynamic processes is the dam break. Since this test case is thoroughly discussed in the literature (e.g., Dumbser 2011; Kao & Chang 2012) but not addressed using the simple multi-fluid SPH algorithm (Monaghan 2013), the fundamental simulation characteristics are not repeated herein. It is, however, important to demonstrate that the artificial repulsion force cures the numerical instability observed at high density ratios as with air and water (Colagrossi & Landrini 2003) such that the multi-fluid characteristics, e.g., the entrainment of air and the fluid's interaction with the surrounding air, are resolved in detail.

Every 25 time steps a Shepard filter correction is applied to ensure zeroth-order consistency of the kernel interpolation (see, e.g., Gomez-Gesteira et al. 2010). The geometric configurations are chosen as by Colagrossi & Landrini (2003), i.e., the water block's ratio of length to height is 2, whereas the height and width of the simulation domain, which is filled up with air particles, is 3 and 5.37 times the height of the water block. In the multi-phase simulation, the resolution is varied in between 2.104 and 2.105 particles to check convergence. Figure 1 illustrates the pressure and velocity distribution of the water phase for the highest resolution with a sampling distance of in both phases. As depicted in Figure 2, at the time instant after the reflection of the wave at the wall the computed phase interface is in good agreement with the inviscid boundary element method solution by Colagrossi & Landrini (2003). Because a significant amount of air is entrained by the breaking wave, a multi-phase model is required to capture the physics; otherwise the entrained void region instantly collapses. In addition, we observe a more uniform distribution of fluid particles if the air phase is explicitly accounted for. The current implementation ensures that the two fluids remain distinctly separated. In contrast to other multi-phase SPH algorithms as demonstrated by Mokos et al. (2014), no unphysical voids occur in the air phase regardless of the chosen resolution, even if the ratio is increased. Despite the physical importance of using a multi-fluid model, Table 2 illustrates the high computational demand associated with this simulation.

Table 2

Computational characteristics of the single- and multi-phase breaking dam simulation

  Comp. time until t = 1s (min) Restrictive speed of sound (m s−1Air particles 
Single-phase 8.3 20  
Multi-phase 1,055 257.5 72,500 
  Comp. time until t = 1s (min) Restrictive speed of sound (m s−1Air particles 
Single-phase 8.3 20  
Multi-phase 1,055 257.5 72,500 
Figure 1

Pressure and velocity distribution of the water phase in the multi-phase dam break simulation at and .

Figure 1

Pressure and velocity distribution of the water phase in the multi-phase dam break simulation at and .

Figure 2

Air entrainment in the reflecting wave of a breaking dam. The computed phase interface is compared to an inviscid boundary element method solution (solid line, Colagrossi & Landrini 2003).

Figure 2

Air entrainment in the reflecting wave of a breaking dam. The computed phase interface is compared to an inviscid boundary element method solution (solid line, Colagrossi & Landrini 2003).

Rising air bubble

As an example of forced aeration occurring in environmental engineering, we study the rise of an air bubble in a stagnant liquid tank. Specified by particles with a uniform sampling distance in both phases, a circular bubble of radius is placed at the bottom of a water-filled tank with a width and height of and , respectively. The liquid's speed of sound is set to . This configuration yields a Reynolds number such that a turbulence model is not required. The Bond number confirms that buoyancy forces predominate surface tension effects and thus a corresponding model is not included in the simulation.

Bubble contours

Figure 3 depicts the particle distribution, pressure and velocity field of the rising bubble at . The bubble contour (black solid line) is compared to the SPH solution by Grenier et al. (2009). We analysed the sensitivity of the deformation of the bubble contours on the artificial repulsion strengthness parameter which is recommended to be set to (Meister & Rauch 2014). In consequence, the stiffness of the phase interface is relaxed and hence the bubble is deformed appropriately, i.e., it attains a horseshoe shape which eventually splits into three segments. The optimized algorithm's evolution of the bubble contours as a function of the dimensionless time is plotted in Figure 4. Owing to the symmetry of the problem, only the left half-plane is shown. The current method's results agree with other SPH (Colagrossi & Landrini 2003; Grenier et al. 2009) and level set method data (Sussman et al. 1994; Wang et al. 2009). Minor deviations like the underestimation of the necking at the lower end of the horseshoe shape are identified.

Figure 3

Pressure and velocity distribution of the air bubble rising in a water tank at . The bubble contour is compared to the reference solution by Grenier et al. (2009).

Figure 3

Pressure and velocity distribution of the air bubble rising in a water tank at . The bubble contour is compared to the reference solution by Grenier et al. (2009).

Figure 4

Evolution of the contour of an air bubble rising in water.

Figure 4

Evolution of the contour of an air bubble rising in water.

Bubble coalescence with the surface

Even though the present method accounts for the correct bubble dynamics within the liquid, the description of the air particles' movement fails as the free surface is approached. This is attributed to the weighting of the particle neighbourhood which becomes void if neighbours are missing. This difficulty is addressed in the next section, but we mention at this stage that the coalescence with the free surface is well described if the surface is explicitly filled up with air particles. To demonstrate this the height of the simulation domain is reduced so that the rising bubble approaches the surface at . The bubble then begins to interact with adjacent air particles. Single dispersed water particles, which penetrate into the air surface (see Figure 5(a)), are attracted to the liquid by the artificial repulsion forces and thus the different phases remain distinctly separated. The coalescence of the bubble alters the particle alignment at the surface, but a regular distribution without any singular particles or voids is maintained.

Figure 5

The bubble merges with the adjacent air phase. The dispersed water particles in (a) are attracted to the fluid due to the artificial repulsion forces and hence a sharp phase interface is maintained (b).

Figure 5

The bubble merges with the adjacent air phase. The dispersed water particles in (a) are attracted to the fluid due to the artificial repulsion forces and hence a sharp phase interface is maintained (b).

LIMITATIONS OF SPH FOR PRACTICAL APPLICATIONS

In the preceding section, we demonstrated the accuracy of the method in modelling the fundamental processes of air entrainment and rising air bubbles. The application to practical problems, however, is difficult due to the limited description of air inlets and the high computational demand which are summarized in the next section.

Air inlet and treatment of free surfaces

To advance to practical problems the presented single bubble study has to be extended to account for the modelling of continuous air inlets. Since the evolution of a particle's physical parameters is governed by a weighted influence of its neighbourhood, air inlet particles interact in both ways with fluid particles. In consequence, it is not possible to directly insert or delete particles and except for special cases even open boundary conditions cannot be accurately imposed (Ferrand et al. 2014). We extensively studied this problem with a novel inlet algorithm based on the buffer zone technique. If inlet particles are assigned constant velocities and the effective hydrostatic pressure such that the pressure gradient drives the inflow, then stable inflow configurations are attained. It is, however, unfeasible to impose the experimental inflow pressure, because the simulation blows up at an early stage. Adversely the bubble formation and detachment is unphysical. Even though the size of the developing bubbles agrees with experimental data, the detachment is deferred. Large numerical errors associated with this scheme cause a deformation of the symmetric bubble shape and the breakup into dispersed particles.

A physically correct inlet scheme can be realized by driving the inflow by a large enclosed pressurized chamber (Das & Das 2009), which is reviewed by Meister & Rauch (2014) using the present method. Despite the eligibility of this algorithm for academic test cases, this method cannot be applied to practical problems since up to 34,000 gas particles are required for modelling the detachment of a single bubble. With the current resources available, this computational overhead is not manageable except for a few small-scale studies. In addition, this method is not applicable to the modelling of a continuous inflow of air as, e.g., required for simulation of aeration-induced mixing in environmental engineering.

If a continuous inlet of air is attained, it is necessary to treat the degassing of air at the free surface unless it is explicitly described by air particles. In the latter case, as demonstrated, the coalescence of the bubble with the surface is accounted for correctly. In turn, this procedure is computationally expensive and only applicable to a finite simulation domain which is filled up with air particles. If treated as a true free surface, air particles have to be explicitly deleted, either based on a number of neighbours or a surface detection criterion (Ihmsen et al. 2011).

HYBRID SPH MODEL: BUBBLE COLUMN SIMULATION

Despite the recent increase in computing power and the establishment of GPU computing (see, e.g., Herault et al. 2010; Dominguez et al. 2013), the SPH modelling of bubbly air-water flow is currently limited to single bubble studies (Szewc et al. 2013) and the merging of two bubbles if a lower resolution is chosen (Grenier et al. 2013). In consequence, the direct application of conventional multi-phase SPH to the modelling of bubble columns remains out of bounds for all practical problems. To allow for simulations in the field of environmental engineering a novel multi-scale approach is proposed where the bubble characteristics are determined by a detailed multi-fluid simulation (see section ‘Rising air bubble’) and assigned to a separate low resolution model. Building on this approach, we study the homogeneous flow regime of a bubble column with a combination of experiments and the numerical SPH method. Despite the simplified flow regime we should emphasize the complexity of this example since effects of wall-bounded turbulence and the arrangement of aerators significantly affect the flow dynamics (for an introduction see, e.g., Kantarci et al. (2005)).

Benefits of this hybrid SPH model are an increase of both particle sampling distance and time steps such that the simulation of the collective movement of 35 bubbles becomes considerably cheaper than resolving a single bubble in detail. In contrast to the air concentration-based model by Chanson (1996) the present algorithm resolves the flow field, which is a prerequisite for coupling SPH to biochemical processes where the bubble-induced mixing of the biokinetic compounds is important. With the approach taken it is for the first time possible to reach the timescales of biokinetic processes which range from seconds to days. Further, the experimental continuous air inlet configuration can be attained which is an unsolved problem for conventional SPH models. The practicality of the approach is demonstrated for the reference experiment below.

Reference experiment

A reference experiment was conducted to determine the bubble velocities and to locate the homogeneous flow regime. Since the findings are passed on to the numerical simulation, the experimental flow configuration should resemble a two-dimensional model. Hence a sparger pipe with five equally spaced orifices with a diameter of is fixed at the bottom centre cross section of a three-dimensional cylindrical plexiglass column (see Figure 6) with the geometric specifications given in Table 3. This configuration limits wall effects and confines the bubble dynamics to the layer coplanar to the aerator. In the experiment, the aeration is driven by a membrane pump and its rate controlled by a flow-meter. To attain a pseudo-steady bubble configuration without any unsteady structures, i.e., a homogeneous bubbly flow regime, the aeration rate is set to .

Table 3

Geometric specifications of the bubble column

H (cm) D (cm) Aspect ratio H/D ( − ) Orifices ( − ) Orifice diameter (cm) Distance Δxo (cm) 
23.5 29 8.1 10−1  3.25 
H (cm) D (cm) Aspect ratio H/D ( − ) Orifices ( − ) Orifice diameter (cm) Distance Δxo (cm) 
23.5 29 8.1 10−1  3.25 
Figure 6

Setup of the bubble column experiment.

Figure 6

Setup of the bubble column experiment.

The bubble dynamics, which are captured by a high speed camera, are summarized in Table 4. In this flow regime, bubbles of approximately uniform size and velocity are equally spaced throughout the domain and a terminal velocity of is reached. Even though the bubble movement is approximately two-dimensional, the ellipsoidal formed bubbles rotate along their axes with a minimum and maximum width between 1.13 cm and 1.64 cm, respectively.

Table 4

Experimental bubble characteristics

Superficial air velocity (cm s−1Height (cm) Width (cm) Vertical distance (cm) Terminal velocity (cm s−1
4.2 10−3 5 10−1 1.13–1.63 2.86 31.59 ± 0.73 
Superficial air velocity (cm s−1Height (cm) Width (cm) Vertical distance (cm) Terminal velocity (cm s−1
4.2 10−3 5 10−1 1.13–1.63 2.86 31.59 ± 0.73 

Numerical simulation using experimental data

Based on the experimental bubble characteristics, a two-dimensional numerical SPH simulation is performed in the cross section coplanar to the aerator. Despite an otherwise identical geometrical setup as in the reference experiment, three-dimensional bubble characteristics like rotation are not included in the numerical model. Since the experimental bubble distributions and velocities are homogeneous throughout the domain, in the numerical simulation bubbles maintain a solid ellipsoidal shape with a height of and a width of 1.4 cm (see Figure 7) with an assigned velocity from the experiment.

Figure 7

Initial bubble shape in the numerical simulation.

Figure 7

Initial bubble shape in the numerical simulation.

In comparison with the direct computation of the bubble dynamics, this combined study allows for modelling the air-induced mixing in a bubble column with 35 single bubbles at a much lower computational cost than required for resolving a single bubble (see Table 5). First, the high density ratio between air and water is numerically challenging for the SPH method because of the discontinuity in density at the phase interface (Colagrossi & Landrini 2003). Regardless of the choice of the corrective algorithm (e.g., Colagrossi & Landrini 2003; Hu & Adams 2007; Grenier et al. 2009), all solutions require a large number of air particles to characterize a single bubble. Since adaptive resolution is not yet established for multi-phase flow, the high level of detail required for curing the numerical instability is defined on the whole simulation domain even if it is otherwise not required. In contrast, if a kernel influence radius of 2 particles per direction is specified, then with the current method only 66 instead of 4,900 air particles per bubble have to be used.

Table 5

Computational demand required for the simulation of a single bubble in comparison with an approximate bubble column model resolving 35 bubbles

Number of bubbles Comp. time until t = 1s (min) Restrictive speed of sound (cm s−1Particles Particles per bubble 
Single 6,608 29.7 93,400 4,900 
Column with 35 17.7 45 66,300 66 
Number of bubbles Comp. time until t = 1s (min) Restrictive speed of sound (cm s−1Particles Particles per bubble 
Single 6,608 29.7 93,400 4,900 
Column with 35 17.7 45 66,300 66 

Figure 8 indicates that the chosen particle spacing is sufficient for resolving the velocity field of the bubble column and for ensuring a stable long-term simulation. Second, in weakly compressible SPH the speed of sound in air should be 14 times higher than for water (Colagrossi & Landrini 2003) and hence the size of the time steps decreases accordingly. On the other hand, no change of the speed of sound is required if bubbles are modelled as solid air particles. In spite of that, Table 5 shows that the restrictive speed of sound of the bubble column simulation is comparable to the single bubble study since in the former case the bubble velocity is higher.

Figure 8

The computed velocity profile of the inner area of the bubble column filled with water is visualized on top of the snapshot of the experimental bubble characteristics which are captured by a high speed camera. At the bottom of the figure the experimental sparger pipe is clearly visible. In the experiment bubbles are homogeneously distributed throughout the domain which corresponds to the position of the rigid bubbles in the numerical simulation.

Figure 8

The computed velocity profile of the inner area of the bubble column filled with water is visualized on top of the snapshot of the experimental bubble characteristics which are captured by a high speed camera. At the bottom of the figure the experimental sparger pipe is clearly visible. In the experiment bubbles are homogeneously distributed throughout the domain which corresponds to the position of the rigid bubbles in the numerical simulation.

CONCLUSIONS AND FUTURE WORK

This work provides an overview of the SPH modelling of aerated flows. Results are produced with the self-developed code titled SPHASE which abbreviates Smoothed Particle Hydrodynamics Activated Sludge Engine. The software development is geared towards the future application in wastewater treatment but generally stimulates the future establishment of SPH as numerical CFD engine in environmental engineering. Real-time visualization of the computation is an important feature of the code and the corresponding program specifications are provided in this paper.

In multi-phase SPH modelling an instability occurs at the air-water phase interface which is cured by applying a weak artificial repulsion force between particles of different phases. The accuracy of the present method is demonstrated by modelling the fundamental processes of air entrainment and rising air bubbles. We have studied the impact of entrained air on the flow dynamics for the example of a breaking wave. In the rising bubble problem where effects of surface tension are neglected since buoyancy forces are predominant, the bubble contours are in accordance with level set method and reference SPH data. Hence unphysical effects of the artificial repulsion force can be ruled out. The modelling of continuous air inlets and the high computational demand associated with the method are identified as key limitations for establishing SPH in this field. Consequently, the simulation of permanent aeration and the resulting interfacial oxygen transfer remain out of bounds, despite an efficient code implementation and the use of multi-GPUs. Hence, in this paper, we propose a combined experimental and numerical study to model the homogeneous flow regime of a bubble column. Since the experimental bubble characteristics are imposed on the numerical simulation, the computational demand is significantly reduced such that the timescales of oxygen consumption due to biokinetic processes can be attained. In the future, this novel approach is expected to allow for the general application of SPH to environmental engineering problems, e.g., for modelling wastewater hydraulics. For a universal application the method is extended by a bubble coalescence and splitting model. A possible approach is seen in an algorithm developed for the animation of air bubbles (Ihmsen et al. 2011), where the high density ratio between air and water is avoided by separately calculating density and pressure for each phase.

ACKNOWLEDGEMENT

This research is part of a project that is funded by the Austrian Science Fund (FWF): [P26768-N28].

Availability of the developed software

Name of the software: SPHASE – Smoothed Particle Hydrodynamics Activated Sludge Engine

Developers: Gregor Burger, Michael Meister, Unit of Environmental Engineering, University of Innsbruck

Contact address: Michael Meister, Unit of Environmental Engineering, University of Innsbruck, Technikerstrasse 13, 6020 Innsbruck, Austria Michael.Meister@uibk.ac.at

Availability: on application to the author

Year first available: 2014

Hardware requirements: PC with multi-core processor

Software requirements: Linux, Mac OS X, Windows XP/7

Programing language: C ++

Program size: 245K

User interface: command line, graphical user interface

REFERENCES

REFERENCES
Adami
S.
Hu
X.
Adams
N.
2012
A generalized wall boundary condition for smoothed particle hydrodynamics
.
J. Comput. Phys.
231
(
21
),
7057
7075
.
Chanson
H.
1996
Air Bubble Entrainment in Free-Surface Turbulent Shear Flows
.
Academic Press
,
Waltham, MA, USA
.
Chanson
H.
2013
Hydraulics of aerated flows: qui pro quo?
J. Hydraul. Res.
51
(
3
),
223
243
.
Colagrossi
A.
Landrini
M.
2003
Numerical simulation of interfacial flows by smoothed particle hydrodynamics
.
J. Comput. Phys.
191
(
2
),
448
475
.
Dominguez
J.
Crespo
A.
Valdez-Balderas
D.
Rogers
B.
Gómez-Gesteira
M.
2013
New multi-GPU implementation for smoothed particle hydrodynamics on heterogeneous clusters
.
Int. J. Multiphase Flow
184
(
8
),
1848
1860
.
Dumbser
M.
2011
A simple two-phase method for the simulation of complex free surface flows
.
Comput. Methods Appl. Mech. Eng.
200
,
1204
1219
.
Ferrand
M.
Violeau
D.
Mayrhofer
A.
Mahmood
O.
2014
Correct boundary conditions for turbulent SPH
. In:
Advances in Hydroinformatics
(
Gourbesville
P.
Cunge
J.
Caignaert
G.
, eds).
Springer
,
Singapore
, pp.
245
258
.
Gomez-Gesteira
M.
Rogers
B.
Dalrymple
R.
Crespo
A.
2010
State-of-the-art of classical SPH for free-surface flows
.
J. Hydraul. Res.
48
,
6
27
.
Grenier
N.
Antuono
M.
Colagrossi
A.
LeTouzé
D.
Alessandrini
B.
2009
An Hamiltonian interface SPH formulation for multi-fluid and free surface flows
.
J. Comput. Phys.
228
(
22
),
8380
8393
.
Grenier
N.
LeTouzé
D.
Colagrossi
A.
Antuono
M.
Colicchio
G.
2013
Viscous bubbly flows simulation with an interface SPH model
.
Ocean Eng.
69
,
88
102
.
Henze
M.
Grady
C.
Gujer
W.
Marais
G.
Matsuo
T.
1987
Activated Sludge Model No. 1. IAWPRC Scientific and Technical Reports No. 1
.
IAWPRC [now IWA]
,
London, UK
.
Herault
A.
Bilotta
G.
Dalrymple
R.
2010
SPH on GPU with CUDA
.
J. Hydraul. Res.
48
,
74
79
.
Hu
X.
Adams
N.
2007
An incompressible multi-phase SPH method
.
J. Comput. Phys.
213
(
1
),
844
861
.
Hwang
Y.
Huang
H.
Chung
N.
Wang
P.
2013
Particle method of characteristics (PMOC) for unsteady pipe flow
.
J. Hydroinform.
15
(
3
),
780
797
.
Ihmsen
M.
Bader
J.
Akinci
G.
Teschner
M.
2011
Animation of air bubbles with SPH
. In:
Proceedings of the International Conference on Computer Graphics Theory and Applications
,
Algarve. Instituto Politécnico de Setúbal
, pp.
225
234
.
Kantarci
N.
Borak
F.
Ulgen
K.
2005
Bubble column reactors
.
Process. Biochem.
40
(
7
),
2263
2283
.
Meister
M.
Rauch
W.
2014
Simulating rising bubble problems with smoothed particle hydrodynamics
. In:
Proceedings of the International Conference on Mathematical Modelling of Fluid Systems, Bristol. Institute of Mathematics and its Applications
, pp.
1
7
.
Mokos
A.
Rogers
B.
Stansby
P.
2014
A multi-phase particle shifting algorithm for SPH simulations for violent hydrodynamics on a GPU
. In:
Proceedings of the International Conference SPHERIC SPH Workshop
,
Paris. EDF R&D
, pp.
1
8
.
Monaghan
J.
1994
Simulating free surface flows with SPH
.
J. Comput. Phys.
110
,
399
406
.
Monaghan
J.
2005
Smoothed particle hydrodynamics
.
Rep. Prog. Phys.
68
(
8
),
1703
1759
.
Monaghan
J.
2013
A simple SPH algorithm for multi-fluid flow with high density ratios
.
Int. J. Numer. Meth. Fluids
71
(
5
),
537
561
.
Morris
J.
Fox
P.
Zhu
Y.
1997
Modeling low Reynolds number incompressible flows using SPH
.
J. Comput. Phys.
136
,
214
226
.
Sussman
M.
Smereka
P.
Osher
S.
1994
A level set approach for computing solutions to incompressible two-phase flow
.
J. Comput. Phys.
114
(
1
),
146
159
.
Tartakovsky
A.
Meakin
P.
Scheibe
D.
West
R. E.
2007
Simulations of reactive transport and precipitation with smoothed particle hydrodynamics
.
J. Comput. Phys.
222
(
2
),
654
672
.
This is an Open Access article distributed under the terms of the Creative Commons Attribution Licence (CC BY-NC-ND 3.0), which permits copying and redistribution for non-commercial purposes with no derivatives, provided the original work is properly cited (http://creativecommons.org/licenses/by-nc-nd/3.0/).