J. Cent. South Univ. Technol. (2011) 18: 1434-1440
DOI: 10.1007/s11771-011-0858-4
Mathematical model for coupled reactive flow and solute transport during heap bioleaching of copper sulfide
YIN Sheng-hua(尹升华)1, 2, WU Ai-xiang(吴爱祥)1, 2, LI Xi-wen(李希雯)1, WANG Yi-ming(王贻明)1
1. Key Laboratory of Ministry of Education of China for Efficient Mining and Safety of Metal Mines,
University of Science and Technology Beijing, Beijing 100083, China;
2. School of Civil and Environment Engineering, University of Science and Technology Beijing,Beijing 100083, China
? Central South University Press and Springer-Verlag Berlin Heidelberg 2011
Abstract: Based on the momentum and mass conservation equations, a comprehensive model of heap bioleaching process is developed to investigate the interaction between chemical reactions, solution flow, gas flow, and solute transport within the leaching system. The governing equations are solved numerically using the COMSOL Multiphysics software for the coupled reactive flow and solute transport at micro-scale, meso-scale and macro-scale levels. At or near the surface of ore particle, the acid concentration is relatively higher than that in the central area, while the concentration gradient decreases after 72 d of leaching. The flow simulation between ore particles by combining X-ray CT technology shows that the highest velocity in narrow pore reaches 0.375 m/s. The air velocity within the dump shows that the velocity near the top and side surface is relatively high, which leads to the high oxygen concentration in that area. The coupled heat transfer and liquid flow process shows that the solution can act as an effective remover from the heap, dropping the highest temperature from 60 to 38 °C. The reagent transfer coupled with solution flow is also analyzed. The results obtained allow us to obtain a better understanding of the fundamental physical phenomenon of the bioleaching process.
Key words: copper sulphide; heap bioleaching; leaching reaction; solution flow; solute transport
1 Introduction
Bioleaching is found to have important applications in the recovery of copper, gold, silver, and uranium from their ores because of its simple equipment, low investment and operation cost, environmental advantages, and reasonable yields over a period of recirculation [1-2]. Copper obtained in this way accounts for more than 20% of the total world production [3-4]. Heap leaching is a mineral processing technology where large piles of crushed or run-of-mine rocks are leached with various chemical solutions to extract valuable minerals [5]. It is very difficult to evaluate the impact of process design and control systems on the operation of the heap due to the complicated heap structure and physicochemical reactions [6]. However, mathematic models can be used to simulate, optimize and improve the heap leaching operations with respect to lixiviant consumption, overall recovery and efficiency.
There are a number of factors to be accounted for in the heap bioleaching modelling, including the metal extraction, liquid and air flow, solute transport, bacterial action, and heat balance. In some models, the kinetics of ore leaching is described using the shrinking-core model (SCM). The SCM was applied to the leaching process of oxide copper ores because the dissolution reactions are rapid and can be ignored in the rate determining sequence of reaction and transport steps [7]. PAUL et al [8] constructed a computer model for bacterial ferric sulfate leaching of the major sulfides found in porphyry copper deposits. PANTELIS et al [9] presented a model of oxidation of mining heaps containing pyrite, and examined the impact of irrigation rate on the oxidation rate of the heap. To simulate the transport of solutes through the flowing channels and the stagnant pores of an unsaturated heap, three models (mixed side-pore diffusion and profile side-pore diffusion with uniform or distributed pore lengths) are derived by BOUFFARD and DIXON [10]. CASAS et al [11] used the model of Michaelis-Menten to calculate the oxidation rate by bacteria attached onto the ore surface, and assumed that a constant number of bacteria homogeneously distributed through the bed. DIXON [12] carried out an analysis of heat conservation during copper sulphide heap leaching, and postulated a new mode of heap behaviour called evaporative autocatalysis. DE ANDRADE LIMA [13] presented a model to simulate the transient evolution of the dissolved chemical species in the heap and column isothermal leaching processes. CARIAGA et al [14] presented the numerical solution of a two-dimensional partial differential equation model for a two-phase flow in a porous media representing the leaching process. CROSS et al [15] developed a framework for multi-phase flows, gas-liquid-solid, in reactive porous media, capturing several key chemico-physical–bacterial phenomena and their interactions. Unsaturated flow of liquid in a bed of uniform and spherical ore particles was studied numerically and experimentally by SHEIKHZADEH et al [16]. An unsteady and two-dimensional model was developed on the basis of the mass conservation equations of liquid phase in the bed and in the particles. The model equations are solved using a fully implicit finite difference method giving the distribution of the degree of saturation in the particles and in the bed and the vertical velocity of flow in the bed, as well as the effect of periodic infiltration on the above distributions. A two-dimensional dynamic model for bio-leaching of secondary copper minerals from a pile was developed by SIDBORN et al [17]. The set of partial differential equations was solved using commercial software. The governing equations for fully coupled flowing-reaction-deformation behaviour with mass transfer in heap-leaching were developed by WU et al [18]. The model equations were solved using an explicit finite difference method under the conditions of invariable application rate and constant hydraulic head.
The multi-process of heap leaching has been tried individually. However, the physical, chemical, hydraulic, and metallurgical processes are systematically linked. Understanding the interactive and coupled processes within the ore bed is important to enhance the performance of heap leaching from design and operating consideration. Based on the achieved results, the coupled process within the heap is investigated from three scales: 1) micro-scale model for lixiviate transport in particle, 2) meso-scale model of liquid flow between particles by utilizing the X-ray CT technology, and 3) macro-scale models for coupled process of liquid and gas flow, solute and heat transport in heap.
2 Model development
A heap is usually up to several kilometers long, and the edges of the heap are sloped with a flat top. The transversal cut of the heap is assumed to be a 2D uniform porous medium of trapezoidal shape. The length of the up and bottom boundaries are 40 and 80 m, respectively, while the height of the trapezoidal is 40 m. Governing equations for the leaching kinetics, solution and gas flow, solute transport are derived based on the spherical, cubic and rectangular representative elementary volume (REV), respectively.
2.1 Chemical reactions
In acidic solution, chalcopyrite is oxidized mainly by ferric ions according to the following reaction:
(1)
In addition to this reaction, the following reaction occurs in the presence of dissolved oxygen:
(2)
The fractional copper extracted by dissolved oxygen and bacterial can be described by the Michaelis-Menten form:
(3)
where β is the stoichiometric factor, ρo is the density of the ore matrix (kg/m3), G is the copper grade, X is the bacteria population (m-3), Vm is the maximum specific respiration rate of bacteria (kg), CL is the concentration of dissolved oxygen (kg/m3), and KM is the Michaelis- Menten half growth rate constant (kg/m3).
The spherical, porous and uniform REV is assumed to represent the ore particle of radius R forming the ore bed. According to Fick’s law of diffusion, the mass balance for reagent H2SO4 within the spherical REV takes the form of a continuity equation [8]:
(4)
where Deff is the effective diffusivity of reagent within the particle pores (m2/s), C(H2SO4) is the concentration of acid in solution (kg/m3), r is the particle radius (m), εo is the porosity of the particle, R(H2SO4) is the creation rate of sulfuric acid, ke is the first-order kinetics constant including the information regarding the particle and its reaction, which is to be determined experimentally (m3/(kg·s)), and S(Cu) is the content of copper associated with the solid phase.
2.2 Liquid flow
The cubic REV includes only fluid and is consistent with the continuum model normally employed in analyses of fluid flow. The stress at microscopic REV can be completely specified by nine components of the stress tensor τ. Consider an infinitesimal rectangular parallelepiped with faces perpendicular to the coordinate axes. On each face, there is a normal stress and a shear stress, which can be further resolved into two components in the directions of the axes. The diagonal elements τ11, τ22 and τ33 of the stress matrix are the normal stresses, and the off-diagonal elements are the tangential or shear stresses. The cube may be imagined to shrink to a point. The equation of motion for liquid flow within heap is described by the following equation [1]:
(5)
where ρL is the density of the liquid (kg/m3), u is the liquid flow velocity (m/s), p is liquid pressure, ? is liquid viscosity, and i and j are the force directions.
2.3 Air flow
The air flow through the bed is described by Darcy’s Law, the well known porous media fluid flow equation. Similar to the solution flow, the momentum equations of air flow are described by the N-S equations, with a reduced air space for flow based on the porosity and a resistivity to the flow due to the porous media. The steady state flow is described by the equation of continuity and the steady state N-S equations given by [11]:
(6)
where ρa (g/m3), pa (N/m2), ua (m/s), Ba (n/m3) and ?a (kg?s/m) are the air density, air pressure, air flow velocity, air body forces and air flow viscosity, respectively.
2.4 Solute transport
During heap leaching process, the barren solution, mainly containing H+, Fe2+, 
and bacteria, infiltrates the heap from the top, soaking the ore and leaching the metal into solution. The pregnant solution rich in Cu2+ is collected at the base of the heap for subsequent electrowinning. The macroscopic REV of the heap comprises of a solution, gas, and solid particles that make up the porous medium. The size of the macroscopic REV must be chosen sufficiently large so that its porosity is representative of the heap as a whole.
During migration through the porous media (heap), sulfuric acid reacts with the chalcopyrite. Acid consumption is assumed to be proportional to its concentration. The mass balance equations for the transport of H2SO4 can be written as [16]
(7)
where εL is the liquid fraction of ore bed, DL is the diffusion (dispersion) coefficient for the species in liquid, and v is the Darcy’s velocity of liquid (m/s).
2.5 Energy balance
There are several mechanisms by which heat may be transferred in the bed. The most important things are the transport by the downward liquid flow, the gas flow through the bed, including the effect of water evaporation, and the heat generated by the leaching reactions. Local thermal equilibrium between solid, liquid, and gas phases is assumed. Considering these assumptions, the equation for the energy balance of the bed may be written as [12]
(8)
where cp,B is the mean heat capacity of the ore bed, kB is the thermal conductivity of the ore bed, ΔHR denotes the heat of reaction per mineral dissolved, cp,g is the heat capacity of gas, and cp,L is the heat capacity of liquid.
2.6 Boundary and initial conditions
These equations consider the mass balance of solute, momentum balance of solution and air flow. The initial conditions for the model are ambient oxygen concentration and temperature, and zero concentration of sulfide acid, ferric and copper ions. The pore of the bed is filled by air before the irrigation begins.
Neither solution flow nor air flow takes place through the particle surface boundary. For the reason of simplification, only H2SO4 diffusion is considered on this boundary. The initial concentration of is set to 0.1 g/L. A saturated solution flow is assumed in the microscope REM. No air flow takes place, and solute transport in this area is ignored. The process of solution flow, air flow, diffusion of H2SO4, O2 and heat has been simulated within the heap geometry.
For up boundary of the heap:
, 
, 
,
, 
For left and right boundaries of the heap:
, 
, 
,
, 
For bottom boundary of the heap:
, 
, 
,
, 
The parameters needed for these simulations are listed in Table 1. The system of differential equations is solved using the simulation tool COMSOL Multiphysics. When the governing equations and boundary conditions are properly defined, the implementation of the equation system into COMSOL Multiphysics is straightforward.
Table 1 Overview of all model parameters used in these simulations
3 Result and discussion
3.1 Acid transfer process at level of particle scale
Diffusion within the particle is the dominant mechanism during coarse ore bioleaching processes, because the metal deposit inside the ore particle could only be dissolved under the precondition of contacting with the reagent. The dissolved reagents, such as H2SO4, move slowly toward the center of the ore through the pore networks. The oxidation of a spherical particle occurs by the reagent diffusion originating from the outer surface, and the oxidation rate of a particle is therefore highly dependent on the particle size and reactant concentration.
Figure 1 shows the acid concentration inside an ore particle from 18 to 180 d of leaching time. Each line represents an increment of 18 d. The purpose of this figure is to emphasize the variation of acid concentration within ore particle and the long time scale for acid diffusion. Within 180 d of simulation, the acid concentration is very different at various portions of the particle. The highest acid concentration near the particle surface is apparent, and the sulfide acid at the surface is maintained at a relatively constant value, which is equal to that of bulk solution. The oxide, secondary sulfide minerals and alkaline gangue are dissolved rapidly at the beginning of leaching process [19]. At or near the surface of an ore particle, the acid concentration is relatively high due to the development of cracks and fissures as the leaching proceeds. Additionally, the particle surface is the nearest place from the bulk solution. The diffusive constraint is very small, so it shows rapid transfer speed of acid in this region.
Fig.1 Acid concentration distributions within ore particle
From the topmost line, it is clear that the center of the ore has not reached steady state at 180 d. However, the acid concentration at the center part of the particle increases slightly even after 180 d of simulation, as the less-porous ore cores with low permeability prevent the liquid from penetrating. The rate of decrease in the acid concentration slows down considerably from 72 d. A gentle concavity in the acid concentration is observed between 72 and 180 d. This is mostly due to the diffusion limitation in the particle core.
Hence, the simulation results suggest that the minerals at and near the surface of the particle would be leached quickly, while those in the particle core region dissolve relatively slowly due to the retarded diffusion of acid.
3.2 Meso-scale flow between ore particles
The investigation of flow phenomena at the meso-scale level could be an important step in scientific discovery for heap leaching research. Combining the X-ray computed tomography (CT) with robust modeling tools, it is now possible to show the pore scale fluid flow modeled by the Navier-Stokes equation. As shown in Fig.2, the transversal cut image of a small amount of compacted ores is obtained based on the X-ray CT experiment. The CT image is then transferred to DXF files (as shown in Fig.3), which is finally imported into COMSOL Multiphysics.
Fig.2 Cross-sectional image of compacted ores
Fig.3 Particle boundary shown by DXF file of CT image
It is assumed that saturated liquid flow exists within a certain part of region consisting of compacted small particles. Figure 4 shows the solution predicted with a Navier-Stokes analysis for the velocities in the pore spaces among those particles. The figure reveals that the solution applied to the heap surface migrates through the porous matrix of the staked ore in some paths, primarily in the downward direction. Not all liquid is in motion; however, a significant portion remains trapped in pores and crevices between the ore particles. It depends strongly on the distribution of solution flow channels and the relative size of stagnant zones. The velocities are higher in the narrow pores and the highest velocity within the simulation area reaches 0.375 m/s. The fluid velocities tend to decrease in stretches where the cross-sectional area for the flow increases.
3.3 Air flow coupled with O2 transfer within heap
Solution in the heap will not fully saturate the porous space among particles, except potentially near the bottom where poor drainage may cause a hold-up of liquid and consequent full saturation. The residual gas space provides an interconnected pore network for air movement. The air flow is mainly induced by air pressure and temperature gradients. The pressure and temperature gradients in the bed do not change much during the process, and as a consequence, the airflow profile does not change much with time. It is observed from Fig.5 that air flows predominantly from the side of the heap through hot region and then up through the top surface. In doing so, it penetrates, in some places, the cooler interior regions. The closer the position apart from the bed surface (sides and top), the faster the gas flows. The lowest air flow rates occur in the central part of the bed, where air flows mainly in the vertical direction. The highest air flow rates are reached under the bed slope and closer to the bottom.
Fig.4 Velocity field plot showing velocities between ore particles
Fig.5 Air velocity fields in dump
Oxygen, carried by the air flow, can penetrate and the fuel can be oxidized in the cooler parts of the heap. Oxygen is transported largely by convection and diffusion from the atmosphere/heap interface, mainly at the slope sides of the ore bed. The rate of leaching reaction is determined by the amount of oxygen available if the ferrous ions are present. The concentration of oxygen in the ore bed is shown in Fig.6. The oxygen concentration is expressed as a percentage with respect to the atmospheric oxygen concentration. It is indicated that the concentration is large along the slope boundary where oxygen flows into the bed due to the natural convection. A minimum in the oxygen concentration is found at an intermediate height in the central zone of the bed. At this point, the oxygen is practically depleted from the convective air flow. The oxygen concentration increases above this zone, mainly due to direct oxygen diffusion from the atmosphere. From this simulation result, a conclusion could be made that convection is a significant air transport mechanism in ore heap, and diffusive transport dominates over convection on the top part of the heap.
Fig.6 Oxygen content contour in dump
3.4 Coupled heat transfer and liquid flow process
The chemical reactions involved in the exothermic oxidation of sulfide minerals result in the generation of substantial amounts of heat. Under the assumption that gas and liquid flows carry no heat out of the bed, the heat is transported only through conduction, transpiration and radiation at the heap surface and slope. In the central part of the bed, the temperature reaches as high as 60 °C due to the heat accumulation and slow heat transportation.
As shown in Fig.7, water draining at the base of the heap can act as an effective heat remover from the heap. In general, the temperature inside the bed increases with increasing the depth, because the leaching solution accumulates the thermal energy liberated by oxidation reactions as it flows downward. The highest temperature zone is located at the bottom, where a large liquid flow path is found. The highest temperature drops to 38 °C from 60 °C when considering the heat transfer by solution flow. The temperature distribution and average temperature are also altered by the liquid flow. Liquid flow affects the oxidation rate in two ways. One is that the water flow transports heat, and reduces the temperature gradient. The other is that the cooling effect of the water flow can prevent the switch-off of the bacterially catalyzed oxidation process by limiting the temperature rise in parts of the heap.
3.5 Reagent transfer coupled with liquid flow
Figure 8 shows the sulfide acid content at the beginning of leaching. Initially, the acid content in the ore bed is negligible. The sulfide acid is offered into the ore bed continuously by sprinkling at the top of the heap. Acid transfers downwards with the flow of liquid and diffuses into the ore particles, and reacts with the copper mineral. Acid is reduced in this reaction, so the acid content of front flow is relatively low. As the sprinkling goes on, a larger quantity of sulfide acid is transported down by the liquid flow. Eventually, the acid content of the whole heap will reach an almost steady state condition, equal to that of supplied solution.
Fig.7 Temperature profile coupled with liquid flow
Fig.8 Sulfide acid content profile coupled with liquid flow
4 Conclusions
1) Mineral at and near the surface of the particle would be leached quickly due to the high concentration of lixiviant in these regions, while that in the particle core region dissolves relatively slowly due to the retarded diffusion of acid. So, the oxidation of the coarse particle is on the time scale of months, or even several years.
2) It is possible to realize the visualization of velocity field between the ore particles by combining the X-ray CT technology with robust modeling tools. The highest velocity is observed in the narrowest pore throat. The air flow through the porous heap is described by N-S equation. The closer the position apart from the bed surface (sides and top), the faster the gas flows. The lowest air flow rate occurs in the central part of the bed, where air flows mainly in the vertical direction.
3) Liquid flow affects the oxidation rate in two ways. One is that the water flow transports heat, and reduces the temperature gradient. The other is that the cooling effect of the water flow can prevent the switch- off of the bacterially catalyzed oxidation process by limiting the temperature rise in parts of the heap. Initially, acid and other reagents transfer downwards with the flow of liquid, diffuse into the ore particles, and react with the copper mineral. Acid is reduced in this reaction, so the concentration of front flow is relatively low.
References
[1] MOUSAVI S M, JAFARI A, YAGHMAEI S, VOSSOUGHI M, SARKOMAA P. Computer simulation of fluid motion in a porous bed using a volume of fluid method: Application in heap leaching [J]. Minerals Engineering, 2006, 19(10): 1077-1083.
[2] PRADHAN N, NATHSARMA K C, SRINIVASA K, SUKLA L B, MISHRA B K. Heap bioleaching of chalcopyrite: A review [J]. Minerals Engineering, 2007, 21(5): 355-365.
[3] WATLING H R. The bioleaching of sulphide minerals with emphasis on copper sulphides—A review [J]. Hydrometallurgy, 2006, 84(1/2): 81-108.
[4] WU Ai-xiang, YIN Sheng-hua, WANG Hong-jiang, QIN Wen-qing, QIU Guan-zhou. Technological assessment of a mining-waste dump at the Dexing copper mine, China, for possible conversion to an in situ bioleaching operation [J]. Bioresource Technology, 2009, 100(6): 1931-1936.
[5] THIEL R, SMITH M E. State of the practice review of heap leach pad design issues [J]. Geotextiles and Geomembranes, 2004, 22(6): 555-568.
[6] BENNETT C R, MCBRIDE D, CROSS M, GEBHARDT J E, TAYLOR D A. Simulation technology to support base metal ore heap leaching [J]. Mineral Processing and Extractive Metallurgy, 2006, 115(1): 41-48.
[7] BARTLETT R W. Simulation of ore heap leaching using deterministic models [J]. Hydrometallurgy, 1992, 29(3): 231-260.
[8] PAUL B C, SOHN H Y, McCARTER M K. Model for ferric sulfate leaching of copper ores containing a variety of sulfide minerals. Part 1: Modeling uniform size ore fragments [J]. Metallurgical Transactions B, 1992, 23(5): 537-548.
[9] PANTELIS G, RITCHIE A I M, STEPANYANTS Y A. A conceptual model for the description of oxidation and transport processes in sulphidic waste rock dumps [J]. Applied Mathematical Modelling, 2002, 26(7): 751-770.
[10] BOUFFARD S C, DIXON D G. Investigative study into the hydrodynamics of heap leaching processes [J]. Metallurgical and Materials Transactions B: Process Metallurgy and Materials Processing Science, 2001, 32(5): 763-776.
[11] CASAS J M, MARTINEZ J, MORENO L, VARGOS T. Bioleaching model of a copper-sulphide ore bed in heap and dump configurations [J]. Metallurgical and Materials Transactions B: Process Metallurgy and Materials Processing Science, 1998, 29(4): 899-909.
[12] DIXON D G. Analysis of heat conservation during copper sulphide heap leaching [J]. Hydrometallurgy, 2000, 58(1): 27--41.
[13] DE ANDRADE LIMA L R P. A mathematical model for isothermal heap and column leaching [J]. Brazilian Journal of Chemical Engineering, 2004, 21(3): 435-447.
[14] CARIAGA E, CONCHA F, SEPULVEDA M. Flow through porous media with applications to heap leaching of copper ores [J]. Chemical Engineering Journal, 2005, 111(2/3): 151-165.
[15] CROSS M, BENNETT C R, CROFT T N, MCBRIDE D, GEBHARDT J E. Computational modeling of reactive multi-phase flows in porous media: Applications to metals extraction and environmental recovery processes [J]. Minerals Engineering, 2006, 19(10): 1098-1108.
[16] SHEIKHZADEH G A, MEHRABIAN M A, MANSOURI S H, SARRAFI A. Computational modelling of unsaturated flow of liquid in heap leaching-Using the results of column tests to calibrate the model [J]. International Journal of Heat and Mass Transfer, 2005, 48(2): 279-292.
[17] SIDBORN M, CASAS J, MARTINEZ J, MORENO L. Two- dimensional dynamic model of a copper sulphide ore bed [J]. Hydrometallurgy, 2003, 71(1-2): 67-74.
[18] WU Ai-xiang, LIU Jin-zhi, TANG Ling-yan. Simulation of coupled flow in reaction-deformation with mass transfer in heap leaching processes [J]. Applied Mathematics and Mechanics, 2007, 28(3): 327-335.
[19] WU Ai-xiang, YIN Sheng-hua, QIN Wen-qing, LIU Ji-shan, QIU Guan-zhou. The effect of preferential flow on extraction and surface morphology of copper sulphides during heap leaching [J]. Hydrometallurgy, 2009, 95(1/2): 76-81.
(Edited by YANG Bing)
Foundation item: Projects(50934002, 51104011) supported by the National Natural Science Foundation of China; Project(IRT0950) supported by Program for Changjiang Scholars and Innovative Research Team in Chinese University; Project(20100480200) supported by China Postdoctoral Science Foundation
Received date: 2010-08-23; Accepted date: 2011-03-23
Corresponding author: YIN Sheng-hua, PhD; Tel: +86-10-62334680; E-mail: csuysh@yahoo.com.cn
