A Model for Analyzing Lake Water Acidification on a Large Region Scale - Part 1: Model Structure

NOT FOR QUOTATION WITHOUT PERMISSION OF THE AUTHOR

A MODEL FOR

ANALYZING

LAgE WATER ACIDIFICATION ON A LABGE REGIONAX,

SCALE

Juha

K M

Maximilian Posch Lea Kauppi

D e c e m b e r

T985

CP-85-48

Cbllaboratiue P a p e r s r e p o r t work which has not been performed solely

at

the International Institute f o r Applied Systems Analysis and which has received only limited review. Views o r opinions expressed herein do not necessarily represent those of the Insti- tute, its National Member Organizations. o r other organizations supporting the work.

INTERNATIONAL INSIlTUTE FOR APPLIED SYSIZMS ANALYSIS 2361 Laxenburg, Austria

AUTHORS

Juha KBmki and Lea Kauppi a r e both from the Water Research Institute, National Board of Waters, P.O. Box 250, SF-00101 Helsinki, Finland. They were formerly with the International Institute f o r Applied Systems Analysis, Laxenburg, Austria.

Maximilian Posch is with the International Institute f o r Applied Systems Analysis, Laxenburg, Austria.

PREFACE

The IIASA "Acid Rain1' P r o j e c t s t a r t e d in 1983 in o r d e r to provide t h e European decision makers with a tool which c a n b e used t o evaluate policies f o r controlling acid r a i n . This modeling e f f o r t i s p a r t of t h e official cooperation between IIASA and t h e UN Economic Commission of E u r o p e (ECE)

.

The IIASA model c u r r e n t l y contains t h r e e linked compartments: Pollu- tion Generation, Atmospheric P r o c e s s e s and Environmental Impacts. Each of t h e s e compartments c a n b e filled by d i f f e r e n t substitutable submodels.

The submodels c u r r e n t l y available are Energy Pathways and S u l f u r Emis- sions, t h e EMEP Long Range T r a n s p o r t Model. F o r e s t Soil pH and Lake Aci- dity. In addition, two submodels are under development: t h e NO, Emissions submodel and t h e Direct F o r e s t Impacts submodel. The f i r s t version of t h e Lake Acidity submodel was p r e s e n t e d in September 1984 in a UNESCO-IHP Workshop in Uppsala, Sweden. Since t h e n s e v e r a l changes have been imple- mented following t h e advice of e x p e r t s . This p a p e r d e s c r i b e s t h e Lake Aci- dity model s t r u c t u r e

as

i t stands in December 1985.

Leen Hordijk

Acid Rain P r o j e c t L e a d e r

ACKNOWLEDGEMENTS

We wish to thank Professor G . Hornberger, Dr. N . Christophersen, Dr.

J . Cosby, Dr. P . Whitehead and Dr. H . Grip f o r providing encouraging com- ments and valuable advice.

-

^vii

-

ABSTRACT

The International Institute f o r Applied Systems Analysis i s developing a computer model which c a n b e used by decision makers t o evaluate policies f o r controlling t h e impact of acid r a i n in Europe. A s p a r t of t h i s t a s k , a sim- p l e dynamic model h a s been developed f o r describing t h e p r o c e s s e s leading t o acidification of s u r f a c e waters. The simulation model is c o n s t r u c t e d of s e v e r a l modules, e a c h of them providing a n overview of a p a r t i c u l a r a s p e c t of l a k e acidification. The meteorologic module calculates t h e amount of water and deposition e n t e r i n g t h e soil o r t h e l a k e d i r e c t l y e a c h month. The IIASA soil acidity submodel accounts f o r t h e soil solution chemistry. A sim- p l e hydrologic method i s applied f o r simulating t h e routing of i n t e r n a l flows s o t h a t t h e convective flow of ions can be estimated. The lake r e s p o n s e i s calculated according

to

t h e equilibrium reactions of inorganic c a r b o n species. These modules

are

described in t h i s p a p e r . In p a r t 2 t h e applica- tion of t h e model o n a l a r g e regional s c a l e will b e described. Monte Carlo techniques will b e used

to

determine those r a n g e s and combinations of input values t h a t p r o d u c e a n a c c e p t a b l e p r e s e n t day lake acidity distribution, when t h e model is driven by a specified deposition.

TABLE OF

_CONTENTS

Authors Preface

Acknowledgements Abstract

Table of Contents 1. Introduction

2. Strategy f o r Model Application 3. Model Struc!.ure

3 . 1 Meteorologic Module 3 . 2 Hydrologic Module 3.3 Soil Chemistry Module 3.4 Lake Module

4. Model Testing List of Symbols References Appendix

A MODEL FOR ANALYZING

LAKE

WATER ACIDIFICATION ON A LARGE REGIONAL SCALE

PBBT

1: MODEL STBUCTUBE

Juha K h S r i , Maximilian Posch and Lea Kauppi

The harmful e f f e c t s on s u r f a c e waters resulting from acidic deposition h a v e been well documented in various p a r t s of t h e N o r t h e r n Hemisphere.

The c a u s a l relationships leading t o f r e s h w a t e r acidification a r e , however, complex and difficult

to

quantify. Hydrochemical models have provided one way of quantifying and integrating various p r o c e s s e s in t h e e n t i r e catch- ments. Models have been used f o r simulating daily variations of water qual- ity in s t r e a m s , caused by variations in deposition, as well as in catchment hydrology and meteorology (8.g. Christophersen et al. 1982). However, many of t h e s e modeling a p p r o a c h e s have been r e g a r d e d as tools f o r d a t a evaluation r a t h e r t h a n tools f o r predicting long-term acidification of t h e catchments.

Recently t h e need

to

provide estimates of potential f u t u r e impacts of a c i d i c deposition h a s been emphasized. Scientific information c a n a s s i s t in evaluating policies f o r emission c o n t r o l by describing quantitative conse- quences of a l t e r n a t i v e scenarios. T h e r e

are at

p r e s e n t t h r e e basic methods f o r making projection of f u t u r e water chemistry of sensitive a r e a s . The f i r s t i s an empirical a p p r o a c h which allows t h e estimation of f u t u r e steady-

state

chemical composition of l a k e s resulting from changes in loading of s t r o n g a c i d s on t h e basis of o b s e r v e d relationships in p r e s e n t conditions (8. g

.

Henriksen, 1980). The second method utilizes complex, p r o c e s s o r i e n t e d submodels f o r catchment hydrology, canopy chemistry, soil chem- i s t r y as well as f o r s t r e a m and l a k e water quality

to

provide a scientific Link between acidic deposition and l a k e acidification (8.g. Chen e t al. 1983). The t h i r d method defines predictive algorithms t h a t l a r g e l y r e t a i n t h e simplicity

of t h e empirical models but t h a t have mechanistic process oriented explana- tions incorporated in t h e i r s t r u c t u r e ,

to

allow

a

theoretical basis f o r estab- lishing confidence in t h e estimates (Cosby

et

al. 1985).

Simple models can b e applied

as

p a r t of a regionalized model s t r u c t u r e . A t t h e beginning of t h e development of IIASAss lake acidification model, no suitable models

were

available f o r this purpose. Therefore,

a

number of existing process descriptions

were

simplified, modified

to

monthly time s t e p and finally linked together

to

form

a

simple working method f o r t h e evalua- tion of lake acidification. This method will b e used as

a

component of t h e IIASA Regional Acidification Information and Simulation (RAINS) model sys-

tem.

The objective of this model system is

to

provide decision makers with

a

tool, which can assist in evaluating policies f o r controlling t h e impact of acidic deposition in Europe. With t h e lake acidification model, t h e response of sensitive lake

areas to

alternative energy p a t t e r n s and emission reduc- tion

measures

will b e analyzed.

In t h e IIASA acid r a i n study. t h e energy-emission model generates s u l -

f i r

emission scenarios f o r Europe, assuming optional programs f o r energy development and emission control (Alcamo

et al.,

1985). Computed emissions

are

converted into

s u w r

deposition scenurios in t h e atmospheric transformation compartment of t h e model system by using t h e long-range t r a n s p o r t model f o r a i r pollutants developed within t h e EMEP-program (see Eliassen and Saltbones. 1983). Sulfur deposition is finally transformed into a n estimation of a d d s t r e s s (Kauppi

et

al. 1985). which forms t h e basic input data f o r environmental impact submodels. Any long-range t r a n s p o r t model from t h e Atmospheric Processes compartment may be linked

to

all models in t h e Environmental Impact compartment (Figure 1). Presently, t h e EMEP sulfur t r a n s p o r t model forms t h e driving f o r c e of t h e f o r e s t soil aci- dity and lake acidity submodels.

The study includes model development and model application. In this paper, t h e o v e d strategy f o r t h e model application

as

w e l l as t h e present model s t r u c t u r e and some model experiments

are

described. The second p a r t includes t h e model application with model results f o r different lake regions in Europe.

2. STRATEGY

FOR

PODEL APPLICATION

In simulation models of environmental systems, based upon physical, chemical and/or biological mechanisms, t h e model s t r u c t u r e (the theory), model inputs, initial conditions

as

well

as

parameter values all necessarily include uncertainty. On

a

regional scale, t h e uncertainty is even g r e a t e r . It h a s been emphasized

in

several studies t h a t t h e analysis of models should concentrate on identifying ranges of inputs. r a t h e r than on traditional parameter estimation (e.g. Fedra 1983; Hornberger and Cosby 1985). Horn- b e r g e r and Cosby (1985) have successfully investigated t h e ability of

a

sim- ple catchment model of sulfate dynamics,

run

with t h e most probable ranges of inputs,

to

produce t h e d i s t r i b u t i o n of measured, present-day

stream

sulfate concentrations f o r t h e White Oak Run region in Virginia, USA. More- over, t h e same procedure has been applied

to

determine which parameter combinations produced acceptable predictions of soil properties when t h e

(a) COMPARTMENTS

r ^---,

^{I I---}

I ¹

^{r--- I 7}

I

^I

I

pollutant I

1

Atmospheric

I

Environmental I

I

generation

I

processes

I

^impact I

I

I

I

I I I

L-.----J L---A

L-.----l

(bl SUBMODELS

-

Forest soil

acidity

Figure 1. Links in t h e IIASA Acid Rain Model.

Sulfur emissions

model w a s driven by t h e o b s e r v e d

stream

chemistry (Cosby e t al. 1985). The basic principle of t h e i r Regionalized S e n s i t i v i t y Analysis h a s been

to

use Monte Carlo techniques f o r simulating t h e r e s u l t of temporal evolution of a number of individual catchments with varying inputs. H o r n b e r g e r and Cosby (1985) suggest t h a t

" ...

i t is a p p r o p r i a t e and. in a sense, e a s i e r

to

p r e d i c t conditions relating

to

d i s t r i b u t i o m of uncertain systems than i t i s t o p r e d i c t conditions f o r a n individual example of t h e system."

Our s t r a t e g y f o r t h e application of t h e model f o r t h e regional s c e n a r i o analysis h a s two distinct levels. A t t h e f i r s t level t h e catchment model i s a b l e

to

analyze changes o v e r time in t h e chemistry of a lake. The model c a n b e

run

f o r any known system f o r which r e l e v a n t lake. catchment and soil information i s available. Examples of model experiments

at

t h e catchment level a r e shown and discussed in c h a p t e r 4 of t h i s p a p e r . In c o u r s e of t h e regional m o d e l development t h e model i s incorporated into a l a r g e r s t r u c - t u r e which s c a l e s t h e s c e n a r i o s from individual systems up

to a

regional level. W e apply t h e Monte C a r l o p a r a m e t e r estimation p r o c e d u r e

to

t h e regional l a k e acidification assessment in o r d e r

to

model regional l a k e water quality distributions.

The Monte Carlo method i s a trial-and-error p r o c e d u r e f o r t h e solution of t h e inverse problem, i.e. f o r estimating t h e poorly known input values from t h e r e q u i r e d output. The basic s t e p s of t h i s estimation p r o c e d u r e

are as

follows (Fedra. 1983):

acidity

*

EMEP sulfur tnnspon

(1)

For t h e s t r u c t u r e of t h e simulation model, performance c r i t e r i a are formulated describing t h e expected satisfactory behavior of t h e

model.

based on available data.

(2)

To estimate all unknown input values, allowable ranges or probability functions are defined f o r them.

(3)

The Monte Carlo program then randomly samples t h e parameter vectors from t h e allowable ranges, runs t h e simulation

m o d e l

through a selected period and finally tests f o r violations of constraint conditions a f t e r t h e simulation.

(4)

This process i s repeated f o r

a

large number of trials.

In our regional model, t h e Monte Carlo method

is

used to determine t h e com- binations of inputs t h a t produce a n acceptable distribution of output vari- ables, observed in t h e study region. For

all

inputs, ranges are chosen broad enough so t h a t any reasonable value

for

a n input could be

selected

in t h e Monte Carlo simulations. In a n ideal case,

t h e r e

should be

a priori

information on t h e shape of distributions

of all

parameters, initial condi- tions and catchment characteristics. In reality, however. this

is

not t h e case and several inputs have to be

selected

from uniform distributions.

Monte Carlo simulations are then

carried

out by randomly selecting a set of input values

f r o m

within these designated ranges and integrating t h e equations f r o m

¹⁹⁶⁰

on using this particular set

of

values.

A

subset of accepted input values corresponding to t h e actual observed present-day frequency distribution in

¹⁹⁸⁰

in each lake region,

is

obtained. The

model.

now calibrated f o r present-day conditions, provides a tool f o r future pro- jections

of

regional lake water acidity. Assuming t h a t t h e set

of

input values obtained in t h e calibration i s representative

of

real catchments in t h e study region, this ensemble can be used

for

t h e scenario analysis

of

t h e response

of

lake systems to different patterns in acidic deposition.

&

a result, this procedure will give new frequency distributions

for m o d e l

output variables which will s e r v e as impact indicators.

A

critical

aci-

dity level has to be established, and according to this level, t h e number of lakes or t h e lake surface

area threatened in each lake district can be

estimated

for

any y e a r and any energy-emission scenario desired.

Agricultural management practices as

w e l l

as o t h e r sources

of

pollu-

tion in t h e catchment have an overwhelming influence on t h e ionic balance

of

s u r f a c e waters. For example, liming and fertilizing practices cause

eutrophication, high alkalinity and high pH-values in those

surface

waters

receiving agricultural runoff. Therefore

t h e

evaluation of impacts will b e

r e s t r i c t e d to sensitive lake regions t h a t do not receive any significant

anthropogenic input except atmospheric pollutants. The behavior

of

each

sensitive lake-district

is

determined by distributions

of

regional soil, lake

and catchment characteristics

as w e l l as by deposition, precipitation and

a i r temperature.

A l l

information regarding this input

data

i s s t o r e d into a

computerized format. Simulations with

⁷⁰

y e a r s t i m e span

(1960-2030)

are

performed

for

each sensitive lake region in Europe using a monthly time

step.

3.

MODEL STRUCTURE

During t h e development of t h e study o u r modeling philosophy has been

to

use a simplified approach which is warranted for a broad geographical scope. The objective has been

to

retain t h e simplicity of t h e model but still have f e w physically realistic processes incorporated in i t s s t r u c t u r e

to

allow

a

theoretical basis f o r assessing confidence in t h e scenarios. The model consists of f o u r modules t h a t

are

linked together as shown in Figure 2. The processes considered in each module are summarized in Table 1. The meteorologic module regulates t h e input f l o w s of

water

and deposition

to

t h e soil and directly

to

t h e lake. The hydrologic and soil chemistry modules t o g e t h e r determine t h e flow of ions leaching f r o m t h e t e r r e s t r i a l catchment

to

t h e lake. N e w equilibrium concentrations in t h e lake

water are

then com- puted in t h e lake module.

Table 1. P r o c e s s e s considered in t h e IIASA lake acidity model.

P r o c e s s Reference

Meteorology :

Partitioning between snow and r a i n Snow melt

Release of deposition from snowpack Ffydrology :

Evapotranspiration

Percolation from upper t o lower r e s e r v o i r L a t e d f l o w

Soil c h e m i s t r y : Carbonate weathering Silicate weathering Cation exchange

Aluminum equilibrium with gibbsite Lake :

Inorganic carbon equilibrium

Shih

et

al. 1972 Chow 1964

Johannessen and Henriksen 1978

Christophersen e t al. 1984 Chen

et

al. 1982

Chen

et

al. 1982 Ulrich 1983 Ulrich 1983 Ulrich 1983

Christophersen

et

al. 1982 Stumm and Morgan 1981

3.1. Meteorologic Module

The purpose of t h e meteorologic module is

to

determine t h e volume of water and proportion of deposition entering t h e catchment within one time s t e p , T . The division

of

t h e

total

precipitation,

P&,

into r a i n ,

P : ,

and snow,

P : ,

is accomplished by Eqs.1a.b using threshold tempsratwes,

Ts

and

T,

and t h e mean monthly temperature T f (Shih

et

al., 1972, Christophersen

et

al. 1984)

l NPUT

MODEL

1 I

DRIVING VARIABLES

I

STATE VARIABLES

'7 1

ACID STRESS INFILTRATION

HYDROLOGIC

I

I

CONVECTIVE FLOW OF IONS

I

I

LAKE MODULE OUTPUT

Figure 2. The overall s t r u c t u r e

of

t h e IIASA lake acidity model.

Snow accumulates, whereas

all

rain during t h e winter i s assumed to run through t h e snowpack and e n t e r t h e

soil.

Also, t h e melting

of

t h e snowpack, mT,

is

set to be proportional t o t h e mean monthly temperature above t h e threshold temperature Ts, using a melting r a t e coefficient

@

(Chow, 1964;

Chen et

al.

1982)

The snowpack, S P ~ , is obtained by summing t h e individual P,f-values

and subtracting t h e m t v a l u e s , as long as SP? stays above z e r o (Eqs.

3

a

and b; h e r e and t h e r e a f t e r primes

refer

r e f e r to a n intermediate step,

which is used

for

computational purposes only)

SP'

=

spT-l

+

P,T s P T

=

SP'

-

mT

Deposition i s assumed

to

accumulate when snow accumulates, t h e

rate

f o r deposition accumulation being D:. The

same

fraction of total deposi- tion, D&, as of total precipitation i s retained in t h e snowpack each month

as

accumulated deposition, DPT (Eqs. 4a.b)

During t h e snowmelt, t h e

rate

f o r t h e release of deposition from t h e snowpack, D,T, i s assumed

to

b e

t w o

times higher than meltwater (Eqs.5 a and b). The fractionation effect observed during t h e snowmelt (Johannessen and Henriksen. 1978) implies t h a t most of t h e impurities in t h e snowpack

are

found in t h e f i r s t m e l t w a t e r .

LPT

=

DP'

-

^D: (5b)

The deposition entering t h e soil o r t h e l a k e will b e called acid stress, asT, in t h e sequel

3.2. Hydrologic Module

The flowpaths of r a i n and snowmelt water through t h e t e r r e s t r i a l sys-

t e m are

important f a c t o r s in determining t h e susceptibility of lakes

to

aci- dification by atmospheric deposition (Chen

et

al. 1982). To provide

a

method f o r simulating t h e routing of internal flows,

a

simple hydrologic model is applied. A combined version of hydrologic models, Birkenes model and ILWAS model, presented by Christophersen

et

al. (1982) and Chen

et

al.

(1982) respectively, i s used.

The IIASA framework

sets

t h e prerequisite of

a

l a r g e spatial scale.

The ILWAS model is highly mechanistic and contains descriptions of t h e p r o c e s s e s both in t h e canopy and in s e v e r a l soil layers. There is thus r a t h e r little curve-fitting involved. The Birkenes model is very s i t e specific and must b e calibrated against t h e typical f e a t u r e s of

a

given catchment

before i t can b e applied. For t h e IIASA context, t h e simple two-layer s t r u c -

ture

of t h e Birkenes model i s applied. Most of t h e physical descriptions of t h e processes for routing t h e

water

through these

two

l a y e r s and o u t of t h e system

are

simplified

f r o m

t h e ILWAS model.

The t e r r e s t r i a l catchment i s vertically segmented into snowpack and

t w o

soil l a y e r s (A- and B-reservoirs; Figure 3). The

A-reservoir

is defined

to

b e identical with t h e uppermost 0.5

meter

soil l a y e r modeled by t h e soil impact model (Kauppi

et

al. 1985), which i s used

later to

account for soil solution chemistry. Physically, t h e f l o w

f r o m

t h e u p p e r r e s e r v o i r can b e thought of as quickflour, which drains down t h e hillsides as piped f l o w

or

fast throughflow and e n t e r s t h e brooks directly (Christophersen and Wright. 1981). This

water

is mainly in contact with humus and t h e u p p e r mineral layer. The B-reservoir in t h e model provides t h e baseflow, which presumably

comes

largely

f r o m

d e e p e r (> 0.5 m) soil l a y e r s (c.f. Christo- phersen and Wright, 1981).

DEEPER PARTS

...-...

OF SOILS

B

LAKE

Figure 3. Simulating t h e routing of

water

flows

f r o m

precipitation

to

lake discharge;

see

t e x t for explanation (modified

f r o m

Chris- tophersen

et al.

1984).

The basic assumption governing t h e soil hydraulics is t h a t rainfall o r meltwater i n f i l t r a t e s as a whole into t h e A-reservoir (c.f. Christophersen and Wright, 1981). Evapotranspiration,

ET,

i s

set

proportional t o t h e mean monthly t e m p e r a t u r e ,

TT,

above 0 C , using a evapotranspiration coeff i- c i e n t E (c.f. Christophersen

et

al. 1984)

The a c t u a l evapotranspiration

rate

is assumed t o b e equal t o t h e poten- t i a l from t h e A-reservoir, and if A becomes empty, from t h e B-reservoir.

The intermediate

water

balance i s given by Eq.12. which considers t h e

water

fluxes between t h e A-reservoir, t h e atmosphere and t h e snowpack

The percolation of

water

into t h e B-reservoir i s controlled by t h e

max-

imum possible percolation r a t e , Qp('), t h e

water

volume available in t h e A- r e s e r v o i r . Q f ) , and t h e s p a c e l e f t in t h e B-reservoir, Q i 3 ) . Any o n e of t h e s e t h r e e f a c t o r s c a n b e

a

limiting f a c t o r f o r percolation. T h e r e f o r e t h e a c t u a l percolation r a t e , Q;, i s

set

equal

to

t h e minimum of t h e s e t h r e e ,

and

where Ks i s t h e hydraulic conductivity

at

saturation, Bf Bfa t h e soil moisture c o n t e n t

at

field c a p a c i t y in A- and B-Layer, r e s p . and B s J , BSg t h e soil moisture content

at

s a t u r a t i o n in A- and B-layer, r e s p .

Lateral flow.

Q;,

is t h e limiting f a c t o r f o r t h e

rate

with which t h e

water

is discharged from t h e B-reservoir

to

streams and lakes. I t is a func- tion of hydraulic conductivity, K,, surface slope, S , soil moisture content above field capacity, catchment width, W, and t h e t e r r e s t r i a l catchment area, A, (Chen

et

al., 1982)

Quickflow is formed from

t w o

fractions; (1) if t h e soil moisture exceeds t h e saturation value, t h e exceeding volume i s assumed

to

e n t e r t h e brooks directly.

4i1),

and (ii) if t h e soil moisture exceeds t h e field capacity value, a fraction of t h e exceeding volume is discharged from t h e A-reservoir

as Lateral

flow,

4 ; ' ) .

The

total

quickflow

at

time s t e p r

is

then t h e sum of t h e s e two,

The volume of

water

retained in both

reservoirs

is then t h e balance between incoming and outgoing

water

volumes,

As a result, t h e hydrologic module simulates discharges from all

reser-

voirs: snowpack and sail r e s e r v o i r s A and B. The

water

from t h e s e t h r e e r e s e r v o i r s mixes in t h e lake within t h e mixing volume before running out from t h e outlet.

3.3. Soil C h e m i s t r y Module

IIASA's soil acidification model is applied as a component of t h i s model

to

compute t h e ion concentrations of t h e i n t e r n a l flows. Ideal mixing i s assumed in t h e r e s e r v o i r s and equilibrium is assumed

to

b e r e a c h e d a c c o r d - ing

to

computed Hc-concentration. The buffer mechanisms i n c o r p o r a t e d in t h e soil model a r e explained in detail elsewhere (Kauppi

et

al. 1985). In t h i s p a p e r only t h e b a s i c principles of t h e soil model a r e described.

I t i s assumed t h a t soils containing free c a r b o n a t e s ( c a l c a r e o u s soils) always have

a

b u f f e r

rate

high enough

to

neutralize any acid

stress.

In non-calcareous soils, however, neutralization depends on t h e intensity of silicate weathering (silicate buffer r a t e ) . As long

as

this b u f f e r r a t e is l a r g e r t h a n t h e a c i d

stress

t h e r e will b e no change in t h e Hc-concentration of soil o r in t h e quality of drainage. If acid s t r e s s exceeds t h e a c t u a l b u f f e r

rate

of t h e silicates. t h e soil is shifted

to

cation exchange b u f f e r r a n g e . Then t h e c a p a c i t y of t h e cation exchange buffer system, BC&, i s depleted with t h e

rate

of acid s t r e s s , asT, minus t h e buffer

rate

of silicates, bra (Eq.16). A non-linear relationship i s assumed between t h e b a s e s a t m t i o n and t h e soil pH within t h e silicate, cation exchange a n d t h e u p p e r aluminum b u f f e r r a n g e ,

as

long

as

BC& r 0,

at

pH from 5.6

to

4.0 (Eq.17).

The s h a p e of t h e pH

-

b a s e saturation relationship h a s been adopted from r e s u l t s of a n equilibrium model by Reuss (1983).

The assumption t h a t water discharged from t h e soil i s in equilibrium with

a

solid aluminum p h a s e h a s been widely used t o calculate t h e dissolved aluminum c o n c e n t r a t i o n s ( s e e e.g. Christophersen

et

al. 1982). The

s a m e

assumption may t h e n b e applied also t o compute t h e buffering through disso- lution of solid aluminum compounds. Gibbsite (AL (OH)3) i s one mineral o f t e n assumed t o c o n t r o l t h e equillbrium concentration of free aluminum ions and mononuclear hydroxy s p e c i e s in freshwaters. W e assume t h e equilibrium with gibbsite t o c o n t r o l buffering in soil after t h e r e i s no b u f f e r capacity l e f t in t h e cation e x c h a n g e buffer r a n g e , i.e. BC&

=

^0. As precipitation i n f i l t r a t e s into t h e soil and mixes with t h e soil solution, disequilibrium con- c e n t r a t i o n s [H'], a n d [dl3 +], a r e obtained,

where

Vr

is t h e volume of soil solution

at

field capacity and t h e infiltrating

water

volume is assumed

to

^equal

Pq + m T -

^ET. The soil solution volume i s simply defined by

The soil thickness of t h e A layer, z A , is fixed to 50

c m

and t h e volumetric

water

content value

at

field capacity. , i s

estimated

separately f o r e a c h soil t y p e based on t h e grain size distribution in soil.

Aluminum i s dissolved

or

precipitated until t h e gibbsite equilibrium

state

is reached (Eq.26). This p r o c e s s involves a change from disequilibrium con- centrations as defined in Eq.27

Combining Eqs.21 and 22 yields a third o r d e r equation which h a s a single real r o o t

The comparison between t h e

rate

of acid

stress

and t h e buffer

rate

t a k e s place

at

e a c h time s t e p in both r e s e r v o i r s b e f o r e t h e chemical s t a t u s of t h e soil solution i s computed. The water entering t h e

B-reservoir

has t h e quality of t h e soil solution leaving t h e A-reservoir. Acid

stress to

t h e

B-

r e s e r v o i r , Q S ~ , is then

where [H']; is t h e H+-concentration

in

t h e

water

leaving t h e

A-reservoir.

Depending upon t h e acid

stress

t h e r e is e i t h e r a net production of base cations

or

t h e r e i s a n exhaustion of cation exchange capacity. In c a s e t h e deposition

rate

of H+ is

lower

than t h e silicate buffer

rate.

t h e weathering f i r s t has

to

fill up t h e cation exchange complex and a f t e r t h a t a n excess supply of base cations

occurs.

The contribution of t h e soil r e s e r v o i r

to

t h e alkalinity of t h e s u r f a c e

water

i s assumed

to

equal t h e amount of t h e excess base cations (Eq.30).

The monthly leaching of hydrogen ions is simulated on t h e b a s i s of t h e simulated hydrogen ion concentrations and t h e simulated monthly d i s c h a r g e s from both r e s e r v o i r s . Additionally, p a r t of deposition and impurities in t h e meltwater f a l l d i r e c t l y on t h e lake. This s o u r c e of acidity i s simply computed from a c i d

stress

and l a k e

area.

A s a r e s u l t , t h e chemis-

try

module gives t h e quantity of acidity,

M H ,

a n d alkalinity,

M h 3 .

from a l l s o u r c e s to t h e l a k e

3.4. Lake Module

The l a k e module computes t h e time p a t t e r n of

water

quality in t h e lake.

The impact on a q u a t i c life will b e estimated on t h e b a s i s of simple t h r e s h o l d pH-values and aluminum concentrations. These c h a r a c t e r i s t i c s

are

most likely

to

indicate damage

to

fish populations a n d o t h e r aquatic organisms.

The change in l a k e water chemistry will b e p r e d i c t e d by means of t i t r a - tion of t h e b a s e c o n t e n t of t h e l a k e with s t r o n g a c i d originating from t h e atmosphere. The initial conditions

-

t h e preacidification water quality

-

h a s t o b e determined f o r

a

given lake. The

water

quality v a r i a b l e of g r e a t importance i s alkalinity, which e x p r e s s e s t h e

total

buffering c a p a c i t y of t h e l a k e water.

In preacidification conditions t h e only affecting p r o c e s s i s assumed

to

b e t h e weathering of c a r b o n a t e s or silicates. In case t h e soil contains f r e e c a r b o n a t e b e a r i n g minerals, t h e l a k e water c a n b e assumed t o b e v e r y high in alkalinity d u e

to

t h e high

rate

of c a r b o n a t e weathering. F o r s i l i c a t e r o c k s , Ulrich (1983) h a s defined weathering

rates

between 0.2-2.0 keq ha -'m - l y r

-'.

The original b i c a r b o n a t e concentration, of t h e l a k e water i s computed by t h e available information on:

t h e annual weathering

rate

of t h e mineral m a t t e r ( b r ) ;

t h e volume of soil through which t h e incoming

water

d r a i n s (A, (z* + z B ) ) ;

themeanannualrunofftowhichproducedHCO~ismixed(P-E).

The following steady-state b i c a r b o n a t e c o n c e n t r a t i o n in t h e o u t l e t of t h e l a k e may t h e n b e calculated based on t h a t information. The b i c a r b o n a t e concentration obtained i s used as t h e initial alkalinity f o r t h e model r u n s .

In c l e a r w a t e r l a k e s t h e c a r b o n a t e alkalinity c a n b e assumed

to

b e t h e only significant buffering agent. mainly with r e a c t i o n (Eq.35). Reaction (Eq.36) c a n b e neglected s i n c e t h e naturally sensitive s u r f a c e waters con- tain only negligible c o n c e n t r a t i o n s of c a r b o n a t e ions.

Reaction (Eq.35) yields

an

expression f o r t h e equilibrium (Eq.37), where [ H ~ C O ~ ] r e p r e s e n t s t h e sum of [C02] and [H2C03].

Combining this with Henry's Law (Stumm and Morgan, 1981)

one finally g e t s

where

K 1

and

KH _are

thermodynamic equilibrium constants, which depend on temperature.

When t h e drainage water, Q l

+ Qg.

t o g e t h e r with t h e d i r e c t

water

input as r a i n on t h e lake,

82,

e n t e r t h e lake and mix within t h e mixing layer, z o , disequilibrium concentrations (Eqs.40 and 41) result,

During t h e snowrnelt t h e mixing l a y e r is assumed

to

b e t h e topmost

water

layer. The

m e l t w a t e r

i s colder than most of t h e l a k e volume and t h e r e f o r e lighter than t h e 4OC

water at

t h e bottom. In this way t h e episodic spring time alkalinity and pH declines in t h e epilimnion can b e estimated.

The

t w o water

l a y e r s

are

then mixed t o g e t h e r a f t e r t h e r e i s no snow l e f t in t h e catchment. During t h e summer, t h e incoming acidity is mixed with t h e whole lake

water

body.

The b u f f e r r e a c t i o n (Eq.35) continues until a new equilibrium

state

according t o ~ q . 3 7 i s accomplished. Equal amounts of hydrogen and b i c a r - bonate ions are consumed

The new equilibrium concentrations, [H']' and [ H C 0 < I T , c a n b e obtained by solving Eqs.39-41. A second o r d e r equation is obtained, from which t h e positive r o o t f o r b i c a r b o n a t e concentration i s a c c e p t e d . The equilibrium hydrogen ion concentration is t h e n calculated from Eq.42

Finally, t h e equilibrium t o t a l alkalinity is given based on t h e definition by Stumm and Morgan (1981).

where t h e dissolved aluminum concentration is

set

p r o p o r t i o n a l

to

t h e hydrogen ion c o n c e n t r a t l o n according

to

t h e gibbsite equilibrium in Eq.28.

4. YODEL TESTING

The complexity of environmental systems and t h e lack of comprehen- s i v e t h e o r e t i c a l background make i t difficult

to test

any conceptualization of a given system. In

m o s t

applications of mathematical formulations describing physicochemical systems, t h e predictions

are

compared with measurements or samples of t h e system. In case t h e p r o c e s s e s under study are v e r y slow, i t i s p r a c t i c a l l y impossible

to test

t h e conceptualization of t h e system by comparing t h e outcome with measurements. Long time s e r i e s of, e.g., l a k e acidity measurements quite r a r e l y exist.

Recently,

a

need h a s been recognized in s e v e r a l environmental management and planning p r o g r a m s

to

c o n s t r u c t policy o r i e n t e d computer- ized tools

to

b e used in t h e decision making p r o c e s s . These formalized mathematical c o n s t r u c t s usually h a v e t h e b e s t c u r r e n t knowledge incor- p o r a t e d in t h e i r s t r u c t u r e . Such models, viewed

as

scientific t h e o r i e s , h a v e

to

b e t e s t a b l e in o r d e r

to

allow t h e i r use f o r management. A formal a p p r o a c h t o t h e o r y testing, based on numerical simulation and Monte C a r l o methods, h a s been proposed by F e d r a (1983). His analysis allows conclu- sions

to

b e drawn on t h e adequacy of both. t h e t h e o r y and t h e u n c e r t a i n inputs, a n d p r o v i d e s some guidance on how

to

improve a c e r t a i n conceptual- ization, even in t h e p r e s e n c e of a high d e g r e e of uncertainty (Fedra, 1983).

The first tests of the conceptualization of lake water acidification have been performed by simulating long term trends in freshwater acidity of a number of individual Finnish lakes (see K h i i r i e t al., 1985). In this model application, crude estimates f o r parameters and initial conditions were derived f o r 40 lakes in Southern Finland from a variety of sources including soil and geologic maps and the water quality data base of the Finnish National Board of Waters. The aim of the

test

runs w a s

to

evaluate whether the model could correctly distinguish between acidified lakes and lakes where no indications of acidification had been observed, when driven with an assumed historical deposition pattern (see Figure 6). The deposition pat-

tern

from 1960 t o 1980 w a s obtained from the RAINS model, in which the upstream submodels, the energy-emissions submodel and the

EMEP

sulfur long-range transport submodel, supplied the input for the environmental impact submodels (see Alcamo e t al., 1985). The other environmental driving variables, ambient a i r temperature and precipitation,

were

obtained from the thirty year climatic means presented in

MUer

(1982). Example model

runs

are shown f o r two lakes in this paper. Catchment characteristics as w e l l a s values f o r initial conditions and parameters used in the application

are

summarized in Table 2.

These first results suggested that the model was able

to

generate an allowable outcome. For example, a strongly decreasing Lake pH

w a s

estimated f o r lake O r a j W i ,

a

typical acidic lake with an observed summer pH below 5.0 (Figure 4). Moreover, a c o r r e c t pH-level in 1980 was predicted f o r lake Venjkvi, which is a circumneutral lake in South-East Finland (Fig- u r e 5). Therefore, these model results were considered promising f o r further applications and

as

a next step Fedra's (1983) formal approach f o r model testing was utilized.

The logical structure of testing a proposed model (a theory). i.e. com- parison of outcomes with observations,

w a s

retained but

at

the same time an appropriate way of describing uncertain inputs

as

w e l l

as

somewhat uncer- tain expected outcome was included. Instead of giving the inputs and outputs

as

specific values, these two data points

were

extended

to

regions in their respective spaces. In the presence of uncertainty, one has

to

deal with a

set

of vectors instead of one vector in the n-dimensional input and output vec-

tor

spaces (Fedra, 1983).

To

test

a specific model, i t is examined whether, f o r

a

s e t of ranges of initial conditions and parameters, allowable outcomes can be produced. The model has

to

be rejected, if no allowable outcome can be generated from a statistically sufficient number of trials. If

a

simple model version fails

to

give an acceptable behavior over the allowable input ranges, i t can be modified by adding more complex procoss descriptions

to

the model. In our case, however,

w e

wish

to

keep the model

as

simple

as

possible in o r d e r

to

minimize the computational steps required. For the regional application, one objective is

to

retain the simplicity of the model

to

allow its interactive use,

The lake acidification model version described above was subjected

to

the testing procedure. The catchments studied in this application were the same

as

f o r the first model application

to

individual basins. The model was incorporated into a Monte Carlo framework, which randomly sampled a

Table 2 . Catchment c h a r a c t e r i s t i c s , initial conditions and p a r a m e t e r values chosen for t h e t w o example catchments.

Lake Orajimi

A~ ( m S ,

=

2 . 2 - l o 5

A, ( m S ,

=

5 . 6 ~ 1 0 ~

ZA ( m ) zg ( m )

2 , ( m ) ( m )

w

( m )

S ( m / m )

rnctotl

^{(eq m}

-2)

(eq m -2)

a c t , ,

8 (eq m -2) -8 (eq m - 2 )

bra _{(eq m} -3 P r _-') Or

Lake Veniirvi

ZA ( m ) zg ( m )

2 , ( m ) ( m )

w

( m )

S ( m / m )

CECtOt m - 2 ) BecEl (eq m - 2 ) CECtot ^j (eq m - 2 )

BC,&

(eq m - 2 )

bra _{(eq m} -3 P r _-') O f

parameter v e c t o r from t h e assigned allowable r a n g e s (Table 3 ) . The model was

run

with

a

constant deposition pattern through a period of t e n y e a r s

to

allow t h e a r b i t r a r y initial values of t h e

state

variables

to

adjust. Next,

a

simulation of twenty y e a r s was performed assuming a historical deposition p a t t e r n shown in Figure 6 and all information about t h e r u n w a s s t o r e d . This process was r e p e a t e d f o r a l a r g e number of t r i a l s (500 times f o r e a c h lake).

Finally, t h e

set

of r u n s obtained

was

analyzed f o r violations of t h e con- s t r a i n t conditions in t h e

course

of t h e simulations.

The uncertainty in t h e measurements and t h e r e p r e s e n t a t i v e n e s s of t h e sample f o r t h e whole s e a s o n

w a s

considered when assigning t h e constraint r a n g e s on t h e basis of observations. The generalized conditions compiled in Table 4 were formulated.

Figure 4. Simulation of t h e lake-pH of lake Orajiirvi driven by t h e his- torical deposition pattern.

Figure 5. Simulation of t h e lake-pH of lake V e n j h i driven by t h e his- torical deposition pattern.

The model

w a s

expected

to

predict correctly t h e observed long term development of lake

water

acidification together with t h e observed seasonal variation of lake acidity. The reproduction of these

t w o

p a r t s of t h e sys-

t e m s

behavior w a s considered relevant

to

t h e regional application of t h e model. If t h e model s t r u c t u r e could fulfill t h e formulated behavior defini- tion, i t

can

also b e assumed

to

describe satisfactorily t h e development of pH- and aLkaLinity4fstributtons in different Lake r e b o n s in Europe.

Table

3.

Assigned allowable input r a n g e s

for

t h e

t w o

catchments stu- died.

Lake Ven jimi

P a r a m e t e r Min. Max. P a r a m e t e r Min. Max.

K,

( m mo -l)

Z A + Z B ( m )

bra (eq m -3 y r -l) W X P 1 CECtot J

~ G F , B 1 CECtot ,B S ( m / m )

@I

@*

2 , ( m ) z ( m )

K,

( m mo -I)

ZA + ZB (m )

bra (eq m -3 y r -l)

=%rP 1 CECt0t.A mJrg 1 CECtOt 3 S ( m / m )

@f

@s

2 , ( m )

z ( m )

Figure

6.

Sulfur deposition scenarios assumed f o r t h e studied catch- ments; (a) high and (b) low deposition scenario.

Results of t h e Monte-Carlo r u n s indicated t h a t t h e model could fulfill

all of t h e assigned constraint conditions. For l a k e Venjgrvi 18

r u n s out of

500 were

found to produce a n acceptable outcome. For lake OraJZirvi only

5

runs out of

500

gave an acceptable behavior. The model seems t o b e a b l e to

-

²⁰

-

Table

4.

Constraint conditions f o r t h e t w o catchments studied.

generate t h e g e n e d trend in lake acidity, but many

of

the runs had to be rejected, because t h e

m o d e l

gave too

l o w

pH-values or t h e seasonal dynam- ics w a s not predicted accurately enough.

There are t h r e e possible explanations f o r getting

so

few sets of input values to match t h e observed water quality. The f i r s t and m o s t appealing is, t h a t t h e historical deposition p a t t e r n - given by t h e upstream models -

^is

biased;

i.e.

t h e emission estimates f o r t h e

period 1960

to

¹⁹⁸⁰

maybe some- what too high. The second reason could be t h e uncertainty in assigning input ranges, they may be too narrow and badly placed. On t h e o t h e r hand, t h e allowable ranges f o r input values should be constrained as much

as pos-

sible in o r d e r to avoid unrealistic combinations of input values. The third explanation could be t h e

m o d e l

itself: t h e proposed

m o d e l

s t r u c t u r e does not

try

to be a final description of t h e lake acidification phenomena. Some additional

processes

may have to be included

in

t h e course of a f u r t h e r development of t h e model. For example, t h e r e

is

r e c e n t evidence t h a t t h e carbonate s y s t e m

is

not t h e only buffering system operating

in

lakes against acidification. Schindler et

al. (1985)

have reported t h a t

a

significant

p r e

portion of t h e acid neutralization in

lakes

in different parts of t h e world has been accounted

for

by auxiliary buffering (sulfate reduction, denitrifi- cation,

etc.). U p

to now t h e model has been kept simple; nevertheless

it

consists of numerous mathematical descriptions, each of them being a theory by

itself.

Their dynamic nonlinear interactions make

it

difficult to relate a failure in t h e o v e d l performance of t h e model to any of t h e indivi- dual modules used.

The model could, thus, reproduce t h e required behavior, which

is

not to say t h a t

it is

validated now. Rather, one could say that t h e model could not b e falsified yet. Meanwhile one can use t h e model cautiously as

a

tool f o r soenario analysis, representing one interpretation of t h e c u r r e n t understanding

of

long term lake acidification. The constraint conditions had to be formulated

so

broad, t h a t several possible conceptualizations of lake acidification necessarily exist, which are able to reproduce t h e

allow-

able systems behavior. The differences may only arise when t h e models are applied f o r future projections (Fedm,

1983). _{W e}

consider

m o d e l

validation a iterative process, which should

also

provide guidance on how to improve t h e

m o d e l .

Lake Orajhrvi [I]

The mean summer pH in

1965

between

5.1

and

5.7

[2]

The mean summer pH in

1980

between

4.4

and

5.0

[3]

The mean spring pH in

1980

between

4.2

and

4.8

Lake Venj&rvi [I]

The mean summer pH

in 1965

between

6.4

and

7.0

[2]

The mean summer pH in

1980

between

6.3

and

6.9

[3]

The mean spring pH

in

1980

between

5.8

and

6.4

The

set

of inputs, which yielded a c c e p t e d simulations, .i.e. model r u n s fulfilling t h e c o n s t r a i n t conditions, c a n be used

f o r

f u t u r e projections of t h e r e s p o n s e of t h e catchments

to

d i f f e r e n t deposition s c e n a r i o s . Two example s c e n a r i o s were introduced using t h e IIASA R A I N S model (Figure 6;

see Alcamo

et

al., 1985). From 1960

to

1980 t h e s c e n a r i o s assumed identical historical energy-emission trends. From 1980 on t h e s c e n a r i o s d e p a r t e d s o t h a t t h e h i g h deposition s c e n a r i o assumed a n increasing

rate

of e n e r g y use throughout E u r o p e , as defined by t h e ECE ' t r e n d s continued' s c e n a r i o (ECE, 1983), linearly e x t r a p o l a t e d

to

2030. The low deposition s c e n a r i o w a s c o n s t r u c t e d from t h e ECE 'conservation' s c e n a r i o , assuming lower rates of e n e r g y u s e and in addition

to

t h a t , effective measures taken f o r t h e c o n t r o l of s u l f u r emissions.

The model c a n b e used f o r any deposition s c e n a r i o f o r producing a s e r i e s of l a k e simulations with t h e accepted

set

of input values. The varia- bility of t h e individual r u n s can b e i n t e r p r e t e d as t h e u n c e r t a i n t y in t h e model r e s u l t s d u e t o t h e uncertainty in t h e input d a t a . This means t h a t some of t h e initial u n c e r t a i n t y i s taken into account throughout t h e long

term

simulation. The model behavior i s demonstrated in Figures 7 and 8 , in which simulations are continued from 1980 on using only t h e 'accepted' input d a t a . The c o n s t r a i n t s c a n b e looked

at

^as windows through which t h e f u t u r e s c e n a r i o s are f o r c e d t o pass.

Results of t h e model r u n s of l a k e VenjZrvi show a c l e a r i n c r e a s e in time in t h e u n c e r t a i n t y of s c e n a r i o s (Figure 7). However, t h e d e g r e e of uncertainty seems t o depend on t h e type of s c e n a r i o used as w e l l

as

on t h e t y p e of l a k e studied. Running lake OrajZwi with t h e high s c e n a r i o t h e five a c c e p t e d

runs

show v e r y l i t t l e variability, whereas with t h e low s c e n a r i o t h e r e are c l e a r d i f f e r e n c e s in t h e r a t e of r e c o v e r y (Figure 8). Neverthe- less, t h e s c e n a r i o s yield a mean summer-pH quite close to 6.0 in 2030.

F e d r a (1981) h a s concluded, t h a t t h e prediction uncertainty (measured a s t h e coefficient of variation of Monte Carlo outputs) i n c r e a s e s with t h e pred- iction time as well as with t h e amount of change in t h e driving variables.

The r e s u l t s from o u t preliminary analysis on s c e n a r i o uncertainty suggest a similar p a t t e r n . The magnitude of uncertainty may, however, v a r y , depend- ing on t h e model and t h e d a t a used. T h e r e f o r e , besides a more comprehen- s i v e u n c e r t a i n t y analysis, also a b e t t e r understanding of t h e model s t r u c - t u r e is essential.

3.e 9 -

19E0 1970 1980 1990 2@@@ 2@10 202@ 2B3@

YERR

1 . 0

1960 1970 1980 1990 2000 201 0 202@ 2B30

YERR

Figure 7. Lake acidity scenarios for lake Venjiirvi generated by using the s e t of input data combinations. which fulfilled the con- straint conditions (see Table 4). In (a) the high scenario was used and all 35

runs

are displayed: in (b), using the low scenario, the mean and the minimum-maximum envelope i s displayed.

3.0 1

^{I I I I I I}

ⁱ

1969 1970 1980 199rE 2000 20 10 2220 2Q3B

YERR

3 . 0 1

^{I I I I I I}

^I

1960 1970 1980 199D 2000 2010 2020 2@30

YERR

Figure 8. Lake acidity scenarios for lake OraJiirvi generated by using the s e t of input data combinations, which fulfilled the con- straint conditions (see Table 4). In (a) the high scenario and In (b) the low scenario was used.

All

5 runs a r e displayed.

LIST

OF

SYMBOLS

In this list of symbols t h e superscript

^T

- the time step index - ^{i s}

suppressed.

time step index

melting rate coefficient f o r forested area evapotranspiration r a t e coefficient actual soil moisture content

soil moisture content at field capacity soil moisture content at saturation lake area

catchment area acid stress rate

actual

cation exchange capacity silicate buffer rate

total

cation exchange capacity

total

deposition

deposition released from t h e snowpack deposition retained in the snowpack accumulated deposition

monthly evapotranspiration mean annual evapotranspiration

f i r s t dissociation constant f o r carbonic acid constant of Henry's

l a w

hydraulic conductivity at saturation gibbsite equilibrium constant

flow of hydrogen ions to t h e lake flow of bicarbonate ions to t h e lake melting rate of snow

mean

annual

precipitation monthly rainfall

monthly snowfall

total

monthly precipitation

partial pressure of

C02

in lake

w a t e r

lateral flow from t h e B-layer

Qd precipitation onto t h e l a k e

QP percolation flow from A-layer into B-layer Qi1) maximum possible percolation r a t e

QiZ) water volume available in t h e A-layer Qd3) s p a c e l e f t in t h e B - r e s e r v o i r

QQ quickflow from t h e catchment ,:I) surf a c e runoff

QJZ) l a t e r a l flow from t h e A-layer S s u r f a c e slope

SP

snowpac k

T

mean monthly t e m p e r a t u r e

T~

t h r e s h o l d t e m p e r a t u r e above which all precipitation falls as r a i n Ts t h r e s h o l d t e m p e r a t u r e below which all precipitation falls

as

snow

vf

^soil w a t e r volume of field capacity

v~

t o t a l water volume in A-layer

v~

t o t a l water volume in B-layer

W catchment width

Z~ thickness of A-layer z~ thickness of B-layer z mean depth of t h e l a k e

o

mixing l a y e r of t h e l a k e

0,s mixing l a y e r of t h e lake in spring

[ALk] equilibrium total alkalinity in l a k e water [H%O<], initial alkalinity in l a k e water

[H2C0j ] sum of [C02] and [H2C03]

[AL 3+] aluminum ion concentration [H+] hydrogen ion concentration [HCO<] b i c a r b o n a t e ion concentration

Alcamo, J., L.Hordijk, J. Kiimgri, P. Kauppi, M. Posch and E. Runca, 1985.

I n t e g r a t e d analysis of acidification in Europe. J. Environ. Manag., 21:

47-61.

Chen, C.W., J.O. Dean, S.A. Gherini and R.A. Goldstein, 1982. Acid r a i n model

-

hydrologic module

-

J . Env. Eng. ASCE. 1 0 8 No. EE3: 455-472.

Chen, C.W., S.A. Gherini. R.J.M. Hudson and J.D. Dean. 1983. The I n t e g r a t e d Lake-Watershed Acidification Study. V o l . 1: Model p r i n c i p l e s a n d appli- cation p r o c e d u r e . Final R e p o r t , September 1983. Lafayette: T e t r a Tech. Inc. 1 8 6 pp.

Chow, V.T. 1964. Handbook of applied hydrology. McGraw-Hill, New York.

Cosby, B.J., R.F. Wright, G.M. H o r n b e r g e r and J.N. Galloway, 1985. Modelling t h e e f f e c t s of a c i d deposition: Assessment of a lumped p a r a m e t a r model of soil

water

and s t r e a m w a t e r chemistry. Water Resour. Res. 21: 51-63.

Christophersen, N., and R.F. Wright. 1981. Sulfate budget and a m o d e l f o r s u l f a t e concentrations in s t r e a m w a t e r

at

Birkenes a small f o r e s t e d catchment in southernmost Norway. Water Resour. Res. 17: 377-389.

Christophersen, N., H.M. S e i p and R.F. Wright. 1982. A model f o r streamwa-

ter

chemistry

at

Birkenes, Norway. Water Resour. Res. 18: 977-996.

Christophersen. N., L.H. Dymbe. M. Johannessen and H.M. Seip, 1984. A model f o r sulphate in streamwater

at

Storgama, Southern Norway.

Ecol. Modelling 21:35-61.

Economic Commission for Europe (ECE), 1983. An efficient energy future.

Butterworths, London, 259 pp.

Eliassen. A.. and J. Saltbones. 1983. Modelling of long-range transport of sulphur o v e r Europe:

a

two-year model run and some model experi- ments. Atm. Environ. 17: 1457-1473.

Fedra, K., 1981. Estimating model prediction accuracy:

a

stochastic approach

to

ecosystem modelling. In: D.M. Drablos (ed), Progress in ecological engineering and management by mathematical modelling.

Liege: Cebedoc, pp.87-105.

Fedra, K. 1983. Environmental modeling under uncertainty: Monte Carlo simulation. IIASA Research Report, RR-83-28. 78 pp.

Henriksen, A. 1930. Acidification of freshwaters:

a

l a r g e scale titration. In:

D. Drablos and A. Tollan (eds). Proc. of International Conference on Ecological Impact of Acid Deposition. 11-14 March 1980. Sandef jord, Norway. pp .68-74.

Hornberger, G.M. and J.B. Cosby, 1985. Evaluation of

a

model of long-term response of catchments

to

atmospheric deposition of sulfate. Proc. of 7th IFAC Symposium on Identification and System Parameter Estimation, Univ. of York, 3-7 July 1985, York, U.K. Pergamon P r e s s (in press).

Johannessen, M. and A. Henriksen, 1978. Chemistry of snow meltwater:

changes in concentration during melting.

Water

Resour. Res. 14: 615- 619.

K-i, J., M. Posch and L. Kauppi, 1985. Development of

a

model analyzing s u r f a c e

water

acidification on a regional scale: application

to

indivi- dual basins in Southern Finland. In: I. Johansson (ed), Hydrological and hydrogeochemical mechanisms and model approaches

to

t h e acidifica- tion of ecological systems. Nordic Hydrological Programme, NHP Report No. 10, pp.151-170.

Kauppi, P., J. K-i, M. Posch. L. Kauppi and E. Matzner, 1985. Acidifica- tion of forest soils: model development and application f o r analyzing impacts of acidic deposition in Europe. Ecol. Modelling (in press).

Miiller, M.J.. 1982. Selected climatic data for

a

global

set

of standard

sta-

tions for vegetation science. D r . W. Junk Publ., The Hague. 306 pp.

Reuss, J.O. 1983. Implications of the Ca-A1 exchange system for the effect of

acid precipitation of soils. J. Environ. Qual. 12: 591-595.

Schindler, D.W., J.W.M. Rudd, C.A. Kelly and M.A. T u r n e r , 1985. E f f e c t s of acidic p r e c i p i t a t i o n o n

i n

situ alkalinity generation in lakes. P a p e r p r e s e n t e d at t h e Muskoka Conference '85, International Symposium on Acidic P r e c i p i t a t i o n , Sep. 15-20, 1985, Ontario, Canada.

Shih,

G.B.,

R.H. Hawkins a n d M.D. Chambers, 1972. Computer modelling of a coniferous forest watershed. In: Age of Changing P r i o r i t i e s for Land a n d

Water.

Am. Soc. Civil Eng., N e w York.

Stumm. W. a n d J.J. Morgan, 1981. Aquatic chemistry. An Introduction emphasizing chemical equilibria in n a t u r a l

waters.

2nd ed. John Wiley&Sons, N e w York. 780 pp.

Ulrich,

B.

1983. Soil acidity a n d i t s relations

to

acid deposition. In:

B.

Ulrich and S. P a n k r a t h (eds). E f f e c t s of accumulation of a i r pollutants in forest ecosystems. D. Reidel Publ. Co. pp. 127-146.