Skip to content

Conversation

@vineetbansal
Copy link
Collaborator

@vineetbansal vineetbansal commented Jan 12, 2026

This is most likely a bug that we haven't encountered yet, simply because we're getting lucky in the order of the keys (representing species->amount mapping) in the current wrapper.

The tests that should fail in this PR is test_phreeqc.py::test_equilibrate_logC_pH_carbonate_8_3 (among others).

The issue seems to be that for a simple solution like this:

Solution({"CO2(aq)": "0.001 mol/L"}, pH=8.3, volume="1 L", engine="phreeqc")

corresponding to the string:

SOLUTION 0
  temp 25.0
  units mol/kgw
  pH 8.3
  pe 8.5
  redox pe
  water 0.9970480319717386
  C(4) 0.0010029607079434515
SAVE SOLUTION 0
END

In addition to CO2, phreeqc reports a tiny amount of (CO2)2 (whatever that represents, not sure..). The following is the output from phreeqc UI, and is consistent with what the old (and our new) wrappers return.

----------------------------Distribution of species----------------------------

                                               Log       Log       Log    mole V
   Species          Molality    Activity  Molality  Activity     Gamma   cm3/mol

   OH-             2.072e-06   2.019e-06    -5.684    -5.695    -0.011     -4.12
   H+              5.138e-09   5.012e-09    -8.289    -8.300    -0.011      0.00
   H2O             5.551e+01   1.000e+00     1.744    -0.000     0.000     18.07
C(4)          1.003e-03
   HCO3-           9.822e-04   9.575e-04    -3.008    -3.019    -0.011     24.57
   CO2             1.079e-05   1.079e-05    -4.967    -4.967     0.000     34.43
   CO3-2           9.923e-06   8.959e-06    -5.003    -5.048    -0.044     -3.97
   (CO2)2          2.137e-12   2.137e-12   -11.670   -11.670     0.000     68.87
H(0)          3.556e-37
   H2              1.778e-37   1.778e-37   -36.750   -36.750     0.000     28.61
O(0)          2.636e-19
   O2              1.318e-19   1.318e-19   -18.880   -18.880     0.000     30.40

------------------------------Saturation indices-------------------------------

The standardize_formula function standardizes (CO2)2 as CO2, and can possibly end up overwriting the CO2 amount with the (CO2)2 amount (which is tiny), in the following logic in .equilibrate():

        solution.components = FormulaDict({})
        # use the output from PHREEQC to update the Solution composition
        # the .species_moles attribute should return MOLES (not moles per ___)
        for s, mol in self.ppsol.species_moles.items():
            solution.components[s] = mol

We simply get lucky that (CO2)2 is consistenly reported before CO2 by the current wrapper, so we end up with the correct amount of CO2. In the new wrappers, the order of the keys in the dict is reversed (for whatever reason, but most probably because we're building the dictionary in the same order as the output of phreeqc UI indicates), so we do encounter this bug.

The attached snippet simply surfaces the bug without fixing it. I'm not sure of a solution without the background on what it means/entails.

@vineetbansal vineetbansal changed the title snippet in code to surface a possible bug Possible bug in standardize_formula Jan 12, 2026
@codecov
Copy link

codecov bot commented Jan 12, 2026

Codecov Report

✅ All modified and coverable lines are covered by tests.
✅ Project coverage is 84.81%. Comparing base (4d1454c) to head (2a0e8d8).
⚠️ Report is 12 commits behind head on main.

Additional details and impacted files
@@             Coverage Diff             @@
##             main     #309       +/-   ##
===========================================
+ Coverage   74.55%   84.81%   +10.25%     
===========================================
  Files          10       10               
  Lines        1521     1521               
  Branches      265      264        -1     
===========================================
+ Hits         1134     1290      +156     
+ Misses        352      199      -153     
+ Partials       35       32        -3     

☔ View full report in Codecov by Sentry.
📢 Have feedback on the report? Share it here.

🚀 New features to boost your workflow:
  • ❄️ Test Analytics: Detect flaky tests, report on failures, and find test suite problems.

@rkingsbury
Copy link
Member

Ooof, good catch @vineetbansal . I have no idea how to interpret (CO2)2. I suppose for the time being we should add an ad-hoc patch to standardize_formula similar to what was done here (see for example - C2O4[-2]

rkingsbury
rkingsbury previously approved these changes Jan 13, 2026
@rkingsbury rkingsbury changed the title Possible bug in standardize_formula Fix incorrect reduction of dimers/polymers in standardize_formula Jan 13, 2026
@rkingsbury rkingsbury added the fix Bug Fixes label Jan 13, 2026
Comment on lines 740 to 749
# --- Code to demonstrate bug ---
# The standardized formula for (CO2)2 is still CO2, so this amount
# (negligible) ends up overwriting the original CO2 amount.
# We simply get lucky in the loop above that (CO2)2 is encountered
# before CO2 in the old wrapper, but we force surfacing the bug here.
problematic_keys = ("(CO2)2",)
for k in problematic_keys:
if k in self.ppsol.species_moles:
solution.components[k] = self.ppsol.species_moles[k]

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I think this can be removed, now that the fix is in place

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Thanks! I just wanted to see if the tests pass before removing this. Hopefully the fix makes sense.

rkingsbury
rkingsbury previously approved these changes Jan 13, 2026
@rkingsbury rkingsbury merged commit d3e9c63 into main Jan 13, 2026
17 checks passed
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

fix Bug Fixes

Projects

None yet

Development

Successfully merging this pull request may close these issues.

3 participants