PhreeqcUsers Discussion Forum

Conceptual Models => Database selection and modification => Topic started by: SaiP on October 06, 2020, 04:15:49 PM

Title: Turn off saturation index while considering mineral precip/dissolution?
Post by: SaiP on October 06, 2020, 04:15:49 PM
Dear David,

Tools I am using:
phreeqcRM (to solve chemistry) + CFD package (to solve transport equation)

Modelling:
Calcite grain dissolution discussed in paper -
"https://link.springer.com/article/10.1007/s10596-019-09903-x"

Information specified in paper for kinetics -
reaction rate = k*g*c [k = rate const = 0.89 mol/m2s, g = activity coeff = 0.001 m3/mol, c = concentration of H+ = 0.01 mol/l];
As you can notice, there is no saturation index.

My input file:
Code: [Select]
SOLUTION 0 #inject acidic solution of ph=2
        temp 25
        pH 2
        units mmol/kgw
        Cl 10 charge
        C(4) 1e-9
END
SOLUTION 1 #initial solution around calcite grain, similar to influx solution
        temp 25
        pH 2
        units mmol/kgw
        Cl 10 charge
        C(4) 1e-9
END
SOLUTION 2 #within the grain having porosity 0.001
        temp 25
        pH 10 charge
        units mmol/kgw
        EQUILIBRIUM_PHASES
        Calcite 0.
END

and my database manipulation for incorporating the reaction rate mentioned above:
Code: [Select]
Calcite
    -start
1   REM   PARM(1) = specific surface area of Calcite (large grains), [m^2/m^3]
10  si_cc = SR("Calcite")
20  rate = 8.9125e-4*PARM(1)*MOL("H+")
30  mole = rate*time #mol/m3
35  mole = mole #mol/l
40  save mole
    -end
In the above code I use parm1 = 20000 1/m which is the reactive surface area. But i notice, even if I increase parm1 by orders of magnitude the dissolution rate is not improved.

I dont know where I am going wrong but I see no dissolution at all in an hour unlike results published in paper.

So, I decided to use "Thermodynamic eq." with no satuation index. I comment out the calcite phases "-logK" and "-analytic" data. By doing so, my solver just crashes at the 1st time step complaining about convergence issues (below)
Code: [Select]
WARNING: Negative moles in solution 922 for C, -1.425926e-02. Recovering...
WARNING: Negative moles in solution 922 for C, -9.981482e-04. Recovering...
WARNING: Negative moles in solution 922 for C, -6.987039e-05. Recovering...
WARNING: Negative moles in solution 922 for C, -3.171901e-05. Recovering...
WARNING: Maximum iterations exceeded, 100

WARNING: Numerical method failed with this set of convergence parameters.

WARNING: Trying smaller step size, pe step size 10, 5 ...

WARNING: Maximum iterations exceeded, 200

WARNING: Numerical method failed with this set of convergence parameters.

WARNING: Trying reduced tolerance 1e-16 ...

WARNING: Maximum iterations exceeded, 200

WARNING: Numerical method failed with this set of convergence parameters.

WARNING: Trying increased tolerance 1e-14 ...

WARNING: Maximum iterations exceeded, 200

WARNING: Numerical method failed with this set of convergence parameters.

WARNING: Trying diagonal scaling ...

WARNING: Maximum iterations exceeded, 200

WARNING: Numerical method failed with this set of convergence parameters.

WARNING: Trying diagonal scaling and reduced tolerance 1e-16 ...

WARNING: Maximum iterations exceeded, 200

WARNING: Numerical method failed with this set of convergence parameters.

WARNING: Trying scaling pure_phase columns 1e-10 ...

WARNING: Maximum iterations exceeded, 200

WARNING: Numerical method failed with this set of convergence parameters.

WARNING: Trying scaling pure_phase columns and diagonal scale 1e-10 ...

WARNING: Negative moles in solution 669 for C, -1.425926e-02. Recovering...
WARNING: Maximum iterations exceeded, 200

WARNING: Numerical method failed with this set of convergence parameters.

WARNING: Negative moles in solution 669 for C, -9.981482e-04. Recovering...
WARNING: Trying delay removal of equilibrium phases 1 ...

WARNING: Negative moles in solution 669 for C, -6.987039e-05. Recovering...
WARNING: Maximum iterations exceeded, 100

WARNING: Numerical method failed with this set of convergence parameters.

WARNING: Negative moles in solution 669 for C, -3.175318e-05. Recovering...
WARNING: Trying increased scaling 1e-09 ...

WARNING: Maximum iterations exceeded, 200

WARNING: Numerical method failed with this set of convergence parameters.

WARNING: Skipping optimize equations for first 5 iterations ...

WARNING: Maximum iterations exceeded, 100

WARNING: Numerical method failed with this set of convergence parameters.

WARNING: Adding inequality to make concentrations greater than zero.

WARNING: Maximum iterations exceeded, 100

WARNING: Numerical method failed with this set of convergence parameters.

WARNING: Trying reduced tolerance 1e-17 ...

WARNING: Maximum iterations exceeded, 200

WARNING: Numerical method failed with this set of convergence parameters.

WARNING: Trying reduced tolerance 1e-18 ...

WARNING: Maximum iterations exceeded, 200

WARNING: Numerical method failed with this set of convergence parameters.

WARNING: The program has failed to converge to a numerical solution.

The following equations were not satisfied:
ERROR:               A(H2O) Activity of water has not converged.        Residual: 3.140644e-06

ERROR: Numerical method failed on all combinations of convergence parameters, cell/soln/mix 922
Stopping.
ERROR: ERROR:               A(H2O) Activity of water has not converged.         Residual: 3.140644e-06

ERROR: Numerical method failed on all combinations of convergence parameters, cell/soln/mix 922

WARNING: Maximum iterations exceeded, 100

WARNING: Numerical method failed with this set of convergence parameters.

WARNING: Trying smaller step size, pe step size 10, 5 ...

WARNING: Maximum iterations exceeded, 200

WARNING: Numerical method failed with this set of convergence parameters.

WARNING: Trying reduced tolerance 1e-16 ...

WARNING: Maximum iterations exceeded, 200

WARNING: Numerical method failed with this set of convergence parameters.

WARNING: Trying increased tolerance 1e-14 ...

WARNING: Maximum iterations exceeded, 200

WARNING: Numerical method failed with this set of convergence parameters.

WARNING: Trying diagonal scaling ...

WARNING: Maximum iterations exceeded, 200

WARNING: Numerical method failed with this set of convergence parameters.

WARNING: Trying diagonal scaling and reduced tolerance 1e-16 ...

WARNING: Maximum iterations exceeded, 200

WARNING: Numerical method failed with this set of convergence parameters.

WARNING: Trying scaling pure_phase columns 1e-10 ...

WARNING: Maximum iterations exceeded, 200

WARNING: Numerical method failed with this set of convergence parameters.

WARNING: Trying scaling pure_phase columns and diagonal scale 1e-10 ...

WARNING: Maximum iterations exceeded, 200

WARNING: Numerical method failed with this set of convergence parameters.

WARNING: Trying delay removal of equilibrium phases 1 ...

WARNING: Maximum iterations exceeded, 100

WARNING: Numerical method failed with this set of convergence parameters.

WARNING: Trying increased scaling 1e-09 ...

WARNING: Maximum iterations exceeded, 200

WARNING: Numerical method failed with this set of convergence parameters.

WARNING: Skipping optimize equations for first 5 iterations ...

WARNING: Maximum iterations exceeded, 100

WARNING: Numerical method failed with this set of convergence parameters.

WARNING: Adding inequality to make concentrations greater than zero.

WARNING: Maximum iterations exceeded, 100

WARNING: Numerical method failed with this set of convergence parameters.

WARNING: Trying reduced tolerance 1e-17 ...

WARNING: Maximum iterations exceeded, 200

WARNING: Numerical method failed with this set of convergence parameters.

WARNING: Trying reduced tolerance 1e-18 ...

WARNING: Maximum iterations exceeded, 200

WARNING: Numerical method failed with this set of convergence parameters.

WARNING: The program has failed to converge to a numerical solution.

The following equations were not satisfied:
ERROR:               A(H2O) Activity of water has not converged.        Residual: 3.144049e-06

ERROR: Numerical method failed on all combinations of convergence parameters, cell/soln/mix 669
Stopping.
ERROR: ERROR:               A(H2O) Activity of water has not converged.         Residual: 3.144049e-06

ERROR: Numerical method failed on all combinations of convergence parameters, cell/soln/mix 669

ERROR:               A(H2O) Activity of water has not converged.        Residual: 3.144049e-06

ERROR: Numerical method failed on all combinations of convergence parameters, cell/soln/mix 669

ERROR:               A(H2O) Activity of water has not converged.        Residual: 3.140644e-06

ERROR: Numerical method failed on all combinations of convergence parameters, cell/soln/mix 922

ERROR: PhreeqcRM failed.
ERROR: PhreeqcRM::RunCells


terminate called after throwing an instance of 'PhreeqcRMStop'
  what():  Failure in PhreeqcRM

Aborted (core dumped)

I tried to reduce my time step size (upto 1e-7s), but it still crashes stating convergence issues. If I have -logK of calcite upto -3 things are running (i see the grain dissolve faster as expected). Why am i unable to turn off logK (saturation index, by setting it to 0)? Is something wrong in my input file? I just ran the same input data without transport as a batch reaction and it works well with phreeqc (all Calcite is dissolved).

Thanks!!
Title: Re: Turn off saturation index while considering mineral precip/dissolution?
Post by: dlparkhurst on October 06, 2020, 05:23:50 PM
It is hard to know what is going on in your transport simulation; I can only address the functioning of a PHREEQC input file. The following is my interpretation of your calculation. You will find that the results do depend on parameter 1.

Code: [Select]
RATES
Calcite
    -start
1   REM   PARM(1) = specific surface area of Calcite (large grains), [m^2/m^3]
10  si_cc = SR("Calcite")
20  rate = 8.9125e-4*PARM(1)*MOL("H+")
30  mole = rate*time #mol/m3
35  mole = mole #mol/l
40  save mole
    -end
END
SOLUTION
pH 2
Cl 1 charge
END
USE solution 1
KINETICS 1
Calcite
-parm 20000
-time 3600 in 10 steps
USER_GRAPH 1
    -headings               Time TOT(Ca) SI(Calcite)
    -axis_titles            "Time, seconds" "Molality" "Saturation index"
    -initial_solutions      false
    -connect_simulations    true
    -plot_concentration_vs  x
  -start
10 GRAPH_X TOTAL_TIME
20 GRAPH_Y TOT("Ca")
30 GRAPH_SY SI("Calcite")
  -end
    -active                 true
END

Note that the RATE  is not limited by a factor of (1-SR(Calcite)). In the simulation the SI increases to about 2, or 2 orders of magnitude supersaturated in about an hour. Without the SR factor, it is possible that the concentrations of Ca and C could increase unbounded, which could cause a convergence failure.

I am not sure what you did in changing the Calcite PHASES definition, but if you remove log_k and analytic expression, the log K will default to zero.