-1-

Blockkurs Hydrogeochemische Modellierung- Theorie -

Prof. B. Merkel & Dipl.Geol. B. Planer-FriedrichLehrstuhl für Hydrogeologie TU Bergakademie Freiberg

1. Einleitung

Hydrogeochemische Modelle sind numerische Werkzeuge, um wasserchemische Analysen zuinterpretieren und geogene sowie anthropogen beeinflußte aquatische Systeme zu analysieren. AlsEingangsgrößen werden in der Regel thermodynamische Daten und möglichst vollständigewasserchemische Analysen benötigt. Die thermodynamischen Daten stehen - soweit es sich umKomplexbildungskonstanten und Löslichkeitsprodukte handelt - in Form von Standard-Datensätzenfür die jeweiligen Programme zur Verfügung. Die Daten zur Beschreibung oberflächenkontrollierterReaktionen müssen durch eigene experimentelle Befunde ergänzt werden. Im Gegensatz zuGrundwasserströmungs- und Transportmodellen bedürfen hydrogeochemische Modelle an sichkeiner Kalibrierung. Bei Berücksichtigung oberflächenkontrollierter Reaktionen und kinetischkontrollierter Reaktionen müssen aber auch hydrogeochemische Modelle kalibriert werden.

2. Thermodynamische Grundlagen

Um den Ablauf einer chemischen Reaktion beschreiben und vorhersagen zu können, ist es nötig, diean der Reaktion beteiligten Species zu erfassen, gewisse Reaktions-Parameter, wie z.B. Enthalpie,Entropie und Gibbs´sche Energie zu kennen, und die Umgebungsbedingungen (homogene-heterogene Reaktion; offenes-geschlossenes System) zu berücksichtigen.

Unter Species versteht man in der Hydrochemie alle wasserlöslichen organischen und anorganischenVerbindungen, d.h. sowohl freie Kationen und Anionen (Na+, K+, Ca2+, Mg2+, Cl- oder F-) als auchVerbindungen verschiedener Elemente, genannt aquatische Komplexe, wie z.B. OH-, HCO3

-, SO42-,

NO3-, CaSO4

0 oder Komplexe mit organischen Liganden.

Die Enthalpie H ist eine thermodynamische Zustandsfunktion und wird auch als dasthermodynamische Potential einer Substanz bezeichnet. Es ist definiert über H=U + p V, wobei Udie innere Energie, p der Druck und V das Volumen ist. Die Entropie S0 ist ein Maß für denOrdnungszustand eines thermodynamischen Systems bzw.ein Maß für Irreversibilität eines Vorgangs.Sie nimmt mit zunehmender Unordnung höhere Werte an. Die Gibb'sche freie Energie G ist ein

-2-

Maß für die Kraft, eine Reaktion ablaufen zu lassen. Ein positiver Wert bedeutet, daß zusätzlicheEnergie nötig ist, um eine Reaktion ablaufen zu lassen. Ein negativer Wert signalisiert, daß eineReaktion spontan abläuft. Der Zusammenhang zwischen diesen drei Größen ist:

H = G + S0*T oder G = H - S0*T

Homogene Reaktionen spielen sich nur innerhalb einer Phase ab, während heterogeneReaktionen zwischen verschiedenen Phasen ablaufen (Wasser-Gas, Wasser-Gestein, Gas-Gestein). In einem geschlossenen System kann nur ein Energieaustausch mit der Umgebungstattfinden, während in einem offenen System neben Energie auch Materie ausgetauscht werdenkann.

2.1. Massenwirkungsgesetz

Die Beschreibung einer chemischen Reaktion erfolgt entweder über den Ansatz einesthermodynamischen Gleichgewichtes, wenn reversible chemische Reaktionen zu beschreiben sindoder über den Ansatz eines kinetischen Gleichgewichtes für irreversible Reaktionen. Hier soll näherauf den Ansatz des thermodynamischen Gleichgewichtes eingegangen werden, bei dem chemischenReaktionen mit Hilfe des Massenwirkungsgesetzes (MWG) in folgender allgemeiner Formbeschrieben werden:

AB ↔ A + B

In dieser vereinfachten Form werden weder die Wertigkeit noch die Stöchiometrie berücksichtigt.Für eine solche Reaktion kann eine Konstante K experimentell bestimmt werden, die die relativeVerteilung der beteiligten Species beschreibt:

[AB]{B}{A}

K ⋅

=

Die eckigen Klammern in obiger Formel geben die „Aktivitäten“ der jeweiligen Species an. Handeltes sich bei AB um eine Mineralphase, ist zusätzlich per Definition vereinbart, daß die „Aktivität“ desMinerals gleich 1 ist. Damit vereinfacht sich die obige Gleichung zu

K = {A}⋅{B}

Die Konstante K wird in Abhängigkeit der jeweiligen Reaktionen als „Löslichkeitsprodukt“ (Lösung/ Fällung), „Stabilitätskonstante“ (Komplexbildung / Komplex-Auflösung, Redoxreaktionen) oder„Selektivitätskonstante“ (Sorption / Desorption) bezeichnet. Der dekadische Logarithmus derKonstante K (log K) wird als pK-Wert bezeichnet.

-3-

Vertauscht man die Seiten der Edukte und Produkte und schreibt statt AB ↔ A + B(Minerallösung) A + B ↔ AB (Mineralbildung, Fällung) ist die zugehörige Konstante K´= 1/K,bzw. der log K´= -log K.

Sollen stöchiometrische Verhältnisse berücksichtigt werden, ergibt sich die Konstante für dieReaktion aA + bB ↔ cC + dD mit a, b, c, d = mol-Anzahl der Edukte A und B, bzw.Produkte C und D zu:

ba

dc

}B{{A}{D}{C}

K ⋅⋅

=

Was ist nun aber der Unterschied zwischen Aktivität und den im Labor bestimmten Konzentrationeneines Elementes?

2.2. Unterschied zwischen Konzentration und Aktivität

Unter der Aktivität ai wird in diesem Zusammenhang die chemische Aktivität eines Ions mit demIndex i oder eines Komplexes mit dem Index i und nicht die radioaktive Aktivität eines Elementsverstanden. Nur in unendlichen verdünnten Lösungen eines aquatischen Systems ist die Stoff-Konzentration ci gleich der Aktivität ai dieses Stoffes. In allen anderen Fällen gilt:

ai = fi . ci mit fi ≠ 1

Der Aktivitätskoeffizient fi beschreibt den funktionalen Zusammenhang zwischen Konzentration undAktivität. Er beschreibt, inwieweit sich die positiv oder negativ geladenen Ionen durch interionareWechselwirkungen gegenseitig beeinflussen und in ihrer „Aktivität“ behindern oder unterstützen. Dader Aktivitätskoeffizient eine nichtlineare Funktion der Ionenstärke ist, ist auch die Aktivität einenichtlineare Funktion der Konzentration. Im Bereich von Ionenstärken bis zu 1 mol/L ist die Aktivitätgeringer als die Konzentration, der Aktivitätskoeffizient hat somit einen Wert < 1. DieserZusammenhang ist in Abb. 1. bis zu einer Ionenstärke von 0,1 mol/L dargestellt. Ganz deutlich ist derEinfluß der Wertigkeit der Species auf den Aktivitätskoeffizienten zu sehen. Mit Ionenstärken > 1mol/L kann sich das aber ändern: die Aktivität einer Species i kann größer werden als dieKonzentration; damit nimmt der Aktivitätskoeffizient Werte > 1 an.

Die Berechnung der Ionenstärke I erfolgt allgemein gemäß der Gleichung

∑=

=n

iii

zc1

2**5.0 I

ci Konzentration aller Species i = 1 bis n in einem aquatischen System in mol/L

-4-

zi Wertigkeit der Species i

Abb. 1. Zusammenhang zwischen Ionenstärke (bis 0,1 mol/L) und Aktivitätskoeffizient (nachHEM 1985)

Mit Hilfe der Debye-Hückel und ähnlicher Gleichungen ist es möglich, die jeweiligenAktivitätskoeffizienten aus der ermittelten Ionenstärke zu berechnen. Im folgenden sind diegängigsten Varianten und ihr maximaler Gültigskeitsbereich (Ionenstärke) aufgeführt:

-5-

Debye-Hückel-Gleichung:

0.005mol/LIIzA)log(f 2

ii <⋅⋅−=

erweiterte Debye-Hückel-Gleichung:

0.1mol/LIIaB1

IzA)log(f

i

2

ii <

⋅⋅+⋅⋅−

=

„WATEQ“- Debye-Hückel-Gleichung:

1mol/LIIb IaB1

IzA)log(f i

i

2

ii <⋅+

⋅⋅+⋅⋅−

=

Davis-Gleichung:

0.5mol/LI0.3I) - I1

I(zA)log(f 2

ii <+

⋅−=

Güntelberg-Gleichung:

0.1mol/LII1,41

I0.5z)log(f 2

ii <+

−=

A Temperatur abhängiger Parameter (� 0,5)B Temperatur abhängiger Parameter (� 0,33 * 10-8)a,b Species-spezifische Parameter (abhängig vom Ionenradius)

Wie anhand dieser Aufstellung ersichtlich ist, reicht die zulässige Ionenstärke bis maximal 1 mol/L.Einige Autoren vertreten aber die Meinung, daß bereits bei Ionenstärken unter 1 mol/L derGültigkeitsbereich der Ionendissoziationstherie beendet ist. Für größere Ionenstärken kann derAnsatz nach PITZER verwendet werden.

-6-

3. Wechselwirkungen innerhalb der wässrigen Phase

3.1. Komplexbildung

Aquatische Komplexe sind definiert als Verbindung eines positiv geladenen Teils (meist ein Metall;möglich ist aber auch ein H+-Ion oder ein anderer positiv geladener Komplex) und eines negativgeladenen Teils, des sog. Liganden. Liganden können sein freie Anionen (Cl-, F-, Br-), ein negativgeladener anorganischer Komplex (OH-, HCO3

-, SO42-) oder organische Komplexe, wie z.B.

Huminsäuren, EDTA, usw.

Die Ladung eines Komplexes zu kennen, ist wichtig für die Abschätzung seinesKontaminationspotentials, da positiv oder negativ geladene Komplexe effektiv sorbiert oderzumindest retardiert werden können, deutlich weniger aber 0-wertige Komplexe.

3.2. Redoxreaktionen

Der Sauerstoffgehalt eines Wassers bestimmt neben anderen Inhaltsstoffen in erster Linie denReduktions-Oxidations-Zustand eines aquatischen Systems (Redoxzustand):

H2O ↔ H2 + 1/2 O2 (Autoprotolyse des Wassers)

Nach dieser Formel ist jeder Sauerstoffkonzentration (rechnerisch) ein bestimmter Wasserstoffgehaltzugeordnet. Die Konzentration von gasförmigen Stoffen in Lösungen kann auch mit derenPartialdruck beschrieben werden. Aus obiger Gleichung folgt:

2.85

2

2

2

2

2 10O}{H

}{pO*}{pH K −==

Dies bedeutet, daß sauerstoffgesättigtes Wasser mit dem Partialdruck pO2 = 1 bar = 1×105 Pa imGleichgewichtszustand mit einem Wasserstoffpartialdruck pH2 von 10-42,6 bar steht. Umgekehrt stehtwasserstoffgesättigtes (vollkommen reduziertes) Wasser im Gleichgewicht mit Sauerstoff mit demPartialdruck pO2 = 10-85,2 bar. Wasserstoffgesättigtes Wasser steht weiterhin im Gleichgewicht mitWasserstoffionen:

H2 ↔ 2H+ + 2e,

so daß im reduzierten System je Wasserstoffmolekül 2 negative Ladungen entstehen. DieseNegativladungen erzeugen gegenüber einem Bezugssystem ein negatives Potential, das gemessenwerden kann (Elektronenpotential pe-). Definitionsgemäß wird das Bezugspotential durch dieNormalwasserstoffelektrode mit pH2 = 1 bar bei einem pH-Wert 0 dargestellt:

-7-

E0(H+/H2) = 0 mV.

Allgemein wird das Potential durch folgende Formel beschrieben:

{red}{ox}

lnFnTR

EE 0

H ⋅⋅

+=

darin bedeuten:

E0 = Normal-Redox-Spannung eines Systems bei den Aktivitäten Ox. = Red.R = Allgemeine Gaskonstante (8,315 J/K mol)T = Absolute Temperatur (Kelvin)n = Anzahl der freiwerdenden Elektronen (e)F = FARADAY-Konstante (96484 C/mol = J/V mol){Ox.} Aktivität des oxidierten Partners{Red.} Aktivität des reduzierten Partners

Da in thermodynamischen Programmen nicht ohne weiteres mit Redoxpotentialen gerechnet werdenkann, wurde das Konzept des pE - Wertes eingeführt. Der pE-Wert entspricht in Analogie zumpH-Wert dem negativen dekadischen Logarithmus der Aktivität der Elektronen. Es wird also miteiner hypothetischen Aktivität bzw. Konzentrationen von Elektronen gerechnet, die de facto imWasser nicht vorhanden sind. Die Umrechnung zu gemessenen Redoxpotentialen ergibt sich aus:

H

- ETR2.303

F}log{e- pE ⋅

⋅⋅==

Durch Einsetzen von Standard-Werten für F, R und T (für 25 Grad Celsius) ergibt sich

pE � 16.9 EH [EH in V]

Da eine experimentelle Bestimmung des Redoxpotentials EH mit herkömmlichen Platin-KalomelElektroden maximal auf ± 30 mV möglich ist, reicht diese vereinfachte Umrechnung in der Regelaus. Für das System H2/H+ gilt:

}{H}{H

ln F2TR

HH

EE2

2

20

H

+

+ ⋅⋅

+

=

}ln{HFnTR

}ln{HFnTR

0E2

2

H⋅

⋅⋅

−⋅⋅⋅

+= +

-8-

}log{HF2

TR2.303pH

FTR2.303

0E 2H⋅

⋅⋅⋅−⋅⋅⋅−= = -0,0591 pH (25°C)

Eine Änderung um eine pH-Einheit bewirkt somit eine Spannungsänderung um 59,1 mV, eineÄnderung des Wasserstoff- (Sauerstoff-) Partialdrucks um eine Größenordnung bewirkt eineSpannungsänderung um 29,6 (14,8) mV. In Redox-Systemen, bei denen lediglich die Wertigkeitgewechselt wird, ohne daß ein Wasserstoffion beteiligt ist (z.B. Fe2+ ↔ Fe3+), liegt einepH-unabhängige Reaktion vor.

In einem pH-/EH-Diagramm kann für jedes Redox-System die vorherrschende Species dargestelltwerden. Solche Darstellungen nennt man Stabilitätsfelder oder (besser) Prädominanz-Diagramme(Abb. 2. ).

Abb. 2. Prädominanzdiagramm für aquatische Spezies im System U-O2-H2O in destilliertemWasser bei 25°C und 10-5Pa Gesamtdruck; die gestrichelte Linie zeigt die Grenze derUO2 Mineralfällung/ -Lösung bei einer Uran-Gesamtkonzentration von 10-8mol (ausLANGMUIR 1997)

-9-

Das mit einer Kalomel-Platinelektrode, aber auch das mit einer Wasserstoffelektrode gemesseneRedoxpotential ist ein Summenpotential, das Integral aller in einem Wasser vorhandenerRedoxspecies. Theoretisch ist es möglich, die Speciesverteilung von allen redoxsensitiven Elementenanalytisch zu bestimmen und aus den zugehörigen Redox-Einzelpotentialen das Redox-Summenpotential zu berechnen. Wird dies beispielsweise für eine Wasserprobe gemacht, indem dieVerteilung As(3) zu As(5) einerseits und andererseits von U(4) zu U(6) sowie N(5) / N(3) {NO3

-/NO2

-} analytisch bestimmt wird, so wird sich in der Regel für jedes Redoxpaar eine andereRedoxspannung errechnen. Die gewichtete Summe aller berechneten Redoxspannungen wird in etwader gemessenen Redoxspannung im Wasser entsprechen. In der Praxis wird aber die analytischeErfassung aller Redoxpaare aufgrund des großen analytischen Aufwandes selten erfolgen.

4. Wechselwirkungen zwischen wässriger Phase und Mineralphase

4.1. Lösung / Fällung

Sättigungsindex

Für die Berechnung von Minerallösung und -ausfällung nach dem MWG ist der Sättigungsindex SIvon großer Bedeutung. Für die Reaktion

AB ↔ A + B

ist der Sättigungsindex SI definiert als log (IAP / K)

dabei steht IAP für das Ionenaktivitäts-Produkt in wäßriger Phase (im Beispiel [A] * [B]) und K fürdas Löslichkeitsprodukt der zugehörigen Mineralphase [AB]. Der Sättigungsindex SI gibt an, ob eineLösung im Gleichgewicht mit einer festen Phase, oder unter- bzw. übersättigt in Bezug auf die festePhase ist. Hierbei bedeutet ein Wert von z.B. +1 eine 10fache Übersättigung, ein Wert von -2 eine100fache Untersättigung bezüglich einer bestimmten Mineralphase.

In der Praxis kann angenommen werden, daß im Bereich von -0,2 bis +0,2 von einem quasi-Gleichgewicht gesprochen werden kann. Dabei sollte aber nicht vergessen werden, daß der SI-Wertin einer logarithmischen Skala angegeben wird: SI = 0,1 bedeutet, daß ca. 25% mehr Ionenaktivitätim Vergleich zum Löslichkeitsprodukt (umgerechnet auf die jeweilige Temperatur) gemessen wurde;SI = 0,2 ca. 50%, SI = 0,3 ca. 100 % usw. Das Ionenaktivitätsprodukt (IAP) errechnet sich aus derSumme der Aktivitäten der in ihrem stöchiometrischen Verhältnis an der Mineralphase beteiligtenSpecies.

So berechnet sich z.B. der SI für die Mineralphase Calcit aus der Reaktionsgleichung

CaCO3 ↔ Ca2+ + CO32-

Daraus folgt: IAP = [Ca2+] . [CO32-]

-10-

Im Labor wird in der Regel die Konzentration von Calcium als Summe aller Calcium-Speciesbestimmt. Bei einer Berechnung muss man also zwei Dinge berücksichtigen:

1. Einerseits muß aus der Ca-Konzentration die Ca-Aktivität berechnet werden und

2. andererseits muß berücksichtigt werden, daß nicht nur freie Ca2+ Ionen vorliegen, sondernCa auch in Komplexe gebunden wird, die sich mit anorganischen und organischen Ligandenbilden können. Im Beispiel des Calciums können dies u.a. die Komplexe CaOH-, CaCO3

0,CaHCO3

+ sein. Ist im Wasser zudem Sulfat vorhanden können sich auch dieSulfatokomplexe CaHSO4

+ und CaSO40 bilden. Dadurch geht mehr Ca in Lösung als bei

einfacher Betrachtung mit dem MWG berechnet, da dabei die Komplexe nicht erfasstwerden.

Mit hydrogeochemischen Programmen ist es sowohl möglich, eingegebene Konzentrationen direkt inAktivitäten umrechnen zu lassen als auch die unter 2. erwähnte Komplexbildung zu berücksichtigen.Zudem kann man sich Sättigungsindices berechnen lassen und andererseits auch ein Wasser„rechnerisch“ mit einer Mineralphase ins Gleichgewicht oder auch in ein Ungleichgewicht zu setzen,indem ein SI Wert von z.B. 0 (Gleichgewicht) oder +1 (10 fache Übersättigung) vorgegeben wird.

Begrenzende Mineralphasen

Wasser, das mit Mineralphasen in Berührung kommt, kann aus der jeweiligen Mineralphase immernur so viele Ionen lösen, wie sich aus dem Löslichkeitsprodukt ergibt. Dabei kommt es u.a. daraufan, welche Ionen bereits im Wasser gelöst sind und inwieweit diese Ionen untereinander Komplexegebildet haben. Veränderungen der Wassertemperatur und der Gaspartialdrücke (z.B. pCO2)können dazu führen, daß eine Übersättigung bezüglich einer Mineralphase auftritt und dieses Mineraldann ausgefällt wird.

Werden zwei Wässer gemischt (z.B. ein Grundwasser mit einem sauren Grubenwasser), so kanneine völlig neue Situation entstehen. Dies betrifft sowohl die Komplexbildung der im Wasser gelöstenIonen als auch die Mineralphasen. Die Bildung oder Veränderung aquatischer Komplexe ist in derRegel ein sehr schneller Prozeß; demgegenüber ist die Lösung oder Fällung von Mineralen meistensdeutlich langsamer. Im Zuge der Mischung von 2 Wässern kann es dazu kommen, daß dasMischwasser gegenüber einer oder mehrerer Mineralphasen übersättigt ist. Wird beispielsweise einsulfathaltiges Grundwasser mit einem anderen Grundwasser gemischt, das kein Sulfat enthält, aberBarium- und Chloridionen, so wird durch die Mischung dieser beiden Wässer das MineralBariumsulfat als begrenzende Phase wirksam werden, da das Mischwasser sowohl Barium- alsauch Sulfationen enthält. Es wird also so lange Bariumsulfat ausgefällt, bis der Sättigungsindexbezüglich BaSO4 den Wert von Null erreicht hat.

Bei der Interpretation hydrogeochemischer Modellrechnungen und insbesondere berechneterSättigungsindizes (SI) ist aber immer zu prüfen, ob die Mineralphase, für die eine Übersättigungausgerechnet wird, unter den gegebenen Temperatur- und Druckbedingungen tatsächlich gebildet

-11-

wird. Herkömmliche hydrogeochemische Codes (PhreeqC, EQ3/6 etc.) bieten diesbezüglich zumberechneten SI-Wert keine zusätzlichen Informationen. Vor allem muß die Kinetik einerMineralbildung (und Minerallösung) berücksichtigt werden. Das erfordert viel chemisches Wissenund im Zweifelsfall intensives Literaturstudium oder Laborversuche zu einer bestimmtenMineralphase.

Schwieriger wird die Berechnung eines SI, wenn nicht reine Mineralphasen sondernMischmineralphasen (solid solutions) in Betracht gezogen werden. Die Untersuchung natürlicherMinerale zeigt, daß die wenigsten Mineralphasen reine Phasen sind, sondern vielfach insbesondereseltene Elemente in wechselnden Gehalten beinhalten. Klassische Beispiele für solid solutionMinerale sind der Dolomit, das Calcit/Rhodochrosit System oder das System Calcit/Otavit(CdCO3). Im folgenden soll das Konzept der Berechnung von Löslichkeiten gemäß solid solutionMineralen in kurzer Form am Beispiel von Calcit und Strontianit dargestellt werden. Die Löslichkeitbeider reiner Mineralphasen errechnet sich aus:

s3

-2

3

2

Calcit ][CaCO}{CO}{Ca

K⋅

=+

und

s3

-2

3

2

Stronianit ][SrCO

}{CO}{Sr K

⋅=

+

Bei Annahme eines solid solution Minerals, das aus einer Mischung dieser beiden Minerale bestehtergibt sich durch Umformung der beiden Gleichungen:

Calcit

tStrontiani

s3

s3

2

2

KK

][CaCO][SrCO

}{Ca }{Sr

⋅=+

+

Daraus folgt, daß für ein bestimmtes Aktivitätsverhältnis von Sr zu Ca im Wasser sich auch einbestimmtes Verhältnis im solid solution Mineral ergibt. Wenn für das nicht-ideale Verhalten derMineralphase in Analogie zum Aktivitätkoeffizienten der aquatischen Species ein Aktivitäts-Korrekturfaktor fcalcit und fStrontianit eingeführt wird, ergibt sich:

tStrontiani

Calcit

2

2

Calcit

tStrontiani

Calcit

tStrontiani

XX

}{Ca }{Sr

ff

KK

⋅=⋅+

+

Dabei steht X jeweils für den molaren Anteil im solid solution Mineral. Das Verhältnis der beidenAktivitätskoeffizienten der Mineralphasen kann im einfachsten Fall zu einemVerteilungskoeffizienten zusammengefaßt werden. Dieser kann durch semi-empirischeApproximation aus Laborversuchen bestimmt werden. Unter Verwendung der Löslichkeitsproduktevon Calcit und Strontianit und Annahme einer Calcium-Aktivität von 1,6 mmol/L, einemVerteilungskoeffizienten von 0,8, sowie einem Verhältnis von 50:1 im solid solution Mineral ergibtsich somit eine Aktivität für Strontium von

-12-

mol/L102.498.010

101.60.020.810

}{CaX

Xf

fK

K}{Sr

6

48.8

3-9.271-

2

Calcit

tStrontiani

Calcit

tStrontiani

Calcit

tStrontiani2

−

−

++

⋅=⋅

⋅⋅⋅⋅=

=⋅⋅⋅=

Würde man Strontianit als begrenzende Phase ansetzen, so könnte erheblich mehr Strontium(Aktivität ca. 2,4.10-4 mol/L) in Lösung sein als bei Ansatz der solid solution Mineralphase.

Dieses Beispiel zeigt die Tendenz, die sich im Zusammenhang mit solid solution Mineralen ergibt: Esbesteht eine Übersättigung oder ein Gleichgewicht hinsichtlich des solid solution Minerals, aber eineUntersättigung gegenüber beiden reinen Mineralphasen, d.h. das solid solution Mineral wird gebildet,aber keine der beiden reinen Phasen. Die Größe dieses Effektes ist abhängig von den Zahlenwertendes Aktivtätskoeffizienten der solid solution Komponente.

Zwei konzeptionelle Modelle der Bildung von „solid solution“ Minerale müssen unterschiedenwerden:

1. „end member“ Modell (beliebige Mischungen zwischen 2 oder mehreren Mineralen)2. „site mixing“ Modell (Fremdelemente können ein bestimmtes Element nur auf bestimmten

Plätzen der Kristallstrukur ersetzen)

In der Version 7 des Programms EQ3/6 sind beide Modelltypen für solid solution Minerale realisiert.Das Problem der solid solution Minerale ist, daß belastbare Daten nur für wenige Systeme verfügbarsind. Die folgende Tabelle (Tab. 1. ) zeigt die „solid solution“ Minerale, die im Gembochs Datensatzvon EQ3/6 vorhanden sind.

Tab. 1. „Solid Solution“ Minerale aus dem Gembochs Datensatz des Programms EQ3/6

System Mineral-FormelBiotite K(Mg,Fe)3AlSi3O10(OH)2

Carbonate-Calcite (Ca,Mn,Zn,Mg,Fe,Sr)CO3

Chlorite-ss (Fe,Mg)5Al2Si3O10(OH)8

Clinoptilolite-hy-ss (Na,K,Cs,NH4,Ca.5,Sr.5)3.467Al3.45(Fe3+)0.017SiClinoptilolite-ss (Na,K,Cs,NH4,Ca.5,Sr.5)3.467Al3.45(Fe3+)0.017SiEpidote-ss Ca2(Fe,Al)Al2Si3O12(OH)Garnet-ss Ca3(Al,Fe)2Si3O12

Olivine (Fe,Mg)2SiO4

Orthopyroxene (Fe,Mg)SiO3

Plagioclase CaAl2Si2O8-NaAlSi3O8

Sanidine-ss (K,Na)AlSi3O8

Saponite-tri (Ca.5,H,K,Mg.5,Na).33Mg3Al.33Si3.67O10(OH)2

Smectite-di (Na,K,Ca.5,Mg.5).33(Al,Mg,Fe)2(Si,Al)4O10(OH)2

-13-

Für einige Elemente haben begrenzende Mineralphasen (reine Mineralphasen und solid solutionMinerale) ohnehin keine Relevanz, weil die Konzentration dieser Elemente über Sorption anTonmineralen, Eisenhydroxiden oder Huminstoffen limitiert wird.

4.2. Sorption

Unter dem Begriff Sorption werden die Termini Adsorption, Absorption, Chemiesorption,Ionentausch und Oberflächenkomplexierung zusammengefaßt. Unter Absorption wird derdiffusive Ein- und Austrag von Wasserinhaltsstoffen in die poröse Matrix eines Gesteins verstanden.Es ist somit zunächst kein oberflächenkontrollierter Prozeß. Chemiesorption dagegen sind Prozessean Oberflächen unter Beteiligung spezifischer Gitterstellen und/oder funktioneller Gruppen.Ionentausch kann sowohl als Chemiesorption als auch als Adsorption verstanden werden.

Von besonderer Bedeutung im Zusammenhang mit Sorption und Kationentausch sind Tonminerale(Aluminiumsilikate), Metalloxide (überwiegend Eisen- und Manganoxide) und organisches Material.Tonminerale bestehen aus 1-n Schichten von Si-O Tetraedern und 1-n Schichten vonAluminiumhydroxid Oktaedern (Gibbsit). Dabei ersetzt Al sehr oft das Si in der Tetraeder-Schichtund Magnesium das Al in der Oktaeder-Schicht. Am Ende des Verwitterungsprozesses von Eisen-und Manganmineralen steht in der Regel die Bildung von Oxiden. Manganoxide bilden sich meistensin oktaedrischer Anordnung; diese ähnelt weitgehend dem Gibsit. Auch Hämatit (Fe2O3) und Goethit(FeOOH) haben eine dem Gibsit ähnliche oktaedriche Struktur. 70 bis 80 Gew.% der organischenMaterie in der Natur sind nach SCHNITZER (1986) zur Gruppe der Huminstoffe zu rechnen. Diessind kondensierte Polymere von aromatischen und aliphatischen Komponenten, die beim Abbaulebender Zellen von Pflanzen und Tieren durch Mikroorganismen entstehen. Huminstoffe sindhydrophil, von dunkler Farbe und weisen Molekularmassen von weniger Hundert bis zu vielenTausend auf. Sie besitzen unterschiedlichste funktionale Gruppen, die in der Lage, sind mitMetallionen in Wechselwirkung zu treten. Huminstoffe (refraktäre organische Säuren) könnenunterteilt werden in Huminsäuren und Fulvinsäuren. Huminsäuren sind unter alkalischen Bedingungenlöslich und fallen unter sauren Bedingungen aus. Fulvinsäuren sind sowohl unter sauren als auchbasischen Bedingungen löslich.

Die Kationen-Austauschkapazität (KAK) von Tonmineralen liegt zwischen 3 und 160 mol/kg. Dieseextrem hohen Austauschkapazitäten beruhen auf zwei physikalischen Gründen:

- extrem große Oberfläche und- einer elektrischen Ladung der Oberflächen.

Diese elektrische Ladungen können unterschieden werden in

- permanente Ladungen und- variable Ladungen.

-14-

Permanente Oberflächenladungen sind auf die Substitution von Metallen in Kristallgittern(Isomorphie) zurückzuführen. Da in der Regel die Substitution durch Metalle mit geringerer Valenzerfolgt, resultiert daraus ein Gesamtdefizit an positiver Ladung für den Kristall. Um diesauszugleichen, bildet sich ein negatives Potential an der Oberfläche, wodurch positiv geladeneMetallionen sorbiert werden können. Die Oberflächenladungen von Tonmineralen sind überwiegendauf Isomorphie zurückzuführen und somit größtenteils permanenter Natur. Dies gilt aber nicht für alleTonminerale; bei Kaolinit sind es weniger als 50% (BOHN et al. 1979).

Neben der permanenten Oberflächenladung existieren Oberflächenladungen, die vom pH-Wertdes Wassers abhängig sind. Diese entstehen durch Protonierung und Deprotonierung funktionellerGruppen der Oberflächen. Unter sauren Bedingungen werden Protonen an den funktionalenGruppen sorbiert, die zu einer insgesamt positiven Ladung der Oberfläche führen. Das Mineral bzw.Teile des Minerals wirken somit als Anionenaustauscher. Unter hohen pH-Werten verbleiben dieSauerstoffatome an den funktionellen Gruppen deprotoniert und das Mineral bzw. ein Teil desMinerals zeigt insgesamt eine negative Oberflächenladung und kann Kationen sorbieren. Im mittlerenpH-Bereich bewirkt die teilweise Protonierung der funktionellen Gruppen, daß sich die Oberflächeinsgesamt neutral verhält. Der zugehörige pH-Wert wird als PZCpH (Point of Zero Charge) bezeichnet. Abb. 3. zeigt das pH abhängigeSorptionsverhalten von Eisenhydroxid Oberflächen. Das Gesamtpotential der pH abhängigenOberflächenladung ist dabei nicht abhängig von der Ionenstärke des Wassers.

Abb. 3. Schematische Darstellung des pH abhängigen Sorptionsverhaltens von Eisen-hydroxidOberflächen durch Anlagerung von H+ bzw. OH- (nach SPARKS 1986) und desPZCpH.

Natürliche Systeme sind vielfach eine Mischung aus Mineralen mit konstanter und variablerOberflächenladung. Das Verhalten in Bezug auf Anionen- und Kationensorption zeigt Abb. 4.Oberhalb pH 3 nimmt die Anionenaustauschkapazität deutlich ab. Die Kationenaustauschfähigkeitbleibt bis ca. pH 5 konstant, um dann extrem anzusteigen

-15-

Abb. 4. Kationen- und Anionenaustauschverhalten von Mineralen als Funktion des pH-Wertes(nach BOHN et al. 1979)

Mathematische Beschreibung der Sorption

Die einfachste Möglichkeit, den experimentellen Befund einer Wechselwirkung einesWasserinhaltsstoffes zwischen fester und flüssiger Phase mathematisch zu beschreiben, ist einelineare Regressionsgleichung:

C* = Kd CmitC* : Masse der am Mineral sorbierten Substanz (mg/kg)Kd: Verteilungskoeffizient (distribution coefficient)C: Konzentration der Substanz im Wasser (mg/L)

Diese lineare Sorptionsisotherme hat einen Vorteil und zwei gravierende Nachteile. Der Vorteilliegt in der mathematischen Einfachheit und der Möglichkeit, die linerare Regressionsgleichung ineinen Retardationsfaktor rf zu überführen und somit als Korrekturterm in die allgemeineTransportgleichung einzuführen [rf = 1 + Bd/θ Kd mit Bd = Feststoffdichte und θ = Wassergehalt].Der erste Nachteil ist, daß der Zusammenhang linear beschrieben wird und somit die Sorption durchnichts beschränkt ist. Es fehlt also zumindest ein oberes Limit. Der zweite noch gravierendereNachteil liegt in der Nicht-Berücksichtigung weiterer Randbedingungen (pH-Wert, Redoxpotential,Ionenstärke, Konkurrenzreaktion um Bindungsplätze etc.). Die Plumpheit des Modells (lumpedparameter) führt dazu, daß es als Prognosemodell kaum geeignet ist. Ergebnisse von Labor- undFeldversuchen sind damit nicht auf reale Systeme übertragbar.

Mit Hilfe einer FREUNDLICH-Isotherme kann ein exponentieller Zusammenhang zwischensorbierten und gelösten Molekülen beschrieben werden.

-16-

C* = Kd C N

Dabei ist N eine weitere empirische Konstante. Der FREUNDLICH-Isotherme kann dasphysikalische Modell einer multilaminellen Belegung der Feststoffoberfläche zugrunde liegen, wobeidavon ausgegangen wird, daß die Bindungsenergie der elektrostatischen Kräfte mit zunehmendemBelegungsgrad der Oberfläche abnimmt. Mit der Verwendung der FREUNDLICH-Isotherme wirddas Manko der fehlenden Limitierung des linearen Modells behoben, nicht jedoch die Nicht-Berücksichtigung weiterer Randbedingungen. Vorteil ist wiederum die Möglichkeit der Überführungin einen Retardationsfaktor rf.

Die LANGMUIR-Isotherme wurde entwickelt, um die Fähigkeit einer Oberfläche mit begrenzterAnzahl von Sorptionsplätzen zu beschreiben:

C1C

C*

⋅+⋅⋅=

αβα

mitα: Sorptions-Konstante als Funktion der Bindungsenergie (L/mg)β: maximal sorbierbare Masse der Substanz (mg/kg)

Auch die LANGMUIR-Isotherme verwendet ausschließlich nichtdeterministische, empirischeParameter (lumped parameter). Die oben aufgeführten gravierenden Nachteile gelten somitgleichermaßen.

Aus wissenschaftlicher Sicht sind einfache Ansätze im Sinne des Kd-Konzeptes (lineare Isothermebzw. FREUNDLICH oder LANGMUIR) unbefriedigend, da die komplexen Prozesse an denOberflächen nicht mittels empirischer Fittingparameter abgebildet werden können. Die Verwendungdes Kd-Konzeptes für Prognosen ist deshalb nur dann gerechtfertigt, wenn Änderungen derRandbedingungen (pH-Wert, Redox-Potential etc.) nicht zu erwarten sind und Parameterfür einen deterministischen und mechanistischen Ansatz nicht ermittelt werden können.Beschreibung des Ionenaustausches

Der Kationen- oder Anionenaustausch kann in Form chemischer Gleichgewichtsreaktionen gemäßdem Massenwirkungsgesetz formuliert werden und geht somit von einer vollständigen Reversibilitätaus. Der Vorteil dieses Ansatzes ist, daß quasi beliebig viele Species in gegenseitige Wechselwirkungum Bindungsplätze an der Oberfläche eines Minerals treten können. Der Austausch von Natriumgegen Calcium kann z.B. wie folgt geschrieben werden:

Na+ + ½ CaX2 ↔ NaX + ½ Ca2+ daraus folgt

-17-

}{Na][CaX}{Ca[NaX]

K0,5

2

0,52

Na

Ca +

+

⋅⋅

=

Diese Notation wird als GAINES-THOMAS Konvention (1953) bezeichnet. Werden statt derAktivitäten die molaren Konzentrationen benutzt, so entspricht dies der VANSELOW Konvention(1932). GAPON (1933) schlug die folgende Schreibweise des Problems vor:

Na+ + Ca0,5X ↔ NaX + ½ Ca2+ daraus folgt

}{NaX][Ca}{Ca[NaX]

K0.5

0,52

Na

Ca +

+

⋅⋅

=

Die Proportionalitätskonstante KCa/Me (spezifischer Austauschkoeffizient) ist im Gegensatz zu denKomplexbildungskonstanten oder Dissoziationskonstanten nicht nur von Druck, Temperatur und derIonenstärke abhängig. Sie ist vor allem auch eine materialspezifische Konstante der jeweiligenFeststoffphase. Sie ist zudem in gewissem Maße auch abhängig von der verwendeten Konvention.Letztlich ist der Kationenaustausch in seiner Selektivität gegenüber verschiedenen Species und seinerGesamtkapazität eine Funktion der inneren und äußeren Oberflächen der Feststoffphase.

Die mathematische Behandlung des Kationen- oder Anionenaustausches ist mit Hilfethermodynamischer Modelle durchführbar, da diese Reaktionen wie eine normale Reaktion im Sinnedes Massenwirkungsgesetzes geschrieben werden können. Als zusätzliche Randbedingung muß imnumerischen Modell die Gesamtaustauschkapazität der festen Phase berücksichtigt werden.

Ein anderer Weg zur Beschreibung des Kationenaustausches besteht im Ansatz von GOUY-CHAPMAN (diffuse double layer theory). Diese Theorie geht auf die Beobachtung zurück, daßim Bereich der Grenzschicht feste Phase zu flüssiger Phase unabhängig von der OberflächenladungKationen und Anionen in einer erhöhten Konzentration vorliegen. Diese erhöhte Konzentration wirddurch die elektrostatischen Kräfte der Oberflächen verursacht. Diesen Kräften entgegen wirkt dieDiffusion, die dazu führt, daß mit zunehmendem Abstand von der Grenzschicht eine Verdünnungeintritt. Der Zusammenhang kann mathematisch durch die Poisson-Boltzmann-Gleichung beschriebenwerden.

OberflächenkomplexierungZusätzlich zum Ionentausch kann es an bestimmten Feststoffen (Eisen-, Aluminium-, Silicium undManganhydroxide, Huminstoffe) zu Komplexierungsreaktionen an Oberflächen kommen. Solchefunktionalen Gruppen enthalten immer O, N und S. Die Reaktion kann nach STUMM & MORGAN(1981) generell wie folgt geschrieben werden:

{GH} + Mez+ ↔ {GMez-1} + H+

-18-

wobei GH für eine funktionale Gruppe wie z.B (R-COOH)n oder (�AlOH)n steht. DieKomplexbildungsfähigkeit funktionaler Gruppen ist stark vom Säure-Base-Verhalten abhängig undberücksichtigt somit die Veränderung des pH-Wertes in einem aquatischen System. In PhreeqC istdas Konzept von DAVIS & KENT (1990) realisiert. Wahlweise kann in PhreeqC auch das Modellder Oberflächenkomplexierung gemäß „diffuse double-layer“ nach DZOMBAK und MOREL(1990) verwendet werden.

5. Wechselwirkungen zwischen wässriger und gasförmiger Phase

Die Löslichkeit von Gasen im Wasser kann bei bekannter Temperatur und bekanntem Partialdruckmit Hilfe des linearen Henry-Gesetzes beschrieben werden:

mi = KHi * pi

mit mi= Molarität des Gases i [mol/kg]KHi = Henry Konstante des Gases ipi = Partialdruck des Gases i[bar]

Wie aus dem Henry-Gesetz bereits zu ersehen, nimmt die Gaslöslichkeit mit steigendem Partialdruckzu. Mit steigender Temperatur und steigender Gesamtmineralisation dagegen nimmt dieGaslöslichkeit ab.

Streng genommen gilt das Henry-Gesetz nur für Gase, die in Lösung nicht weiter reagieren, wie z.B.

N2, O2, Ar. Bei NH4 und CO2 dagegen ist die Anwendung des Henry-Gesetzes nur bei gleichzeitiger

Berücksichtigung der Folgereaktionen sinnvoll. Obwohl im Falle des Kohlendioxids z.B. nur ca. 1%

der H2CO3 in Abhängigkeit vom pH-Wert in HCO3- und CO3

2- dissoziiert, kommt es durch die sich

anschließenden Komplexbildungen mit Kationen zu einer deutlich höheren Lösung von CO2 in

Wasser. Werden zudem auch noch Mineralphasen (z.B. Calcit) gelöst und dabei Protonen

verbraucht, so sind dies weitere Folgereaktionen, durch die wesentlich mehr CO2 in Lösung geht, als

nach dem HENRY-Gesetz zunächst berechnet wird.

6. Praxis der hydrogeochemischen Modellierung

Hydrogeochemische Modelle sind inzwischen ein Standardwerkzeug in der Bearbeitunghydrogeochemischer und hydrogeologischer Fragestellungen. Eine Auswahl von Programmen ist inTab. 2. aufgelistet.

-19-

Tab. 2. Auswahl hydrogeochemischer Programme und deren zeitliche Entwicklung

Hydrogeochemische Programme

* WATEQ (1974) WATEQF (1976) WATEQ2 (1979 WATEQ3 (1981)WATEQ4F (1987)

* PHREEQE (1982) PHRQPITZ (1988) PHREEQX PHREEQM(1993)PhreeqC (1995, 1997) PHREEQC 2 (2000)

* REDEQL (1972) REDEQL2 (1973) GEOCHEM (1979)

* MINEQL (1976) MINTEQ (1984) HYDRAQL (1988) MINTEQA2 (1991)

* SOLMINEQ (1973) SOLMINEQ.88 (1988/90)

* EQ3/6 (1978 - 2000)

6.1. Genereller Aufbau der Programme

Der Aufbau hydrogeochemischer Modellierungsprogramme besteht aus input, output sowie derRechengrundlage, d.h. einer Datenbank, in der thermodynamische Parameter, wieLöslichkeitskonstanten für bestimmte Reaktionsgleichungen definiert sind.

InputAls Eingabeparameter aller hydrogeochemischen Modellierungen bedarf es einer möglichstvollständigen und exakten wasserchemischen Analyse.

DatenbankenGrundsätzlich ist es bei nahezu allen Programmen möglich, sich eine eigene thermodynamischeDatenbank zu erstellen. Dies ist allerdings ein erheblicher Aufwand und erfordert große Sorgfalt. Ausdiesem Grund wird man in der Regel auf bestehende Datensätze zurückgreifen. Tab. 3. zeigt eineAuswahl von thermodynamischen Datenbanken und die jeweils berücksichtigten Elemente.

Zu beachten ist beim Vergleich von Löslichkeitsprodukten aus verschiedenen Datenbanken, obdieselben Reaktionsgleichungen zugrunde liegen. Unterschiedliche Reaktionsgleichungen für dieBildung desselben Minerals liefern unterschiedliche Löslichkeitsprodukte. Löslichkeitsprodukte undKomplexbildungskonstanten müssen daher - soweit sie der Literatur entnommen werden – stetseindeutig mit einer Reaktionsgleichung verknüpft sein.

Tab. 3. Auswahl thermodynamischer Datenbanken und der darin berücksichtigten Elemente

-20-

database NEA PhreeqC WATEQ4F CHEM-VAL6

GEM-BOCHS.V2EQ8

HATCHES NAGRANTB

MINEQL-PSI

MINTEQA2

form forprogram

PHREEQEQ 3/6

PhreeqC WATEQ4FPhreeqC

PHREEQEdBASE

EQ3/6 PHREEQE PHREEQE

PHREEQE

MINTEQA2PhreeqC

Ag + + + +Al + + + + + + + + +Am + + + + + +Ar +As + + + + + + + +Au +B + + + + + + + + +Ba + + + + + + + + +Be + +Br + + + + + + + + +C + + + + + + + + +Ca + + + + + + + + +Cd + + + + + +Ce +Cl + + + + + + + + +Cm +Co +Cr + +Cs + + + +Cu + + + + + +Dy +Er +Eu +F + + + + + + + + +Fe + + + + + + + + +Ga +Gd +H + + + + + + + + +Ha +He +Hg + + +Ho +I + + + + + + + +In +K + + + + + + + + +Kr +La +Li + + + + + + + + +Lu +Mg + + + + + + + + +Mn + + + + + + + +Mo + + + +N + + + + + + + + +Na + + + + + + + + +Nb +Nd + +Ne +Ni + + + + + + +Np + + + + +O + + + + + + + + +P + + + + + + + + +Pa +Pb + + + + + +Pd + + + +Pm +Pr + +

-21-

Pu + + + + +Ra + +Rb + + + +Re +Rn +Ru +S + + + + + + + + +Sb + +Sc +Se + + + + + + +Si + + + + + + + + +Sm + +Sn + + + + +Sr + + + + + + + + +T b + +Tc + + + + +T h + + + + +Ti +Tl + +Tm + +U + + + + + + + +V + +W +Xe +Y +Yb + +Zn + + + + +Zr + + + +

Konsistenz thermodynamischer Daten

Thermodynamische Daten werden durch Laborversuche unter definierten Randbedingungen(Temperatur und Ionenstärke) ermittelt. Dabei werden häufig relativ hohe Ionenstärken (0,1 oder 1molare Lösungen) verwendet. Für die Rückrechnung der Komplexbildungs-Konstanten oder auchder Löslichkeitsprodukte auf eine Ionenstärke von Null werden die gleichen Verfahren wie zurBerechnung der Aktivitäten aus gemessenen Konzentrationen verwendet (z.B. erweiterte Debye-Hückel Gleichung). Da die Gültigkeit der Ionendissoziationstheorie spätestens bei 1-molarenLösungen endet, bewegen sich solche Experimente zum Teil in einem Bereich, der nicht mehr odernur eingeschränkt durch die Ionendissoziationstheorie abgedeckt ist.

Werden Löslichkeitsprodukte und Komplexbildungskonstanten der Literatur entnommen, so werdenauf diese Weise Daten zusammengetragen, die unter unterschiedlichen experimentellenRandbedingungen und unterschiedlichen Berechnungsmethoden zur Ermittlung von Konstanten beiIonenstärke Null gewonnen wurden. Damit ergibt sich zwangsläufig eine Inkonsistenz der Daten. Esstellt sich somit die Frage, ob es besser ist, mit einem in sich konsistenten und damit möglicherweisekleineren und unvollständigen Datensatz zu arbeiten oder aber einen sehr umfangreichen undhinsichtlich der Species quasi vollständigen aber möglicherweise inkonsistenten Datensatz zuverwenden. Diese Frage könnte nur von Fall zu Fall durch eine Verifikation mit unabhängigenMeßdaten (Benchmarks) entschieden werden. Dies ist aber leider nur selten realisierbar, da dieanalytische Bestimmung der Species in einem natürlichen Wasser sehr aufwendig ist.

-22-

OutputIm Output werden die berechneten Ergebnisse angezeigt. Standardmäßig sind das z.B. für dasProgramm PhreeqC2 die Wiederholung der eingelesenen Datenbank und des input-files, dieLösungs-Zusammensetzung (einzelne Elemente in mol/L), die Beschreibung der Lösung (pH, pe,Ionenstärke, Temperatur, Ladungsgleichgewicht), die Aufführung der Redoxpaare (falls vorhanden),die Speziesverteilung und die Sättigungsindizes der Eingabelösung. Erst danach folgt unter „Beginningof reaction calculations“ die Ausgabe der Aufgaben-spezifischen Ergebnisse (veränderte Lösung).

6.2. Vergleich der Programme

Die gebräuchlichsten hydrogeochemischen Modelle sind WATEQ4F, MINTEQA2, PhreeqC2 undEQ3/6. Mit WATEQ4F können lediglich Speziationsberechnungen durchgeführt werden. DasAnwendungsspektrum von PhreeqC2 und EQ3/6 geht weit darüber hinaus. Typische Anwendungensind:

• Speziierung/Komplexierung• Gleich/Ungleichgewichtseinstellung bezüglich Mineralen oder Gaspartialdrücken (Sättigungsindex)• Mischen verschiedener Wässer• Auswirkungen von Temperaturänderungen simulieren• stöchiometrische Reaktionen (z.B. Titrationen)• Verdunstungseinflüsse berechnen• Reaktionen mit festen, flüssigen und gasförmigen Phasen (offenes/geschlossenes System)• Sorption (Kationenaustausch, Oberflächenkomplexierung)• Genese von Wässern durch inverse Modellierung ermitteln (reaction path finding)• Reaktiver StofftransportDas Programm PHREEQE, Vorläufer des heutigen PhreeqC2, wurde 1982 (PARKHURST et al.1982) in der Programmiersprache FORTRAN geschrieben und später durch eine Version mitPitzer-Gleichungen für Ionenstärken > 1 mol/L ergänzt (PHRPITZ). PHREEQX und PHREEQM(APPELO & POSTEMA 1994) sind Versionen, die den ursprünglichen Code in bestimmten Dingenerweiterten. 1995 wurde das Programm in der Programmiersprache C neu geschrieben (PhreeqC:PARKHURST 1995). Dabei wurde der Gleichungslöser überarbeitet und eine Reihe vonzusätzlichen Optionen vorgesehen. Ferner wurde ein einfaches 1-dimensionales Transportmodulimplementiert, das vergleichbar ist zu dem in PHREEQM (APPELO & POSTEMA 1994).PHRKIN ist ein Zusatzmodul zu PHREEQE zur Modellierung kinetisch gesteuerter Reaktionen.

PhreeqC2 (PARKHURST & APPELO 1999) ist die aktuellste Version, die u.a. folgendezusätzliche Simulationen ermöglicht:

• Dispersion oder Diffusion• kinetische Reaktionen mit benutzerdefinierten Umsatzraten• Bildung idealer oder nichtidealer solid solution minerals

-23-

• Veränderung der Zahl an Austauscherplätzen mit Lösung oder Fällung von Reaktanten• Einbeziehung von Isotopenbilanzen in Berechnungen inverser Modelierung

Zudem ist es möglich, den Daten-output benutzerdefiniert zu verkürzen und in einem spreadsheetkompatiblem Dateiformat zu exportieren. Zur Programmierung für Benutzer-spezifischerFragestellungen ist ein BASIC Code im Programm implementiert.

EQ3/6 besteht im Gegensatz zu PhreeqC2 aus zwei Programmen: EQ3 ist ein reinerSpeciierungscode, dessen Ergebnisse für weitergehende Fragestellungen innerhalb EQ6weiterverarbeitet werden. Insgesamt birgt EQ 3/6 insbesondere im Vergleich zum neuen PhreeqC2kaum mehr Vorteile. Vor allem das sehr viel voluminösere Datenformat für die thermodynamischenGrunddaten und die Inputfiles ist als Nachteil zu betrachten. Tab. 4. und Tab. 5. zeigen denVergleich des Datenformats in EQ 3/6 und PhreeqC2 exemplarisch für die Definition des MineralsRutherfordine (UO2CO3).

Tab. 4. Auszug aus der NEA - Datenbasis für EQ3/6; Definition des Minerals Rutherfordine (fettmarkiert sind die Elemente, die in ähnlicher Form auch im PhreeqC Datensatz stehen)

UO2CO3

date last revised = 02-jul-1993

keys = solid

V0PrTr = 0.000 cm**3/mol (source =

* mwt = 330.03690 g/mol

3 chemical elements =

1.0000 C 5.0000 O 1.0000 U

4 species in data0 reaction

-1.0000 UO2CO3 -1.0000 H+

1.0000 HCO3- 1.0000 UO2++

* log k grid (0-25-60-100/150-200-250-300 C) =

-3.8431 -4.1434 -4.4954 -4.7855

-5.0616 -5.2771 500.0000 500.0000

* Extrapolation algorithm: constant enthalpy approximation

* gflag = 1 (reported delG0f used)

* basic source = 92gre/fug

* delG0f = -1563.046 kj/mol

* delH0f = -1689.647 kj/mol

* S0PrTr = 144.200 j/(mol*K)

Tab. 5. Auszug aus der WATEQ4F - Datenbasis für PhreeqC; Definition des MineralsRutherfordine

Rutherfordine 606

UO2CO3 = UO2+2 + CO3-2

log_k -14.450

-24-

delta_h -1.440 kcal

Der Vorteil des PhreeqC Datenformats liegt eindeutig in seiner Kompaktheit und der Tatsache, daßdie Reaktionsgleichung in der Syntax chemischer Formeln geschrieben ist.

Wichtig ist zu beachten, daß EQ 3/6 und PhreeqC2 für die Definition des Minerals Rutherfordinezwei unterschiedliche Reaktionsgleichungen verwenden. Während PhreeqC2 die Reaktionsgleichung

UO2CO3 = UO22+ + CO3

2-

verwendet, ist es in EQ3/6 die Gleichung

UO2CO3 + H+ = HCO3- + UO2

2+

Deshalb kann auch das Löslichkeitsprodukt, wie bereits einleitend zum Thema „Datenbanken“erwähnt, nicht gleich sein kann.

Im folgenden wird gezeigt, wie eine einfache Aufgabe mit PhreeqC und EQ3/6 formuliert wird. Dazusoll die Lösung des Minerals Rutherfordine in einem Wasser mit jeweils 1 mmol/L Natrium undChlorid und geringen Sulfatkonzentrationen (0,0001 mmol/L) unter oxidierenden Bedingungen (pe =14) bei 25°C und einem CO2-Partialdruck von 0,00033 bar (Atmosphärendruck) simuliert werden.Auch hier wird deutlich, daß die Definition eines Problems in PhreeqC2 leichter und damit schnellerdurchgeführt werden kann.

Tab. 6. Beispiel für INPUT-File EQ3/6 (Lösung von Rutherfordine als Funktion des CO2-Partialsdrucks)

EQ3NR input file name= co3aqui.3iDescription= "Uranium Carbonate solution"Version level= 7.2

endit. Tempc= 2.50000E+01 rho= 1.00000E+00 tdspkg= 0.00000E+00 tdspl= 0.00000E+00 fep= 0.00000E+00 uredox= tolbt= 0.00000E+00 toldl= 0.00000E+00 tolsat= 0.00000E+00 itermx= 0* 1 2 3 4 5 6 7 8 9 10 iopt1-10= 0 0 0 0 0 0 0 0 0 0 iopg1-10= 0 0 0 0 0 0 0 0 0 0 iopr1-10= 0 0 0 0 0 0 0 0 0 0 iopr11-20= 0 0 0 0 0 0 0 0 0 0 iodb1-10= 0 0 0 0 0 0 0 0 0 0 uebal= H+ nxmod= 0data file master species= Na+ switch with species= jflag= 0 csp= 1.00000E-03data file master species= UO2++ switch with species= jflag= 19 csp= 0. Mineral= UO2CO3data file master species= HCO3-

-25-

switch with species= jflag= 21 csp= -3.481 gas= CO2(g)data file master species= SO4-- switch with species= jflag= 0 csp= 1.00000E-10data file master species= Cl- switch with species= jflag= 0 csp= 1.00000E-03data file master species= H+ switch with species= jflag= 16 csp= -7.0000data file master species= O2(aq) switch with species= jflag= -1 csp= 0.data file master species= H2(aq) switch with species= jflag= -1 csp= 0.data file master species= HS- switch with species= jflag= -1 csp= 0.endit.

Tab. 7. Beispiel für PhreeqC-INPUT-File (Lösung von Rutherfordine als Funktion des CO2

Partialdrucks)

TITLE Beispiel Lösung Rutherfordine als Funktion des CO2 Partialdruckes

SOLUTION 1 Wasser mit 1 mmol/L Na und Clunits mmol/kgwtemp 25pH 7pe 14Na 1S(6) 1E-7Cl 1

EQUILIBRIUM_PHASES 1 CO2(g) -3.481Rutherfordine 0

END

Der input-Aufwand für PhreeqC2 hat sich weiter vereinfacht durch eine neue Windows-Oberfläche(Abb. 5. ), die es über einen key-word Index rechts im Programmfenster ermöglicht, im input-file dienotwendigen Programm-Befehle und Syntaxformulierungen per Mausklick schnell zu aktivieren undim output-file gezielt nach den gewünschten Ergebnissen zu suchen. Die gewählte Datenbank kann ineinem extra Fenster gelesen, nicht aber verändert werden. Zusätzlich gibt es die Fenster „GRID“(Daten im spreadsheet-Format) und „CHART“ zur graphischen Darstellung der Ergebnisse.

Abb. 5. Einfaches Beispiel einer chem. Analyse in PhreeqC2 for Windows key wordIndex

Key word SOLUTIONfür Eingabe der Chem.Analyse

möglicheSubbefehleunterSOLUTION

chem. Analyse

Berechnung starten

-26-

Im Blockkurs Hydrogeochemische Modellierung wird daher mit dem Programm PhreeqC2 forWindows gearbeitet. Da das Programm vom USGS (United States Geological Survey) entwickeltwurde, ist es zudem public domain und kann vom Internet heruntergeladen werden:

Homepage des USGS: http://wwwbrr.cr.usgs.gov/projects/GWC_coupled/phreeqc/index.html

• Programm: PhreeqC2.2 mit Windows-Oberfläche: “PHREEQC for Windows (1.0M)”(http://www.geo.vu.nl/users/posv/phreeqc.html)

• Manual: „ PDF Documentation (1.1M)“(ftp://brrcrftp.cr.usgs.gov/geochem/unix/phreeqc/manual.pdf)

Für weitere Informationen / Daten zum Blockkurs:homepage Prof. B. Merkel (http://www.geo.tu-freiberg.de/~merkel/vorlesung/)

7. Literatur

APPELLO C.A.J. & POSTEMA D. (1994): Geochemistry, groundwater and pollution. Rotterdam,Balkema

BALL J.W. & NORDSTROM D.K. (1991): User´s Manual for WATEQ4F. US GeologicalSurvey Open-File Report 91-183

BOHN H.L., McNEAK B.L. & O´CONNOR G.A. (1979): Soil Chemistry. Wiley-InterscienceBOHN et al. (1983)

DAVIS J. & KENT D.B. (1990): Surface complexation modeling in aqueous geochemistry. In:Hochella M.F. & White A.F (eds):Mineral-Water Interface Geochemistry; Washington D.C.,Mineralogical Society of America, Reviews in Mineralogy, Vol.23, 5

DZOMBAK D.A & MOREL F.M.M. (1990): Surface complexation modeling - Hydrous ferricoxide, New York, John Wiley

GAPON E.N. (1933) Theory of exchange adsorption. J.Gen. Chem. (USSR) 3.

GAINES G.L. & THOMAS H.C. (1953): Adsorption studies on clay minerals. A formulation of thethermodynamics of exchange adsorption. J. Chem. Phys. 21

GLYNN P.D. (1990) : Modeling solid solution reactions in low-temperature aqueous systems. In:MELCHIOR D.C. & BASSET R.L. (eds). Chemical modeling of aqueous systems II. ACS Sympl.Ser 416, 74-86

GLYNN P.D. & REARDON E.J. (1990): Solid-solution aqueous-solution equililibria:thermodynamic theory and representation. Am.J.Sci. 290, 164-201

-27-

HEM J.D. (1985): Study and interpretation of the chemical characteristics of natural waters.-U.S.Geol.Surv.Water-Supply Paper 2254, 3.Aufl.

LANGMUIR D. (1997): Aqueous environmental geochemistry.-Prentice Hall; New Jersey.

PARKHURST, D.L., THORSTENSON D.C. & PLUMMER L.N. (1982): Phreeqe – a computerprogram for geochemical calculations.- U.S.Geol.Surv. Water-Resour.Invest.80-96:210.

PARKHURST, D.L. (1995): User´s guide to PhreecQ – a computer program for speciation,reaction-path, advective-transport, and inverse geochemical calculations.- U.S.Geol.Surv. Water-Resour.Invest.Rept.95-4227.

PARKHURST, D.L., AND APPELO, C.A.J. (1999): User's guide to PHREEQC (Version 2)--acomputer program for speciation, batch-reaction, one-dimensional transport, and inversegeochemical calculations: U.S. Geological Survey Water-Resources Investigations Report 99-4259,312 p.

PERASON F.J., BERNER U. & HUMMEL W. (1992): Nagra thermochemical Data Base. II.Supplemental Data 05/92. Technical Report 91-18

SCHNITZER M. (1986): Binding of humic substances by soil mineral colloids. In: Interactions ofsoil Minerals With Natural Organics and Microbes. In: HUANG P.M. & SCHNITZER M. (Eds).Soil. Sci. Soc. Am. Publ. No. 17, Madison, WI

SPARKS D.L. (1986): Soil Physical Chemistry. CRC Press Inc., Boca Raton, FL

STUMM W. & MORGAN J.J. (1981): Aquatic Chemistry.- Wiley & Sons, New York.

VANSELOW A.P.(1932): Equilibria of the base-exchange reactions of bentonites, permutites, soilcolloids and zeolites. Soil Sci. 33