Conceptual Models > Kinetics and rate controlling factors
High pCO2 values usage in CO2 capture script
(1/1)
Erdinc:
Hello PHREEQC gurus,
I have been working on script to calculate the dissolution of minerals and CO2 capture through weathering of slag mineralization. I have run my simulations from 0.1% to up to 99% CO2. from 0.1 to around 6 % makes some sense but anything around 50% and 99% shows erratic readings which I believe could be PHREEQC sensitivity. Would you please let me know your thoughts? Many thanks.
Best Regards,
Erdinc
--- Code: ---TITLE DATABASE C:\Users\ErdincCosgun\Documents\Projects\Databases\carbfix_2022.dat
TITLE TESTING KINETICS OF SPECIFIC SOLIDS
TITLE 100 CELLS, 1CM HIGH
TITLE FLOW-THROUGH TEST, 3 LITRE
TITLE === FULL ANTIMAGNETIC SUITE ===
TITLE KARDEMIR PROCESS WATER MEAN VALUES
TITLE REDOX CONTROLLED BY IRON SPECIES
TITLE CHEMICAL SPECIATION
TITLE PHASES FROM Carbfix2022
TITLE CALCITE DISSOLVES AS KINETIC PHASE
TITLE CALCITE PRECIPITATION AS EQUILIBRIUM PHASE
TITLE CO2 INJECTED FROM DAC 50% TO SOLUTION, log(0.5):-0.301
TITLE Li ADDED IN TRACE TO KARDEMIR SUPPLIED WATER
PHASES # Minerals, Directly taken from Carbfix2022.dat PHASES.
Calcite_diss # Calcite that dissolves only. Alphanumeric name of phase; no spaces are allowed.
CaCO3 + H+ = Ca+2 + HCO3- # Dissolution reaction
log_k 1.8487 # Log K at 25 ?C for the reaction
-delta_H -25.7149 kJ/mol # Enthalpy of reaction at 25 ?C for the reaction
# deltafH -288.552 kcal/mol
-analytic -1.4978e2 -4.8370e-2 4.8974e3 6.0458e1 7.6464e1 # Identifier for coefficients for an analytical expression for the temperature dependence of log K.
# Range 0-350
-Vm 36.934 # molar_volume, the molecular weight divided by the density of the solid at 25 ?C. Ex: Molecular Weight/Mass (g/mole)=100.09, Density (g/ml or g/cm3)= 2.71 therfore, Vm= 100.09 / 2.71 = 36.934
# Extrapol supcrt92
# Ref HDN+78
Calcite_ppt # Calcite that percipitates only
CaCO3 + H+ = Ca+2 + HCO3-
log_k 1.8487
-delta_H -25.7149 kJ/mol
# deltafH -288.552 kcal/mol
-analytic -1.4978e2 -4.8370e-2 4.8974e3 6.0458e1 7.6464e1
# Range 0-350
-Vm 36.934 # molar_volume, the molecular weight divided by the density of the solid at 25 ?C. Ex: Molecular Weight/Mass (g/mole)=100.09, Density (g/ml or g/cm3)= 2.71 therfore, Vm= 100.09 / 2.71 = 36.934
# Extrapol supcrt92
# Ref HDN+78
Brucite
Mg(OH)2 + 2 H+ = Mg+2 + 2 H2O
log_k 16.2980
-delta_H -111.34 kJ/mol
# deltafH -221.39 kcal/mol
-analytic -1.0280e2 -1.9759e-2 9.0180e3 3.8282e1 1.4075e2
# Range 0-350
-Vm 24.63 # molar_volume, the molecular weight divided by the density of the solid at 25 ?C. Ex: Molecular Weight/Mass (g/mole)=58.32, Density (g/ml or g/cm3)= 2.39 therfore, Vm= 58.32 / 2.39 = 24.40
# Extrapol supcrt92
# Ref HDN+78
Larnite
Ca2SiO4 + 4 H+ = SiO2 + 2 Ca+2 + 2 H2O
log_k 38.4665
-delta_H -227.061 kJ/mol
# deltafH -551.74 kcal/mol
-analytic 2.6900e1 -2.1833e-3 1.0900e4 -9.5257 -7.2537e4
# Range 0-300
-Vm 51.6 # HDN+78 # molar_volume, the molecular weight divided by the density of the solid at 25 ?C. Ex: Molecular Weight/Mass (g/mole)=172.24, Density (g/ml or g/cm3)= 3.28 therfore, Vm= 172.24 / 3.28 = 52.512
# Extrapol Cp integration
# Ref 82sar/bar
Portlandite
Ca(OH)2 + 2 H+ = Ca+2 + 2 H2O
log_k 22.5552
-delta_H -128.686 kJ/mol
# deltafH -986.074 kJ/mol
-analytic -8.3848e1 -1.8373e-2 9.3154e3 3.2584e1 1.4538e2
# Range 0-300
-Vm 33.056 # thermo.com.V8.R6+.tdat # molar_volume, the molecular weight divided by the density of the solid at 25 ?C. Ex: Molecular Weight/Mass (g/mole)=74.09, Density (g/ml or g/cm3)= 2.23 therfore, Vm= 74.09 / 2.23 = 33.224
# Extrapol Cp integration
# Ref RHF79
# All other solids are defined in the chosen database
SOLUTION 1 Kardemir Process/grid water # Directly taken from KARDEMIR Steel Manufacturer Plant. Define the composition of an aqueous solution. Speciate an analysis, calculate saturation indices.
temp 25 # temperature in degrees Celcius, default = 25 C, Kardemir 25C
pH 8.25 # +-0.1,default pH = 7, Kardemir 8.5
pe 4 # pe = -log(electron activity), default = 4.
redox pe # reduction-oxidation reactions
units ppm # same as mg/l, default units ppm, mmol/kg water
density 1 # density, default = 1 kg/L, Kardemir 1
Ca 160 # +-0.1,total dissolved Ca, Kardemir 140-180
#Mg 3.4 # +-0.1,total dissolved Mg # Speciate that can be in grid water. No presence of Mg in Kardemir water
#N(5) 0.84 # 0.84% ppm or mg/l of Nitrogen # Speciate that can be in grid water. No presence of N(5) in Kardemir water
N(-3) 0.35 # Kardemir NH3 0.2-0.5
#Na 22 # 11-33 mg/l,total dissolved Na # Speciate that can be in grid water. No presence of Na in Kardemir water
#K 0.82 # total dissolved K # Speciate that can be in grid water. No presence of K in Kardemir water
P 0.22 # 0.3% ppm or mg/l of Phosforous, Kardemir P2O5 0.15-0.30
Fe(2) 0.25 # total dissolved Fe,generally <0.01, Kardemir 0.2-0.3, assumed all to be as ferrous
Fe(3) 1e-9 # ferric iron assumed as present in negligible quantities
Alkalinity 240 # 30-400ppm reasonable range alkalinity for domestic drinking water,low alkalinity 0.4- 0.8, Kardemir 210-270 CaCo3
C(4) 4.620e-03 mol/kgs # equivalent of Alkalinity 240 and pH 8.25
Cl 85 charge # total dissolved Cl, Kardemir 50-120. Charged balanced.
S(6) 70 # ppm or mg/l of Sulfate, total dissolved Sulfate, Kardemir 60-80
#Li 6.941E-03 # Lithium (1umol/L) is a tracer for Process Water
#Si 2.80855 # mg/L = 0.1 mmol/L
-water 1 # kg water, default = 1 kg
SAVE SOLUTION 1
END
SOLUTION 2 Median Black Sea rainfall concentrations from Alagha and Tunkel (2003). CO2(atm), O2(atm) added. Chloride charge balance
temp 25
pH 5.15
pe 4
redox pe
units ppm
density 1
Al 3.37
Ca 1.8
Cl 1.01 charge # adapt Cl to obtain charge balance
C(4) 1 CO2(g) -3.373 # equilibrium with atmospheric CO2
Fe 2.82
K 1.72
Mg 0.88
Mn 0.07
N(-3) 0.56
N(5) 1.5
Na 2
S(6) 2.52
O(0) 0.01 O2(g) -0.678 # equilibrium with atmospheric O2
-water 1 # kg
END
EQUILIBRIUM_PHASES 1-100
# Equilibrate SOLUTION 1 water with calcite precipitation
# Brucite 0 4.733 dissolve_only # 4.733 (mol_mineral/L_pores)calculated. PHREEQC output assigns this value (4.733 )to each cell rather than the total depth of the core therfore dividing the moles of mineral by the number of cells is required, 4.733/100=0.04733 cells will give mols /cell.See supporting document
# Calcite_diss 0 9.49 dissolve_only # Slag calcite dissolves, 9.49 (mol_mineral/L_pores)calculated. PHREEQC output assigns this value (9.49 )to each cell rather than the total depth of the core therfore dividing the moles of mineral by the number of cells is required, 9.49/100=0.0949 cells will give mols /cell.See supporting document
Calcite_ppt 0 0 # Switch on calcite as a phase that can precipitate within days
# Larnite 0 9.87 dissolve_only # 9.87 (mol_mineral/L_pores)calculated. PHREEQC output assigns this value (9.87 )to each cell rather than the total depth of the core therfore dividing the moles of mineral by the number of cells is required, 9.87/100=0.0987 cells will give mols /cell.See supporting document
# Portlandite 0 14.43 dissolve_only # 14.43 (mol_mineral/L_pores)calculated. PHREEQC output assigns this value (14.43 )to each cell rather than the total depth of the core therfore dividing the moles of mineral by the number of cells is required, 14.43/100=0.1443 cells will give mols /cell.See supporting document
END
MIX 1 Combine rainwater with (process/grid water + DAC)
1 0.88 # 1-(Rainfall_rate/(Rainfall_rate+Procwtr_rate), use fraction 0.88 of solution 1
2 0.12 # 1-(Procwtr_rate/(Rainfall_rate+Procwtr_rate), use fraction 0.12 of solution 2
USE equilibrium_phases 1
SAVE SOLUTION 0
END
SOLUTION 1-101 Simulated mixture of (process water, DAC) and rainwater.
temp 25
pH 6.963
pe 13.64
redox pe
units mol/kgw
density 1
-water 1 # kg
C(4) 5.06E-03
Ca 3.52E-03
Cl 2.12E-03
Fe 3.94E-06
K 1.07E-06
Mg 1.93E-06
N 2.79E-05
Na 1.07E-05
P 6.26E-06
S(6) 6.45E-04
END
REACTION 0
# Remove / Evaporate water in steps. All solutions defined by SOLUTION input are scaled to have exactly 1 kg (approximately 55.5 mol) of water, unless -water identifier is used.
H2O -1
# Evaporate 11.1012 moles(= 0.2 * 55.506), 20% evaporation from pond in 4 step. # 11.1012 mole in 4 step. # The first step removes(-) 2.775e+00 moles of reaction (units are ?moles?) to the initial solution; the second step removes 5.551e+00 moles of reaction to the initial solution; the third 8.326e+00 moles; and the fourth 1.110e+01 moles; each reaction step begins with the same initial solution and removes only the amount of reaction specified.
11.1012 moles in 4 steps
EQUILIBRIUM_PHASES 101 DAC CO2 capture to 50% # BOTTOM CELL IS THE EFFLUENT FROM THE COLUMN, AS IT EXPERIENCES CONTACT WITH DAC-AMENDED AIR
CO2(g) -0.301 10 # Calculated CO2 SI= -0.301, Amount of 10 mol available
O2(g) -0.679 10 # Calculated O2 SI= -0.679, Amount of 10 mol available
# N2(g) -0.107 10
Calcite_ppt 0 0
SELECTED_OUTPUT 0
-file 230828_1125_Antimag_1cm_101cell_DAC_FlowThru_12pw-flush_dt-dx_Cal-KinEqu_CV5_Si_Init-solns20March2025V0_50pc500Days.tsv
-ionic_strength true
-distance true
-totals Ca C(4) Si
-molalities Ca+2 CaCO3 CO3-2 HCO3- CO2
-saturation_indices Calcite Larnite Brucite Portlandite CO2(g)
-equilibrium_phases Calcite Calcite_ppt CO2(g) O2(g)
-kinetic_reactants Calcite_diss Larnite Brucite Portlandite
END
USER_PRINT
10 p = SYS("equi", count , name$ , type$ , moles )
20 FOR i = 1 to count
30 dv = PHASE_VM(name$(i))*EQUI_DELTA(name$(i))
40 delta_v = delta_v + dv
50 PRINT name$(i), dv, " cm^3"
60 NEXT i
70 PRINT "Total change in mineral volume: ", delta_v, " cm^3"
80 END
RATES
Calcite_diss # calcite dissolving only
-start # Identifier marks the beginning of a Basic program by which the moles of reaction for a time subinterval are calculated.
5 REM Formula = CaCO3
10 REM Rate law modified from Zhang et al. 2019
20 REM Rate law from Marty et al 2015
30 REM unit should be J, m2, mol,kgw-1 and second-1
40 REM parm(1) is surface area in the unit of m2/kgw
50 REM calculation of surface area can be found in the note
60 REM M is current moles of minerals. M0 is the initial moles of minerals
70 REM parm(2) is a correction factor
80 REM Affinity ON
90 sr_mineral = SR("Calcite") # Calculate calcite saturation ratio
95 si_mineral = SI("Calcite") # Calculate calcite saturation index
100 moles = 0
110 REM rate=0 if mineral supersaturated in solution
120 IF(M <= 0 OR sr_mineral > 1) THEN GOTO 320
130 a0 = PARM(1)
132 REM Volume is 1 litre solution
135 v = 1
140 temp = 1/TK
145 REM acid solution parameters
150 k_acid = 0
160 eapp_acid = 0
170 n_acid = 0
175 REM neutral solution parameters
180 k_neut = 6.59E+04
190 eapp_neut = 66000
195 REM base solution parameters
200 k_base = 1.04E+09
210 eapp_base = 67000
220 n_base = 1.6
225 REM Rate Law in TST form
230 hplus = ACT("H+")
240 r = 8.314
250 r_acid = k_acid*EXP((-eapp_acid/r)*temp)*(hplus^n_acid)
260 r_neut = k_neut*EXP((-eapp_neut/r)*temp)
270 REM Base in terms of bicarbonate
280 r_base = k_base*EXP((-eapp_base/r)*temp)*(ACT("HCO3-")^n_base)
290 r_all = r_acid+r_neut+r_base
300 rate = (a0/v)*(M/M0)^0.67*r_all*(1-sr_mineral)*PARM(2) # Calculate overall rate, factor
310 moles = rate*TIME # Calculate moles of reaction over time interval given by TIME. Note that the multiplication of the rate by TIME must be present in one of the Basic lines.
320 SAVE moles # Return moles of reaction for time subinterval with ?SAVE?. A SAVE statement must always be present in a rate program.
-end # Identifier marks the end of a Basic program by which the number of moles of a reaction for a time subinterval is calculated.
Brucite
-start
5 REM Formula = Mg(OH)2
10 REM Rate parameters from Zhang et al. 2019 database
20 REM Rate law from Palandri and Kharaka (2004)
25 REM experimental condition range T=25-75C, pH=1-5
30 REM unit should be J, m2, mol,kgw-1 and second-1
40 REM parm(1) is surface area in the unit of m2/kgw
50 REM calculation of surface area can be found in the note
60 REM M is current moles of minerals. M0 is the initial moles of minerals
70 REM parm(2) is a correction factor
80 REM Affinity ON
90 sr_mineral = SR("Brucite")
95 si_mineral = SI("Brucite")
100 moles = 0
110 REM rate=0 if mineral supersaturated in solution
120 IF(M <= 0 OR sr_mineral > 1) THEN GOTO 380
130 a0 = PARM(1)
140 REM Volume is 1 litre solution
150 v = 1
160 temp = 1/TK
170 REM acid solution parameters
180 k_acid = 4.00E+05
190 eapp_acid = 59000
200 n_acid = 0.500
210 REM neutral solution parameters
220 k_neut = 1.30E-01
230 eapp_neut = 42000
240 REM base solution parameters
250 k_base = 0
260 eapp_base = 0
270 n_base = 0
280 REM Rate Law in TST form
290 hplus = ACT("H+")
300 r = 8.314
310 r_acid = k_acid*EXP((-eapp_acid/r)*temp)*(hplus^n_acid)
320 r_neut = k_neut*EXP((-eapp_neut/r)*temp)
340 r_base = k_base*EXP((-eapp_base/r)*temp)*(hplus^n_base)
350 r_all = r_acid+r_neut+r_base
360 rate = (a0/v)*(M/M0)^0.67*r_all*(1-sr_mineral)*PARM(2)
370 moles = rate*TIME
380 SAVE moles
-end
Larnite
-start
5 REM Formula = Ca2SiO4
10 REM Rate parameters from Carbfix_2022 database
20 REM Rate law from Palandri and Kharaka (2004)
30 REM unit should be J, m2, mol,kgw-1 and second-1
40 REM parm(1) is surface area in the unit of m2/kgw
50 REM calculation of surface area can be found in the note
60 REM M is current moles of minerals. M0 is the initial moles of minerals
70 REM parm(2) is a correction factor
80 REM Affinity ON
90 sr_mineral = SR("Larnite")
95 si_mineral = SI("Larnite")
100 moles = 0
110 REM rate=0 if mineral supersaturated in solution
120 IF(M <= 0 OR sr_mineral > 1) THEN GOTO 380
130 a0 = PARM(1)
140 REM Volume is 1 litre solution
150 v = 1
160 temp = 1/TK
170 REM acid solution parameters
180 k_acid = 5.25e8
190 eapp_acid = 70400
200 n_acid = 0.44
210 REM neutral solution parameters
220 k_neut = 0
230 eapp_neut = 0
240 REM base solution parameters
250 k_base = 8.25e5
260 eapp_base = 60900
270 n_base = 0.22
280 REM Rate Law in TST form
290 hplus = ACT("H+")
300 r = 8.314
310 r_acid = k_acid*EXP((-eapp_acid/r)*temp)*(hplus^n_acid)
320 r_neut = k_neut*EXP((-eapp_neut/r)*temp)
330 IF(k_base = 0) THEN GOTO 350
340 r_base = k_base*EXP((-eapp_base/r)*temp)*(hplus^n_base)
350 r_all = r_acid+r_neut+r_base
360 rate = (a0/v)*(M/M0)^0.67*r_all*(1-sr_mineral)*PARM(2)
370 moles = rate*TIME
375 REM PRINT "SI(Larnite) = " si_mineral ", Molar adjustment = " moles
380 SAVE moles
-end
Portlandite
-start
5 REM Formula = Ca(OH)2
10 REM Rate law modified from Zhang et al. 2019
20 REM Rate law from marty et al 2015
25 REM experimental condition range T=25-80C, pH=5-7
30 REM units should be J, m2, mol,kgw-1 and second-1
40 REM parm(1) is surface area in the unit of m2/kgw
50 REM calculation of surface area can be found in the note
60 REM M is current moles of minerals. M0 is the initial moles of minerals
70 REM parm(2) is a correction factor
80 REM Affinity ON
90 sr_mineral = SR("Portlandite")
95 si_mineral = SI("Portlandite")
100 moles = 0
110 REM rate=0 if mineral supersaturated in solution
120 IF(M <= 0 OR sr_mineral > 1) THEN GOTO 320
130 a0 = PARM(1)
132 REM Volume is 1 litre solution
135 v = 1
140 temp = 1/TK
145 REM acid solution parameters
150 k_acid = 1.10E+10
160 eapp_acid = 75000
170 n_acid = 0.600
175 REM neutral solution parameters
180 k_neut = 3.04E+05
190 eapp_neut = 75000
195 REM base solution parameters
200 k_base = 0
210 eapp_base = 0
220 n_base = 0
225 REM Rate Law in TST form
230 hplus = ACT("H+")
240 r = 8.314
250 r_acid = k_acid*EXP((-eapp_acid/r)*temp)*(hplus^n_acid)
260 r_neut = k_neut*EXP((-eapp_neut/r)*temp)
270 REM Base in terms of bicarbonate
280 r_base = k_base*EXP((-eapp_base/r)*temp)*(ACT("OH-")^n_base)
290 r_all = r_acid+r_neut+r_base
300 rate = (a0/v)*(M/M0)^0.67*r_all*(1-sr_mineral)*PARM(2)
310 moles = rate*TIME
320 SAVE moles
-end
kin_HCO3
-start
10 REM
20 kin_ca=1*KIN_DELTA("Portlandite")+2*KIN_DELTA("Larnite")+0*KIN_DELTA("Brucite")+(0.5)*KIN_DELTA("Calcite_diss")
30 kin_mg=0*KIN_DELTA("Portlandite")+0*KIN_DELTA("Larnite")+1*KIN_DELTA("Brucite")+0*KIN_DELTA("Calcite_diss")
60 moles = -2*(kin_ca+kin_mg)
70 print kin_ca, kin_mg, moles
320 SAVE moles
-end
KINETICS 1-100
Portlandite
-formula Ca(OH)2 1
-m 0.1443 # 14.43/100 cells = 0.1443 cells will give mols /cell. Initial moles of reactant. Mass conc port vs pores/ Molecular-mass_port
-m0 0.1443 # Default is equal to initial moles
-parms 2.303 1 # Bulk area * Area% Port
-tol 1e-08 # Default is 1?10-8 mol
Larnite
-formula Ca2SiO4 1
-m 0.0987 # 9.87/100 cells = 0.0987 cells will give mols /cell
-m0 0.0987
-parms 2.499 1
-tol 1e-08
Calcite_diss
-formula CaCO3 1
-m 0.0949 # 9.49/100 cells = 0.0949 cells will give mols /cell
-m0 0.0949
-parms 1.682 1
-tol 1e-08
Brucite
-formula Mg(OH)2 1
-m 0.04733 # 4.733/100 cells = 0.04733 cells will give mols /cell
-m0 0.04733
-parms 0.556 1
-tol 1e-08
-steps 1
-step_divide 0.1
-runge_kutta 6
-bad_step_max 500
-cvode true
-cvode_steps 200
-cvode_order 5
END
USE SOLUTION 0
USE REACTION 0
USE equilibrium_phases 101
USE KINETICS 1
SAVE equilibrium_phases 1-10
SAVE solution 1-10
TRANSPORT
-cells 101 # 100 cells slag heap, last cell closed pond (1 litre)
# -shifts 303 # 3 passes of pore water through the slag heap
-shifts 1245 # 12.326 passes of pore water through the slag heap
# -shifts 101 # 1 passes of pore water through the slag heap
-lengths 0.01 # 1cm cells
-time_step 34722.2 # seconds retention time in cells - Antimagnetic = 0.401 days. Shift (101) x Timestep(34722.2) = 3506942.2 sec = 3506942.2 /60/ 60/ 24 = The total time of the simulation (40.5896 days).Shift (303) x Timestep(34722.2) = 10520826.6 sec = 10520826.6/60/ 60/ 24 = The total time of the simulation (121.768 days).
# Shift (1245) x Timestep(34722.2) = 43200000sec = 43200000 /60/ 60/ 24 = The total time of the simulation (500 days)
-print_cells 1, 100-101
# -print_cells 10-11
-print_frequency 1245 # Print only at end
# -print_frequency 33 # Print only at end
USER_PUNCH 1 # This keyword data block is used to define Basic programs that print user-defined quantities to the selected-output file. Any Basic ?PUNCH? statement will write to the selected-output file.
-headings kin_hco3 Pore_vol # Cum-C-heading
-start
7 REM IF(CELL_NO < 101) THEN GOTO 100
10 IF (TOTAL_TIME = 0) THEN PUT(0,1)
20 kin_ca=1*KIN_DELTA("Portlandite")+2*KIN_DELTA("Larnite")+0*KIN_DELTA("Brucite")+(0.5)*KIN_DELTA("Calcite_diss")
30 kin_mg=0*KIN_DELTA("Portlandite")+0*KIN_DELTA("Larnite")+1*KIN_DELTA("Brucite")+0*KIN_DELTA("Calcite_diss")
60 kin_hco3 = -2*(kin_ca+kin_mg)
65 PUT(GET(1)+kin_hco3,1)
67 cum_c = GET(1)
68 PUNCH kin_hco3
#70 PUNCH cum_c
75 REM punch pore volumes, no_cells = 100
77 no_cells = 100
80 PUNCH (STEP_NO + .5) / no_cells
100 REM END
-end
SELECTED_OUTPUT 1
-file 230828_1125_Antimag_1cm_101cell_DAC_FlowThru_12pw-flush_dt-dx_Cal-KinEqu_CV5_Si_adv-react20March2025V0_50pc500Days.tsv
-ionic_strength true
-distance true
-totals Li Ca C(4) Si # Lithium is a tracer for Process Water
-molalities Ca+2 CaCO3 CO3-2 HCO3- CO2
-saturation_indices Calcite Larnite Brucite Portlandite CO2(g)
-equilibrium_phases Calcite Calcite_ppt CO2(g) O2(g)
-kinetic_reactants Calcite_diss Larnite Brucite Portlandite kin_HCO3
END
--- End code ---
dlparkhurst:
I don't understand what you are simulating. I have simplified your script and only simulated 20 cells.
You have evaporated 20% of the water in the first 10 cells. However, that evaporated water advects through the column and will be gone after 20 shifts (for my simulation), after which, there is no effect of evaporated water. Calcite is always in equilibrium.
You will have explain what you mean by "erratic readings". If possible, make a batch reaction where a change in CO2 concentration demonstrates the erratic readings. Otherwise, please use a column with fewer cells that runs quickly and explain what results are unexpected.
The following shows the "bubble" of evaporated water moving down the column after 5 shifts and the equilibrium with calcite.
--- Code: ---TITLE DATABASE C:\Users\ErdincCosgun\Documents\Projects\Databases\carbfix_2022.dat
TITLE TESTING KINETICS OF SPECIFIC SOLIDS
TITLE 100 CELLS, 1CM HIGH
TITLE FLOW-THROUGH TEST, 3 LITRE
TITLE === FULL ANTIMAGNETIC SUITE ===
TITLE KARDEMIR PROCESS WATER MEAN VALUES
TITLE REDOX CONTROLLED BY IRON SPECIES
TITLE CHEMICAL SPECIATION
TITLE PHASES FROM Carbfix2022
TITLE CALCITE DISSOLVES AS KINETIC PHASE
TITLE CALCITE PRECIPITATION AS EQUILIBRIUM PHASE
TITLE CO2 INJECTED FROM DAC 50% TO SOLUTION, log(0.5):-0.301
TITLE Li ADDED IN TRACE TO KARDEMIR SUPPLIED WATER
RATES
Calcite_diss # calcite dissolving only
-start # Identifier marks the beginning of a Basic program by which the moles of reaction for a time subinterval are calculated.
5 REM Formula = CaCO3
10 REM Rate law modified from Zhang et al. 2019
20 REM Rate law from Marty et al 2015
30 REM unit should be J, m2, mol,kgw-1 and second-1
40 REM parm(1) is surface area in the unit of m2/kgw
50 REM calculation of surface area can be found in the note
60 REM M is current moles of minerals. M0 is the initial moles of minerals
70 REM parm(2) is a correction factor
80 REM Affinity ON
90 sr_mineral = SR("Calcite") # Calculate calcite saturation ratio
95 si_mineral = SI("Calcite") # Calculate calcite saturation index
100 moles = 0
110 REM rate=0 if mineral supersaturated in solution
120 IF(M <= 0 OR sr_mineral > 1) THEN GOTO 320
130 a0 = PARM(1)
132 REM Volume is 1 litre solution
135 v = 1
140 temp = 1/TK
145 REM acid solution parameters
150 k_acid = 0
160 eapp_acid = 0
170 n_acid = 0
175 REM neutral solution parameters
180 k_neut = 6.59E+04
190 eapp_neut = 66000
195 REM base solution parameters
200 k_base = 1.04E+09
210 eapp_base = 67000
220 n_base = 1.6
225 REM Rate Law in TST form
230 hplus = ACT("H+")
240 r = 8.314
250 r_acid = k_acid*EXP((-eapp_acid/r)*temp)*(hplus^n_acid)
260 r_neut = k_neut*EXP((-eapp_neut/r)*temp)
270 REM Base in terms of bicarbonate
280 r_base = k_base*EXP((-eapp_base/r)*temp)*(ACT("HCO3-")^n_base)
290 r_all = r_acid+r_neut+r_base
300 rate = (a0/v)*(M/M0)^0.67*r_all*(1-sr_mineral)*PARM(2) # Calculate overall rate, factor
310 moles = rate*TIME # Calculate moles of reaction over time interval given by TIME. Note that the multiplication of the rate by TIME must be present in one of the Basic lines.
320 SAVE moles # Return moles of reaction for time subinterval with ?SAVE?. A SAVE statement must always be present in a rate program.
-end # Identifier marks the end of a Basic program by which the number of moles of a reaction for a time subinterval is calculated.
Brucite
-start
5 REM Formula = Mg(OH)2
10 REM Rate parameters from Zhang et al. 2019 database
20 REM Rate law from Palandri and Kharaka (2004)
25 REM experimental condition range T=25-75C, pH=1-5
30 REM unit should be J, m2, mol,kgw-1 and second-1
40 REM parm(1) is surface area in the unit of m2/kgw
50 REM calculation of surface area can be found in the note
60 REM M is current moles of minerals. M0 is the initial moles of minerals
70 REM parm(2) is a correction factor
80 REM Affinity ON
90 sr_mineral = SR("Brucite")
95 si_mineral = SI("Brucite")
100 moles = 0
110 REM rate=0 if mineral supersaturated in solution
120 IF(M <= 0 OR sr_mineral > 1) THEN GOTO 380
130 a0 = PARM(1)
140 REM Volume is 1 litre solution
150 v = 1
160 temp = 1/TK
170 REM acid solution parameters
180 k_acid = 4.00E+05
190 eapp_acid = 59000
200 n_acid = 0.500
210 REM neutral solution parameters
220 k_neut = 1.30E-01
230 eapp_neut = 42000
240 REM base solution parameters
250 k_base = 0
260 eapp_base = 0
270 n_base = 0
280 REM Rate Law in TST form
290 hplus = ACT("H+")
300 r = 8.314
310 r_acid = k_acid*EXP((-eapp_acid/r)*temp)*(hplus^n_acid)
320 r_neut = k_neut*EXP((-eapp_neut/r)*temp)
340 r_base = k_base*EXP((-eapp_base/r)*temp)*(hplus^n_base)
350 r_all = r_acid+r_neut+r_base
360 rate = (a0/v)*(M/M0)^0.67*r_all*(1-sr_mineral)*PARM(2)
370 moles = rate*TIME
380 SAVE moles
-end
Larnite
-start
5 REM Formula = Ca2SiO4
10 REM Rate parameters from Carbfix_2022 database
20 REM Rate law from Palandri and Kharaka (2004)
30 REM unit should be J, m2, mol,kgw-1 and second-1
40 REM parm(1) is surface area in the unit of m2/kgw
50 REM calculation of surface area can be found in the note
60 REM M is current moles of minerals. M0 is the initial moles of minerals
70 REM parm(2) is a correction factor
80 REM Affinity ON
90 sr_mineral = SR("Larnite")
95 si_mineral = SI("Larnite")
100 moles = 0
110 REM rate=0 if mineral supersaturated in solution
120 IF(M <= 0 OR sr_mineral > 1) THEN GOTO 380
130 a0 = PARM(1)
140 REM Volume is 1 litre solution
150 v = 1
160 temp = 1/TK
170 REM acid solution parameters
180 k_acid = 5.25e8
190 eapp_acid = 70400
200 n_acid = 0.44
210 REM neutral solution parameters
220 k_neut = 0
230 eapp_neut = 0
240 REM base solution parameters
250 k_base = 8.25e5
260 eapp_base = 60900
270 n_base = 0.22
280 REM Rate Law in TST form
290 hplus = ACT("H+")
300 r = 8.314
310 r_acid = k_acid*EXP((-eapp_acid/r)*temp)*(hplus^n_acid)
320 r_neut = k_neut*EXP((-eapp_neut/r)*temp)
330 IF(k_base = 0) THEN GOTO 350
340 r_base = k_base*EXP((-eapp_base/r)*temp)*(hplus^n_base)
350 r_all = r_acid+r_neut+r_base
360 rate = (a0/v)*(M/M0)^0.67*r_all*(1-sr_mineral)*PARM(2)
370 moles = rate*TIME
375 REM PRINT "SI(Larnite) = " si_mineral ", Molar adjustment = " moles
380 SAVE moles
-end
Portlandite
-start
5 REM Formula = Ca(OH)2
10 REM Rate law modified from Zhang et al. 2019
20 REM Rate law from marty et al 2015
25 REM experimental condition range T=25-80C, pH=5-7
30 REM units should be J, m2, mol,kgw-1 and second-1
40 REM parm(1) is surface area in the unit of m2/kgw
50 REM calculation of surface area can be found in the note
60 REM M is current moles of minerals. M0 is the initial moles of minerals
70 REM parm(2) is a correction factor
80 REM Affinity ON
90 sr_mineral = SR("Portlandite")
95 si_mineral = SI("Portlandite")
100 moles = 0
110 REM rate=0 if mineral supersaturated in solution
120 IF(M <= 0 OR sr_mineral > 1) THEN GOTO 320
130 a0 = PARM(1)
132 REM Volume is 1 litre solution
135 v = 1
140 temp = 1/TK
145 REM acid solution parameters
150 k_acid = 1.10E+10
160 eapp_acid = 75000
170 n_acid = 0.600
175 REM neutral solution parameters
180 k_neut = 3.04E+05
190 eapp_neut = 75000
195 REM base solution parameters
200 k_base = 0
210 eapp_base = 0
220 n_base = 0
225 REM Rate Law in TST form
230 hplus = ACT("H+")
240 r = 8.314
250 r_acid = k_acid*EXP((-eapp_acid/r)*temp)*(hplus^n_acid)
260 r_neut = k_neut*EXP((-eapp_neut/r)*temp)
270 REM Base in terms of bicarbonate
280 r_base = k_base*EXP((-eapp_base/r)*temp)*(ACT("OH-")^n_base)
290 r_all = r_acid+r_neut+r_base
300 rate = (a0/v)*(M/M0)^0.67*r_all*(1-sr_mineral)*PARM(2)
310 moles = rate*TIME
320 SAVE moles
-end
kin_HCO3
-start
10 REM
20 kin_ca=1*KIN_DELTA("Portlandite")+2*KIN_DELTA("Larnite")+0*KIN_DELTA("Brucite")+(0.5)*KIN_DELTA("Calcite_diss")
30 kin_mg=0*KIN_DELTA("Portlandite")+0*KIN_DELTA("Larnite")+1*KIN_DELTA("Brucite")+0*KIN_DELTA("Calcite_diss")
60 moles = -2*(kin_ca+kin_mg)
70 print kin_ca, kin_mg, moles
320 SAVE moles
-end
PHASES # Minerals, Directly taken from Carbfix2022.dat PHASES.
Calcite_diss # Calcite that dissolves only. Alphanumeric name of phase; no spaces are allowed.
CaCO3 + H+ = Ca+2 + HCO3- # Dissolution reaction
log_k 1.8487 # Log K at 25 ?C for the reaction
-delta_H -25.7149 kJ/mol # Enthalpy of reaction at 25 ?C for the reaction
# deltafH -288.552 kcal/mol
-analytic -1.4978e2 -4.8370e-2 4.8974e3 6.0458e1 7.6464e1 # Identifier for coefficients for an analytical expression for the temperature dependence of log K.
# Range 0-350
-Vm 36.934 # molar_volume, the molecular weight divided by the density of the solid at 25 ?C. Ex: Molecular Weight/Mass (g/mole)=100.09, Density (g/ml or g/cm3)= 2.71 therfore, Vm= 100.09 / 2.71 = 36.934
# Extrapol supcrt92
# Ref HDN+78
Calcite_ppt # Calcite that percipitates only
CaCO3 + H+ = Ca+2 + HCO3-
log_k 1.8487
-delta_H -25.7149 kJ/mol
# deltafH -288.552 kcal/mol
-analytic -1.4978e2 -4.8370e-2 4.8974e3 6.0458e1 7.6464e1
# Range 0-350
-Vm 36.934 # molar_volume, the molecular weight divided by the density of the solid at 25 ?C. Ex: Molecular Weight/Mass (g/mole)=100.09, Density (g/ml or g/cm3)= 2.71 therfore, Vm= 100.09 / 2.71 = 36.934
# Extrapol supcrt92
# Ref HDN+78
Brucite
Mg(OH)2 + 2 H+ = Mg+2 + 2 H2O
log_k 16.2980
-delta_H -111.34 kJ/mol
# deltafH -221.39 kcal/mol
-analytic -1.0280e2 -1.9759e-2 9.0180e3 3.8282e1 1.4075e2
# Range 0-350
-Vm 24.63 # molar_volume, the molecular weight divided by the density of the solid at 25 ?C. Ex: Molecular Weight/Mass (g/mole)=58.32, Density (g/ml or g/cm3)= 2.39 therfore, Vm= 58.32 / 2.39 = 24.40
# Extrapol supcrt92
# Ref HDN+78
Larnite
Ca2SiO4 + 4 H+ = SiO2 + 2 Ca+2 + 2 H2O
log_k 38.4665
-delta_H -227.061 kJ/mol
# deltafH -551.74 kcal/mol
-analytic 2.6900e1 -2.1833e-3 1.0900e4 -9.5257 -7.2537e4
# Range 0-300
-Vm 51.6 # HDN+78 # molar_volume, the molecular weight divided by the density of the solid at 25 ?C. Ex: Molecular Weight/Mass (g/mole)=172.24, Density (g/ml or g/cm3)= 3.28 therfore, Vm= 172.24 / 3.28 = 52.512
# Extrapol Cp integration
# Ref 82sar/bar
Portlandite
Ca(OH)2 + 2 H+ = Ca+2 + 2 H2O
log_k 22.5552
-delta_H -128.686 kJ/mol
# deltafH -986.074 kJ/mol
-analytic -8.3848e1 -1.8373e-2 9.3154e3 3.2584e1 1.4538e2
# Range 0-300
-Vm 33.056 # thermo.com.V8.R6+.tdat # molar_volume, the molecular weight divided by the density of the solid at 25 ?C. Ex: Molecular Weight/Mass (g/mole)=74.09, Density (g/ml or g/cm3)= 2.23 therfore, Vm= 74.09 / 2.23 = 33.224
# Extrapol Cp integration
# Ref RHF79
# All other solids are defined in the chosen database
SOLUTION 1 Kardemir Process/grid water # Directly taken from KARDEMIR Steel Manufacturer Plant. Define the composition of an aqueous solution. Speciate an analysis, calculate saturation indices.
temp 25 # temperature in degrees Celcius, default = 25 C, Kardemir 25C
pH 8.25 # +-0.1,default pH = 7, Kardemir 8.5
pe 4 # pe = -log(electron activity), default = 4.
redox pe # reduction-oxidation reactions
units ppm # same as mg/l, default units ppm, mmol/kg water
density 1 # density, default = 1 kg/L, Kardemir 1
Ca 160 # +-0.1,total dissolved Ca, Kardemir 140-180
#Mg 3.4 # +-0.1,total dissolved Mg # Speciate that can be in grid water. No presence of Mg in Kardemir water
#N(5) 0.84 # 0.84% ppm or mg/l of Nitrogen # Speciate that can be in grid water. No presence of N(5) in Kardemir water
N(-3) 0.35 # Kardemir NH3 0.2-0.5
#Na 22 # 11-33 mg/l,total dissolved Na # Speciate that can be in grid water. No presence of Na in Kardemir water
#K 0.82 # total dissolved K # Speciate that can be in grid water. No presence of K in Kardemir water
P 0.22 # 0.3% ppm or mg/l of Phosforous, Kardemir P2O5 0.15-0.30
Fe(2) 0.25 # total dissolved Fe,generally <0.01, Kardemir 0.2-0.3, assumed all to be as ferrous
Fe(3) 1e-9 # ferric iron assumed as present in negligible quantities
Alkalinity 240 # 30-400ppm reasonable range alkalinity for domestic drinking water,low alkalinity 0.4- 0.8, Kardemir 210-270 CaCo3
C(4) 4.620e-03 mol/kgs # equivalent of Alkalinity 240 and pH 8.25
Cl 85 charge # total dissolved Cl, Kardemir 50-120. Charged balanced.
S(6) 70 # ppm or mg/l of Sulfate, total dissolved Sulfate, Kardemir 60-80
#Li 6.941E-03 # Lithium (1umol/L) is a tracer for Process Water
#Si 2.80855 # mg/L = 0.1 mmol/L
-water 1 # kg water, default = 1 kg
END
SOLUTION 2 Median Black Sea rainfall concentrations from Alagha and Tunkel (2003). CO2(atm), O2(atm) added. Chloride charge balance
temp 25
pH 5.15
pe 4
redox pe
units ppm
density 1
Al 3.37
Ca 1.8
Cl 1.01 charge # adapt Cl to obtain charge balance
C(4) 1 CO2(g) -3.373 # equilibrium with atmospheric CO2
Fe 2.82
K 1.72
Mg 0.88
Mn 0.07
N(-3) 0.56
N(5) 1.5
Na 2
S(6) 2.52
O(0) 0.01 O2(g) -0.678 # equilibrium with atmospheric O2
-water 1 # kg
END
EQUILIBRIUM_PHASES 1-100
# Equilibrate SOLUTION 1 water with calcite precipitation
# Brucite 0 4.733 dissolve_only # 4.733 (mol_mineral/L_pores)calculated. PHREEQC output assigns this value (4.733 )to each cell rather than the total depth of the core therfore dividing the moles of mineral by the number of cells is required, 4.733/100=0.04733 cells will give mols /cell.See supporting document
# Calcite_diss 0 9.49 dissolve_only # Slag calcite dissolves, 9.49 (mol_mineral/L_pores)calculated. PHREEQC output assigns this value (9.49 )to each cell rather than the total depth of the core therfore dividing the moles of mineral by the number of cells is required, 9.49/100=0.0949 cells will give mols /cell.See supporting document
Calcite_ppt 0 0 # Switch on calcite as a phase that can precipitate within days
# Larnite 0 9.87 dissolve_only # 9.87 (mol_mineral/L_pores)calculated. PHREEQC output assigns this value (9.87 )to each cell rather than the total depth of the core therfore dividing the moles of mineral by the number of cells is required, 9.87/100=0.0987 cells will give mols /cell.See supporting document
# Portlandite 0 14.43 dissolve_only # 14.43 (mol_mineral/L_pores)calculated. PHREEQC output assigns this value (14.43 )to each cell rather than the total depth of the core therfore dividing the moles of mineral by the number of cells is required, 14.43/100=0.1443 cells will give mols /cell.See supporting document
END
MIX 1 Combine rainwater with (process/grid water + DAC)
1 0.88 # 1-(Rainfall_rate/(Rainfall_rate+Procwtr_rate), use fraction 0.88 of solution 1
2 0.12 # 1-(Procwtr_rate/(Rainfall_rate+Procwtr_rate), use fraction 0.12 of solution 2
USE equilibrium_phases 1
SAVE SOLUTION 0
END
SOLUTION 1-101 Simulated mixture of (process water, DAC) and rainwater.
temp 25
pH 6.963
pe 13.64
redox pe
units mol/kgw
density 1
-water 1 # kg
C(4) 5.06E-03
Ca 3.52E-03
Cl 2.12E-03
Fe 3.94E-06
K 1.07E-06
Mg 1.93E-06
N 2.79E-05
Na 1.07E-05
P 6.26E-06
S(6) 6.45E-04
END
REACTION 0
# Remove / Evaporate water in steps. All solutions defined by SOLUTION input are scaled to have exactly 1 kg (approximately 55.5 mol) of water, unless -water identifier is used.
H2O -1
# Evaporate 11.1012 moles(= 0.2 * 55.506), 20% evaporation from pond in 4 step. # 11.1012 mole in 4 step. # The first step removes(-) 2.775e+00 moles of reaction (units are ?moles?) to the initial solution; the second step removes 5.551e+00 moles of reaction to the initial solution; the third 8.326e+00 moles; and the fourth 1.110e+01 moles; each reaction step begins with the same initial solution and removes only the amount of reaction specified.
11.1012 moles #in 4 steps
EQUILIBRIUM_PHASES 101 DAC CO2 capture to 50% # BOTTOM CELL IS THE EFFLUENT FROM THE COLUMN, AS IT EXPERIENCES CONTACT WITH DAC-AMENDED AIR
CO2(g) -0.301 10 # Calculated CO2 SI= -0.301, Amount of 10 mol available
O2(g) -0.679 10 # Calculated O2 SI= -0.679, Amount of 10 mol available
# N2(g) -0.107 10
Calcite_ppt 0 0
END
KINETICS 1-100
Portlandite
-formula Ca(OH)2 1
-m 0.1443 # 14.43/100 cells = 0.1443 cells will give mols /cell. Initial moles of reactant. Mass conc port vs pores/ Molecular-mass_port
-m0 0.1443 # Default is equal to initial moles
-parms 2.303 1 # Bulk area * Area% Port
-tol 1e-08 # Default is 1?10-8 mol
Larnite
-formula Ca2SiO4 1
-m 0.0987 # 9.87/100 cells = 0.0987 cells will give mols /cell
-m0 0.0987
-parms 2.499 1
-tol 1e-08
Calcite_diss
-formula CaCO3 1
-m 0.0949 # 9.49/100 cells = 0.0949 cells will give mols /cell
-m0 0.0949
-parms 1.682 1
-tol 1e-08
Brucite
-formula Mg(OH)2 1
-m 0.04733 # 4.733/100 cells = 0.04733 cells will give mols /cell
-m0 0.04733
-parms 0.556 1
-tol 1e-08
-steps 1
-step_divide 0.1
-runge_kutta 6
-bad_step_max 500
-cvode true
-cvode_steps 200
-cvode_order 5
END
USE SOLUTION 0
USE REACTION 0
USE equilibrium_phases 101
USE KINETICS 1
SAVE equilibrium_phases 1-10
SAVE solution 1-10
END
TRANSPORT
-cells 20 # 101 # 100 cells slag heap, last cell closed pond (1 litre)
-shifts 5 # 1245 # 12.326 passes of pore water through the slag heap
-lengths 0.01 # 1cm cells
-time_step 34722.2 # seconds retention time in cells - Antimagnetic = 0.401 days. Shift (101) x Timestep(34722.2) = 3506942.2 sec = 3506942.2 /60/ 60/ 24 = The total time of the simulation (40.5896 days).Shift (303) x Timestep(34722.2) = 10520826.6 sec = 10520826.6/60/ 60/ 24 = The total time of the simulation (121.768 days).
-print_cells 0-20
-punch_cells 1-20
-print_frequency 5 # 1245 # Print only at end
-punch_frequency 5
USER_GRAPH 1
-headings dist Portlandite Larnite Calcite Brucite H2O_kg
-axis_titles "Distance" "SI" "Water, kg"
-initial_solutions false
-connect_simulations true
-plot_concentration_vs x
-start
10 graph_x DIST
20 graph_y SI("Portlandite"), SI("Larnite"), SI("Calcite_diss"), SI("Brucite")
30 graph_sy TOT("water")
100 END
-end
END
--- End code ---
Erdinc:
Thanks David. I have calculated the shift by number of days I want to simulate (500 days in my case) to understand the dissolution of the mineral kinetics. I have run the simulations from different percentages of pCO2 values and from 0.1% up to 6% there is good indication on how much of the minerals dissolve during the life of the 500 days. When I switch the pCO2 value 7% and above up to 99%, I see the simulations stops completing all the way to 500 days and therefore same for minerals dissolutions and CO2 capture potentials are not known. I was thinking the effect of accurate thermodynamic solubility data, the non-ideal behavior of CO2 fugacity, and potential limitations in the way PHREEQC handles gas-water interactions would be a reason for the high pCO2 values. Is this something you come across before? Unfortunately, I do not have options to attach what I am simulating in regard to graphs but this is what I can explain to the best of my knowledge.
Many Thanks for the support
dlparkhurst:
Okay, I ran your (very) long simulation with llnl.dat and a version of carbfix.dat. Both ran to completion.
What error am I supposed to see?
Will a 20 cell model with no mixruns demonstrate the problem? The current run time is too long to debug usefully.
I'll send you a message and you can send me a script or graph that shows me what you are talking about. I don't see anything wrong.
Navigation
[0] Message Index
Go to full version