Simple example of kinetic dissolution of calcite with Phreeqc and Aquifer Lab - Tutorial

This is an example of calcite dissolution with kinetics where a simplified rate and initial molality are defined. The tutorial covers creation of an input file for Phreeqc, database selection, Phreeqc simulation in Aquifer App and analysis of processed output data on a Jupyter notebook. Finally the molality of Calcium and pH are plotted with time.

Tutorial

Phreeqc input file

TITLE Kinetic dissolution of calcite

SOLUTION 1
    temp 10
    pH   6
    C    1

EQUILIBRIUM_PHASES 1
    CO2(g)  -1.5

KINETICS 1
    Calcit2
    formula CaCO3
    -m0 1
    -steps 30000 in 20 steps

RATES
    Calcit2
    -start
    10 rate=10^-6.91 - 10^-1.52*(tot("Ca"))^2
    20 save 1 * rate * time
    -end

END

Code

import pandas as pd
import numpy as np
import matplotlib.pyplot as plt

#solution composition 
compDf = pd.read_csv('solutionComposition.csv', index_col=[0,1,2])

compDf.head()

			Element	Molality	Moles
Simulation	Type	Number
1	initial	1	C	0.001000	0.001000
	batch	1	C	0.002326	0.002326
		1	Ca	0.000184	0.000184
		2	C	0.002688	0.002688
		2	Ca	0.000365	0.000365

#caDf = compDf.loc[(compDf.index.get_level_values('Type') == 'batch') & (compDf.Element == 'Ca')]
caDf = compDf.loc[compDf.Element == 'Ca']
caDf.head()

			Element	Molality	Moles
Simulation	Type	Number
1	batch	1	Ca	0.000184	0.000184
		2	Ca	0.000365	0.000365
		3	Ca	0.000540	0.000540
		4	Ca	0.000707	0.000707
		5	Ca	0.000863	0.000863

#solution description for ph values
descDf = pd.read_csv('solutionDescription.csv', index_col=[0,1,2])
descDf.head()

			Parameter	Value
Simulation	Type	Number
1	initial	1	pH	6.00000
		1	pe	4.00000
		1	Specific Conductance (µS/cm, 10°C)	8.00000
		1	Density (g/cm³)	0.99972
		1	Volume (L)	1.00033

#phDf = descDf.loc[(descDf.index.get_level_values('Type') == 'batch') & (descDf.Parameter == 'pH')]
phDf = descDf.loc[descDf.Parameter == 'pH']
phDf.head()

			Parameter	Value
Simulation	Type	Number
1	initial	1	pH	6.000
	batch	1	pH	6.017
		2	pH	6.210
		3	pH	6.338
		4	pH	6.432

#define time list
timeList = np.linspace(0,30000,21)
timeList

array([    0.,  1500.,  3000.,  4500.,  6000.,  7500.,  9000., 10500.,
       12000., 13500., 15000., 16500., 18000., 19500., 21000., 22500.,
       24000., 25500., 27000., 28500., 30000.])

#plot Ca molality
plt.plot(timeList[1:],caDf.Molality*1000)
plt.ylabel('Ca (mMoles)')
plt.xlabel('Time Seconds')

Text(0.5, 0, 'Time Seconds')

#plot pH
plt.plot(timeList,phDf.Value)
plt.xlabel('pH')
plt.xlabel('Time (Seconds)')

Text(0.5, 0, 'Time (Seconds)')

Input data

You can download the input files from this link.