USER_GRAPH RETURNING ZEROS WITH KINETIC SIMULATIONS

[swhit]

I am trying to simulate a kinetic beaker reaction where I add Calcite to a beaker of acidic solution, and I want to plot the results using USER_GRAPH.  However, for some reason the USER_GRAPH function is returning zeros for some quantities.

When I export the data to excel using SELECTED_OUTPUT, the values are non-zero.

I would prefer to be able to use USER_GRAPH, as it allows for a more rapid model development process since I don't have to go into excel (or python, or whatever) after each model execution to view the results.

Any thoughts as to why USER_GRAPH is returning zero values?


CODE BELOW:

SOLUTION 1# Analytical Results for Pond 3, from Table 20 of EE/CA
pH 2.5 # Default: 7.
pe 4.0 #11.59 # Default: 4., 11.59 from tailings beaker test
redox pe
-units mg/L # Default: mmol/kgw
Al 100 mg/L
Ba 1.7 mg/L
Cd 1.7 mg/L
Ca 500.0 mg/L
Cu 8.9 mg/L
Fe 3100 mg/L
Mg 97 mg/L
Mn 1200 mg/L
Zn 610 mg/L
Na 0.0 mg/L
Cl 5.8
F 12.0 mg/L
S(6) 15000 charge  # concentration is adapted to charge balance
C(4) 0.6 as HCO3 # Concentration is in mg of HCO3 = 0.6 / 0.61 = 0.98 mmol/L
Si 1.0 mg/L

EQUILIBRIUM_PHASES 1 # cell number or a range of cells. Default: 1.
#atmospheric gases
CO2(g) -3.5
O2(g) -0.67

#soil phases (update to match XRD)
#Calcite 0 0 # name (must be defined in PHASES), Saturation Index (Default: 0), initial amount (moles, default: 10 moles)
Dolomite 0 0 #diss # Dolomite can dissolve only
Quartz 0 0 #precipitate # Quartz can precipitate only if uncommented
Chalcedony 0 0

#secondary phases
Fe(OH)3(a) 0 0
Gypsum 0 0

SAVE SOLUTION 1

END

SOLUTION 2 #pure water
pH 7
pe 4
temp 25

USE EQUILIBRIUM_PHASES 1

SAVE SOLUTION 2

END

MIX 1 #based on previous experimental procedures
1 .35
2 .65

SAVE SOLUTION 3

END

USE SOLUTION 3
USE EQUILIBRIUM_PHASES 1

KINETICS 1
Calcite
   -m       0
   -m0      1
   -parms   12500.0      30
   -tol     1.e-8
-steps 30000 in 100 steps #seconds
INCREMENTAL_REACTIONS true


RATES
Calcite
-start
1   rem M = current number of moles of calcite
2   rem M0 = number of moles of calcite initially present
3   rem PARM(1) = A/V, cm^2/L
4   rem PARM(2) = exponent for M/M0
10  si_cc = SI("Calcite")
20  if (M <= 0 and si_cc < 0) then goto 200
30    k1 = 10^(0.198 - 444.0 / TK )
40    k2 = 10^(2.84 - 2177.0 / TK)
50    if TC <= 25 then k3 = 10^(-5.86 - 317.0 / TK )
60    if TC > 25 then k3  = 10^(-1.1 - 1737.0 / TK )
70    t = 1
80    if M0 > 0 then t = M/M0
90    if t = 0 then t = 1
100   area = PARM(1) * (t)^PARM(2)
110   rf = k1*ACT("H+")+k2*ACT("CO2")+k3*ACT("H2O")
120   rem 1e-3 converts mmol to mol
130   rate = area * 1e-3 * rf * (1 - 10^(2/3*si_cc))
140   moles = rate * TIME
200 SAVE moles
-end

REACTION 3 # cell number or a range of cells. Default: 1.
Calcite 1 # Reactant, followed by  stoichiometric coefficient (useful for multiple react.)
2e-1 in 100 steps # Adds x amount in x equal steps

USER_GRAPH 1
    -headings Mol_KCalcite #Mol_KCalcite_Delta Mol_ECalcite Mol_Ferrihydrite #why are moles returning as zero?
    -axis_titles "Time (s)" "Moles"
    -initial_solutions true
    #-axis_scale y_axis 0 1e-1 1e-2
  -start
  10 PLOT_XY SIM_TIME, KIN("Calcite"), color = Blue, line_w = 0, \
            symbol = Circle, sybol_size = 7
  -end
USER_GRAPH 2
    -headings Mol_Calcium
    -axis_titles "Time (s)" "Moles"
    -initial_solutions true
    #-axis_scale y_axis 0 1e-1 1e-2
  -start
  10 PLOT_XY SIM_TIME, TOT("Ca"), color = Red, line_w = 0, \
      symbol = XCross
  -end
USER_GRAPH 3
    -headings SI_Calcite
    -axis_titles "Time (s)" "SI_Calcite"
    -initial_solutions true
    #-axis_scale y_axis -2 2 1e-1
  -start
  10 PLOT_XY SIM_TIME, SI("Calcite"), color = Blue, line_w = 0, \
            symbol = Circle, sybol_size = 7
  -end
USER_GRAPH 4
    -headings pH
    -axis_titles "Time (s)" "pH"
    -initial_solutions true
    #-axis_scale y_axis -2 2 1e-1
  -start
  10 PLOT_XY SIM_TIME, -LA("H+"), color = Blue, line_w = 0, \ #why is pH returning as 0?
            symbol = Circle, sybol_size = 7
  -end


SELECTED_OUTPUT
-file D:\GoogleDrive\GradSchool\Research\Caselton\PHREEQC\BatchRxnSulfideTailingsAlluvium.csv # file name
-selected_out true
-user_punch true
-high_precision false
-reset false
-pH true
-pe true
-water true
-ionic_strength true
-time true
-solution true
-equilibrium_phases Calcite
-kinetic_reactants Calcite
-saturation_indices Calcite

[dlparkhurst]

You misspelled symbol_size.

Also, you probably do not want M = 0 and M0 = 1 while using the multiplication factor of M/M0 in RATES.

[swhit]

Doh!

I guess I didn't realize that it wouldn't produce an error message so the behavior confused me.

Thanks for the note about M/M0 -- I hadn't gotten around to setting that up properly yet since I struggling with plotting.