- Email: [email protected]
A model for the leaching of major inorganics from retorted Rundle oil shale
War. Res. Vol. 27, No. 2, pp. 277--291,1993 Printed in Great Britain.All rightsreserved
0043-1354193 $5.00+0.00 Copyright © 1993 Pergamon Press Ltd
A MODEL FOR THE LEACHING OF MAJOR INORGANICS FROM RETORTED RUNDLE OIL SHALE A. A. KROL, P. R. F. BELL@ and P. F. GREENFIELD') Department of Chemical Engineering, The University of Queensland, St Lucia, Qld 4072, Australia (First received August 1991; accepted in revised form July 1992) Abstract--A model is developed to describe the leaching of major soluble inorganics (the sulphate and chloride salts of calcium, magnesium, sodium and potassium) from processed oil shale sourced from the Rundle resource, Queensland, Australia. The model accounts for rapid dissolution of readily soluble salts, the formation of non-charged ion pairs in solution, cation exchange and slow dissolution/precipitation of calcium sulphate (anhydrite/gypsum). I.~,achate results from an experiment involving irrigation of an initially dry column of processed shale are used to validate the model. A piston model is used to account for unsaturated water flow, and this is expressed in the water quality model by a moving bottom boundary condition. Tracer experiments show that a good basis for the leaching model is the classic convection-dispersion model. A moving node finite element method is used to obtain a numerical solution to this leaching model. The model provides a good description of the leaching of each inorganic species. Deviations between model predictions and experimental results are greatest at slower flowrates and for monovalent cations. It is concluded that the main chemical mechanisms controlling leaching of major inorganics from processed shale have been identified. Key worda--Rundle, oil shale, retorted shale, leaching, chemical equilibria, modelling, finite element, unsaturated flow, salt transport
INTRODUCTION The exploitation of oil shale resources will generate a large quantity of solid waste consisting o f processed (or spent) shale, low grade unprocessed shale and overburden. Retorting operations change the physical and chemical characteristics o f the shale and increase the potential mobility of salts (Slettevold et al., 1978; Stollenwerk, 1980). Lcachates are saline and generally alkaline in nature. There are two broad environmental concerns associated with such leachate: • revegetation of dumps must not be handicapped by concentration of salts in the surface zone due to an evaporation-driven upward flux; and • leachate containing toxic species must not enter ground or surface water in quantifies likely to pose a hazard. The Rundle oil shale resource is located on the central Queensland coast. Spent shale used in this study was part of a bulk sample extracted from the resource in the late 1970s which was retorted in Germany using the Lurgi-Ruhrgas process. Leachate studies on this material have been in progress since that time. The earliest studies concentrated on establishing the characteristics of leachate generated by irrigation of laboratory columns (Bell et aL, 1982), These were supplemented by studies designed to determine the mechanisms controlling leaching of
inorganic constituents (Bell et al., 1986; Krol et al. 1986, 1988). The work described in this paper involves development of a model to predict the quality of leachate generated from Rundle spent shale. Previous water quality models for oil shale leachate have been developed for material extracted from the Green River resource in the US. Rundle shale is softer and more similar to a claystone than US shale (Rowley and Brown, 1981). Most of these previous models have concentrated on mass transfer aspects of the leaching process rather than the chemistry involved (e.g. Kuo et al., 1979; Ramirez et al., 1983). They have also tended to neglect leaching under conditions of unsaturated water flow. An exception is the model proposed by Margheim (1975) which considered three solute phases: • a solution phase of Ca, Mg and Na salts; • an exchangeable phase consisting of Ca 2+, Mg 2+ and N a + adsorbed to the spent shale; and • a soluble salt, CaSO4. However, Margheim's model did not include any correcfion for the effect of ionic activity coefficients and its success in fitting experimental data was largely due to curvefitting of model parameters. In our studies it was decided a priori that the developed model should as far as possible be calibrated by independent experiments designed to determine which mechanisms should be included in the 277
A. A. KROLet al.
278
model and provide model parameters. The number of parameters obtained by curvefitting of column leachate data was to be minimized and such data were in the main to be used for model validation. The developed model was thus a test of whether the chemical mechanisms controlling the leaching of major inorganic constituents (calcium, magnesium, sodium, potassium, chloride and sulphate) from the shale had been correctly identified. A simple physical situation was chosen for the validation experiment. This involved irrigation of an initially dry column of spent shale with a constant flux of distilled water. The physics of the column experiment and previous experimental studies (Bell et al., 1986; Krol et al., 1986, 1988) showed that it was necessary to consider
inclusion of at least the following mechanisms in the model of this experiment: (1) Infiltration of water into a dry column followed by steady state flow. (2) Hydrodynamic dispersion and diffusion, and possibly mass transfer between "immobile" and "mobile" zones of water. (3) Rapid dissolution of most salts and slower dissolution of CaSO4. (4) Ion exchange of cations. (5) Formation of the ion pairs CaSO ° and MgSO °. (6) The effect of changes in ionic strength on activity coefficients and thus on the reactions (3)-(5) above.
A
constant flow we steel plate ~brane :t to hollow
s ~ t l o n In c ~ t ~ t sis via small holes
iy E gamma ray 8ource
bl mtforTn
ked with shale
wring temdometew
i x sheet rio
KL
leachata take off
i
Fig. I. One metre leachate column.
J
to vacumn controNer
279
Model for leaching of major inorganics
Experiments D E F
Superficial velocity, v (cm/h) 0.406 0.694 0.148
Table i. Column operating conditions Mean dry bulk Average pore Water content, 0 densities, Ps volume (PV) (cm3/cm3) (g/cm3) (litres) 0.585 0.888 8.44 0.616 0.887 8.96 0.551 0.881 8.00
The colum experiments which provided the model validation data are described below. Tracer studies were carried out on these columns in order to select and parameterize the hydrodynamic model on which the overall leaching model was based. The experimental section is followed by a description of the development of the model, a summary of the numerical solution scheme employed and a comparison between model predictions and experimental results. EXPERIMENTAL
Tracer volume, fraction of PV 0.684 0.792 0.823
Leachate volume collected(PV) (litres) 79.3 (9.4) 78.1 (8.7) 88.5 (11.1)
D, E and F) were carried out under different hydraulic conditions. The operating conditions of the columns and the tracer studies are summarized in Table 1. COLUMN HYDRODYNAMICS
Water flow during infiltration The water content profiles during infiltration into the dry columns were all very steep. The profiles measured on column E are presented in Fig. 2. It was necessary to develop a model for such water flow in order to provide spatial water content and velocity profiles for use by the leaching model. The usual approach to such a situation is to use an algebraic or numerical solution for the coupling of the unsaturated form of Darcy's equation with the continuity equation (e.g. Van Genuchten, 1980). This method was not used here because a situation involving
The leaching experiments were carried out in a constant temperature room (25 + 2°C) using a I m long perspex column of 12cm square internal cross section (Fig. I). The column was repacked for each experiment with airdried spent shale to a consistent density of 890 + 20 kg/m 3. Distilled, deionized water was added at a constant rate to the top of the column by a peristaltic pump. An even distribution of water was achieved by placing an 8/~m steep wetting fronts is difficult to solve numerically. cellulose acetate membrane between the top water reservoir F u r t h e r , it h a s been established t h a t with the excepand the shale. The bubbling pressure of this membrane was tion o f chloride the range o f flow h y d r o d y n a m i c s such that the water reservoir remained full throughout the experiment and water was preferentially supplied to any less represented b y experiments D, E a n d F has little saturated parts of the shale surface, i.e. where the capillary effect o n the leaching b e h a v i o u r o f m a j o r inorganics pressure was highest. In all experiments the water was ( K r o l et al., 1988). Hence, a simple " p i s t o n " w a t e r supplied at a rate below the saturated hydraulic conductivity flow m o d e l was used to a p p r o x i m a t e the w a t e r conof the shale (0.8cm/h), and thus the flow through the tent profiles by a m o v i n g step function. Such a column remained unsaturated. The progress of the wetting front was followed using m o d e l was used by H i l d e b r a n d a n d H i m m e l b l a u a y-ray attenuation system (a modified GEC 026 nuclear (1977) to predict the t r a n s p o r t o f nitrate d u r i n g density gauge containing a 780 MBq mCs source). When the infiltration into a sand column. T h e required i n p u t wetting front reached the base of the column, leachate was d a t a is the e q u i l i b r i u m water c o n t e n t , 0~, for the t o p extracted by a vacuum equivalent to the capillary pressure b o u n d a r y infiltration rate, vi. T h e position o f the measured by the tensiometer at the top of the column. This procedure achieved even capillary pressure and water conwetting f r o n t at a n y time, t, for a c o l u m n o f length tent profiles down the length of the column. The volume, L, c a n t h e n be predicted using the equations given in pH and conductivity of the discharged leachate were monT a b l e 2. T h e application o f the p i s t o n flow m o d e l to itored and samples were analysed for calcium, magnesium, c o l u m n e x p e r i m e n t E is s h o w n in Fig. 2. T h e model sodium, potassium, chloride and sulphate (Krol et al., 1988). predicts the p o s i t i o n o f the wetting f r o n t to within Once the column hydrodynamics had reached a steady state (i.e. the flowrates into and out of the column were equal and the water content profile was even and conTable 2. Equations describing the piston model for water infiltration stant) they were investigated using a tritium tracer. Tritium into a column of shale is relatively non-adsorbing, does not change water flow properties and does not interfere with the leaching mechVi anisms of those inorganic components which were being (i.e. before the wetting front reaches the bottom of the column) monitored. The tracer was injected into the column as a finite pulse wetting front position, L~ ffi tv~ with a tritium activity of approx. 107 dps/l. During tracer Or experiments the discharged leachate was directed to an column profiles, 0 ~ z ~ L~: v ffi v~, 0 = Or automatic fraction collector which collected samples at set L,01L lation counter and a toluene-based scintillant as the countvi ing medium. The count rate was corrected on the basis of (i.e. after the wetting front has reached the bottom of the column) external standard ratios and tritium standards using the method of Fox (1976). L~ ffi L Further details of all experimental procedures are proO ~ z ~ L : v f f i v l , O=Oj vided by Krol (1985). Three leaching experiments (labelled
280
A.A. I ~ o L et al.
0.7] --I
0.6l
II 8 3 hours "E "E
O'S l
l I
I
0.4
!
@
I I I I
0.3 O
0
" ~
(first leachate)
0.2
5 4 hours
a '~,30.5 hours ~,~
,1 1.3 hours
J
KEY ~, ~, experiment
L
-----
/
! I I
0.1
piston model
I 0
O
10
20
30
40
50
60
7O
80
90
100
Distance from Top of Column Z (cm)
Fig. 2. Application of piston water flow model to column E.
5 cm, an accuracy which was considered adequate for our purposes. In the column experiments a maximum of approx. 2% of the total shale weight was leached. The impact of such leaching on porosity and permeability, and hence on water flow characteristics, is regarded as insignificant given the level of accuracy of the overall modelling procedures. Tracer studies Tracer breakthrough results were interpreted using CFITIM, a non-linear least squares fitting programme developed by Van G-enuchten (1981). Two models were fitted to the data, a convectiondispersion model (model A) and a physical nonequilibrium model (model B). Model B assumes the presence of "mobile" and "immobile" zones of water with diffusive mass transfer occurring between the zones, as described in detail by Van G-enuchten and Wierenga (1976). Both models assume that any adsorption of the tracer is rapid relative to the water flow rate. The one-dimensional variation of solution and solid concentrations of a species (c and s, respectively)
with time (t) and distance (z) are described by the following equation: 00c Os ~2c Oc ~ + pa ~ -- K ~zz-- V-~z
(!)
K is the dispersion/diffusion coefficient. By making use of the linear Freundlich isotherm (s = kc) and suitable substitutions, model A can be expressed in the dimensionless form used by CFITIM (Van Genuchten, 1981): R ~c = I ~2c
0x
B 0X2
~c
(2)
0X
The dimensionless variables are defined as follows: vL ; B --- -~-
C = ~ oc .,
p 0k , "
vt
Z
x = Z ;
R=I+
(3)
~ = o-2"
The dimensionless form of model B can be similarly developed (Van G-cnuchten, 1981; Krol, 1985). Results obtained from fitting model A to tracer breakthrough data are given in Table 3. Figure 3
Table 3. Fit of convection--dispersion model (A) to tracer data Model parameters fitted by CFITIM
Calculated parameters
Experiment code
B
R
T
K (cm2/h)
Tracer recovery (%)
D E F
272+30 228+22 232+37
1.011 +0.004 1.003+0.005 1.004+0.007
0.661 +0.007 0.791 +0.006 0.828+0.010
0.15+0.02 0.30+0.03 0.066+0.011
97 100 101
Note: + errors are 95% confidence limits.
Model for leaching of major inorganics shows the model fit to the tracer breakthrough data from column E. Model A provided a good fit to the data in all cases. The fitted variables were the Brenner number, B, the retardation factor, R, and the input pulse size. Prediction of the latter effectively represents numerical integration of the breakthrough data, and comparison with the known input pulse size enabled the tracer recoveries (97-101%) to be calculated. As expected the tritium tracer showed negligible adsorption (R values close to 1). The calculated dispersion coefficients (K) show a positive correlation with irrigation rate, i.e. E > D > F. Model B also showed a good fit to the experimental data, indicating that zones of "immobile" water may be present within the column. Krol (1985) presents arguments which conclude that model B may in fact represent a better approximation of the hydrodynamics in the columns than model A. However, it was the convection-dispersion model (A) which was coupled with the piston water flow model (Table 2) to provide the hydrodynamic basis for the leaching model. Model A was selected because it requires less input data and it makes inclusion of a sophisticated chemistry model much simpler. MECHANISMS INCLUDED IN THE CHEMISTRY MODEL
281
the form of anhydrite dissolve more slowly at a rate which depends on solution concentrations. The batch experiments also showed that most readily soluble components are present on, or close to, the readily accessible surface of spent shale particles. It is assumed in the column leaching model that the rapidly dissolving components (Na, Mg, K, Cl and that Ca and SO4 which is not present as anhydrite) are instantaneously dissolved at the wetting front. This assumption is affected in the model by considering dissolution as a finite dispersive flux across the moving wetting front into the wetted region. The size of this flux equates with the quantity of "component i" which must dissolve at the wetting front, and determines the concentration gradient of i at the wetting front boundary (z --Lw):
Oci [
"~z ,,=z,,
= vPBS~
~
The value of K in this equation is assumed to be that determined in the tracer studies (Table 3). Dissolution of Ca and SO4 from slowly dissolving anhydrite is characterized by a kinetic model used previously for dissolution of gypsum (e.g. Keisling et al., 1978; Glas et al., 1979): OS16
p~-~- : Sk, O(cl ± c,s).
Dissolution of salts Batch experiments have shown that, when extraparticle mass transfer effects are minimized, approx. 90% of leachable Mg, Na, K and CI dissolves in the first hour (Krol et al., 1988). Ca and SO4 present in
Key Input tracer pulse output tracer, experiment output tracer, model
EXPERIMENT water content 0 . 6 1 6 om3/©m 3 Darcy velocity 0.694 cm/hr
1.0
........
~ 0.8
I
1
2"
(5)
The parameter k, is the reaction rate constant for dissolution of anhydrite and, as included in equation (5), implicitly includes surface area. Subscripted concentrations can be identified by reference to Table 4.
1.2
--•
(4)
"\
i0.6 L)
0.2
oJ0 Pore Volumes of Leachata Eluted since start of Tracer Injection
Fig. 3. Fit of convection-dispersion model to tracer breakthrough curve for column E.
A. A. K.ROLet al.
282
Table 4. Symbolsused to representindividualspecies Solution "Solid exchangeable" phase phase Solid "leachable"phase (mequiv/l) (mequiv/g) (mequiv/g) cl ct s, (as CaCI,)
Species Ca2+ Mg2+
c2
Na +
c3
K+ CISO2CaSO° MgSO°
c( cs c6 c,6 c26
c2
$16 (as CaSO4) S2
c3
s3
--
__
c( ----
The term (c~-cts) represents the "saturation deficit" of the solution with respect to anhydrite. The concentration of Ca t+ which would result in the solution becoming saturated with respect to anhydrite, C~s, is calculated using the solubility expression: K~'6 cts = m
s( ss s6 (exceptthat as CaSO() s,6 (as CaSOD
adsorbed species. The Gapon model used in this previous work is applied in the column leaching model: c0.5~, ~r
= - -
(6)
C6
~K = c3cl
c,~ c~-~
~K = ~ . K~6 =
.
YlY6
Speciation in solution The presence in solution of the uncharged ion pairs (CaSO ° and MgSO °) increases the solubility o f anhydrite and affects the ion exchange of Ca t+ and Mg 2+. The ion pair concentrations are calculated using their equilibrium relationships: C16
C! C6
(8)
C26
K~'6 = c2 c6"
(9)
The superscript "*" on the equilibrium constant identifies it as a value corrected for the effect of activity coefficients: Kt6 YiY6
K~'6 = - )'16
~i~26y2~/6
K * , ffi ~
(14)
(7)
The parameter S in equation (5) is a switch used in the numerical solution. It is set on (S = 1) or off (S --0) for a given time at each spatial position in the column, depending on whether anhydrite is present or not. Further details are provided by Krol (1985).
K*, = --
(13)
c?.~
where K ~ is the solubility constant for anhydrite, Kt6, adjusted by the solution activity coefficients of Ca t+ , ?l and SO~-, Y6: Kt6
(12)
c~°~
(10)
(11)
Equations (12)-(14) are written in terms of solution concentrations rather than activities. It was decided not to incorporate activity coefficients into the equations as the resultant improvement in accuracy is minimal when compared with the noted variability (standard deviations o f up to 80%) in the mean measured values of the selectivity coefficients (Krol et ai., 1986). However, a mean correction factor of 20% was applied to the measured values of ~K and ~K in order to account for the average effect of the difference in salt concentration between the experiments used to measure selectivity coefficients and the column leachate studies. These selectivity coefficients for monovalent-divalent cation are those most sensitive to the effect of activitys coefficients. The cation exchange capacity (CEC) was assumed to remain constant at the value measured for the spent shale on completion of the column leaching experiments (0.255 mequiv/g; Krol et al., 1988): CEC = c~ + c2 + c3 + c4.
(15)
Activity coefficients The activity coefficients of charged species were calculated using the extended Debye--Huckel equation. F o r uncharged species the following empirical relationship was used (Wigley, 1977):
T26
Cation exchange Krol et al. (1986) demonstrated that cation exchange on the shale is significant and that equilibrium is rapidly attained between the solution and
lOgl0 Yi= 0.1L
(16)
Summary of chemistry model The incorporation o f the above mechanisms into the leaching model results in the total concentration
Model for leaching of major inorganics of any component at a wetted part of the column being represented as follows: Ca ----Ca2+ + CaSO°
+ Ca2+-shale + CaSO4(s)
Mg = Mg2+ + MgSO°
+ Mg 2+-shale
Na = Na +
+ Na+-shale
K=K
+
283
assumption [equation (4)]. When the column is fully wetted (i.e./.~ = L) the wetting front boundary condition is replaced by the usual zero concentration gradient: ~(~C5 z
+ K +-shale
z'L
=
O.
(18)
The appropriate top boundary condition for irrigation with distilled, deionized water is the constant flux condition:
CI = CISO, = SO~- + CaSO ° + MgSO°
+ CaSO4(s)
In solution, transportable
az-
Associated with solid, not transportable
Rapidly dissolving solid species are not associated with the solid because such species are assumed to dissolve instantaneously on contact with the wetting front. THE COLUMN LEACHING M O D E L FOR C H L O R I D E
Chloride is the only component which is assumed both to dissolve rapidly and to be independent of the other ieachate components in the chemistry model. The leachate model for CI was thus developed separately from the general leachate model. It provides a check on the validity of three assumptions: (1) the applicability of the piston water flow model to water infiltration; (2) the applicability of the convection-dispersion model for column hydrodynamics; (3) the instantaneous dissolution of rapidly soluble salts at the wetting front.
Development of the chloride model
- -
(17)
el has no chemical reaction interactions with other species and thus equation (17) is also the transport model for el. The region of applicability of equation (17) (0 < z
Model predictions and discussion The prediction of the chloride leaching model is compared with discharged leachate concentrations for column E in Fig. 4(A). At this flow velocity (0.694 cm]h), the model predictions agreed very well with the experimental results. However, the goodness of fit deteriorated for the slower and more unsaturated flow conditions of columns F [(Fig. 4(B)] and D (Krol, 1985). In fact the profiles predicted by the model were relatively insensitive to flow
Table 5. Distribution of leachable solids on dry shale
Component
Total* concentration (mequiv/g)
Rapidly dissolving udts (mequiv/g)
Anhydrite (mequiv/g)
Ca Mg Na K CI SO.
0.3496 0.0589 0.0397 0.0147 0.0125 0.1954
0 0.0089 0.0097 0.0057 0.0125 0.0108
0.1846 ----0.1846
*grol (1985). w a ~/2--0
(19)
The model was solved using a moving node Galerkin finite element method based on that proposed by O'Neill and Lynch (1980). This method minimizes problems of numerical dispersion which would otherwise be very pronounced, particularly when the wetting front initially starts to move into the column of dry shale and the effective Brenner number (B) based on the wetted column length (vLw/K) is at a minimum. Integration in time was by a Crank-Nicbolson finite difference scheme. The finite element nodes were positioned at each time step by placing the last node at the wetting front and evenly spacing the remaining nodes over the wetted region. No attempt was made to optimize the position of the nodes to provide a more economical and accurate solution to the model because the full leachate model, incorporating all the major inorganic components, required a more evenly spaced mesh. A detailed outline of the development and solution of the nodal equations is given by Krol (1985). The chloride leaching model was solved for all columns (i.e. D, E and F). Required model input data are included in Tables 1 (o, Pn and v), 3 (K) and 5 (ss).
For CI, the convection--dispersion model [equation (1)] can be simplified to:
00c5 02c5 0c5 Ot = K~z2 v O'zz"
VCs)]=.0--0.
Exchange sites (estimated)* (mequiv/g) 0.165 0.05 0.03 0.01
284
A.A. KIgOLet al.
conditions, whereas the experimental results showed greater dispersion under slow flow conditions. A possible reason for this behaviour is a break-
down in the instantaneous dissolution assumption due to: • a greater proportion of immobile water being present as the flow becomes more unsaturated, with dissolution being much slower in immobile regions; and • the wetting front becoming less sharp, resulting in CI dissolution occurring over a finite bandwidth rather than at a sharp front.
Despite the short-comings of the chloride leaching model, it was used as the basis for modelling the leaching of other major inorganics. This was justified because Cl was rapidly leached from the columns during the passage of approximately one third of a pore volume of leachate, whilst the other inorganics remained present for many pore volumes (Krol e t al., 1988). Thus the shape of the wetting front and the column hydrodynamics were not the major factors determining the nature of these other leaching profiles, and their simulation using a simple hydrodynamic model with minimum data requirements was reasonable.
400 350 300 I 28(]
KEY
J
leo ~s~p~/mentnl I
200
i
15C
(A)
100
4oo~ ~'°I i aoo~
i
x.y
I. =--,-, II
1.°I
~I°°F"'\. O o
I I[
,-, I I
~
Eluted Leechate (lflre)
Fig. 4. Comparison of predicted and observed chloride results. (A) Column E, v = 0.694 cm/h and (B) column F, v = 0.148 cm/h.
Model for leaching of major inorganics THE COLUMN LEACHING MODEL FOR
ALL MAJOR INORGANICS
Model development
The transport equations for each of the other components (Ca, Mg, Na, K and SO4) cannot be solved independently, as was done for CI, because these equations are coupled through cation exchange, solution speciation and anhydrite dissolution reactions. The transport equations for these components arc:
Calcium:
oOC,+ Pa~t + o-Td Ocl~ + p.-Osl~ ~ Ot
02ct
= x ~
Ocl
- ~~
O"cl6
Oct6
+ r-~F~ - ~ - F ;
(20)
Magnesium: 0 0c2
O~ . ,, 0 % 02c~ 0c2 . 02e26 0c26 = K ~Tz2 -- V ~z + K-~z2 - v w
(21)
Sodium.
0 0c3
O~
O=c3
Oe3
O~
02c,
0c4
+ PB'~ =
(22)
Potassium:
00c,
+ p,-~ =
(23)
The concentration of sulphate is obtained by charge balance: (76 = C1 + (:72+ (73 + C4 - (75.
(24)
The appropriate top and bottom boundary conditions for equations (20)-(23) are similar to those for the chloride model. Equations (17) and (20)-(24) are non-linear, stiff and interacting. Expressions for the interacting chemistry, i.e. anhydrite dissolution [equation (5)], ion pair formation [equations (8) and (9)] and cation exchange [equations (12)-(14)], were substituted into the cation transport equations (20)-(23). The charge balance relationship [equation (24)] was used to eliminate SO~- to give a set of four equations. This equation set was then solved in a similar manner to the chloride model. However, because of the non-linear nature of the nodal equations, an iterative procedure was used to solve the problem. Full details are given by Krol (1985). Model predictions and discussion
The developed model was used to simulate the leaching of column E. The only model coefficient varied was the anhydrite dissolution rate constant, k,. Mass balance errors due to the numerical technique
285
were generally less than 1% of the total leachable concentration. A reduction in the time step or spatial grid spacing had negligible effect on the predicted results. The figures referred to below present total concentrations of Ca, Mg or SO4 (i.e. the free ion plus that portion existing as the appropriate ion pair). Sulphate. Figure 5 compares the predicted and experimental SO4 concentrations in the leachate. The SO4 profile is characterized by four regions: (I) a rapid initial fall; (II) a rise to a peak after about 2.4 I. (0.3 pore volumes of leachate); (III) a gradual drop in concentration until 50-601. (6.3-7.7 pore volumes) of leachate have been discharged; and (IV) another rapid fall. Some oscillation in the numerical results was experienced for region IV due to a numerical stability problem when the remaining anhydrite only occupies the last few nodes in the column. Regions I-IV are common to all of the figures and thus discussion of the mechanisms affecting the leaching of one component in a particular region is also relevant to the leaching of other components in that region. Figure 5 shows that model predictions of SO4 concentration are comparable with experimental results over the whole profile. In region I the C! and other component concentrations (Figs 4, 6, 7 and 8) are all high and this would normally result in enhanced anhydrite dissolution clue to the high ionic strength of the solution. However, the concentration of Ca in this region is also high due to Ca 2÷ being displaced from exchange sites as the wetting front approaches the bottom of the column (see below for discussion on this exchange phenomenon), resulting in gypsum being precipitated in the wetting front region and a net drop in SO4 concentration. This effect is delayed due to slow precipitation kinetics [equation (5)] and the initial leachate achieves supersaturation with respect to gypsum. The level of supersaturation reduces with larger values of the dissolution/precipitation rate constant, k, [see Fig. 5(B)]. In region II, the SO4 concentration rises because more anhydrite dissolves as Ca is readsorbed by the exchange sites and its concentration in solution falls (Fig. 6). The gradual fall of SO4 concentration in the region III is due to a drop in anhydrite (and gypsum) solubility as ionic strength falls in this region (see Figs 7 and 8) whilst the Ca concentration remains relatively constant due to ion exchange effects (Fig. 6). In both regions II and III, there is relatively little effect ofk. values in the range 0.1-1 h -~ because there is sufficient length of shale containing anhydrite to ensure that the discharged leachate remains saturated with respect to anhydrite. An effect becomes evident in region IV when the leachate becomes unsaturated with respect to anhydrite. In region IV, the remaining anhydrite/gypsum is washed out of the column giving a predicted convex curve in the S04 profile. The experimental concave curve suggests that the final leaching of CaSO4 salts is not well modelled
286
A.A. ]¢J~OL et al.
by equation (5). This is possibly due to the changing wetted area of the salts as leaching proceeds, and final control of CaSO4 leaching by intra-particle diffusion. Comparison of model and experimental results indicates that a k, value of the order of 0.1 h-m is reasonable. This is in line with the 0.05-2 h - ~ range reported for gypsum by Glas et al. (1979). Calcium. Predicted total Ca concentrations in the discharged leachate are shown in Fig. 6. The profiles are characterized by a rapid initial fall followed by an almost plateau level of concentration in region III and an eventual rapid fall as CaSO4 is washed out.
More details of regions I and II are shown in Fig. 6(B). The initial Ca level is high because Ca is displaced by other cations from dry cation exchange sites at the wetting front. The Ca profile in these regions is broadened by dispersive effects, and to a lesser extent, by slow precipitation of displaced Ca to form gypsum. As for Cl, the model fits better in regions I and II at higher pore water velocities and the fit deteriorates as the pore water velocity falls (Krol, 1985). Region III is characterized by almost level (slightly rising) Ca concentrations. In this region, the source of Ca in solution is the dissolution of anhydrite and
8O 7O
i °° q
50
(A)
40 •
IV
3O
0
\
w 10
%
I, 10
I 20
I 30
I 40
60
!1"I
•
60
70
Volume of Eluted L e a e h a t e (litre)
80
•
* J .........
J
e_
70 60
O
• ,/;/
o
-4,,2"
i .° 40
U
III (iS)
30 0 o~ 20
M
10
O0
I 1
I 2
I 8
i 4
6
Volume of Eluted L m m h a t e (lltre)
Fig. 5. Comparison of predictod and observed sulphatc results, column E. (A) Whole profile and (B) initial profile ( , ka=O.1 h-I; - ~ , ka = 0.25 h-n; - - - , k. = l.Oh-n).
Model for leaching of major inorganics
287
the general trend of the experimental data but proportionaUy larger deviations exist than did for Ca and SO4. These deviations are attributed to the use of constant cation exchange selectivity coefficients in the model. Values measured experimentally (Krol et al., 1986) showed variation with solution composition. The high initial concentration peak is due to the dissolution of readily soluble Mg salts at the wetting front. The spread of this peak over regions I and II is due to dispersion and cation exchange. As with Ca and CI, the fit at low pore water velocities is not as good as that at high pore velocities (Krol, 1985).
gypsum. However, the profile of Ca does not follow that of SO 2- (Fig. 5) because ieachate within the column is now passing over exchange sites deficient in Ca. Ca therefore exchanges with adsorbed Mg, Na and K, displacing them into solution (Figs 7 and 8). The combination of these two mechanisms serves to buffer the calcium concentration at an approximately constant value. Magnesium. Figure 7 compares predicted total Mg concentrations with experimental results. The profile is characterized by a rapid initial fall in regions I and II followed by a gradually reducing fall in region III to a zero level in region IV. Model predictions follow
1
120
80
(AI
b
0
10
20 30 40 50 V o l u m e o f Eluted L o a c h a t e (Utre)
60
70
140
120
~
I0© (U)
8C
lU
u I
u
40
°e
•
o
•
o
•
o
•
o
,
20 0 0
rl 2
I 4
I 6
I 8
10
V o l u m e o f Eluted I . e ~ h a t o (Into)
Fig. 6. Comparison of predicted and observed calcium results, column E. (A) Whole profile and (B) initial profile( , k, f O . l h-J; - - m , k, ffi O.25 h-I; - - - , k, ffi l.O h-I).
288
A.A. KRoL et al. 100
80
11 I=
60
0 m g
io
III
IV
o
0
D
2O
•
Q 0
10
20
30
40
oo
L 50
-O~) -
-
-70
Volume of Eluted L e a c h a t e (lltre)
Fig. 7. Comparison of predicted and observed magnesium results, column E. In region III, the concentration of Mg in solution is controlled by cation exchange (including the effect of the existence of the uncharged MgSO~4 ion pair). Mg is stripped out of the cation exchange sites by Ca which remains at a high level (Fig. 6). Once the Ca concentration drops away in region IV the predicted concentrations of Mg fall to low levels. Sodium and potassium. Figure 8 compares the predicted Na concentrations with experimental results.
The shape of the profile is controlled by similar mechanisms to those controlling the leaching of Mg. However, the deviation between model predictions and experimental data is larger. This is attributed to the fact that the experimentally measured values of the monovalent-divalent cation exchange selectivity coefficient, 3K, showed a much larger variation than those of the divalent-divalent cation exchange selectivity coefficient, ~K, and a marked dependence on
100
8O
i
80 "
-0
lU
10
:0
SO
IV
40
SO
- -e6
-
Volume of Eluted L e a o h a t e (lltre)
Fig. 8. Comparison of predicted and observed sodium results, colunm E.
-
-70
Model for leaching of major inorganics Na concentration (Krol et al., 1986). Modelling the effect of a varying ] K value might result in a better fit to the experimental data. The concentrations of K in leachate are much lower than Na (Krol et al., 1988) and the model predictions compared with experimental data in a simple manner to those for Na. Description of solute transport mechanisms as the wetting front moves down through the dry spent shale column. Figure 9 shows how the predicted concentrations of cations and anions at the wetting front vary as the wetting front infiltrates a column of spent shale. The predicted cation and anion concentration
289
profiles within the column at the time when the wetting front reaches the base of the column are given in Fig. 10. The four regions (I-IV) marked on Fig. 10 roughly correspond to those marked on the discharged ieachate concentration plot (i.e. Figs 5-8). The mechanisms of leaching during infiltration into the column are now discussed with reference to Figs 9 and 10. When the distilled water infiltrates the shale, the readily soluble Mg, Ha and K salts dissolve instantaneously at the wetting front. Ca and SO4 concentrations rise because Ca is stripped out of the cation exchange sites and anhydrite dissolves.
"° ! 120
o
,i 4o o
0
--
I ,, I I I 20 40 60 80 Dlstsnoo of Wetting Front from Top of Column (cm)
100
300
~
240
i
180
i
120
•
!
CI
...................
0
0
i 20
S04
I 40
I 60
I 80
i(
Distance of Wetting Front from Top of Cohmm (ore)
Fig. 9. Concentrations at the wetting front. (A) Cations and (B) anions (k, =0.25 h-l).
290
A.A. KROLet
Initially, therefore, cation and anion concentrations at the front rise (Fig. 9). However, the solution at the wetting front eventually becomes saturated with respect to anhydrite and the amount of anhydrite dissolving at the front falls despite the increasing ionic strength. Ca desorbing from exchange sites then causes precipitation of gypsum (the hydrated form of anhydrite) at the wetting front because it introduces excess Ca into solution and the solubility product of gypsum (or anhydrite in the present model) is exceeded. Thus, SO, concentrations at the front start to fall (Fig. 9) and the solid anhydrite concen-
o2.
tration at some point in the column exceeds that originally present on the dry spent shale. Because the precipitation rate of gypsum is slow, it is possible for leachate to become supersaturated with respect to anhydrite]gypsum. Na, Mg and to a lesser extent, K are desorbed from the exchange sites into solution behind the wetting front [regions III and IV of Fig. 10(A)]. Thus, despite the fact that these cations do not exist as slowly dissolving salts, they have more dispersed solution concentration profiles than C1 [Fig. 10(A) and (B)].
160 i (A) 120
80
IV
.//
•
/J, '
m
/
4o
O~)
20
40
60
80
100
Distance from Top of Column (cm)
300 (a)
i
200
L IV
III
100
0
0
l
•
I 20
I I ~...t.~¥ 40 60 80 Distance from Top of Column (cm)
TI'I 100
Fig. 10. Concentrations within the column at wetting front breakthrough. (A) Cations and (B) anions (ko =ffi0.25 h-').
Model for leaching of major inorganics CONCLUSIONS The chemical mechanisms controlling the leaching of Ca, Mg, Na, K, CI and S04 from spent shale have been identified. A model has been developed which incorporates expressions for these mechanisms and predicts leachate compositions up to the time when CaSO4 salts are totally dissolved. Uncertainties in the characterization of anhydrite dissolution and cation exchange are reflected by the presence of some deviations between model predictions and experimental results. Problems with the numerical solution technique include small errors in component mass balances and oscillating numerical solutions. Further work in progress aims to develop the current model for field application and address these problems. Acknowledgements--The financial assistance provided under the National Energy Research Development and Demonstration Program and the provision of spent shale samples by the Esso Rundle Group are gratefully acknowledged. Thanks are also expressed to Dr M. T. Van Genuchten for supplyinga copy of his curvefitting program, CFITIM. REFERENCES
Bell P. R. F., Geronimos G. L., Dobson K. and Greenfield P. F. (1982) Characteristics of leachates from Rundle shale processed in a Lurgi pilot retort. Envir. Technol. Lett. 3, 433-440. Bell P. R. F., Krol A. A. and GreenfieldP. F. (1986) Factors controlling the leaching of major and minor constituents from processed oil shale. War. Res. 20, 741-750. Fox B. W. (1976) Techniques of Sample Preparation for Liquid Scintillation Counting. North-Holland, Amsterdam. Glas T. K., Klute A. and McWhorter D. B. (1979) Dissolution and transport of gypsum in soils: theory. Soil Sci. Soc. Am. J. 43, 265-268. Hildebrand M. A. and Himmeiblau D. M. (1977) Transport of nitrate ion in unsteady, unsaturated flow in porous media. A.I. Ch. JI 23, 326-335. Keisling T. C., Ran P. S. C. and Jessup R. E. (1978) Pertinent criteria for describing the dissolution of
291
gypsum beds in flowing water. Soil Sci. Soc. Am. J. 42, 234-236. Krol A. A. (1985) Development of a model for leaching of retorted Rundle oil shale. Ph.D. thesis, Department of Chemical Engineering, The University of Queensland. Krol A. A., Bell P. R. F., Greenfield P. F. and Dunstan M. J. (1986) Ion exchange properties of retorted Rundle oil shale. Wat. Res. 10, 1299-1306. Krol A. A., Bell P. R. F., Greenfield P. F. and Dunstan M. J. (1988) Batch and column studies of the leaching of major inorganics from spent Rundle oil shale. Envir. Technol. Left. 9, 1073-1088. Kuo T. M. C., Park W. C. and Lindemanis A. (1979) Inorganic leaching of spent shale from modified in situ processing. In Proceedings of the 12th Oil Shale Symposium, Colorado, pp. 81-93. Margheim G. A. (1975) Water pollution from spent shale. Ph.D. thesis, Colorado State University. O'Neill K. and Lynch D. R. (1980) Effective and highly accurate solution of diffusion and convection--diffusion problems using moving, deforming coordinates. In Proceedings of the 3rd International Conference on Finite Elements in Water Resources, Mississippi, pp. 3.67-3.76. Ramirez W. F., Feerer J. L. and Morelli P. T. (1983) Mechanisms of leachate generation and transport in retorted oil shale. AIChE Meeting (Fall). Pap. 65D. Rowley P. G. and Brown T. (1981) Australian oil shale development. In Proceedings of Oil Shale, the Environmental Challenges II, Colorado, pp. 1-28. Slettevold C. A., Biermann A. H. and Burnham A. K. (1978) A surface area and pore volume study of retorted oil shale. Lawrence Livermore Lab. Rep. UCRL-52610. Stollenwerk K. G. (1980) Geochemistry of leachate from retorted and unretorted Colorado oil shale. Ph.D. thesis, University of Colorado. Van Genuchten M. T. (1980) A comparison of numerical solutions of the one dimensional unsaturated-saturated flow and mass transport equations. In Proceedings of the 3rd International Conference on Finite Elements in Water Resources, Mississippi, pp. 3.49-3.66. Van Genuchten M. T. (1981) Non equilibrium transport parameters from miscible displacement experiments. US Salinity Laboratory, Research Report No. 119. Van Genuchten M. T. and Wierenga P. J. (1976) Mass transfer studies in sorbing porous media: I. Analytical solutions. Soil Sci. Soc. Am. J. 40, 473-480. Wigiey T. M. L. (1977) WATSPEC--a computer program for determining the equilibrium speciation of aqueous solutions. British Geomorphological Research Group, Tech. Bull. No. 20.
0043-1354193 $5.00+0.00 Copyright © 1993 Pergamon Press Ltd
A MODEL FOR THE LEACHING OF MAJOR INORGANICS FROM RETORTED RUNDLE OIL SHALE A. A. KROL, P. R. F. BELL@ and P. F. GREENFIELD') Department of Chemical Engineering, The University of Queensland, St Lucia, Qld 4072, Australia (First received August 1991; accepted in revised form July 1992) Abstract--A model is developed to describe the leaching of major soluble inorganics (the sulphate and chloride salts of calcium, magnesium, sodium and potassium) from processed oil shale sourced from the Rundle resource, Queensland, Australia. The model accounts for rapid dissolution of readily soluble salts, the formation of non-charged ion pairs in solution, cation exchange and slow dissolution/precipitation of calcium sulphate (anhydrite/gypsum). I.~,achate results from an experiment involving irrigation of an initially dry column of processed shale are used to validate the model. A piston model is used to account for unsaturated water flow, and this is expressed in the water quality model by a moving bottom boundary condition. Tracer experiments show that a good basis for the leaching model is the classic convection-dispersion model. A moving node finite element method is used to obtain a numerical solution to this leaching model. The model provides a good description of the leaching of each inorganic species. Deviations between model predictions and experimental results are greatest at slower flowrates and for monovalent cations. It is concluded that the main chemical mechanisms controlling leaching of major inorganics from processed shale have been identified. Key worda--Rundle, oil shale, retorted shale, leaching, chemical equilibria, modelling, finite element, unsaturated flow, salt transport
INTRODUCTION The exploitation of oil shale resources will generate a large quantity of solid waste consisting o f processed (or spent) shale, low grade unprocessed shale and overburden. Retorting operations change the physical and chemical characteristics o f the shale and increase the potential mobility of salts (Slettevold et al., 1978; Stollenwerk, 1980). Lcachates are saline and generally alkaline in nature. There are two broad environmental concerns associated with such leachate: • revegetation of dumps must not be handicapped by concentration of salts in the surface zone due to an evaporation-driven upward flux; and • leachate containing toxic species must not enter ground or surface water in quantifies likely to pose a hazard. The Rundle oil shale resource is located on the central Queensland coast. Spent shale used in this study was part of a bulk sample extracted from the resource in the late 1970s which was retorted in Germany using the Lurgi-Ruhrgas process. Leachate studies on this material have been in progress since that time. The earliest studies concentrated on establishing the characteristics of leachate generated by irrigation of laboratory columns (Bell et aL, 1982), These were supplemented by studies designed to determine the mechanisms controlling leaching of
inorganic constituents (Bell et al., 1986; Krol et al. 1986, 1988). The work described in this paper involves development of a model to predict the quality of leachate generated from Rundle spent shale. Previous water quality models for oil shale leachate have been developed for material extracted from the Green River resource in the US. Rundle shale is softer and more similar to a claystone than US shale (Rowley and Brown, 1981). Most of these previous models have concentrated on mass transfer aspects of the leaching process rather than the chemistry involved (e.g. Kuo et al., 1979; Ramirez et al., 1983). They have also tended to neglect leaching under conditions of unsaturated water flow. An exception is the model proposed by Margheim (1975) which considered three solute phases: • a solution phase of Ca, Mg and Na salts; • an exchangeable phase consisting of Ca 2+, Mg 2+ and N a + adsorbed to the spent shale; and • a soluble salt, CaSO4. However, Margheim's model did not include any correcfion for the effect of ionic activity coefficients and its success in fitting experimental data was largely due to curvefitting of model parameters. In our studies it was decided a priori that the developed model should as far as possible be calibrated by independent experiments designed to determine which mechanisms should be included in the 277
A. A. KROLet al.
278
model and provide model parameters. The number of parameters obtained by curvefitting of column leachate data was to be minimized and such data were in the main to be used for model validation. The developed model was thus a test of whether the chemical mechanisms controlling the leaching of major inorganic constituents (calcium, magnesium, sodium, potassium, chloride and sulphate) from the shale had been correctly identified. A simple physical situation was chosen for the validation experiment. This involved irrigation of an initially dry column of spent shale with a constant flux of distilled water. The physics of the column experiment and previous experimental studies (Bell et al., 1986; Krol et al., 1986, 1988) showed that it was necessary to consider
inclusion of at least the following mechanisms in the model of this experiment: (1) Infiltration of water into a dry column followed by steady state flow. (2) Hydrodynamic dispersion and diffusion, and possibly mass transfer between "immobile" and "mobile" zones of water. (3) Rapid dissolution of most salts and slower dissolution of CaSO4. (4) Ion exchange of cations. (5) Formation of the ion pairs CaSO ° and MgSO °. (6) The effect of changes in ionic strength on activity coefficients and thus on the reactions (3)-(5) above.
A
constant flow we steel plate ~brane :t to hollow
s ~ t l o n In c ~ t ~ t sis via small holes
iy E gamma ray 8ource
bl mtforTn
ked with shale
wring temdometew
i x sheet rio
KL
leachata take off
i
Fig. I. One metre leachate column.
J
to vacumn controNer
279
Model for leaching of major inorganics
Experiments D E F
Superficial velocity, v (cm/h) 0.406 0.694 0.148
Table i. Column operating conditions Mean dry bulk Average pore Water content, 0 densities, Ps volume (PV) (cm3/cm3) (g/cm3) (litres) 0.585 0.888 8.44 0.616 0.887 8.96 0.551 0.881 8.00
The colum experiments which provided the model validation data are described below. Tracer studies were carried out on these columns in order to select and parameterize the hydrodynamic model on which the overall leaching model was based. The experimental section is followed by a description of the development of the model, a summary of the numerical solution scheme employed and a comparison between model predictions and experimental results. EXPERIMENTAL
Tracer volume, fraction of PV 0.684 0.792 0.823
Leachate volume collected(PV) (litres) 79.3 (9.4) 78.1 (8.7) 88.5 (11.1)
D, E and F) were carried out under different hydraulic conditions. The operating conditions of the columns and the tracer studies are summarized in Table 1. COLUMN HYDRODYNAMICS
Water flow during infiltration The water content profiles during infiltration into the dry columns were all very steep. The profiles measured on column E are presented in Fig. 2. It was necessary to develop a model for such water flow in order to provide spatial water content and velocity profiles for use by the leaching model. The usual approach to such a situation is to use an algebraic or numerical solution for the coupling of the unsaturated form of Darcy's equation with the continuity equation (e.g. Van Genuchten, 1980). This method was not used here because a situation involving
The leaching experiments were carried out in a constant temperature room (25 + 2°C) using a I m long perspex column of 12cm square internal cross section (Fig. I). The column was repacked for each experiment with airdried spent shale to a consistent density of 890 + 20 kg/m 3. Distilled, deionized water was added at a constant rate to the top of the column by a peristaltic pump. An even distribution of water was achieved by placing an 8/~m steep wetting fronts is difficult to solve numerically. cellulose acetate membrane between the top water reservoir F u r t h e r , it h a s been established t h a t with the excepand the shale. The bubbling pressure of this membrane was tion o f chloride the range o f flow h y d r o d y n a m i c s such that the water reservoir remained full throughout the experiment and water was preferentially supplied to any less represented b y experiments D, E a n d F has little saturated parts of the shale surface, i.e. where the capillary effect o n the leaching b e h a v i o u r o f m a j o r inorganics pressure was highest. In all experiments the water was ( K r o l et al., 1988). Hence, a simple " p i s t o n " w a t e r supplied at a rate below the saturated hydraulic conductivity flow m o d e l was used to a p p r o x i m a t e the w a t e r conof the shale (0.8cm/h), and thus the flow through the tent profiles by a m o v i n g step function. Such a column remained unsaturated. The progress of the wetting front was followed using m o d e l was used by H i l d e b r a n d a n d H i m m e l b l a u a y-ray attenuation system (a modified GEC 026 nuclear (1977) to predict the t r a n s p o r t o f nitrate d u r i n g density gauge containing a 780 MBq mCs source). When the infiltration into a sand column. T h e required i n p u t wetting front reached the base of the column, leachate was d a t a is the e q u i l i b r i u m water c o n t e n t , 0~, for the t o p extracted by a vacuum equivalent to the capillary pressure b o u n d a r y infiltration rate, vi. T h e position o f the measured by the tensiometer at the top of the column. This procedure achieved even capillary pressure and water conwetting f r o n t at a n y time, t, for a c o l u m n o f length tent profiles down the length of the column. The volume, L, c a n t h e n be predicted using the equations given in pH and conductivity of the discharged leachate were monT a b l e 2. T h e application o f the p i s t o n flow m o d e l to itored and samples were analysed for calcium, magnesium, c o l u m n e x p e r i m e n t E is s h o w n in Fig. 2. T h e model sodium, potassium, chloride and sulphate (Krol et al., 1988). predicts the p o s i t i o n o f the wetting f r o n t to within Once the column hydrodynamics had reached a steady state (i.e. the flowrates into and out of the column were equal and the water content profile was even and conTable 2. Equations describing the piston model for water infiltration stant) they were investigated using a tritium tracer. Tritium into a column of shale is relatively non-adsorbing, does not change water flow properties and does not interfere with the leaching mechVi anisms of those inorganic components which were being (i.e. before the wetting front reaches the bottom of the column) monitored. The tracer was injected into the column as a finite pulse wetting front position, L~ ffi tv~ with a tritium activity of approx. 107 dps/l. During tracer Or experiments the discharged leachate was directed to an column profiles, 0 ~ z ~ L~: v ffi v~, 0 = Or automatic fraction collector which collected samples at set L,
280
A.A. I ~ o L et al.
0.7] --I
0.6l
II 8 3 hours "E "E
O'S l
l I
I
0.4
!
@
I I I I
0.3 O
0
" ~
(first leachate)
0.2
5 4 hours
a '~,30.5 hours ~,~
,1 1.3 hours
J
KEY ~, ~, experiment
L
-----
/
! I I
0.1
piston model
I 0
O
10
20
30
40
50
60
7O
80
90
100
Distance from Top of Column Z (cm)
Fig. 2. Application of piston water flow model to column E.
5 cm, an accuracy which was considered adequate for our purposes. In the column experiments a maximum of approx. 2% of the total shale weight was leached. The impact of such leaching on porosity and permeability, and hence on water flow characteristics, is regarded as insignificant given the level of accuracy of the overall modelling procedures. Tracer studies Tracer breakthrough results were interpreted using CFITIM, a non-linear least squares fitting programme developed by Van G-enuchten (1981). Two models were fitted to the data, a convectiondispersion model (model A) and a physical nonequilibrium model (model B). Model B assumes the presence of "mobile" and "immobile" zones of water with diffusive mass transfer occurring between the zones, as described in detail by Van G-enuchten and Wierenga (1976). Both models assume that any adsorption of the tracer is rapid relative to the water flow rate. The one-dimensional variation of solution and solid concentrations of a species (c and s, respectively)
with time (t) and distance (z) are described by the following equation: 00c Os ~2c Oc ~ + pa ~ -- K ~zz-- V-~z
(!)
K is the dispersion/diffusion coefficient. By making use of the linear Freundlich isotherm (s = kc) and suitable substitutions, model A can be expressed in the dimensionless form used by CFITIM (Van Genuchten, 1981): R ~c = I ~2c
0x
B 0X2
~c
(2)
0X
The dimensionless variables are defined as follows: vL ; B --- -~-
C = ~ oc .,
p 0k , "
vt
Z
x = Z ;
R=I+
(3)
~ = o-2"
The dimensionless form of model B can be similarly developed (Van G-cnuchten, 1981; Krol, 1985). Results obtained from fitting model A to tracer breakthrough data are given in Table 3. Figure 3
Table 3. Fit of convection--dispersion model (A) to tracer data Model parameters fitted by CFITIM
Calculated parameters
Experiment code
B
R
T
K (cm2/h)
Tracer recovery (%)
D E F
272+30 228+22 232+37
1.011 +0.004 1.003+0.005 1.004+0.007
0.661 +0.007 0.791 +0.006 0.828+0.010
0.15+0.02 0.30+0.03 0.066+0.011
97 100 101
Note: + errors are 95% confidence limits.
Model for leaching of major inorganics shows the model fit to the tracer breakthrough data from column E. Model A provided a good fit to the data in all cases. The fitted variables were the Brenner number, B, the retardation factor, R, and the input pulse size. Prediction of the latter effectively represents numerical integration of the breakthrough data, and comparison with the known input pulse size enabled the tracer recoveries (97-101%) to be calculated. As expected the tritium tracer showed negligible adsorption (R values close to 1). The calculated dispersion coefficients (K) show a positive correlation with irrigation rate, i.e. E > D > F. Model B also showed a good fit to the experimental data, indicating that zones of "immobile" water may be present within the column. Krol (1985) presents arguments which conclude that model B may in fact represent a better approximation of the hydrodynamics in the columns than model A. However, it was the convection-dispersion model (A) which was coupled with the piston water flow model (Table 2) to provide the hydrodynamic basis for the leaching model. Model A was selected because it requires less input data and it makes inclusion of a sophisticated chemistry model much simpler. MECHANISMS INCLUDED IN THE CHEMISTRY MODEL
281
the form of anhydrite dissolve more slowly at a rate which depends on solution concentrations. The batch experiments also showed that most readily soluble components are present on, or close to, the readily accessible surface of spent shale particles. It is assumed in the column leaching model that the rapidly dissolving components (Na, Mg, K, Cl and that Ca and SO4 which is not present as anhydrite) are instantaneously dissolved at the wetting front. This assumption is affected in the model by considering dissolution as a finite dispersive flux across the moving wetting front into the wetted region. The size of this flux equates with the quantity of "component i" which must dissolve at the wetting front, and determines the concentration gradient of i at the wetting front boundary (z --Lw):
Oci [
"~z ,,=z,,
= vPBS~
~
The value of K in this equation is assumed to be that determined in the tracer studies (Table 3). Dissolution of Ca and SO4 from slowly dissolving anhydrite is characterized by a kinetic model used previously for dissolution of gypsum (e.g. Keisling et al., 1978; Glas et al., 1979): OS16
p~-~- : Sk, O(cl ± c,s).
Dissolution of salts Batch experiments have shown that, when extraparticle mass transfer effects are minimized, approx. 90% of leachable Mg, Na, K and CI dissolves in the first hour (Krol et al., 1988). Ca and SO4 present in
Key Input tracer pulse output tracer, experiment output tracer, model
EXPERIMENT water content 0 . 6 1 6 om3/©m 3 Darcy velocity 0.694 cm/hr
1.0
........
~ 0.8
I
1
2"
(5)
The parameter k, is the reaction rate constant for dissolution of anhydrite and, as included in equation (5), implicitly includes surface area. Subscripted concentrations can be identified by reference to Table 4.
1.2
--•
(4)
"\
i0.6 L)
0.2
oJ0 Pore Volumes of Leachata Eluted since start of Tracer Injection
Fig. 3. Fit of convection-dispersion model to tracer breakthrough curve for column E.
A. A. K.ROLet al.
282
Table 4. Symbolsused to representindividualspecies Solution "Solid exchangeable" phase phase Solid "leachable"phase (mequiv/l) (mequiv/g) (mequiv/g) cl ct s, (as CaCI,)
Species Ca2+ Mg2+
c2
Na +
c3
K+ CISO2CaSO° MgSO°
c( cs c6 c,6 c26
c2
$16 (as CaSO4) S2
c3
s3
--
__
c( ----
The term (c~-cts) represents the "saturation deficit" of the solution with respect to anhydrite. The concentration of Ca t+ which would result in the solution becoming saturated with respect to anhydrite, C~s, is calculated using the solubility expression: K~'6 cts = m
s( ss s6 (exceptthat as CaSO() s,6 (as CaSOD
adsorbed species. The Gapon model used in this previous work is applied in the column leaching model: c0.5~, ~r
= - -
(6)
C6
~K = c3cl
c,~ c~-~
~K = ~ . K~6 =
.
YlY6
Speciation in solution The presence in solution of the uncharged ion pairs (CaSO ° and MgSO °) increases the solubility o f anhydrite and affects the ion exchange of Ca t+ and Mg 2+. The ion pair concentrations are calculated using their equilibrium relationships: C16
C! C6
(8)
C26
K~'6 = c2 c6"
(9)
The superscript "*" on the equilibrium constant identifies it as a value corrected for the effect of activity coefficients: Kt6 YiY6
K~'6 = - )'16
~i~26y2~/6
K * , ffi ~
(14)
(7)
The parameter S in equation (5) is a switch used in the numerical solution. It is set on (S = 1) or off (S --0) for a given time at each spatial position in the column, depending on whether anhydrite is present or not. Further details are provided by Krol (1985).
K*, = --
(13)
c?.~
where K ~ is the solubility constant for anhydrite, Kt6, adjusted by the solution activity coefficients of Ca t+ , ?l and SO~-, Y6: Kt6
(12)
c~°~
(10)
(11)
Equations (12)-(14) are written in terms of solution concentrations rather than activities. It was decided not to incorporate activity coefficients into the equations as the resultant improvement in accuracy is minimal when compared with the noted variability (standard deviations o f up to 80%) in the mean measured values of the selectivity coefficients (Krol et ai., 1986). However, a mean correction factor of 20% was applied to the measured values of ~K and ~K in order to account for the average effect of the difference in salt concentration between the experiments used to measure selectivity coefficients and the column leachate studies. These selectivity coefficients for monovalent-divalent cation are those most sensitive to the effect of activitys coefficients. The cation exchange capacity (CEC) was assumed to remain constant at the value measured for the spent shale on completion of the column leaching experiments (0.255 mequiv/g; Krol et al., 1988): CEC = c~ + c2 + c3 + c4.
(15)
Activity coefficients The activity coefficients of charged species were calculated using the extended Debye--Huckel equation. F o r uncharged species the following empirical relationship was used (Wigley, 1977):
T26
Cation exchange Krol et al. (1986) demonstrated that cation exchange on the shale is significant and that equilibrium is rapidly attained between the solution and
lOgl0 Yi= 0.1L
(16)
Summary of chemistry model The incorporation o f the above mechanisms into the leaching model results in the total concentration
Model for leaching of major inorganics of any component at a wetted part of the column being represented as follows: Ca ----Ca2+ + CaSO°
+ Ca2+-shale + CaSO4(s)
Mg = Mg2+ + MgSO°
+ Mg 2+-shale
Na = Na +
+ Na+-shale
K=K
+
283
assumption [equation (4)]. When the column is fully wetted (i.e./.~ = L) the wetting front boundary condition is replaced by the usual zero concentration gradient: ~(~C5 z
+ K +-shale
z'L
=
O.
(18)
The appropriate top boundary condition for irrigation with distilled, deionized water is the constant flux condition:
CI = CISO, = SO~- + CaSO ° + MgSO°
+ CaSO4(s)
In solution, transportable
az-
Associated with solid, not transportable
Rapidly dissolving solid species are not associated with the solid because such species are assumed to dissolve instantaneously on contact with the wetting front. THE COLUMN LEACHING M O D E L FOR C H L O R I D E
Chloride is the only component which is assumed both to dissolve rapidly and to be independent of the other ieachate components in the chemistry model. The leachate model for CI was thus developed separately from the general leachate model. It provides a check on the validity of three assumptions: (1) the applicability of the piston water flow model to water infiltration; (2) the applicability of the convection-dispersion model for column hydrodynamics; (3) the instantaneous dissolution of rapidly soluble salts at the wetting front.
Development of the chloride model
- -
(17)
el has no chemical reaction interactions with other species and thus equation (17) is also the transport model for el. The region of applicability of equation (17) (0 < z
Model predictions and discussion The prediction of the chloride leaching model is compared with discharged leachate concentrations for column E in Fig. 4(A). At this flow velocity (0.694 cm]h), the model predictions agreed very well with the experimental results. However, the goodness of fit deteriorated for the slower and more unsaturated flow conditions of columns F [(Fig. 4(B)] and D (Krol, 1985). In fact the profiles predicted by the model were relatively insensitive to flow
Table 5. Distribution of leachable solids on dry shale
Component
Total* concentration (mequiv/g)
Rapidly dissolving udts (mequiv/g)
Anhydrite (mequiv/g)
Ca Mg Na K CI SO.
0.3496 0.0589 0.0397 0.0147 0.0125 0.1954
0 0.0089 0.0097 0.0057 0.0125 0.0108
0.1846 ----0.1846
*grol (1985). w a ~/2--0
(19)
The model was solved using a moving node Galerkin finite element method based on that proposed by O'Neill and Lynch (1980). This method minimizes problems of numerical dispersion which would otherwise be very pronounced, particularly when the wetting front initially starts to move into the column of dry shale and the effective Brenner number (B) based on the wetted column length (vLw/K) is at a minimum. Integration in time was by a Crank-Nicbolson finite difference scheme. The finite element nodes were positioned at each time step by placing the last node at the wetting front and evenly spacing the remaining nodes over the wetted region. No attempt was made to optimize the position of the nodes to provide a more economical and accurate solution to the model because the full leachate model, incorporating all the major inorganic components, required a more evenly spaced mesh. A detailed outline of the development and solution of the nodal equations is given by Krol (1985). The chloride leaching model was solved for all columns (i.e. D, E and F). Required model input data are included in Tables 1 (o, Pn and v), 3 (K) and 5 (ss).
For CI, the convection--dispersion model [equation (1)] can be simplified to:
00c5 02c5 0c5 Ot = K~z2 v O'zz"
VCs)]=.0--0.
Exchange sites (estimated)* (mequiv/g) 0.165 0.05 0.03 0.01
284
A.A. KIgOLet al.
conditions, whereas the experimental results showed greater dispersion under slow flow conditions. A possible reason for this behaviour is a break-
down in the instantaneous dissolution assumption due to: • a greater proportion of immobile water being present as the flow becomes more unsaturated, with dissolution being much slower in immobile regions; and • the wetting front becoming less sharp, resulting in CI dissolution occurring over a finite bandwidth rather than at a sharp front.
Despite the short-comings of the chloride leaching model, it was used as the basis for modelling the leaching of other major inorganics. This was justified because Cl was rapidly leached from the columns during the passage of approximately one third of a pore volume of leachate, whilst the other inorganics remained present for many pore volumes (Krol e t al., 1988). Thus the shape of the wetting front and the column hydrodynamics were not the major factors determining the nature of these other leaching profiles, and their simulation using a simple hydrodynamic model with minimum data requirements was reasonable.
400 350 300 I 28(]
KEY
J
leo ~s~p~/mentnl I
200
i
15C
(A)
100
4oo~ ~'°I i aoo~
i
x.y
I. =--,-, II
1.°I
~I°°F"'\. O o
I I[
,-, I I
~
Eluted Leechate (lflre)
Fig. 4. Comparison of predicted and observed chloride results. (A) Column E, v = 0.694 cm/h and (B) column F, v = 0.148 cm/h.
Model for leaching of major inorganics THE COLUMN LEACHING MODEL FOR
ALL MAJOR INORGANICS
Model development
The transport equations for each of the other components (Ca, Mg, Na, K and SO4) cannot be solved independently, as was done for CI, because these equations are coupled through cation exchange, solution speciation and anhydrite dissolution reactions. The transport equations for these components arc:
Calcium:
oOC,+ Pa~t + o-Td Ocl~ + p.-Osl~ ~ Ot
02ct
= x ~
Ocl
- ~~
O"cl6
Oct6
+ r-~F~ - ~ - F ;
(20)
Magnesium: 0 0c2
O~ . ,, 0 % 02c~ 0c2 . 02e26 0c26 = K ~Tz2 -- V ~z + K-~z2 - v w
(21)
Sodium.
0 0c3
O~
O=c3
Oe3
O~
02c,
0c4
+ PB'~ =
(22)
Potassium:
00c,
+ p,-~ =
(23)
The concentration of sulphate is obtained by charge balance: (76 = C1 + (:72+ (73 + C4 - (75.
(24)
The appropriate top and bottom boundary conditions for equations (20)-(23) are similar to those for the chloride model. Equations (17) and (20)-(24) are non-linear, stiff and interacting. Expressions for the interacting chemistry, i.e. anhydrite dissolution [equation (5)], ion pair formation [equations (8) and (9)] and cation exchange [equations (12)-(14)], were substituted into the cation transport equations (20)-(23). The charge balance relationship [equation (24)] was used to eliminate SO~- to give a set of four equations. This equation set was then solved in a similar manner to the chloride model. However, because of the non-linear nature of the nodal equations, an iterative procedure was used to solve the problem. Full details are given by Krol (1985). Model predictions and discussion
The developed model was used to simulate the leaching of column E. The only model coefficient varied was the anhydrite dissolution rate constant, k,. Mass balance errors due to the numerical technique
285
were generally less than 1% of the total leachable concentration. A reduction in the time step or spatial grid spacing had negligible effect on the predicted results. The figures referred to below present total concentrations of Ca, Mg or SO4 (i.e. the free ion plus that portion existing as the appropriate ion pair). Sulphate. Figure 5 compares the predicted and experimental SO4 concentrations in the leachate. The SO4 profile is characterized by four regions: (I) a rapid initial fall; (II) a rise to a peak after about 2.4 I. (0.3 pore volumes of leachate); (III) a gradual drop in concentration until 50-601. (6.3-7.7 pore volumes) of leachate have been discharged; and (IV) another rapid fall. Some oscillation in the numerical results was experienced for region IV due to a numerical stability problem when the remaining anhydrite only occupies the last few nodes in the column. Regions I-IV are common to all of the figures and thus discussion of the mechanisms affecting the leaching of one component in a particular region is also relevant to the leaching of other components in that region. Figure 5 shows that model predictions of SO4 concentration are comparable with experimental results over the whole profile. In region I the C! and other component concentrations (Figs 4, 6, 7 and 8) are all high and this would normally result in enhanced anhydrite dissolution clue to the high ionic strength of the solution. However, the concentration of Ca in this region is also high due to Ca 2÷ being displaced from exchange sites as the wetting front approaches the bottom of the column (see below for discussion on this exchange phenomenon), resulting in gypsum being precipitated in the wetting front region and a net drop in SO4 concentration. This effect is delayed due to slow precipitation kinetics [equation (5)] and the initial leachate achieves supersaturation with respect to gypsum. The level of supersaturation reduces with larger values of the dissolution/precipitation rate constant, k, [see Fig. 5(B)]. In region II, the SO4 concentration rises because more anhydrite dissolves as Ca is readsorbed by the exchange sites and its concentration in solution falls (Fig. 6). The gradual fall of SO4 concentration in the region III is due to a drop in anhydrite (and gypsum) solubility as ionic strength falls in this region (see Figs 7 and 8) whilst the Ca concentration remains relatively constant due to ion exchange effects (Fig. 6). In both regions II and III, there is relatively little effect ofk. values in the range 0.1-1 h -~ because there is sufficient length of shale containing anhydrite to ensure that the discharged leachate remains saturated with respect to anhydrite. An effect becomes evident in region IV when the leachate becomes unsaturated with respect to anhydrite. In region IV, the remaining anhydrite/gypsum is washed out of the column giving a predicted convex curve in the S04 profile. The experimental concave curve suggests that the final leaching of CaSO4 salts is not well modelled
286
A.A. ]¢J~OL et al.
by equation (5). This is possibly due to the changing wetted area of the salts as leaching proceeds, and final control of CaSO4 leaching by intra-particle diffusion. Comparison of model and experimental results indicates that a k, value of the order of 0.1 h-m is reasonable. This is in line with the 0.05-2 h - ~ range reported for gypsum by Glas et al. (1979). Calcium. Predicted total Ca concentrations in the discharged leachate are shown in Fig. 6. The profiles are characterized by a rapid initial fall followed by an almost plateau level of concentration in region III and an eventual rapid fall as CaSO4 is washed out.
More details of regions I and II are shown in Fig. 6(B). The initial Ca level is high because Ca is displaced by other cations from dry cation exchange sites at the wetting front. The Ca profile in these regions is broadened by dispersive effects, and to a lesser extent, by slow precipitation of displaced Ca to form gypsum. As for Cl, the model fits better in regions I and II at higher pore water velocities and the fit deteriorates as the pore water velocity falls (Krol, 1985). Region III is characterized by almost level (slightly rising) Ca concentrations. In this region, the source of Ca in solution is the dissolution of anhydrite and
8O 7O
i °° q
50
(A)
40 •
IV
3O
0
\
w 10
%
I, 10
I 20
I 30
I 40
60
!1"I
•
60
70
Volume of Eluted L e a e h a t e (litre)
80
•
* J .........
J
e_
70 60
O
• ,/;/
o
-4,,2"
i .° 40
U
III (iS)
30 0 o~ 20
M
10
O0
I 1
I 2
I 8
i 4
6
Volume of Eluted L m m h a t e (lltre)
Fig. 5. Comparison of predictod and observed sulphatc results, column E. (A) Whole profile and (B) initial profile ( , ka=O.1 h-I; - ~ , ka = 0.25 h-n; - - - , k. = l.Oh-n).
Model for leaching of major inorganics
287
the general trend of the experimental data but proportionaUy larger deviations exist than did for Ca and SO4. These deviations are attributed to the use of constant cation exchange selectivity coefficients in the model. Values measured experimentally (Krol et al., 1986) showed variation with solution composition. The high initial concentration peak is due to the dissolution of readily soluble Mg salts at the wetting front. The spread of this peak over regions I and II is due to dispersion and cation exchange. As with Ca and CI, the fit at low pore water velocities is not as good as that at high pore velocities (Krol, 1985).
gypsum. However, the profile of Ca does not follow that of SO 2- (Fig. 5) because ieachate within the column is now passing over exchange sites deficient in Ca. Ca therefore exchanges with adsorbed Mg, Na and K, displacing them into solution (Figs 7 and 8). The combination of these two mechanisms serves to buffer the calcium concentration at an approximately constant value. Magnesium. Figure 7 compares predicted total Mg concentrations with experimental results. The profile is characterized by a rapid initial fall in regions I and II followed by a gradually reducing fall in region III to a zero level in region IV. Model predictions follow
1
120
80
(AI
b
0
10
20 30 40 50 V o l u m e o f Eluted L o a c h a t e (Utre)
60
70
140
120
~
I0© (U)
8C
lU
u I
u
40
°e
•
o
•
o
•
o
•
o
,
20 0 0
rl 2
I 4
I 6
I 8
10
V o l u m e o f Eluted I . e ~ h a t o (Into)
Fig. 6. Comparison of predicted and observed calcium results, column E. (A) Whole profile and (B) initial profile( , k, f O . l h-J; - - m , k, ffi O.25 h-I; - - - , k, ffi l.O h-I).
288
A.A. KRoL et al. 100
80
11 I=
60
0 m g
io
III
IV
o
0
D
2O
•
Q 0
10
20
30
40
oo
L 50
-O~) -
-
-70
Volume of Eluted L e a c h a t e (lltre)
Fig. 7. Comparison of predicted and observed magnesium results, column E. In region III, the concentration of Mg in solution is controlled by cation exchange (including the effect of the existence of the uncharged MgSO~4 ion pair). Mg is stripped out of the cation exchange sites by Ca which remains at a high level (Fig. 6). Once the Ca concentration drops away in region IV the predicted concentrations of Mg fall to low levels. Sodium and potassium. Figure 8 compares the predicted Na concentrations with experimental results.
The shape of the profile is controlled by similar mechanisms to those controlling the leaching of Mg. However, the deviation between model predictions and experimental data is larger. This is attributed to the fact that the experimentally measured values of the monovalent-divalent cation exchange selectivity coefficient, 3K, showed a much larger variation than those of the divalent-divalent cation exchange selectivity coefficient, ~K, and a marked dependence on
100
8O
i
80 "
-0
lU
10
:0
SO
IV
40
SO
- -e6
-
Volume of Eluted L e a o h a t e (lltre)
Fig. 8. Comparison of predicted and observed sodium results, colunm E.
-
-70
Model for leaching of major inorganics Na concentration (Krol et al., 1986). Modelling the effect of a varying ] K value might result in a better fit to the experimental data. The concentrations of K in leachate are much lower than Na (Krol et al., 1988) and the model predictions compared with experimental data in a simple manner to those for Na. Description of solute transport mechanisms as the wetting front moves down through the dry spent shale column. Figure 9 shows how the predicted concentrations of cations and anions at the wetting front vary as the wetting front infiltrates a column of spent shale. The predicted cation and anion concentration
289
profiles within the column at the time when the wetting front reaches the base of the column are given in Fig. 10. The four regions (I-IV) marked on Fig. 10 roughly correspond to those marked on the discharged ieachate concentration plot (i.e. Figs 5-8). The mechanisms of leaching during infiltration into the column are now discussed with reference to Figs 9 and 10. When the distilled water infiltrates the shale, the readily soluble Mg, Ha and K salts dissolve instantaneously at the wetting front. Ca and SO4 concentrations rise because Ca is stripped out of the cation exchange sites and anhydrite dissolves.
"° ! 120
o
,i 4o o
0
--
I ,, I I I 20 40 60 80 Dlstsnoo of Wetting Front from Top of Column (cm)
100
300
~
240
i
180
i
120
•
!
CI
...................
0
0
i 20
S04
I 40
I 60
I 80
i(
Distance of Wetting Front from Top of Cohmm (ore)
Fig. 9. Concentrations at the wetting front. (A) Cations and (B) anions (k, =0.25 h-l).
290
A.A. KROLet
Initially, therefore, cation and anion concentrations at the front rise (Fig. 9). However, the solution at the wetting front eventually becomes saturated with respect to anhydrite and the amount of anhydrite dissolving at the front falls despite the increasing ionic strength. Ca desorbing from exchange sites then causes precipitation of gypsum (the hydrated form of anhydrite) at the wetting front because it introduces excess Ca into solution and the solubility product of gypsum (or anhydrite in the present model) is exceeded. Thus, SO, concentrations at the front start to fall (Fig. 9) and the solid anhydrite concen-
o2.
tration at some point in the column exceeds that originally present on the dry spent shale. Because the precipitation rate of gypsum is slow, it is possible for leachate to become supersaturated with respect to anhydrite]gypsum. Na, Mg and to a lesser extent, K are desorbed from the exchange sites into solution behind the wetting front [regions III and IV of Fig. 10(A)]. Thus, despite the fact that these cations do not exist as slowly dissolving salts, they have more dispersed solution concentration profiles than C1 [Fig. 10(A) and (B)].
160 i (A) 120
80
IV
.//
•
/J, '
m
/
4o
O~)
20
40
60
80
100
Distance from Top of Column (cm)
300 (a)
i
200
L IV
III
100
0
0
l
•
I 20
I I ~...t.~¥ 40 60 80 Distance from Top of Column (cm)
TI'I 100
Fig. 10. Concentrations within the column at wetting front breakthrough. (A) Cations and (B) anions (ko =ffi0.25 h-').
Model for leaching of major inorganics CONCLUSIONS The chemical mechanisms controlling the leaching of Ca, Mg, Na, K, CI and S04 from spent shale have been identified. A model has been developed which incorporates expressions for these mechanisms and predicts leachate compositions up to the time when CaSO4 salts are totally dissolved. Uncertainties in the characterization of anhydrite dissolution and cation exchange are reflected by the presence of some deviations between model predictions and experimental results. Problems with the numerical solution technique include small errors in component mass balances and oscillating numerical solutions. Further work in progress aims to develop the current model for field application and address these problems. Acknowledgements--The financial assistance provided under the National Energy Research Development and Demonstration Program and the provision of spent shale samples by the Esso Rundle Group are gratefully acknowledged. Thanks are also expressed to Dr M. T. Van Genuchten for supplyinga copy of his curvefitting program, CFITIM. REFERENCES
Bell P. R. F., Geronimos G. L., Dobson K. and Greenfield P. F. (1982) Characteristics of leachates from Rundle shale processed in a Lurgi pilot retort. Envir. Technol. Lett. 3, 433-440. Bell P. R. F., Krol A. A. and GreenfieldP. F. (1986) Factors controlling the leaching of major and minor constituents from processed oil shale. War. Res. 20, 741-750. Fox B. W. (1976) Techniques of Sample Preparation for Liquid Scintillation Counting. North-Holland, Amsterdam. Glas T. K., Klute A. and McWhorter D. B. (1979) Dissolution and transport of gypsum in soils: theory. Soil Sci. Soc. Am. J. 43, 265-268. Hildebrand M. A. and Himmeiblau D. M. (1977) Transport of nitrate ion in unsteady, unsaturated flow in porous media. A.I. Ch. JI 23, 326-335. Keisling T. C., Ran P. S. C. and Jessup R. E. (1978) Pertinent criteria for describing the dissolution of
291
gypsum beds in flowing water. Soil Sci. Soc. Am. J. 42, 234-236. Krol A. A. (1985) Development of a model for leaching of retorted Rundle oil shale. Ph.D. thesis, Department of Chemical Engineering, The University of Queensland. Krol A. A., Bell P. R. F., Greenfield P. F. and Dunstan M. J. (1986) Ion exchange properties of retorted Rundle oil shale. Wat. Res. 10, 1299-1306. Krol A. A., Bell P. R. F., Greenfield P. F. and Dunstan M. J. (1988) Batch and column studies of the leaching of major inorganics from spent Rundle oil shale. Envir. Technol. Left. 9, 1073-1088. Kuo T. M. C., Park W. C. and Lindemanis A. (1979) Inorganic leaching of spent shale from modified in situ processing. In Proceedings of the 12th Oil Shale Symposium, Colorado, pp. 81-93. Margheim G. A. (1975) Water pollution from spent shale. Ph.D. thesis, Colorado State University. O'Neill K. and Lynch D. R. (1980) Effective and highly accurate solution of diffusion and convection--diffusion problems using moving, deforming coordinates. In Proceedings of the 3rd International Conference on Finite Elements in Water Resources, Mississippi, pp. 3.67-3.76. Ramirez W. F., Feerer J. L. and Morelli P. T. (1983) Mechanisms of leachate generation and transport in retorted oil shale. AIChE Meeting (Fall). Pap. 65D. Rowley P. G. and Brown T. (1981) Australian oil shale development. In Proceedings of Oil Shale, the Environmental Challenges II, Colorado, pp. 1-28. Slettevold C. A., Biermann A. H. and Burnham A. K. (1978) A surface area and pore volume study of retorted oil shale. Lawrence Livermore Lab. Rep. UCRL-52610. Stollenwerk K. G. (1980) Geochemistry of leachate from retorted and unretorted Colorado oil shale. Ph.D. thesis, University of Colorado. Van Genuchten M. T. (1980) A comparison of numerical solutions of the one dimensional unsaturated-saturated flow and mass transport equations. In Proceedings of the 3rd International Conference on Finite Elements in Water Resources, Mississippi, pp. 3.49-3.66. Van Genuchten M. T. (1981) Non equilibrium transport parameters from miscible displacement experiments. US Salinity Laboratory, Research Report No. 119. Van Genuchten M. T. and Wierenga P. J. (1976) Mass transfer studies in sorbing porous media: I. Analytical solutions. Soil Sci. Soc. Am. J. 40, 473-480. Wigiey T. M. L. (1977) WATSPEC--a computer program for determining the equilibrium speciation of aqueous solutions. British Geomorphological Research Group, Tech. Bull. No. 20.