Table of Contents
ISRN Materials Science
Volume 2014, Article ID 454691, 11 pages
Research Article

Gas Bubbles Expansion and Physical Dependences in Aluminum Electrolysis Cell: From Micro- to Macroscales Using Lattice Boltzmann Method

NSERC/Alcoa Industrial Research Chair MACE3 and Aluminium Research Centre (REGAL), Laval University, Quebec, QC, Canada G1V 0A6

Received 29 October 2013; Accepted 5 January 2014; Published 23 February 2014

Academic Editors: C. Carbonaro and A. O. Neto

Copyright © 2014 Mouhamadou Diop et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.


This paper illustrates the results obtained from two-dimensional numerical simulations of multiple gas bubbles growing under buoyancy and electromagnetic forces in a quiescent incompressible fluid. A lattice Boltzmann method for two-phase immiscible fluids with large density difference is proposed. The difficulty in the treatment of large density difference is resolved by using nine-velocity particles. The method can be applied to simulate fluid with the density ratio up to 1000. To show the efficiency of the method, we apply the method to the simulation of bubbles formation, growth, coalescence, and flows. The effects of the density ratio and the initial bubbles configuration on the flow field induced by growing bubbles and on the evolution of bubbles shape during their coalescence are investigated. The interdependencies between gas bubbles and gas rate dissolved in fluid are also simulated.

1. Introduction

A steadily increasing computational power, new developments, and improvement of numerical methods allow numerical simulation to model more and more physical phenomena. In this sense, the in situ gas generated by alumina particles and carbon engenders the formation of gas bubbles in aluminum cells. These gas bubbles are produced by chemical reaction of alumina with the carbon, which consists of complex process of electrochemical reduction of alumina in the cells. The continuously released gas of the alumina reaction generates nuclei underneath the surface of carbon anode and detaches to dissolve in the cryolite which subsequently grows to bubbles. The bubbles arrange into a cellular structure in the cryolite, which is preserved by a rapid cyclic flow of gas.

A major aim of the simulations presented here is to gain a better fundamental physical understanding of gas bubbles because of the complexity of gas bubbles- liquid flow in chemical system [1, 2], since it determines largely the dynamic of the cryolite. The detailed understanding may help to reduce the energy consumption and to improve the present energy efficiency by controlling the dynamics of bubbles in cells.

In contrast with actual aluminum electrolysis cell models, which focus rather on specific aspects and for the most part on gas bubbles system, a comprehensive and general gas bubbles model is needed to be developed here which includes the wealth of phenomena occurring in aluminum electrolysis cell process. The model will elucidate the influence of gas bubbles as described in [3, 4] and processing parameters such as the cryolite viscosity, the surface tension, the diffusion coefficient of the propellant in the cryolite, the ambient pressure, the heating rate, and the cell geometry. This information is then helpful for optimizing the gas bubbles structure. However, as an imitation of the real world and due to the limited computational speed, some physical processes must be omitted or simplified in the physical model. The comprehensiveness and generality of the physical model demand for a numerical solution. The choice of the numerical method is made in favor of the lattice Boltzmann model (LBM). LBM explores in particular complex flow problems to which the bubbles clearly belong with its large, dynamically moving gas/liquid interfaces [5, 6]. Advantages of the method are the computational locality and parallelism, the relative ease of implementation, the possibility for the incorporation of additional features, and the strong physical background as a kinetic scheme. These advantages allowed developing a complete simulation tool within the scope of this work [712]. The main part of this work is concerned with the development and adaptation of the new method called LBM to specify the aluminum electrolysis process. As a main goal the complex gas bubbles evolution in aluminum electrolysis cells is studied with an efficient method that is able to overcome difficulties in the treatment of complex interfaces [1319]. For capturing the interface, novel-free surface boundary conditions are developed. In consistency with the hydrodynamical problem, the LBM is utilized for description of the diffusive processes of the propellant. Furthermore, additional model features can be incorporated within the framework of the LBM to tackle a multiphysical system of millions of degrees of freedoms.

2. Physics of Gas Bubbles

This section gives a detailed introduction to the technological aspects and the physical phenomena of the gas bubbles process by in situ gas generation. In particular, the focus relies upon the anodic chemical reaction of carbon.

The properties of the alumina, particle size and distribution, atomization atmosphere, purity, and heat treatment are crucial for the aluminum electrolysis process. The alumina and the cryolite assumed are mixed homogeneously. The subsequent nucleation method needs to assure a high densification in order to suppress gas loss by percolation at the beginning of the process. Optionally, the compacted material can be cut into near net-shape forms to study the phenomena of channelization and penetration of fluids in fracture.

The carbon anode considered as compacted material is constituted of one or several pieces of the compacted precursor. In the smelter there are flows from cell to cell in aluminum busbars while, in each cell, the current flow (I) is set to 25 KA and run downwards through the cell composed of 8 anodes, cryolite, molten metalpad, and a cathode carbon block.

In a matter of minutes after the beginning of the electrolysis process, gas bubbles nucleate and expand underneath the anodes in the cryolite and the process is interrupted by removing the gas out of the cell in order to stabilise the system. The final product consists of a closed cell with a cyclic formation of gas bubbles generated from the periodical reaction of alumina in the cryolite from feeders and surrounded by a fully dense skin called the top crust. The density in the multiphase fluid flows ranges from 0.4 to 0.8 g/cm3, which corresponds to a relative density of 15–30%.

Focusing on the actual electrolysis process, five expansion stages, dependent on the temperature, are identified in Figure 1: thermal expansion of the precursor (0), bubble nucleation under surface of anodes (I), expansion in the semisolid and liquid range (II), expansion in the liquid temperature range (III), and finally the gas bubbles collapse (IV).

Figure 1: Schematic expansion of gas divided into five regimes. Phase (0) depicts the thermal expansion of the precursor, phase (I) denotes the bubbles nucleation and expansion under surface of anodes, the expansion in the semisolid and liquid range is in phase (II), phase (III) is the expansion in the liquid temperature range, and finally phase (IV) is the gas bubbles collapse.

The beginning of the gas formation, regime 0–II in Figure 1, is marked by decomposition of the alumina, melting the precursor and nucleation. The decomposition of about 6 wt% up to 960°C releases approximately a gas volume which exceeds 10 times the precursor volume. The decomposition begins while the alumina matrix is still solid (at about 950°C). Partial pressures in the order of the matrix yield strength may arise and tear open the powder compact precursor.

The crack-like openings in the carbon anode are problematic because the CO2 can percolate through these channels to the ambience being lost for the multiphase system. In addition to these cracks, the small pores over the anode serve as nuclei when the alumina melts. In the melt, heterogeneous nucleation may occur at nucleation sites such as impurities or propellant particles.

The major gas bubbles expansion takes place while the cryolite is liquid, regime III-IV in Figure 1. The strong increase in expansion relies on the substantial oxide carbon release during the considered temperature interval. The oxide carbon dissolves and diffuses into numerous nuclei, which consequently expand. At a certain point, the gas bubbles impinge upon each other. They must give up their favorable shape. Thin fluid films separate the bubbles. The variations of gas bubble sizes and thus interior pressure as well as preceding film ruptures lead to bent and corrugated cell walls.

The stabilization of the gas bubbles films is essential. Without stable cell walls, bubbles would immediately undergo coalescence when they contact each other. By coalescence with the atmosphere enforced by buoyancy, the bubbles and consequently massive gas portions from alumina are lost from the system. A typical gas bubbles structure cannot evolve. The gas bubbles may owe stable films to the presence of additives such as oxides. Since the alumina possesses oxide, that oxide naturally resides as a network in the precursor. The repulsive force is termed disjoining pressure in the following. Additionally, network fragments make thin films less penetrable to the melt so that the cryolite remains in the films. These effects resemble a drastically enlarged viscosity. This assumption is supported by the observed shape stability of a once melted alumina precursor.

During the gas bubbles expansion, film ruptures, growth coalescence, capillary, gravitational drainage, and electromagnetic, interbubble diffusion, gas bubbles flow and gas bubbles collapse occur. They are sketched in the following.

Film ruptures are inevitable. One observes a minimal film thickness in gas bubbles of the order of 50 μm. This leads to the assumption that a film rupture event occurs when this critical minimal film thickness dc is exceeded. A rupture event causes the collapse of the entire film. Surface tension forces withdraw the film remnants towards the surrounding structure. The structure time , starting from the initial defect in the film to its complete disappearance, is a matter of seconds. One cause for the film thinning is simply the expansion of the bubbles. If gas bubbles pore inflates, its cell wall thickness decreases due to mass conservation. At a certain point, is reached leading to rupture and coalescence. Thereby, the gas bubbles system coarsens, that is, the mean pore diameters increase. Considering a representative cell of the gas bubbles, the ratio of pore volume to pore shell is approximately proportional to .

Gas bubbles drainage is the flow of liquid through and out of the gas system under the influence of external and capillary forces. In the gravitational drainage, the cryolite in the plateau borders is drawn downwards by graviton so that the cryolite accumulates downwards at the cryolite layer. In the capillary drainage, the different pressures in cell films and plateau borders originating from the surface tension cause the cryolite to withdraw from the films. Naturally, these processes will reduce the film thickness in the gas bubbles system which leads to cell wall ruptures. Estimating the times scale of gas bubbles expansion, it is found to be very short compared to usual bubble in fluid. Only a drastically enlarged viscosity and in particular the static disjoining pressure in the films hinder the gravitational and capillary drainage to cause gas bubbles collapse. When the stabilization is significant, the drainage may even become negligible.

Distinct gas bubble sizes and consequently distinct bubbles may lead to interbubble diffusion of oxide carbon via the cell walls. Tiny bubbles suffer from gas loss, while large bubbles earn oxide carbon, similar to Ostwald ripening. A study of the diffusion flow of oxide carbon across cell borders demonstrates the large time scales, compared to the bubbly times, of the gas transport in cryolite. Therefore, the interbubble diffusion is negligible in the cryolite bubbly evolution. During the evolution, gas bubbles are permanently exposed to mechanical forces due to the expansion itself or spatial restriction by the cell design. Deformation of pores, rearrangement or destruction of cell films, and gas bubbles flow are the consequences. Topological rearrangements and film ruptures occur often in an avalanche-like manner. In the final stage of bubbles expansion, large bubbles are present compared to the gas bubbles system size. Gas bubbles collapse at the cryolite surface leads now to severe oxide carbon loss, in the system so that the gas system decays; that is, the collapse starts.

3. Physical Model

Reviewing gas bubbles in aluminum electrolysis cell in Figure 2, almost no theoretical models can be found. Considering the bubbly evolution in liquid, it becomes apparent that specific aspects of gas bubbles are commonly studied in isolation. In order to include many physical phenomena of the electrolysis process, a comprehensive model is constructed based on fundamental physical equations. In order to obtain the full information of gas bubbles structure, it is required to directly simulate gas bubbles structure in aluminum electrolysis process. Nevertheless, the model will naturally contain simplifications because it is a realistic imitation of the real industrial problem and because of inaccessibility in cells.

Figure 2: Vertical cut short of an aluminum electrolysis cell in which the gas expands from liquids. The points , , and are the references points for the evolution of gas in the cell.

The gas bubbles are stochastically distributed and formed simultaneously at the starting point. Dissolved alumina particles reacting with the carbon anode diffuse into the nuclei over the surface of carbon anodes. The decomposition and dissolution of alumina are modeled by a homogeneous volume source in the cryolite domain . The nuclei grow over the undersurface of anodes to bubbles. Consequently, the surrounding cryolite is displaced. Interactions between bubbles are accounted for at first over the anodes by a short-range disjoining pressure. Cell walls ruptures occur when a minimal critical film thickness is attained. As a simplification, the model uses a static electromagnetic force term in the expression of Navier-Stokes equation. Furthermore, the starting point of the model indicates already that the melting of the alumina precursor and an explicit description of the chemical decomposition of the propellant are neglected. Similarly, the final bubbly structure is achieved by simply stopping the simulation at certain point of time. In this model the temperature remains constant at 960°C in space and time. Consequently, the material parameters are constants. A temperature and heating rate enters the model only implicitly by the homogeneous volume source term of the dissolved oxide carbon.

The solute oxide carbon in the cryolite obeys the diffusion equation in : where   is the volume source term, is the gas diffusion, is the concentration field of oxide carbon, and denotes the time.

The liquid domain contains islands, which are the nuclei or gas bubbles. At the gas-liquid interface , Sievert’s law holds as follows: where the pressure is distinct for each bubble and is Sievert’s constant. Furthermore, the gradient at the interface is parallel to the unit normal direction of the interface: The dynamics of cryolite flow caused predominately by gas bubbles expansion is described by the famous incompressible Navier-Stokes equation in : Since the viscous stresses of the gas are negligible, free surface boundary conditions are applicable at : where the subscripts and denote the normal and tangential component with regard to . The gas pressures of bubbles , the surface tension pressure , and the disjoining pressure enter the model via these boundary conditions. In the bubbles, the ideal gas law is considered: where   is the gas pressure in bubble given by the ideal gas equation, denotes the ideal gas constant and denotes the amount of gas, denotes the temperature, and denotes the volume of bubble . The change of bubble volume per time is given by the normal velocity of the bubble-liquid interface: The gas exchange between cryolite and bubbles is given by Frick’s law so that the following change of the amount of substance per time applies: A special type of bubble is the atmosphere. It has a constant pressure at all times. Its oxide carbon partial pressure, entering Sieverts’ law, can be chosen differently from .

A heuristic approach condenses the film stabilization in a single function, the disjoining pressure . is a repulsion experienced by bubble in the presence of another bubble : The distance from bubble to bubble is obtained by moving along the normal direction. When the neighboring bubble is too distant; that is, exceeds , no interactions of the bubbles occur. At present, the functional form of is chosen arbitrarily.

When the film thickness falls below a critical threshold , the cell wall is subjected to rupture. These coalescence criteria are extracted in agreement with numerical experimental data. To conclude, desired features, such as rheology or drainage, evolve automatically in the model because of its generality. Furthermore, processing parameters, such as viscosity or ambient pressure, are incorporated naturally and hence accessible to investigations.

3.1. Numerical Model

The liquid domain in which gas bubbles develop contains two domains: the gas and the liquids . This domain can be open or closed as shown in Figure 2. One notes the entire walls of the domain. We will seek to establish a system of equations verified by the velocity field and the pressure in all the domain , that is, in as in . The velocity field and the pressure are defined in the liquid part in the same manner they are also in the gaseous part, but this time taking into account that . In the mixture gas liquid, one applies the conservation laws. In addition, one applies a condition of sticking contact on walls of the entire domain .

The last point is disputable and constitutes a first approximation. This leads to the following formulation.

Find the velocity , , and the pressure , such as : where is the gaseous volume fraction in the mixture within the cryolite before the rise of gas bubbles over the surface of the liquid. is a field, depending on the point where we are in the gaseous domain and time. The viscosity depends on , as well as the gas-liquid system. By neglecting the volume mass of gas () before the volume mass of the liquid (about ), on can write The cryolite is assumed to be incompressible, so that The micro-macro-coupling is carried out then through the expansion rate of the gas that remains to be determined. In the microscopic case, the expansion of gas bubbles is governed by the difference of pressures between the internal and the ambient, and the latter constitutes a macroscopic parameter of the problem, exterior of our microscopic model. Thus, it follows that one uses where denotes the internal pressure of the gas and denotes the expansion rate.

This expansion rate is explicitly dependent on time because it implies a gas creation by chemical reaction from the alumina. The time appears explicitly in the expression of . Finally, the system in velocity pressure to solve writes as follows.

Find the velocity , and the pressure , such as : After the germination of gas bubbles and their detachment from the undersurface of the anodes, the gas bubbles still expand in the gas-liquid mixture and are assumed to be spherical. They rise to the surface and start to reduce their interaction between them. At the surface of the liquid they occupy a small percentage of the total volume, to a structure of mixture gas and liquid with deformed gas bubbles using more than 85% of the total volume.

3.2. Micro-Macro-Coupling

The micro-macro-coupling in the aluminum electrolysis cell is carried out by closing (15) via the expansion rate (see (14)). One proposes to determine this expansion rate from simulations using a microscopic model. This model takes into account the expansion of a representative volume during a moment of the macroscopic flow. This model takes as input variables the initial configuration of gas as well as its initial pressure, and it intervenes as a variable of macroscopic order, , in the ambient pressure at a considered point of the macroscopic flow. The expansion rate is expressed first when the interactions of gas bubbles are not taken into account. In this study, a volume representative in expansion, of volume , contains a gaseous part of volume and a liquid part of a constant over time. The gas rate is defined by . One recovers the following relation between and : It is assumed that the gas in the cryolite is composed of spherical bubbles of volume , radius , and the pressure . Whereas the cryolite is regarded as Newtonian and of consistency . The pressure is assumed to be constant in the liquid and equal to the external pressure that represents the constraint applied on the wall by the external medium. Thus, gas bubbles grow by difference of pressure, and their interactions are neglected here. Consequently, one finds that So that an algebraic transformation will give Without the interactions between bubbles, one applies over each bubble the expression of the velocity (), pressure (), and . Finally, it gives The above equation now is simplified because the bubbles are taken with the same pressure in the entire volume of the domain . One obtains The differential of (16) and its combination with (20), in function of pressure’s difference writes Finally, one expresses the pressure of gas in function of . is supposed to be homogeneous in the gaseous domain, by using the ideal gas law for each bubble and by summing over all bubbles, are obtained that . Using (16), one can write that Thus, by neglecting the interactions between bubbles, and by assuming a homogeneity of the pressure, the gas rate expansion in the representative volume gives the expression in function of and . Now, following, a macroscopic approach, the expansion rate of the gas is a field. One obtains then To solve numerically (23), its domain of validity must be extended in all the domain , in a manner to retrieve its second term null in the liquid domain .

4. Numerical Implementation

Reviewing the physical model from a numerical point of view, the most challenging aspect is the parallel implementation of the large, complex, and dynamically moving interfaces between the cryolite and the gas bubble. Among the large variety of numerical methods, the lattice Boltzmann method (LBM) [2023] is a promising candidate for numerical investigations in aluminum electrolysis cell. Suitable free surface boundary conditions need to be newly developed for the lattice Boltzmann method. In addition, the other features of the physical model, for example, the diffusion equation and the disjoining pressure, need to be included into the LBM.

At first, the focus remains on the fluid dynamics. Since the beginning of the last century, it is known that the kinetic Boltzmann equation can approximate the Navier-Stokes equation. By numerically solving the kinetic Boltzmann equation, the lattice Boltzmann model supplies an approximate solution for the Navier-Stokes. The LBM is consequently a kinetic scheme, based on particle dynamics. In detail, the space is divided into a 2D regular lattice. On each lattice site, single-particle distribution functions represent the number of particles at lattice site and time in state .

The distinct states are defined by a set of nine discrete velocities: In each time step , the distribution functions are advected along their corresponding velocity and subsequently they undergo a collision operation locally at each time site. In the simplest case, the collision step consists of a relaxation of the function towards an equilibrium function . Both steps are graphically displayed in Figure 8 and combined in the lattice Boltzmann equation (LBE) [24]: The macroscopic variables and are obtained by the zeroth and first moments: The collision operator is constructed to fulfill mass, momentum, and energy conservation at each lattice site and to ensure the correct symmetries of the Navier-Stokes: and is the sound speed squared.

The force term includes the gravity g and electromagnetic force [23, 25, 26]: By Chapman-Enskog analysis, the LBM is proved to approximate the Navier-Stokes in the incompressibility, , and for small Reynolds numbers. The Chapman-Enskog analysis starts with a separation of the time scales and expansion of the distribution functions: where is a small number. The following constraints are imposed on the expanded distribution functions: Further, the left hand side of the LBE is expanded by a second order Taylor series and the resulting terms in (25) are substituted by the above series expansions. An equation is derived containing various orders of epsilon. Evaluating the first two moments of this equation results in the hydrodynamic equations. The first order of epsilon provides the Euler equations. Including the second order of epsilon, the NSE is derived: where .

The presented LBM is applied in the fluid region [2729]. In the gas bubbles system, gas regions and the outer cell designs are present as well. In the regular lattice, one may define the distinct regions by distinct cells: Fluid , bubble , interface , and wall cells . The LBM is executed on the fluid cells. In the interface cells, the free surface boundary conditions in (5) and need to be implemented in (6). After the advection step and before the collision takes place, the interface cells obtain known distribution functions from liquid regions and unknown ones from the gas regions. The missing populations are reconstructed by a first order Chapman-Enskog approximation; see (30), where in all terms of are neglected: In , (27), in (33). The macroscopic quantities , , and are unknown. The free surface boundary conditions, (5) and (6), are utilised to eliminate the unknowns , , and . Therefore, the scalar products are expressed in the local coordinate system at the interfacial consisting of the normal and the vector tangential as depicted in Figure 3. Furthermore, the quadratic velocity terms in are neglected in order to obtain linear system of equations. Then, Here, , . is the pressure at the interface of the considered gas bubble including surface tension and disjoining pressure.

Figure 3: Interface cell after advection step with unknowns function (dotted line) and known function (full line).

After reconstruction the usual collision operation takes place. Yet, the interface advection remains to be discussed. For this reason, interface cells contain additional information which indicates the amount of fluid inside each cell. This liquid volume fraction is updated according to the in- and outflow of the interface cell. For more accuracy, the difference between the incoming and the outgoing distribution functions is calculated on each link and weighted by factor which depends on the position of the interface:

With is the area fraction obtained through , , and is related to via . The total flux in the interface cell is obtained by summing over all .

Empty or full interface cells are detected when exceeds its limits, . They are changed to gas or fluid cells, respectively. When cells are transformed, the closeness of the interface needs to be guaranteed; that is, unknown distribution functions must be advected into interface cells. This demands sometimes transformations of cells in the close neighborhood of the initially transformed interface cell. Unknown distribution functions in these cells are obtained by linear extrapolations. The simplest wall boundary conditions are utilized in the wall cells, which represent the cell design.

Bubble cells surrounded by interface cells build connected and closed areas of cells which represent the bubbles. The volume and amount of gas in the bubbles need to be determined in each time step. Initially, the bubble volume is known. Then, it is sufficient to calculate the change of volume per time step. It is given by the change of the volume fractions in the interface cells of the bubble : Similarly, only the oxide carbon diffusion across the interface is recorded in order to obtain the amount of substance in the bubbles. In the lattice Boltzmann formalism, the flux is given by the exchange of distribution functions across the interface which is approximated by the cell borders between fluid and interface cells.

5. Numerical Validation

These numerical implementations need to be validated. For this, comparisons between numerical and analytical results give information on the validity range of the parameters and the accuracy of this scheme. This test will merely focus on newly developed implementation issues and on single relaxation time, whereas the interface advection is tested by streaming a gas bubble in uniform flow field. Distortions of the round bubble shape and the small deviations from analytically predetermined way of the bubble center appear after the advection. Overall, the advection appears to be very promising results.

The same holds for the curvature algorithm. The results of template sphere methods are better than those of standard finite difference methods (FDM) often applied to volume of fluid methods (VOF).

The Laplace law has been used for several parameters by considering bubble in liquid flow. In Figure 4, the numerical result depicts a small deviation from the theoretical expectation marked by the dotted lines. During these investigations, small spurious currents occur. For the tested variables comprised between , , and for the fluid, are used all over the computations.

Figure 4: Numerical tests of Laplace law. The dashed lines represent the theoretical equation.

They remain small overall the experiment. Dynamically, the free surface boundary conditions are investigated by standing capillary and external forces waves such as gravity and the electromagnetic force. These waves are damped by the viscosity of the liquid. Frequency and damping rate obtained numerically are in good agreement with theoretical expectations. Good agreement is also depicted by a numerical investigation on the free surface boundary condition as shown in (5) at the interface of a rising gas bubble. In addition, the gas bubble shapes and their rising velocity are observed. The later seems to decrease very fast below expectations with a shrinking effect, when the surface tension is much larger than the effects of the external forces impact small or negligible. In conclusion, one can summarise that the numerical experiment exhibits good agreements with theoretical expectations in regard to parameter ranges of our numerical tests.

6. Results and Discussion

Having provided and tested the physical model and its numerical implementation, the simulations of gas bubbles expansion in aluminum cell process are finally performed. The tests carried out show that the consistency of the liquid intervenes as multiplicative factor in the expression of the viscosity in the system bubbles liquid . In this study, one considers a system of bubbles from 84 to 1500, that is, having different densities of gas bubbles. The choice of the structure of bubbles does not influence the viscosity of bubbles liquid and helps to keep the integrality of bubbles in the computational domain during the bubbles expansion. Thus, the effects of edges are restricted. In Figure 5, one plotted the values of the ratio in function of gas rate, computed during the expansion of bubbles. For a number of bubbles over 84, we have uniquely represented the obtained values before the effect of edges intervenes (rate <70%), whose description of bubbles-liquid structure is good. Comparing the curve of Figure 5 with the right 1-W plotted on the same figure, one observes a nonlinear dependence of the viscosity to the gaseous fraction . One distinguishes three different areas in the curves of Figure 5. On one hand, for a gas rate inferior to 20%, the more the density of bubbles is important, the smaller the bubbles are. They offer less resistance to shearing. This explains the fact that the equivalent gas rate is high. As long as the interaction between bubbles is negligible, the gas fraction remains small.

Figure 5: Function obtained by numerical simulation in function of the number of gas bubbles.

Finally, when the gas fraction tends towards 1, the viscosity also tends asymptotically towards zero. Thus, we have put in evidence, on one hand, the nonlinear dependence of the viscosity of mixture and, on the other hand, the dependence of the viscosity on the bubbles density. One remarks that the density of bubbles is consequent when the viscosity tends to an asymptotic curve.

6.1. Dependence of Viscosity on Bubbles-Liquid Mixture

There are several numerical simulations of bubbles-liquid structure by varying the bubbles density. The initial gas fraction is the same for each simulation, , as well, and the number of bubbles considered is important, more the sizes of whose, given at same gas rate are small. One used in this study an ordinated structure (bubbles dispersed regularly on the grid) and bubbles having different sizes and same initial pressure in a manner to conserve the integrality of the structure in the computational domain, on another hand to satisfy the hypothesis of the model. Because the interaction among bubbles of same pressure is absent, the gas rate model predicts an accurate profile for 84, 520, 789, 1000 and 1500 bubbles, comparing with the analytical expression of gas rate. The results show that the gas rate and the expansion rate are strongly influenced by the density of bubbles. The former will decrease with the increase of latter. This reflects the interactions between bubbles and the constraints exerted by the ones to the others when the expansion of 84 bubbles is considered. As long as the gas fraction does not exceed 45%, its evolution remains similar to that obtained without interaction. The divergence shown between the two fractions beyond 50% explains the fact that the gas goes out of the computational volume causing a decrease of the amount of liquid present in this volume. Consequently, beyond a certain limit, the expansion rate of gas fraction is overestimated compared to the reality. This approach is applicable for each simulation.

Figure 6 shows that the more the number of considered bubbles is important, the more these effects of edges are negligible, compared to whole bubble remaining inside the system. We choose to express the interaction of bubbles in the expression of gas rate by introducing simply in (23), a parameter determined from the curves in Figure 6. Its value is in function of bubbles density (Figure 7) and then depends also on the radius of these bubbles. The curves corresponding to (23) related to the parameter are plotted in Figure 7. Finally, a micro-macro-coupling is used to express the expansion rate of the mixture bubbles liquid in which the liquid is Newtonian and of consistency .

Figure 6: Evolution of gas rate over simulations of volume expansion of bubbles liquid with 84 to 1500 bubbles.
Figure 7: Evolution of gas rate during expansion of bubbles-liquid structure. The points correspond to the obtained values by computation of the value of parameter. The curves are the approximations estimated from (23).
Figure 8: Evolution of gas fraction, viscosity, velocity, and of the pressure at different points in the cell during the gas expansion.
6.2. Free Expansion of Gas Bubbles in Cell

The computational domain is regarded as an open system; that is, the gaseous matter can freely go out of the boundary part located at a certain height on which one applies a constraint null such as () as indicated in Figure 8 depicting the aluminium electrolysis cell. We have chosen, on Figure 8, the evolution of the gas fraction, viscosity, and of the pressure at different points A, B, C, and D. The pressures lower than the previous one in Figure 8(d), particularly at , due to the outgoing face on which ; the velocity at the end of the expansion is itself higher than the previous one in Figure 8(c).

The micro-macro-coupling used in this study is based on the fact that the simulation of the expansion of the representative volume needs a macroscopic variable, that is, the ambient pressure. The expansion happens through the fluctuations of pressure in the mixture bubbles-liquid. Thus, heterogeneities of gas fraction and then the size of bubble (or the density of bubbles) appear in the mixture. The expansion of the mixture depends on its structure, via the viscosity and the volume mass.

7. Conclusion

In this study, a new scheme has been developed to simulate the expansion of gas bubbles in aluminum electrolysis cell in order to describe the evolution of the viscosity and the volume mass of a representative volume of bubbles-liquid mixture during the expansion of the gas. It is found that the dependence of these two parameters is on the bubbles density in the mixture. The influence of gas bubbles interactions on these parameters has been revealed. From these results, we have carried out a micro-macro-coupling to express the evolution of the mixture volume mass during expansion in function of pressures fluctuation in the mixture. This approach helps to take into account the electrolysis process and also the influence of the mixture structure during the simulations of expansion. These simulations provide the variations of gas fraction in the cell. The simulations also show that, close to the cooled edges, the gas expansion is slower and even stopped. Thus, a fine surface skin is created on the cryolite.

Conflict of Interests

The authors declare that there is no conflict of interests regarding the publication of this paper.


The authors would like to acknowledge the financial support of the Natural Sciences and Engineering Research Council of Canada and Alcoa. A part of the research presented in this paper was financed by the Fonds de Recherche du Quebec-Nature et Technologies (FRQ-NT) by the intermediary of the Aluminium Research Center-REGAL. The authors wish also to thank Dr. Donald Ziegler from the Alcoa Canada Primary Metals for the scientific discussions.


  1. R. Keunings, “Advances in the computer modeling of the flow of polymeric liquids,” Computational Fluid Dynamics Journal, vol. 9, pp. 449–458, 2001. View at Google Scholar
  2. L. Lefebvre and R. Keunings, “Numerical simulation of chemically-reacting polymer flows,” in Proceedings of the 1st European Computational Fluid Dynamics Conference, vol. 2, pp. 1133–1138, Brussels, Belgium, September 1992.
  3. A. Arefmanesh and S. G. Advani, “Diffusion-induced growth of a gas bubble in a viscoelastic fluid,” Rheologica Acta, vol. 30, no. 3, pp. 274–283, 1991. View at Publisher · View at Google Scholar · View at Scopus
  4. D. F. Baldwin, C. B. Park, and N. P. Suh, “Microcellular sheet extrusion system process design models for shaping and cell growth control,” Polymer Engineering & Science, vol. 38, no. 4, pp. 674–688, 1998. View at Publisher · View at Google Scholar · View at Scopus
  5. R. Benzi, S. Succi, and M. Vergassola, “The lattice Boltzmann equation: theory and applications,” Physics Report, vol. 222, no. 3, pp. 145–197, 1992. View at Publisher · View at Google Scholar · View at Scopus
  6. P. L. Bhatnagar, E. P. Gross, and M. Krook, “A model for collision processes in gases. I. Small amplitude processes in charged and neutral one-component systems,” Physical Review, vol. 94, no. 3, pp. 511–525, 1954. View at Publisher · View at Google Scholar · View at Scopus
  7. J. H. Han and C. D. Han, “Bubble nucleation in polymeric liquids. I. Bubble nucleation in concentrated polymer solutions,” Journal of Polymer Science B, vol. 28, no. 5, pp. 711–741, 1990. View at Publisher · View at Google Scholar · View at Scopus
  8. J. H. Han and C. D. Han, “Bubble nucleation in polymeric liquids. II. Theoretical considerations,” Journal of Polymer Science B, vol. 28, no. 5, pp. 743–761, 1990. View at Publisher · View at Google Scholar · View at Scopus
  9. S. Hutzler, D. Weaire, and F. Bolton, “The effects of Plateau borders in the two-dimensional soap froth. III. Further results,” Philosophical Magazine B, vol. 71, no. 3, pp. 277–289, 1995. View at Publisher · View at Google Scholar · View at Scopus
  10. C. Isenberg, The Science of Films and Soap Bubbles, Dover, New York, NY, USA, 1992.
  11. D. Weaire and R. Phelan, “A counter-example to Kelvin's conjecture on minimal surfaces,” Philosophical Magazine Letters, vol. 69, no. 2, pp. 107–110, 1994. View at Publisher · View at Google Scholar
  12. J. U. Brackbill, D. B. Kothe, and C. Zemach, “A continuum method for modeling surface tension,” Journal of Computational Physics, vol. 100, no. 2, pp. 335–354, 1992. View at Publisher · View at Google Scholar · View at Scopus
  13. N. S. Ramesh, D. H. Rasmussen, and G. A. Campbell, “Numerical and experimental studies of bubble growth during the microcellular foaming process,” Polymer Engineering & Science, vol. 31, no. 23, pp. 1657–1664, 1991. View at Publisher · View at Google Scholar
  14. L. E. Scriven, “On the dynamics of phase growth,” Chemical Engineering Science, vol. 10, no. 1-2, pp. 1–13, 1959. View at Publisher · View at Google Scholar · View at Scopus
  15. D. L. Weaire and S. Hutzler, The Physics of Foams, Oxford University Press, New York, NY, USA, 1999.
  16. R. Clift, J. R. Grace, and M. Weber, Bubbles, Drops and Particles, Academic Press, New York, NY, USA, 2002.
  17. F. Bolton and D. Weaire, “The effects of Plateau borders in the two-dimensional soap froth I. Decoration lemma and diffusion theorem,” Philosophical Magazine B, vol. 63, no. 4, pp. 795–809, 1991. View at Publisher · View at Google Scholar · View at Scopus
  18. J. M. Buick, Lattice Boltzmann methods in interfacial wave modelling [Ph.D. thesis], University of Edinburgh, Edinburgh, Scotland, 1997.
  19. J. W. Bullard, E. J. Garboczi, W. C. Carter, and E. R. Fuller Jr., “Numerical methods for computing interfacial mean curvature,” Computational Materials Science, vol. 4, no. 2, pp. 103–116, 1995. View at Publisher · View at Google Scholar · View at Scopus
  20. N. Cao, S. Chen, S. Jin, and D. Martinez, “Physical symmetry and lattice symmetry in the lattice Boltzmann method,” Physical Review E, vol. 55, no. 1, pp. R21–R28, 1997. View at Publisher · View at Google Scholar · View at Scopus
  21. X. He and L.-S. Luo, “A priori derivation of the lattice Boltzmann equation,” Physical Review E, vol. 55, no. 6, pp. R6333–R6336, 1997. View at Publisher · View at Google Scholar · View at Scopus
  22. X. He and L.-S. Luo, “Theory of the lattice Boltzmann method: from the Boltzmann equation to the lattice Boltzmann equation,” Physical Review E, vol. 56, no. 6, pp. 6811–6817, 1997. View at Publisher · View at Google Scholar · View at Scopus
  23. Z. Guo, C. Zheng, and B. Shi, “Discrete lattice effects on the forcing term in the lattice Boltzmann method,” Physical Review E, vol. 65, no. 4, Article ID 046308, 6 pages, 2002. View at Publisher · View at Google Scholar · View at Scopus
  24. J. M. Buick and C. A. Created, “Gravity in a lattice Boltzmann model,” Physical Review E, vol. 61, no. 5, pp. 5307–5320, 2000. View at Publisher · View at Google Scholar · View at Scopus
  25. S. Chen, D. Martínez, and R. Mei, “On boundary conditions in lattice Boltzmann methods,” Physics of Fluids, vol. 8, no. 9, pp. 2527–2536, 1996. View at Publisher · View at Google Scholar · View at Scopus
  26. M. R. Swift, E. Orlandini, W. R. Osborn, and J. M. Yeomans, “Lattice Boltzmann simulations of liquid-gas and binary fluid systems,” Physical Review E, vol. 54, no. 5, pp. 5041–5052, 1996. View at Publisher · View at Google Scholar · View at Scopus
  27. A. K. Gunstensen, D. H. Rothman, S. Zaleski, and G. Zanetti, “Lattice Boltzmann model of immiscible fluids,” Physical Review A, vol. 43, no. 8, pp. 4320–4327, 1991. View at Publisher · View at Google Scholar · View at Scopus
  28. D. Grunau, S. Chen, and K. Eggert, “A lattice Boltzmann model for multiphase fluid flows,” Physics of Fluids A, vol. 5, no. 10, pp. 2557–2562, 1993. View at Publisher · View at Google Scholar · View at Scopus
  29. M. Do-Quang, E. Aurell, and M. Vergassola, “An inventory of lattice Boltzmann models of multiphase flows,” Tech. Rep., Parallel and Scientific Computing Institute, Department of Scientific Computing, Uppsala University, Uppsala, Sweden, 2000. View at Google Scholar