Get different selected output value with same script in IPHREEQC and PHREEQC

[Yongqiang]

Hi PHREEQC users and developers,

I am using IPHREEQC and PHREEQC to calculation dissolution. However, I got different results with same input. Please see the following input code:

IPHREEQC

Code: [Select]

# -*- coding: utf-8 -*-
"""
Created on Tue Mar  7 09:29:09 2023

@author: che458
"""

from win32com.client import Dispatch

def result(input_str):
    db_path = 'LLNL.DAT'#'carbfix.dat'
    dbase = Dispatch('IPhreeqcCOM.Object')
    dbase.LoadDatabase(db_path)
    dbase.RunString(input_str)
    return dbase.GetSelectedOutputArray()

if __name__ == '__main__':
    a =     """
DATABASE c:\phreeqc\database\LLNL.DAT

SOLUTION 1
EQUILIBRIUM_PHASES 1
CO2(g) 3 100
SAVE SOLUTION 1
END

SELECTED_OUTPUT
    -file  selected.out
    -reset false
USER_PUNCH
    -headings pH Mol(CaCO3) Mol(CO3-2) Mol(Ca) TotalC Al An
    -start
    10 PUNCH -LA("H+")#PR_P("CO2(g)")
    20 PUNCH MOL("CaCO3")*1000
    30 PUNCH MOL("CO3-2")*1000
    40 PUNCH MOL("Ca+2")*1000
    50 PUNCH TOT("C")
    60 PUNCH -EQUI_DELTA("Albite")*1000
    70 PUNCH -EQUI_DELTA("Anorthite")*1000
    -end
END



USE SOLUTION 1
EQUILIBRIUM_PHASES 2
Albite 0 24.9
Anorthite 0 24.5
SAVE SOLUTION 2
END
 

"""
    tem = result(a)     
    print(tem)

PHREEQC:

Code: [Select]

DATABASE c:\phreeqc\database\LLNL.DAT

SOLUTION 1
EQUILIBRIUM_PHASES 1
CO2(g) 3 100
SAVE SOLUTION 1
END

SELECTED_OUTPUT
    -file  selected.out
    -reset false
USER_PUNCH
    -headings pH Mol(CaCO3) Mol(CO3-2) Mol(Ca) TotalC Al An
    -start
    10 PUNCH -LA("H+")#PR_P("CO2(g)")
    20 PUNCH MOL("CaCO3")*1000
    30 PUNCH MOL("CO3-2")*1000
    40 PUNCH MOL("Ca+2")*1000
    50 PUNCH TOT("C")
    60 PUNCH -EQUI_DELTA("Albite")*1000
    70 PUNCH -EQUI_DELTA("Anorthite")*1000
    -end
END



USE SOLUTION 1
EQUILIBRIUM_PHASES 2
Albite 0 24.9
Anorthite 0 24.5
SAVE SOLUTION 2
END

Look forward to your insights.

Regards,
Michael

[Yongqiang]

The amount of albite dissolution is different though other outputs are same. 

[dlparkhurst]

Are you sure you are using the same database? Or perhaps you are using different versions of PHREEQC between the com and the batch version.

Add the log K of albite to your output file.

Code: [Select]

LK_PHASE("Albite")

Otherwise, print the output file and look for differences.

[Yongqiang]

Hi David,

Thanks for your help.

I used your suggested method to output the logK of albite. Both of them outputted 2.6574 for albite dissolution.

I also checked the thermodynamic database. Both of them used LLNL.DAT.

The IPHREEQC and batch PHREEQC versions are 3.7.3-15968. I tried to output the file with the COM version in python. However, it seems that no output file is available for COM version in Python.

I tried to lower the pressure for CO2(g). Same results can be obtained at lower pressure, such as -3 and -1 for Equilibrium_phase command. However, deviations emerge at higher pressure.

Best Regards,
Michael

[dlparkhurst]

Can you give me an example of which numbers differ? and by how much?

[Yongqiang]

Hi David,

Thanks for your reply. The following are the numbers:
IPHREEQC-COM version:

Code: [Select]

(('pH', 'Mol(CaCO3)', 'Mol(CO3-2)', 'Mol(Ca)', 'TotalC', 'Al', 'An', 'LogK(Al)'),
(4.794098585416058, 0.08333842715243284, 0.011002936982880721, 112.39594014981397, 35.108246694244414, -8.071054935498978e-08, 333.8527245630374, 2.6574779847717167))The -EQUI_DELTA("Albite")*1000 is -8.071054935498978e-08.

PHREEQC batch version:

Code: [Select]

          pH   Mol(CaCO3)   Mol(CO3-2)      Mol(Ca)       TotalC           Al           An     Logk(Al)
  4.7941e+00   8.3338e-02   1.1003e-02   1.1240e+02   3.5108e+01   1.9279e-08   3.3385e+02   2.6575e+00
The -EQUI_DELTA("Albite")*1000 is 1.9279e-08.

We can see other values are the same. Only the data for albite is different.

Best Regards,
Michael

[dlparkhurst]

The default convergence criterion is 1e-8, so relative differences in this range are possible.

Add the following before the first END and see if your results become more similar:

Code: [Select]

KNOBS
-conv 1e-10


[Yongqiang]

Hi David,

Yes, the COM version gave a same value as the batch version. Thanks.

So the COM version is different with the batch version, correct?

Cheers,
Michael

[dlparkhurst]

No, it is exactly the same code, but there can be different rules for floating point operations between the com and the batch version. Every processor and environment can have different rounding conventions, different numbers of extra bits in operations, and other details that can lead to slightly variable results. If there are still differences after lowering the tolerance to between 10e-10 and 10-12, then there may be a problem. 

[Yongqiang]

Totally got you! That makes sense.

Many thanks, David.