Hi, professor,
when I learn this case (Kinetic Dissolution of Quartz), I don't know how to get the value of quartz_dissolved.
When I add the command "print(quartz_dissolved)" in the "ratefun()" function, I can get a lot of value (just one time_step). However, this value is not the same with the value I got from the customer software.
I hope to get your help.
best wish!!
%pylab inline
import phreeqpython
from scipy.integrate import odeint
pp = phreeqpython.PhreeqPython('phreeqc.dat')
def ratefun(sol, quartz_dissolved, m0, A0, V):
m = m0-quartz_dissolved
rate = (A0/V)*(m/m0)*0.6710*-13.7(1-sol.sr("Quartz"))
print(quartz_dissolved) # I try this method, but it is wrong.
return rate * 1e3
solution1 = pp.add_solution({})
t = np.array([])
y = []
year = 365243600
for time, sol in solution1.kinetics('SiO2',
rate_function=ratefun,
time=np.linspace(0,5*year, 2),
m0=158.5,
args=(23.13,0.16)):
t = np.append(t, time)
y.append(sol.total_element('Si', units='mmol'))
Hi, professor, when I learn this case (Kinetic Dissolution of Quartz), I don't know how to get the value of quartz_dissolved. When I add the command "print(quartz_dissolved)" in the "ratefun()" function, I can get a lot of value (just one time_step). However, this value is not the same with the value I got from the customer software. I hope to get your help. best wish!!
%pylab inline import phreeqpython from scipy.integrate import odeint pp = phreeqpython.PhreeqPython('phreeqc.dat') def ratefun(sol, quartz_dissolved, m0, A0, V): m = m0-quartz_dissolved rate = (A0/V)*(m/m0)*0.6710*-13.7(1-sol.sr("Quartz")) print(quartz_dissolved) # I try this method, but it is wrong. return rate * 1e3 solution1 = pp.add_solution({})
t = np.array([]) y = []
year = 365243600
for time, sol in solution1.kinetics('SiO2', rate_function=ratefun, time=np.linspace(0,5*year, 2), m0=158.5, args=(23.13,0.16)): t = np.append(t, time) y.append(sol.total_element('Si', units='mmol'))