I have an assignment to plot two consecutive first order reactions, then find the maximum concentration for reaction B.
I've managed to plot a graph of the three functions, i'm struggling with
finding the maximum value. My teacher told me to use optimize.fmin()
(I figure he wants me to create a negative of the function for reaction B,
and find the minimum for that function, which should be the maximum.),
only trouble is it's not working!
here's what I have so far, I have tried other value other than 0.75 for the second argument in the optimize.fmin()
function.
Where am I going wrong here ? I see the error is saying it is expecting an array but getting a sequence ? is this in relation to the t=linspace(0,tmax,20)
line
of code where I create 20 evenly spaced points of to a total runtime of 20 minutes for the experiment
%pylab inline
from matplotlib import *
from scipy import *
k1 = 0.15
k2 = 0.10
A0 = 2
tmax = 21
t = linspace(0,tmax,20)
e1= e**(-k1*t)
e2= e**(-k2*t)
def conc_A(t):
A = A0 * e1
return A
def conc_B(t):
B = A0 *(k1 / (k2-k1)) * (e1 - e2)
return B
def conc_C(t):
C = (A0/ (k2-k1)) * (k2 * ((1 - e1 ))* - (k1 *(1-e2)))
return C
pylab.plot(t,conc_A(t),label ='[A]')
plot(t,conc_B(t),label= '[B]')
plot(t,conc_C(t),label= '[C]')
pylab.legend(loc='upper right')
plt.xlabel("Time (minutes)" )
plt.ylabel("Concentration Mol $Dm^{-3}$")
plt.title("Rates of reaction of two consecutive first order reactions")
def neg_B(t):
return -conc_B(t)
optimize.fmin(neg_B,0.75)
The error I get is
---------------------------------------------------------------------------
ValueError Traceback (most recent call last)
<ipython-input-88-f481d0e274f9> in <module>()
48
49
---> 50 optimize.fmin(neg_B,0.75)
51
52
C:\Users\gsandle1\AppData\Local\Continuum\Anaconda2\lib\site-packages\scipy\optimize\optimize.py in fmin(func, x0, args, xtol, ftol, maxiter, maxfun, full_output, disp, retall, callback, initial_simplex)
391 'initial_simplex': initial_simplex}
392
--> 393 res = _minimize_neldermead(func, x0, args, callback=callback, **opts)
394 if full_output:
395 retlist = res['x'], res['fun'], res['nit'], res['nfev'], res['status']
C:\Users\gsandle1\AppData\Local\Continuum\Anaconda2\lib\site-packages\scipy\optimize\optimize.py in _minimize_neldermead(func, x0, args, callback, maxiter, maxfev, disp, return_all, initial_simplex, xatol, fatol, **unknown_options)
515
516 for k in range(N + 1):
--> 517 fsim[k] = func(sim[k])
518
519 ind = numpy.argsort(fsim)
ValueError: setting an array element with a sequence.
Taking a lot more points in t
, e.g. 2000, you may find the time of the maximum by finding the index of the array at which the numerical maximum occurs.
t = np.linspace(0,tmax,2000)
#... rest of code
print t[np.argmax(conc_B(t))] # prints 8.11
If you want to use optimize.fmin
, I would suggest you first read the documentation. It states that the first argument needs to be a function. So you need to provide a function that should be minimized.
import numpy as np
import scipy.optimize as optimize
k1 = 0.15
k2 = 0.10
A0 = 2.
tmax = 21
e1= lambda t: np.exp(-k1*t)
e2= lambda t: np.exp(-k2*t)
conc_A = lambda t: A0 * e1(t)
conc_B= lambda t: A0 *(k1 / (k2-k1)) * (e1(t) - e2(t))
conc_C = lambda t: (A0/ (k2-k1)) * (k2 * ((1. - e1(t) ))* -(k1 *(1.-e2(t))))
print optimize.fmin(lambda t: -conc_B(t),0.75)
# Optimization terminated successfully.
# Current function value: -0.888889
# Iterations: 23
# Function evaluations: 46
# [ 8.10930176]