Dynamic simulation of the Dead Sea

Advances in Water Resources 25 (2002) 263–277 www.elsevier.com/locate/advwatres

Dynamic simulation of the Dead Sea B.N. Asmar

*,1,

Peter Ergenzinger

Berlin Environmental Research Group, Fachbereich Geowissenshaften, Freie Universit€at, Berlin, Germany Received 12 January 2001; received in revised form 14 November 2001; accepted 5 December 2001

Abstract The Dead Sea is a hypersaline terminal lake located in the Rift Valley between Jordan and Israel. In this work a generalised mathematical model describing the behaviour of the Dead Sea has been developed. The model established the condition of the Sea by evaluating a series of ordinary differential equations describing mass balances on the water and major chemical species in the Sea. The Sea was modelled as a two-layer system. The model was validated by comparing its predictions to measured level records. The results obtained highlighted the importance of detailed evaporation modelling, showed the necessity to model the Sea as a two-layer system, validated the usage of average distribution data to estimate the flowrates of rivers, and justified ignoring diffusion effects in further modelling. The model predicted that in the case of continuing current conditions, the level will continue to decline, at a decelerating rate, because the area and evaporation rate are both decreasing. Under these conditions, the model shows that the salinity of both layers will continue to increase, and that seasonal stratification and seasonal crystallisation of gypsum and aragonite will continue. Ó 2002 Published by Elsevier Science Ltd. Keywords: Dead Sea; Mathematical modelling; Dynamic simulation; Crystallisation; Salinity

1. Introduction The Dead Sea (DS) is a closed lake situated at the deepest part of the Jordan Rift Valley, between Moab (Kerak) Mountains in the east and the Judean Desert in the west. Its catchment area covers approximately 40 000 km2 in six countries. See Fig. 1. The Sea’s eastern shore belongs to Jordan, and its western shore is currently under Israeli control. (The north-western shore is part of the West Bank, where ongoing negotiations between Israel and the Palestinian National Authority will determine its final status). The ‘‘classical’’ DS consisted of two basins (northern and southern) divided by the Lisan Peninsula. In 1966, the Sea covered an area of 940 km2 (76% in the northern basin). Its length was 76 km, with an average width of 14 km. The total volume of the water in the Sea was estimated to be 142 km3 (only 0.5% in the southern basin). The level of the Sea was 399 m below mean sea level (bmsl). The maximum depth was 749 m bmsl, whereas * Corresponding author. Tel.: +44-115-9513053; fax: +44-1159514181. E-mail address: [email protected] (B.N. Asmar). 1 Current address: School of Chemical, Environmental and Mining Engineering, University of Nottingham, University Park, Nottingham NG7 2RD, UK.

the average depth of the southern basin was only 10 m. In 1994, the Sea consisted only of a shrinking northern basin, whereas the southern basin has been converted entirely into industrial evaporation pans. The DS area was reduced to 683 km2 and its level was 406 m bmsl. Fig. 2 shows the DS in 1966 and 1994 illustrating the difference in the area. The climate in the DS area is highly arid, with a mean annual rainfall of around 65 mm. Average January temperatures range between 14 and 17 °C, and the August mean temperature is 34 °C, with a recorded maximum of 51 °C. The main feature of the DS is its salinity, which was 332 g/l (the highest in the world) resulting in an average density of 1.23 g/l. The DS had been in a stratified state for centuries. This was due to the density difference between a dense deep layer and a more dilute upper layer. A transient zone existed between the two layers, and was characterised by a steep density gradient [1]. The continuous drop in the Sea level caused by less water inflow into the Sea resulted in an increase in the salinity of the surface layer, which meant an increased water density. In 1979, these factors caused a complete overturn. The two-layer system disappeared and was replaced by one homogeneous layer with temperature, salinity and density higher than those of the former deep water [2].

0309-1708/02/$ - see front matter Ó 2002 Published by Elsevier Science Ltd. PII: S 0 3 0 9 - 1 7 0 8 ( 0 1 ) 0 0 0 6 3 - X

264

B.N. Asmar, P. Ergenzinger / Advances in Water Resources 25 (2002) 263–277

Fig. 2. The Dead Sea in 1966 and 1994.

Fig. 1. The Dead Sea location and catchment area.

The decline in water flow into the Sea and the massive utilisation of the Sea’s water for industrial purposes are still causing a continuous drop in the level of the Sea and a shrinkage in its area. This problem prompted numerous studies into its causes and effects. The result was a series of proposed solutions of how to halt the decline. However, most studies dealt with specific topics, and made no attempt to analyse their results within a broader frame treating the Sea as an integrated system, for example looking into the salinity, solubility of the salts and evaporation. Very few models were found in literature, which looked into simulating the behaviour of the DS as a whole, e.g. [3–5]. Mero and Simon in two papers presented a detailed model of the Sea’s behaviour. They took into account the existence of two layers in the Sea, and developed a relationship between the salinity and the density of the Seawater. However, instead of using direct measurements as data, they estimated all inputs based on the rainfall in Jerusalem and Nazareth, 20 and 100 km from the DS, respectively. They also used a simplified approach to estimate evaporation, and included the southern basin of the Sea in

their calculation. These factors affected the accuracy of their model significantly. On the other hand, the model developed by Harza was over-simplified. It ignored the two-layer system in the DS, and did not include the chemical’s material balances. Consequently, their estimates of an annual average of 90 cm decline in the Sea level were noticeably high when compared to the recorded average of 55 cm in the last 20 years. In this work, a dynamic model has been developed to predict the behaviour of the DS. The mathematical model of the Sea is presented next. The model, although simple and based on basic mass balance concepts, is detailed enough to achieve satisfactory accuracy. The implementation of the simulation model and its validation is presented later, and is followed by examples of the results to illustrate the model capabilities. The paper concludes with a summary of the main results.

2. Mathematical model The approach presented can be used as a general framework to model saline water bodies. The state of any system is determined by solving water and chemicals mass balances, where the detailed chemicals budget is used to enhance the estimation of evaporation and hence improve the overall accuracy of the model. For each system localised relationships are used to describe

B.N. Asmar, P. Ergenzinger / Advances in Water Resources 25 (2002) 263–277

the details, for example the volume-area-level relationships. For the DS, the mathematical model describing the Sea’s behaviour consists of a set of algebraic and a set of differential equations. The state of the Sea is established by solving a set of ordinary differential equations formed by writing a set of mass balance equations for the water and eight major chemicals in the Sea. The principle of energy conservation is accounted for implicitly through evaporation calculations. In addition to that, several state equations are used that correlate various Sea properties. The model development is divided into two parts: the water mass balance, and the mass balances over other chemicals. The overall state of the Sea is then determined combining both parts together. 2.1. Water mass balance Several water budget calculations for the DS were performed earlier [3,4,6,7], both as full studies or as approximate calculations required as complementary parts of other studies. Almost all these studies assumed the evaporation rate to be constant by ignoring the effects of salinity changes on evaporation, consequently the results showed inaccuracies in the water budget calculations. The model derived here tries to rectify this by accounting for the evaporation dynamics. The water budget calculation is treated here as a simple mass balance problem, in which the principle of conservation of mass in a closed system is applied. The general mass balance equation is Accumulation ¼ Input  Output þ Generation  Consumption:

ð1Þ

Some general assumptions are made in order to simplify the derivation of the mathematical model. These assumptions are: 1. The DS water body is considered to consist of two layers, which are characterised by having different salinity and density. The depth of each layer is not considered fixed or constant, but is determined solely by the results of the mass balances in the model. This is significant as it allows the non-existence of the top layer in some situations, and differs from other models, which assume a fixed top layer depth (e.g. [4]). 2. The time step for the entire calculation is one month. This choice is inevitable since this model uses measured data, which is limited and in most cases is reported only on monthly basis. However, in some cases only yearly data were found. Using a daily step will introduce large uncertainties as the input data have to be estimated. 3. The density of the water streams is considered to be that of fresh water [1 g=cm3 ] for the inflows for which

265

salinity is less than [10 g/l]. This can be justified as the change in density is negligible below that salinity [8]. 4. No water is generated or consumed by means of chemical reactions in the DS. Thus the third and the fourth terms on the right hand side of Eq. (1) are omitted. 5. Each stream is considered chemically homogeneous and its chemical and physical properties are fixed at the specific time step. 6. The calculations are performed only on the northern basin of the DS as the southern part is totally transformed into salt pans and controlled by industry (since 1982 the two basins are effectively separated). This is justified, as the plans to restore the DS to its previous level exclude the southern basin where dams and dikes are to be constructed to prevent the incoming water from reaching the southern basin. In this model the water state is determined by two mass balance equations, one per layer. Each equation follows the general form of Eq. (1). The inflows into the DS consist of both natural and man-made inflows. Fig. 1 shows the catchment area of the DS and identifies its main water resources. The Jordan River was the main water source supplying the DS. Historically, it accounted to 75% of all supply to the Sea. But with the increasing demand on water in the region, most of its water is now diverted. Several other small rivers flow into the Sea from the east, and are (from north to south) Al-U’zaymi, Zerqa-Ma’in-Zara, Mujib, Ibn Hammad, Kerak, and Hasa. The latter flows to the southern basin of the DS. From the west there is no permanent inflows into the Sea. Other natural inflows consist of the rainfall directly into the Sea, the storm water run off in the dry valleys, the surface and submerged springs inflow, and the infiltration from the surrounding aquifers. Manmade inflow is basically the end-brine from industry. For centuries, the only output from the DS was evaporation. But since 1950s a new output exists, in the form of the industry intake water pumped from the Sea to be processed to produce chemicals. Fig. 3 shows a block diagram of the water mass balance over the DS

Fig. 3. Water mass balance over the Dead Sea.

266

B.N. Asmar, P. Ergenzinger / Advances in Water Resources 25 (2002) 263–277

including all input and output streams. The surface water is assumed to flow into the top layer, as well as the direct rainfall over the Sea. The infiltration is assumed to flow into the bottom layer. Industries around the Sea extract water from depths of up to 40 m, thus the industrial pumping is assumed to be always from the bottom layer, as well as the returning end-brines that are considered sinking directly to the bottom layer. Evaporation is included as an output in both layers. The calculation method used assumes that evaporation will start at the top layer, and proceeds to the bottom layer only if the top layer is completely evaporated. In the top layer, the water mass balance is expressed as dW1 ¼ dt

X

X

FWMRV þ FWMSP þ FWMR X þ FWMSV  E1 ;

ð2Þ P

where W1 is the mass of water in the top layer, FWMRV is Pthe sum of the rivers water mass inflow into the Sea, P FWMSP is the sum of the springs water mass inflow, FWMSV is the sum of the storm run-off in the side valleys water mass inflow, FWMR is the water mass rainfall over the Sea, and E1 is the evaporation from the top layer. In the bottom layer, the water mass balance is expressed as dW2 ¼ FWMIN þ FWMEB  FWMP  E2 ; dt

FVRV;i ¼ FYRV;i  xj;RV;i ;

where FVRV;i is the river flowrate [m3 /month], FYRV;i is the average yearly flowrate [m3 =y], and xj;RV;i is the monthly fraction of the flowrate, i is the river index and j is the month index. The water mass flow [kg/month] for each river FWMRV;i is obtained by solving a simple mass balance. The volumetric flowrate FVRV;i is multiplied by the dif3 ference between the river’s density ½kg=m  qRV;i and 3 salinity ½kg=m  SRV;i FWMRV;i; ¼ FVRV;i ðqRV;i  SRV;i Þ:

ð5Þ

The monthly fractions of inflow values for the Jordan River were obtained by analysing statistically the available monthly data over a period of 30 years (see [7,9]). Similar analysis was also conducted on the data obtained for the Mujib River [9]. Figs. 4 and 5 show the monthly flowrate fractions for Jordan River and Mujib River, respectively. For the other four eastern rivers considered, no data is available concerning their monthly flowrates, therefore the same flow fractions calculated for the Mujib River are considered to apply for them as they are located in an area with similar meteorological conditions.

ð3Þ

where W2 is the mass of water in the bottom layer, FWMIN is the infiltration water mass inflow into the DS, FWMEB is the end-brine water mass inflow, FWMP is the water mass of the industry intake, and E2 is the evaporation from the bottom layer. 2.1.1. Rivers flowrate Due to insufficient data of measured discharges for the rivers in the area, and the fact that the accuracy of the measured data is poor [4,9], there is a need to estimate the natural rivers inflow. The recent reduction of these inflows, due to the intensive water consumption in the area and the need to use all the available water sources, mean that using simple average values for the inflows is insufficient. The method assumes that the fraction of each monthly flow relative to the total flow is constant and is independent of the flow quantity. These fractions are calculated by obtaining the dimensionless ratios of flowrate in each month to the total flowrate. The mean values of these ratios are calculated over a number of years, so as to obtain dimensionless quantities that can best describe the behaviour of the river independent from fluctuations. The general formula to describe the rivers inflow into the DS is expressed as

ð4Þ

Fig. 4. Monthly flowrate fractions for the Jordan River.

Fig. 5. Monthly flowrate fractions for the Mujib River.

B.N. Asmar, P. Ergenzinger / Advances in Water Resources 25 (2002) 263–277

The average yearly flowrates for all rivers can be found in Appendix A. 2.1.2. Rainfall The method used in estimating the flowrate monthly fractions of the rivers is again used here to calculate monthly fractions for the rainfall data. The data for the monthly rainfall were obtained [10] for Deir-Alla in the Jordan Valley some 50 km north of the DS (1965–1995), and at Safi at the southern tip of the southern basin of the DS (1975–1995). The assumption above was tested by analysing statistically the available monthly data for a time period of 21 years (1975–1995) from records for both stations. Fig. 6 shows the rainfall monthly fractions in both stations. In this model, the fractions obtained for Safi are used. The rainfall water mass flowrate equation is expressed as FWMR ¼ ðqR  SR Þ  FYR  xj;R  A;

ð6Þ

where FWMR is the monthly water mass flowrate of rainfall into the Sea [kg/month], FYR is the average annual rainfall [m/y], xj;R is the monthly fraction of the rainfall, qR is the density of the rainfall water [kg=m3 ], SR is the salinity of the rainfall water [kg=m3 or g/l], and A is the DS surface area [m2 ]. 2.1.3. Storm run-off flowrate It was assumed that the storm run-off water inflow is related to the rainfall. This was suggested previously by Arad and Michaeli [11], who considered the flowrate in the valleys surrounding the DS to be less than 10% of the total rainfall in the catchment area of each valley. They stated that in the southern parts of the area the flow might even be discarded. El-Naqa [12] used a similar concept in the Mujib River area. Based on the data available for Mujib River flood values [12], and Tsin Valley that reaches the southern basin of the DS [6] a value of 2% of the total rainfall of a valley’s catchment area is assumed to flow as storm run-off water into the Sea.

Fig. 6. Monthly rainfall fractions for Deir Alla and Safi Stations.

267

The storm run-off water mass flowrate equation is expressed as FWMSV ¼ ðqR  SR Þ  FYSV  xj;R ;

ð7Þ

where FWMSV is the monthly water mass flowrate of storm run-off water into the DS [kg/month], FYSV is the yearly volumetric flowrate of storm run-off water [m3 =y], xj;R is the monthly fraction of the rainfall, qR is the 3 density of the rainfall water [kg=m ], and SR is the sa3 linity of the rainfall water [kg=m or g/l]. The water characteristics and the distribution fractions for the rainfall are used here as the water flow is assumed to be generated directly from the rainfall. The estimates for the catchment areas and average rainfall of the dry valleys are found in Appendix A. 2.1.4. Springs flowrate Springs flowrate data are largely unavailable in the DS area. Several estimations for the flow of the springs were given in literature. The estimates varied considerably. For example, Stiller and Chung [13] estimated the measured springs around the DS as 30 106 m3 =y whereas Chan and Chung [14] estimated the springs flow as 150 106 m3 =y. The flowrates used for each spring are found in Appendix A. For each spring, the water mass flowrate can be calculated as FWMSP;i ¼ ðqSP;i  SSP;i Þ  FYSP;i  xj;RV;Mujib ;

ð8Þ

where FWMSP;i is the monthly water mass flowrate of each spring into the DS [kg/month], FYSP;i is the yearly volumetric flowrate of each spring [m3 =y; qSP;i is the 3 density of each spring’s water [kg=m ], and SSP;i is the 3 salinity of each spring’s water [kg=m or g/l]. For the calculations in the model, the flowrate monthly distribution is assumed to follow the distribution of Mujib River, hence the monthly fraction of Mujib River xj;RV;Mujib is used in Eq. (8). 2.1.5. Infiltration flowrate Several estimates in the literature appear to evaluate the water inflow into the DS from the surrounding aquifer, which usually includes the submerged springs. These values vary considerably because all of them are based on calculations. For example, Bentor [6] estimated the flow at 9:85 106 m3 =y, and Salameh [15] estimated the infiltration value from the eastern side only at 100 106 m3 =y. Steinhorn [2] estimated the flow of the submerged springs which included the seepage at 21–45 106 m3 =y. This estimation is based on a comprehensive heat balance, which depended heavily on measured values of all the variables. The accuracy of his estimate is considered more reliable than the other estimates reported, and therefore is used in this model.

268

B.N. Asmar, P. Ergenzinger / Advances in Water Resources 25 (2002) 263–277

The infiltration water mass flowrate is expressed as FWMIN ¼ ðqIN  SIN Þ  FYIN

mth ; 365:25

ð9Þ

where FWMIN is the monthly infiltration water mass flowrate [kg/month], FYIN is the yearly infiltration volumetric flowrate [m3 =y], qIN is the density of infiltra3

tion’s water [kg=m ], and SIN is the salinity of infiltration’s water [kg=m3 or g/l]. For the calculations in the model, the infiltration flowrate is assumed to be distributed equally on a daily basis (i.e. a constant daily flowrate), and then corrected for each month. This introduces the final term in Eq. (9), where mth is the number of days in each month. 2.1.6. Evaporation from the DS It is well known that evaporation is less from saline than from fresh water surfaces when all other conditions, apart from salinity are the same. Consequently, the general evaporation models cannot be applied directly to calculate evaporation from the DS, but have to be modified to satisfy the lake conditions. The water budget method was used for a long time as the sole estimation method for evaporation in the Jordan valley [7,16] with estimations around 1.6 m/y for the DS. Considering its limitations, the need to determine evaporation rate independently is necessary to overcome the inaccuracies in this method. Recently several studies were carried out to estimate the evaporation rate of the DS, since this data is an essential part of many other studies concerning the future development projects in the DS area. See [17–22]. However, none of the above studies took into account all the factors that are involved in the evaporation process in the DS, and none tried to combine all these effects in terms of one general formula. Here, evaporation is calculated using a modified Penman formula developed by Asmar and Ergenzinger [23] that includes all the parameters. Full details of the derivation are given in their paper. Their results using this approach varied within 3% of the results obtained by Neuman [24], Calder and Neal [17] and Oroud [22]. The generalised formula is expressed as E ¼ f ðS; u; Ta ; hÞ;

ð10Þ 2

where E is the evaporated water flux [kg=m s], S is the water salinity [kg=m3 ], u is the wind speed [m/s], Ta is the air temperature [°C] and h is the air humidity expressed in terms of air vapour pressure [mbar]. 2.1.7. Industrial intake and end-brine rejection flowrate The industries at the southern tip of the Sea pump considerable amount of the DS water to their evaporation pans (approximate ratio of 1.3:1 between Israel and Jordan). After extracting the commercial salts, the companies return roughly 40% of the original intake in

Fig. 7. Monthly flowrate fractions for water pumped from the Dead Sea and end-brine rejected into the Sea.

form of end-brine [25]. The intake reached 424:396 106 m3 =y in 1997, with a planned intake of 531 106 m3 =y by the year 2000. The monthly fractions for the intake and the endbrine rejection data are calculated based on the assumption that dimensionless fractions represent the true distribution throughout the year, and are shown in Fig. 7. The method is similar to that used in estimating the rivers and the rainfall flow. Here, the monthly intake and rejection data were obtained from APC [26]. For the intake flow pumped from the DS, the density and salinity are the same as that of the DS. The water mass flowrate equation is expressed as FWMP ¼ ðqDS  SDS Þ  FYP  xj;P ;

ð11Þ

where FWMP is the monthly brine water mass flowrate [kg/y], FYP is the yearly brine volumetric flowrate [m3 =y; qDS is the density of the DS water [kg=m3 ; SDS is the salinity of the DS water [kg=m3 or g/l], and xj;P is the monthly pumping flowrate fraction. For the end-brine water flowing into the DS, the water mass flowrate equation is expressed as FWMEB ¼ ðqEB  SEB Þ  FYEB  xj;EB ;

ð12Þ

where FWMEB is the monthly end-brine water mass flowrate [kg/y], FYEB is the yearly end-brine volumetric flowrate [m3 =y; qEB is the density of end-brine’s water 3 3 [kg=m ; SEB is the salinity of end-brine’s water [kg=m or g/l], and xj;EB is the monthly end-brine flowrate fraction. 2.2. Chemicals mass balances For centuries chemicals accumulated in the DS. The source of these chemicals was a subject of several studies, e.g. [6,7], who concluded that the saline springs were the major contributor of the salts despite their insignificant water contribution. In this model, a series of mass balances over the main chemical species are required in addition to the water mass balance to establish the

B.N. Asmar, P. Ergenzinger / Advances in Water Resources 25 (2002) 263–277

overall state of the DS at any point in time. The chemical inputs to the Sea consist of the chemicals dissolved in the water inflows into the Sea, which were discussed earlier. Until the 1950s the chemicals output from the DS was negligible. However, the chemical industries are now a major chemical output stream from the Sea, where millions of tons of salt are extracted each year. To solve the mass balances the general form of Eq. (1) is also applied here. The following general assumptions are made to simplify the problem. These assumptions hold for all the chemical species involved. 1. The mass balances are performed on ionic bases. 2. The mass balances in the Sea are performed on the following ions: Naþ1 , Mgþ2 , Caþ2 , Kþ1 , Cl1 , Br1 , 1 SO2 4 , HCO3 , which exist in significant quantities. 3. The trace materials mass balances are ignored 4. No reactions in which any of the major chemicals is generated or consumed take place. Thus, the third and fourth terms in Eq. (1) are equal to zero. 5. The change in the concentration of salts is considered one-dimensional (i.e. constant in the horizontal x–y plane and variable only vertically), which means that at a specific depth the concentration is the same all over the Sea. 6. As in the water budget, the sources that flow into the southern basin were excluded, although some of them are known to be very saline such as Sodom Springs. Two mass balance equations for each chemical, one per layer, are evaluated to establish its state in the Sea. Each equation follows the general form of Eq. (1). Fig. 8 is a block diagram illustrating the chemicals mass balances showing the input and output streams. In the top layer, the mass balance of chemical i is expressed as X dCh1;i X ¼ MChRV;i þ MChSP;i þ MChR;i dt X þ MChSV;i  Chrys1;i þ Diffi ; ð13Þ

Fig. 8. Chemicals’ mass balances over the Dead Sea.

269

where P Ch1;i is the total mass of chemical i in the top layer, MChRV;i is the P sum of the rivers chemical i mass inflow into the Sea, MCh PSP;i is the sum of the springs chemical i mass inflow, MChSV;i is the sum of the storm run-off in the side valleys chemical i mass inflow, MChR:i is the chemical i mass in rainfall over the Sea, Chrys1;i is the mass of chemicals precipitated in the specific time period [kg/month], and Diffi is the mass of chemical i transferred by diffusion from or into this layer. In the bottom layer, the water mass balance of chemical i is expressed as dCh2;i ¼ MChIN;i þ MChEB;i  MChP;i  Chrys2;i dt  Diffi ;

ð14Þ

where Ch2;i is the total mass of chemical i in the bottom layer, MChIN;i is the infiltration chemical i mass inflow, MChEB;i is the end-brine chemical i mass inflow, MChP;i is chemical i mass extracted for industry, Chrys2;i is the mass of chemicals precipitated in the specific time period [kg/month], and Diffi is the mass of chemical i transferred by diffusion from or into this layer. This term is the same as in Eq. (13), but with an opposite sign. A general formula that describes each chemical input and output into and from the DS is expressed as MChi ¼ Ci  FV ;

ð15Þ

where MChi is the mass of component i in each input or output [kg/month], Ci is the concentration of the chemical i in the stream [kg=m3 ], and FV is the monthly volumetric flowrate of each stream [m3 / month]. The chemical composition of the natural inflows into the DS is considered uniform throughout the year. However, the unavailability of the data is a major problem in this study, as the data for some sources are missing completely, and for many others they are either incomplete or are known only for a period of time less than several months, hence are not representative for the specific source. Nevertheless such data were used whenever they were available, and inevitably some data had to be estimated. As a result, the following additional assumptions were made: 1. The chemical composition of the storm run-off water is the same as the rain. This is a simplification as the floods in the valleys carry a lot of mud and salts with it, but the small amount that flows and the lack of data leaves no choice. 2. The chemical composition of the bottom springs is similar to the surrounding aquifer composition. The chemical compositions of all the DS water inflows (including the end-brines) used in this model are found in Appendix A.

270

B.N. Asmar, P. Ergenzinger / Advances in Water Resources 25 (2002) 263–277

2.3. Salts precipitation and crystallisation Crystallisation of salts in the DS was observed for a long time in the form of what was called ‘‘whitening’’ periods of the lake. In simple terms, each salt starts to precipitate when the solution is saturated with the specific salt. In the DS the process is not that simple because the other salts in the solution interact due to ion-pairing effect [27], which interferes with the crystallisation process. Several studies were performed on the precipitation of chemicals in the DS. Sass and Ben-Yaakov [28] studied the aragonite precipitation, Katz et al. [27] and Krumgalz and Millero [29] studied the gypsum. Several studies concerning halite precipitation were also conducted [27,29–33], and Gavrieli et al. [32] studied the carnallite. The precipitation of aragonite, gypsum, halite and carnallite is discussed below as their saturation degrees are either reached or close to. The studies conducted found that in the DS, aragonite (CaCO3 ) is the first material precipitated, followed by gypsum (CaSO4  2H2 O), then halite (NaCl). The carnallite (KCl  MgCl2  6H2 O) is next to crystallise, however, the latter is only expected when the brine 3 reaches a specific gravity of 1.3 g=cm [34]. Based on the above studies, this model assumes that the ions are crystallised in the following forms: 1. All Naþ1 as NaCl. 2. All Kþ1 as KCl in the combined form carnallite KCl  MgCl2  6H2 O. 3. All Br1 as MgBr2 . 4. All HCO1 3 as CaCO3 . 5. All SO2 4 as CaSO4  2H2 O. 6. Remaining Mgþ2 as MgCl2 in the combined form carnallite KCl  MgCl2  6H2 O. 7. Remaining Caþ2 as CaCl2 . To quantify the amount of salts precipitated, a standard solubility calculation over each salt is conducted [35]. In general terms, in a saturated solution of the salt Mm Nn , an equilibrium exists between the salt and its ions as Mm Nn  mM þn þ nN m :

ð16Þ

The solubility product constant of the salt Ksp is equal to the product of the concentrations of the ions ½M, ½N  produced in a saturated solution expressed in molality terms [mol Chemical/kg water], each raised to a power equal to its coefficient in the balanced equilibrium equation n

m

Ksp ¼ ½M ½N  :

ð17Þ

A saturated solution only exists when the ion product is exactly equal to Ksp. So when the ion products exceeds Ksp, a supersaturated solution exists and precipitation occurs. The degree of saturation of the salt DSt is defined as the ratio of the ion product

at any point to the ion product at saturation conditions (i.e. Ksp). If Dst is larger than one then precipitation occurs. n

DSt ¼

m

½M ½N  : Ksp

ð18Þ

However, when more than one salt in a solution are common ions, the solubility of these salts in the solution becomes less than that in pure water due to what is called the common ion effect. It can easily be understood on the basis of Le Ch^atelier’s Principle, where the existence of a common ion increases the concentration of that ion driving the equilibrium in Eq. (16) to the left. So to maintain the value of Ksp, the concentration of the other ion is reduced causing the reduction in solubility. Suppose the cation M þn is the common ion in a solution, then the amount of the salt Mm Nn precipitated is determined by the concentration of the anion N m . The saturation concentration of ½N s is calculated as  1=m Ksp ½N s ¼ : ð19Þ ½Mn The total amount (per kg water) of the salt crystals precipitated MPMN can then be calculated as MPMN ¼ ð½N   ½N s ÞMWMN ;

ð20Þ

where MWMN is the molecular weight of the salt Mm Nn . The amount of the anion precipitated is calculated as: MPN ¼

MPMN MWN : MWMN

ð21Þ

And the amount of the cation precipitated is calculated as MPM ¼

MPMN MWM : MWMN

ð22Þ

Aragonite: The crystals of aragonite were observed in 1943 and 1959 [34]. At that time the surface of the lake became milky white and the transparency was considerably reduced. Sass and Ben-Yaakov [28] assumed that the water of the DS is saturated with respect to aragonite. They calculated that the solubility of aragonite Kspa is 6:7 104 , and considered that it is independent of temperature. If the calculated degree of saturation DSta is larger than unity, precipitation occurs. Here HCO1 is the 3 determining factor as Caþ2 is the common ion and is always in excess in the DS. Gypsum: The DS is saturated with respect to gypsum, and gypsum crystals were observed in the Sea and as sediments [33]. The solubility constant of gypsum Kspg was found to be dependent on the ionic strength I [27,29]. An empirical equation to evaluate it was suggested by Katz et al. [27] Kspg ¼ 103  ð5:498  0:3709 IÞ:

ð23Þ

B.N. Asmar, P. Ergenzinger / Advances in Water Resources 25 (2002) 263–277

The ionic strength is calculated as [36] X I ¼ 0:5 ðCmi  zi2i Þ;

ð24Þ

where zii is the ion charge, and Cmi is the ion’s molality. If the calculated degree of saturation of gypsum DStg is larger than unity, precipitation occurs. Here SO2 4 is the determining factor as Caþ2 is again the common ion and is always in excess in the DS. Halite: The halite crystals on the lake surface were observed at the central part of the DS far from any effects of near shore perturbations. At that time the den3 sity of the lake was 1:23348 g=cm . Also the salt sediments found in the shallow part in the south are due to mixing of the DS water with highly concentrated endbrines from potash processing. The salt reefs are now seen in the Lisan Strait area after the level of the Sea dropped [30]. Anati [33] studied the salt precipitation in the DS, he ignored both aragonite and gypsum as their values were too small. He calculated that 77% of the halite precipitated at the periphery of the lake, mostly at the southern tip, and only 23% precipitated in the interior of the lake. For halite, the degree of saturation of halite DSth is calculated as [27]  þ1  1  2 Cl c NaCl Na DSth ¼ ; ð25Þ Ksph where Ksph is the solubility constant of NaCl, and c NaCl is the mean activity coefficient. The solubility constant depends on temperature. From the values reported in Gavrieli et al. [32] for the constant at temperatures between 25 and 50 °C, the following relationship is derived using linear regression Ksph ¼ 0:1531Tw þ 34:6842

ðr ¼ 0:987Þ;

ð26Þ

where Tw is the water temperature in the DS. Formulae to calculate the mean activity coefficient were suggested by Gavrieli et al. [32] and Katz et al. [27]. Both formulae lead to very similar results. The formula suggested by the latter authors is used in this model X ln c NaCl ¼ ln c NaClð0Þ þ ðaNaClðiÞ  Ii Þ; ð27Þ i

where c NaClð0Þ is the mean activity coefficient of pure NaCl solution, at ionic strength I. The mean activity coefficient of pure NaCl is calculated as follows ln c NaClð0Þ ¼ 0:0022I 2 þ 0:07512  0:3873; 6 < I < 13;

ð28Þ

where aNaClðiÞ is the interactive coefficient between NaCl and the ion i. For the DS, Gavrieli et al. [32] suggested formulae to calculate these coefficients at different temperatures. In the DS the temperature is mostly below 35 °C, so average formulae between the temperature 25 and

271

35 °C are used for aNaCl–MgCl2 and aNaCl–CaCl2 . A constant value of )0.0134 is used for aNaCl–KCl . aNaCl–MgCl2 ¼ 0:00272I  0:010885;

ð29Þ

aNaCl–CaCl2 ¼ 0:001785I  0:012505:

ð30Þ

Gavrieli et al. [32] found that the DS is slightly oversaturated with respect to halite, and thus the value of DSth of 1.06 is taken as the basis after which crystallisation occurs. This agrees with Beyth [37] who observed that the DS is slightly supersaturated with respect 3 to halite, and it precipitated at a density of 1.24 g=cm 3 compared to the expected value of 1.15–1.21 g=cm . Here Naþ1 is the determining factor as Cl1 is the common ion, and is always in excess in the DS. Carnallite: Under the current conditions of the DS carnallite precipitation does not take place. Gavrieli et al. [32] calculated the conditions needed for carnallite to reach its saturation point. They concluded that this point will be reached when the molal concentration of all ions except Naþ and Cl increase by a factor of 1.67. They estimated that this point will be reached in 20/ 3 years when the density of the DS reaches 1.3 g=cm . In this model, the precipitation of carnallite is not included. 2.3.1. Diffusion In the DS, diffusion can occur in two ways: 1. The diffusion from the sediments into the water bulk of the DS, which adds chemicals to the DS. 2. The diffusion between the surface layer and the bulk layer of the Sea water, which only changes the concentration in the water layers. Several studies concerning the first case were conducted. Nishri [38] used a modified Fick’s Law to calculate diffusion from the sediments in the DS. Stiller et al. [39] examined the molecular diffusion coefficient equation suggested by Lerman in 1978. They found that the molecular diffusion in the DS is hindered by tortuosity of the sediments where porosity is considered to be a measure of it. In this case, the time the salts need to diffuse from the bottom of the Sea layer is extremely long, and as the Sea basin is covered with a layer of halite and other precipitating crystals the driving force for mass transfer is reduced to a minimum as the Sea becomes saltier and even saturated resulting in no driving force for diffusion. On the other hand, diffusion of other chemicals existing below the halite layer is even slower, because these chemicals need first to diffuse through a solid phase, which is extremely slow. Therefore in this model the diffusion from the Seabed is ignored. Several studies were also conducted concerning the diffusion within the DS water layers. Larman in 1969 [34] found that mixing by diffusion is extremely slow and it would take thousands of years to achieve mixing across the pycnoline layer. Ganor and Katz [31] studied the

272

B.N. Asmar, P. Ergenzinger / Advances in Water Resources 25 (2002) 263–277

effect of double diffusion between the surface layer and the bulk of the DS. They used the double-diffusive instability test suggested by Turner [40], and found that for the DS the double diffusive instability was ruled out. Asmar [8] proposed a detailed diffusion model to describe the phenomena in the DS. The model was based on Fick’s Law, with the assumption that the Sea follows the case of diffusion of each chemical component through non-diffusing water, and ignoring the multicomponent interaction effects. He found that mass transfer by diffusion between the two layers in the DS is negligible and accounts for less than 0.1% of the total annual mass transfer for each chemical. Therefore, in this model diffusion is neglected from the chemicals mass balance equations. 2.4. Final state combined calculations The final state of the Sea is established by combining the water mass balance and the chemicals mass balances. The results of these calculations give the total mass and the salinity value in each layer, which enables the calculation of the level and the area of the DS. 2.4.1. Salinity calculation The salinity in this model is defined as the total mass 3 of salts per volume of DS solution [kg=m ]. In each layer, the salinity is calculated as Si ¼

Chi qi ; ðChi þ Wi Þ

ð31Þ

where Ch is the total mass of salts in the layer, W is the water mass, and q is the density. However, since the density of the DS is a function of salinity as shown in Eq. (32) obtained from Mero and Simon [3], an iterative solution is required to evaluate the salinity of each layer. q ¼ 1:0  105 ðTw  4:0Þ

1:865

þ 0:0777 expð0:00325SÞ :

ð32Þ

3. Simulation model The simulation model was written using a FORTRAN-based commercial package Advanced Continuous Simulation Language (ACSL). The model consists of 18 differential equations evaluating the mass balances of water and the major eight chemicals in both layers in the Sea. These equations are solved numerically using Runge–Kutta Second-Order algorithms with an integration step of 0.2 month. The time step used in formulating the model equations is one month. This choice is inevitable due to data availability limitations. A builtin implicit solver is used for the iterative calculation to evaluate the salinity in both layers. This proved to be the critical point of the program and its main problem, because it needs a fairly close initial condition to converge. The model was written using double precision calculations. Since the level is the main concern of the model, it has been chosen as the variable used to evaluate the numerical error in the calculations, especially when determining the integration step size. The model is validated by comparing its performance with measured DS levels. The year 1989 was taken as the starting point for comparison since a set of consistent data for the DS level and salinity within an acceptable time span is available. Fig. 9 shows a comparison between the measured levels of the DS [26] and the predicted levels calculated using the simulation model for the period 1989–1997 using both the average annual data, and a corrected set of data based on the actual rainfall in the validation period. As seen from Fig. 9, the predicted levels of the simulation model compare well with the measured level, so that satisfactory agreement is achieved. The differences are due to inevitable errors in estimating the evaporation, infiltration, ground water and storm run-off water. The weather conditions are not stable in the region and droughts (late 1980s) and rainy seasons (1992) are

2.4.2. Level, area and volume calculations The volume of the DS is obtained by summing the volumes of both layers of the DS. The volume of each layer is found simply dividing its mass by its density. Eq. (33) describes the volume calculation, where V is the volume of the DS [m3 ]. V ¼

ðW1 þ Ch1 Þ ðW2 þ Ch2 Þ þ : q1 q2

ð33Þ

To calculate the area and the level of the DS, two relationships between the two variables and the volume of the DS are fitted using data obtained from Harza [5], where A is the area in [km2 ], and L is the level in [m]. A ¼ 1874:11 þ 63:686V  0:5887V 2 þ 0:0019V 3 ;

ð34Þ

L ¼ 695:1334 þ 2:9865V  0:0057V 2 :

ð35Þ

Fig. 9. Comparison between monthly measured and predicted Dead Sea levels between 1989 and 1996.

B.N. Asmar, P. Ergenzinger / Advances in Water Resources 25 (2002) 263–277

273

recurring regularly in the region. Using average monthly values for estimating some parameters (e.g. flow distributions) proved to be acceptable as shown in Fig. 9, where the monthly averages compensated for the effects of the droughts and the rainy seasons. As a result the two models came to a fairly close agreement after the period of 8 years.

4. Results and discussion The simulation model is used to predict the behaviour of the DS for short periods. Using the model to predict the behaviour for long periods can only be done with caution, as the parameters used may change significantly and thus give misleading results. The results presented in Figs. 10–16 cover a period of 100 years between 1989 and 2088 assuming that the current conditions continue. The industrial intake is increased at the year 2000 by the projected amounts planned by the industries (25%). 4.1. Level and area Assuming the continuation of the current conditions, both the level and the area of the DS will continue to recede. This is shown in Figs. 10 and 11. This is easily explained as the water output from the Sea exceeds the input, and hence the net water balance results in a negative accumulation. A closer look at the two figures, however, shows that the rate of reduction in both level and area is decreasing. This is clearer in the area case where the rate of reduction is slower. The deceleration is explained as a negative feedback process. As the water in the DS drops, the salinity increases, leading to lower evaporation, and hence lower water output. 4.2. Evaporation In this model evaporation can be expressed in several forms. Fig. 12 shows evaporation expressed in two

Fig. 10. Predicted Dead Sea level between 1989 and 2088.

Fig. 11. Predicted Dead Sea area between 1989 and 2088.

forms. In Fig. 12(a) total evaporation is shown in kg/ month. Expressed in this form, evaporation can directly be included into mass balance calculations. However, in most evaporation literature, evaporation is expressed in terms of m/y or mm/y. Fig. 12(b) shows total evaporation in terms of m/month. In this form, evaporation is expressed per unit area. As shown in Fig. 12, evaporation varies considerably between winter and summer months, with considerably higher values in the summer. The salinity of the top layer decreases during winter months reaching its lowest values at the beginning of the spring. This coincides with increasing temperatures as spring starts. Both factors result in increasing evaporation. However as evaporation increases salinity also increases leading to lower evaporation later in the summer despite higher temperatures. Fig. 13 shows monthly evaporation rate [kg/month] and salinity values where the interaction between the two parameters is evident. Monthly evaporation rates were reported by Newman [24] and by Oroud [22]. The results here agree with their observation in terms of evaporation distribution. But the quantitative values differ slightly as they report results for different salinity conditions. In the longer term, as the salinity in the Sea increases, evaporation decreases. However, in Fig. 12(b) evaporation rate (expressed in height units) appears to be constant or slightly increasing. This is misleading, as the apparent increase is due to the way evaporation is expressed. The amount of evaporation in kg water decreases with increasing salinity, and at the same time the area is decreasing. The rate of area reduction is higher than that of evaporation reduction. When evaporation is expressed in height units, the amount of evaporation in kg is divided by the area, therefore the result is a false impression of rate increase. The evaporation calculation in this model assumes its occurrence on the top layer, and only when the layer is totally evaporated, evaporation starts from the bottom layer. This assumption produced results close to the values reported. It is clear that the low values found by

274

B.N. Asmar, P. Ergenzinger / Advances in Water Resources 25 (2002) 263–277

(a)

(b) Fig. 12. Predicted evaporation in the Dead Sea between 1989 and 2088: (a) mass rate; (b) (m/month).

(a)

(b)

Fig. 13. Example of predicted salinity in a 12-month period (plotted with total evaporation from the Dead Sea): (a) top layer; (b) bottom layer.

Fig. 14. Predicted salinity in the Dead Sea between 1989 and 2088.

Stanhill [20] were as a result of ignoring the existence of the thin top layer in his model. Using this method also does not require using correction factors in evaporation calculations as Harza [5] did.

Fig. 15. Example of a predicted stratification in the Dead Sea in a 12month period.

4.3. Salinity and density In this model the density of the DS is calculated as a function of the relevant salinity as shown in Eq. (32),

B.N. Asmar, P. Ergenzinger / Advances in Water Resources 25 (2002) 263–277

275

upon the assumption that the two-layer system exists as long as the salinity is lower in the top layer. When salinity is equal there is only one layer. Any additional salts input during this condition is corrected numerically by a flow of salts to the bottom layer. This leads to the fixed salinity condition seen between July and November in Fig. 13(a). As can be seen from Fig. 15, this model predicts seasonal stratification, and that a condition of two layers exists in the Sea in winter, spring and most summer months. This result agrees with the behaviour of the DS reported by Gavrieli et al. [32], who noted that normally seasonal stratification takes place in the DS, which ends usually when overturn takes place in the autumn. Fig. 16. Example of predicted crystallisation of chemicals in the Dead Sea in a 12-month period.

therefore the behaviour of density is tied to and can be expressed through the behaviour of salinity. The salinity in the top layer oscillates throughout the year between its lowest values in spring to highest values late summer when there is only one layer (see the following section). As shown in Fig. 14 the changes in the salinity of the top layer are significant, and ignoring them results in large errors in estimating the evaporation. Using only the salinity in the bottom layer underestimates evaporation. Fig. 13 illustrates the salinity profile in both layers of the DS. In the top layer the salinity values oscillate but with increasing both the minimum and the maximum peaks. This behaviour is strongly linked to both evaporation as more salts accumulation leads to less evaporation, and to the salinity of the bottom layer which when increases means that top layer will exist for a longer period until its salinity reaches the salinity of bottom layer to cause overturn. In the bottom layer, the situation is simpler; the salinity is always increasing as the quantity of water in the bottom layer is always decreasing. However, a closer look at Fig. 13(b) reveals a slight decrease in salinity in spring time, this is due to the fact that no endbrine flows into the Sea during winter times and there is no industry intake also. It should be noted here, see Fig. 13(a), that the condition of a one-layer system can be seen when the salinity remains constant while the evaporation rate decreases. Comparisons of the model results with published data were close. For example Vengosh and Rosenthal [41] reported a salinity of 330.6 g/l, the model values were 331.5–331.9 g/l. Also Encyclopaedia Britannica [42] reported a value of 332 g/l, the model predictions were 332–332.4 g/l.

4.5. Crystallisation and salts precipitation The degree of saturation in each layer is used to detect the occurrence of crystallisation. As explained earlier it is assumed that crystallisation occurs when the water is saturated with aragonite (DSta ¼ 1) and gypsum (DStg ¼ 1), and slightly oversaturated with halite (DSth ¼ 1:06). In the case analysed, seasonal aragonite and gypsum precipitation occurs during late summer and autumn months as shown in Fig. 16, where the disappearance of the top layer results in higher local salinity producing saturation conditions. This is due to the fact that the input streams contain relatively higher quantities of carbonate and sulphate which result in higher local concentrations. On the other hand the halite degree on saturation is far below 1.06 and thus its crystallisation does not happen. The situation is different in the bottom layer where all three chemicals are far below saturation points, and thus no crystallisation is expected. These results confirm what was noticed in the DS of whitening phenomenon, and gypsum crystallisation reported by Sass and Ben-Yaakov [28], Katz et al. [27], Steinhorn [30] and Krumgalz and Mellero [29]. But they contradict the observations of Anati [33] concerning halite crystallisation (as this model shows), that halite is not close to saturation conditions. This can be explained as follows, in this model it is assumed that at a certain level concentration is equal throughout the Sea, which does not account for local differences in concentration at specific areas such as close to the mouth of River Jordan. Anati [33] noticed that most halite crystallisation occurs at the southern tip of the DS close to the endbrine flow. In this area local salinity is much higher and can cause saturation conditions that lead to this local crystallisation.

4.4. Number of layers

5. Conclusions

The salinity conditions in the Sea are used as an indication of the number of layers. Calculations are based

1. A mathematical model that describes the dynamic behaviour of the DS has been developed based on

276

B.N. Asmar, P. Ergenzinger / Advances in Water Resources 25 (2002) 263–277

simple mass balances on water and the major chemical species in the Sea. The model has been implemented with an ACSL code and validated by comparing its performance with recorded DS data. The model is highly versatile and can be modified easily to reflect specific needs, or for further development. The model can predict the future behaviour of the DS based on the current conditions and under Table 1 Flowrates of the inputs into the Dead Sea Source

all sorts of disturbances including a proposed Dead Sea–Red Sea Canal. The lack of credible records for the flowrates of rivers, springs and valleys, and the unavailability of other physical and chemical data in the region is the main obstacle that hinders accurate predictions. However, to overcome this, the method using average distribution values to describe the flowrates of the rivers and the rainfall proved to be successful as the results obtained using the model were close to the measured values. An accurate model to describe evaporation behaviour is essential to achieve accurate results. Ignoring the existence of a thin top layer leads to large errors in the estimation of evaporation. This can be easily seen in the results reported by Stanhill [20], which underestimated evaporation by almost 30%. Assuming that current conditions continue to exist in the DS region, the model predicts a continuous decline in the DS level. The decline rate is decelerating as the area and evaporation rate are both decreasing. Under these conditions, the model shows that the salinity of both layers will continue to increase, and seasonal crystallisation of gypsum and aragonite will continue. The model predicts current seasonal stratification conditions. It predicts that this situation will eventually turn into permanent stratification condition, as the continuous decline in area and evaporation will reach a point at which a permanent thin top layer persists. The time this situation is reached cannot be predicted accurately as the conditions in the area,

2.

3.

Flowrate 3

Rivers [m =y] Jordan Mujib Zerqa-Ma’in-Zara al-U’zaymi Ibn Hammad Kerak Rainfall [mm/y] Storm run-off [m3 =y] Springs [m3 =y] W. Haditha Ain Thra’a Ein Fashkha/Zuqim – Main Ein Fashkha – Minor Qaneh-Samar/Ghuweir-Tureiba Ein Gedi Yesha Infiltration [m3 =y] Pumping [m3 =y] End-brine [m3 =y]

400 106 42:4 106 25 106 3:252 106 4:75 106 2:86 106 77:19 8:961 106

4.

8:205 106 1:389 106 39:447 106 12:623 106 13:149 106 3 106 3 106 33 106 424:396 106 169:758 106

5.

Table 2 Chemical composition of the inputs into the Dead Sea Source Rivers [mg/l] Jordan Mujib Zerqa-Ma’in-Zara al-U’zaymi Ibn Hammad Kerak

Naþ

Caþ2

Kþ

765 316 428 734 124 96

51 40 56 96 10 8

Mgþ2

330 106 139 238 86 67

276 108 29 50 39 30

Cl

Br

1630 573 709 1215 253 197

2 0.3 7 12 5 4

SO2 4 452 279 198 339 119 93

HCO 3 436 140 301 516 264 205

Salinity 4032 1562.3 1867 3200 900 700

Rainfall [lg/l] 6.493 Springs [mg/l] W. Haditha Ain Thra’a Ein Fashkha/Zuqim – Main Ein Fashkha – Minor Qaneh-Samar/Ghuweir-Tureiba Ein Gedi Yesha

0.411

4.569

124 39 868 6841 771 30 27 800

10 8.6 96 1083 120 – 3500

86 76 343 2789 365 – 10 300

19 800

–

56 056

2.479 39 14.6 381 5448 707

14.004

–

11.430

–

39.386

– 18 700

253 88 2820 32 012 3848 45 11 753

5 0 99 714 84 – 2206

119 85 71 186 91 – 1060

264 183 260 263 293 – 142

900 494.2 4938 49 336 6280 75 181 238

8141

66 125

–

1971

152

152 245

–

–

–

Groundwater [mg/l] End-brine [g/l] 3.849

4.208

34.244

87.867

326.516

456.684

B.N. Asmar, P. Ergenzinger / Advances in Water Resources 25 (2002) 263–277

especially industrial production and hence water intake, are changing. Acknowledgements The authors thank Dr. Paul Langston for his critical reading and editing of the manuscript. They also thank Prof. Marc Parlange and the anonymous reviewers for their useful comments and suggestions. Appendix A. Model parameters See Tables 1 and 2.

References [1] Carmi I, Gat JR, Stiller M. Tritium in the Dead Sea. Earth Planet Sci Lett 1984;71:377–89. [2] Steinhorn I. The disappearance of the long-term meromictic stratification of the Dead Sea. Limnol Oceanogr 1985;30:451–72. [3] Mero F, Simon E. A daily simulation model for evaluation future Dead Sea levels. In: Diskin M, editor. Scientific basis for water resources management (Proceedings of the Jerusalem Symposium, September 1985), vol. 153. IAHS Publication; 1985. p. 265–76. [4] Simon E, Mero F. The simulation of 121-year series of natural Dead Sea inflows. In: Diskin M, editor. Scientific basis for water resources management (Proceedings of the Jerusalem Symposium, September 1985), vol. 153. 1985. p. 381–92. [5] Harza, Jordan Rift Valley integrated development study. Interim reports, 1996. [6] Bentor YK. Some geochemical aspects of the Dead Sea and the question of its age. Geochim Cosmochim Acta 1961;25:239–60. [7] Abed AM. Geology of the Dead Sea: waters, salts and evolution. Amman, Jordan: Dar Al-Arqam; 1985 [in Arabic]. [8] Asmar BN. Dynamic behaviour of the Dead Sea water body and the effects of the Dead Sea–Red Sea Canal. Ph.D. Dissertation, Freie Universit€ at Berlin, 1998. [9] Central Water Authority (CWA) – HKJ, Technical Report No. 33, Amman, Jordan, 1964. [10] Meterological Department (MD) – HKJ, Monthly data, Amman, Jordan, 1995. [11] Arad A, Michaeli A. Hydrogeological investigations in the western catchment of the Dead Sea. Isr J Earth Sci 1967;16:181–96. [12] El-Naqa A. Hydrological and hydrogeological characteristics of Wadi el Mujib catchment area, Jordan. Environ Geol 1993;22:257–71. [13] Stiller M, Chung YC. Radium in the Dead Sea: a possible tracer for the duration of meromixis. Limnol Oceanogr 1984;29:574–86. [14] Chan LH, Chung Y. Barium and radium in the Dead Sea. Earth Planet Sci Lett 1987;85:41–53. [15] Salameh E. Country relevant water resources. Dirasat – Nat Resour 1983;10:95–108. [16] Stanhill G. Evaporation from Lake Tiberias: an estimate by the combined water balance-mass transfer approach. Isr J Earth Sci 1969;18:101–8. [17] Calder IR, Neal C. Evaporation from saline lakes: a combination equation approach. Hydrol Sci – Journal-des Sciences Hydrologiques 1984;29:89–97.

277

[18] Salhorta AM, Adams EE, Harleman DRF. Effect of salinity and ionic composition on evaporation analysis of Dead Sea evaporation ponds. Water Resour Res 1985;21:1336–44. [19] Salhorta AM, Adams EE, Harleman DRF. The alpha, beta, gamma of evaporation from saline water bodies. Water Resour Res 1987;23:1769–74. [20] Stanhill G. Changes in the rate of evaporation from the Dead Sea. Int J Climatol 1994;14:465–71. [21] Oroud IM. Evaluation of saturation vapour pressure over hypersaline water bodies at the southern edge of the Dead Sea – Jordan. Sol Energy 1994;53:497–503. [22] Oroud IM. Effect of salinity upon evaporation from pans and shallow lakes near the Dead Sea. Theoret Appl Climatol 1995;52:231–40. [23] Asmar BN, Ergenzinger P. Estimation of evaporation from the Dead sea. Hydrol Process 1999;13:2743–50. [24] Neumann J. Tentative energy and water balances for the Dead Sea. Bull Res Counc Isr 1958;7G:137–63. [25] Arab Potash Company Ltd., The Arab Potash Company into the next century with confidence, Amman, Jordan, 1996. [26] Arab Potash Company Ltd., Personal communications, 1997. [27] Katz A, Stranisky A, Taitel-Goldman N, Beyth M. Solubilities of gypsum and halite in the Dead Sea and its mixtures with seawater. Limnol Oceanogr 1981;26:709–16. [28] Sass E, Ben-Yaakov S. The carbonate system in hypersaline solutions: Dead Sea brines. Mar Chem 1977;5:183–99. [29] Krumgalz BS, Millero FJ. Physico-chemical study of the Dead Sea waters, III. On gypsum saturation in Dead Sea waters and their mixtures with Mediterranean Sea water. Mar Chem 1983; 13:127–39. [30] Steinhorn I. In situ precipitation at the Dead Sea. Limnol Oceanogr 1983;28:580–3. [31] Ganor J, Katz A. The geochemical evolution of halite structures in hypersaline lakes: the Dead Sea, Israel. Limnol Oceanogr 1989;34:1214–23. [32] Gavrieli I, Starinsky A, Bein A. The solubility of halite as a function of temperature in the highly saline Dead Sea brine systems. Limnol Oceanogr 1989;34:1224–34. [33] Anati DA. How much salt precipitates from the brines of a hypersaline lake? The Dead Sea as a case study. Geochim Cosmochim Acta 1993;57:2191–6. [34] Steinhorn I, Gat JR. The Dead Sea. Sci Am 1983;249: 84–91. [35] Brady JE, Humiston GE. General chemistry: principles and structures. 4th ed. New York, USA: Wiley; 1986. [36] Skoog DA, West DM, Holler FJ. Fundamentals of analytical chemistry. 7th ed. Fort Worth, TX: Saunders College Publishing; 1996. [37] Beyth M. Recent evolution and present stage of Dead Sea brines. In: Nissenbaum A, editor. Hypersaline brines and evaporatic environments (Proceedings of the Bat Sheva Seminar on Saline Lakes and Natural Brines). Dev. Sedimentol. 1980;28: 155–66. [38] Nishri A. The geochemistry of manganese in the Dead Sea. Earth Planet Sci Lett 1984;71:415–26. [39] Stiller M, Kaushansky P, Carmi I. Recent climatic changes recorded by the salinity of pore waters in the Dead Sea sediments. Hydrobiologia 1983;103:75–9. [40] Turner JS. Buoyancy effects in fluids. Cambridge, UK: Cambridge University Press; 1973. [41] Vengosh A, Rosenthal E. Saline groundwater in Israel: its bearing on the water crisis in the country. J Hydrol 1994;156: 389–430. [42] Encyclopaedia Britannica – Online Service, 1995. Available from: http://www.eb.com.