
zeros (( n, 1 )) d_depnt_d_indep = - rAdash / FA0 d_depnt_d_indep = - alpha / ( 2 * y ) * ( 1 + eps * X ) return d_depnt_d_indep # The "driver" that will integrate the ODE(s): # - # Start by specifying the integrator: # use ``vode`` with "backward differentiation formula" r = integrate.

Import numpy as np from scipy import integrate from matplotlib.pylab import ( plot, grid, xlabel, ylabel, show, legend ) def PBR_model ( indep, depnt ): """ Dynamic balance for the packed bed reactor (PBR) demo problem class 05C indep: the independent ODE variable, such as time or length depnt: a vector of dependent variables X = depnt = the conversion y = depnt = the pressure ratio = P/P_0 = y Returns d(depnt)/d(indep) = a vector of ODEs """ # Assign some variables for convenience of notation X = depnt y = depnt # Constants FA0 = 0.1362 # kdash = 0.0074 # alpha = 0.0367 # eps = - 0.15 #Algebraic equations rAdash = - kdash * ( 1 - X ) / ( 1 + eps * X ) * y # Output from this ODE function must be a COLUMN vector, with n rows n = len ( depnt ) d_depnt_d_indep = np. Check the video out for a detailed explanation of the code. So rather use one of the other options below, if possible.Įnter the following Polymath code (as demonstrated in class on 11 February. Unfortunately this is only a 15 day version, which is pretty useless for this course, since you will require the code for at least the next 30 to 45 days. The tutorial below is similar, but uses a reactor design example.Ĭonsider the pair of differential equations (covered in class on 11 February):ĭownload and install Polymath from the CD/DVD included with your course textbook (Windows only sorry Mac and Linux versions are not available). Here is a tutorial a wrote a few years ago when I taught 3E4 in 2010.

If you are not familiar with this topic, it is your responsibility to quickly catch up, because we will use this intensively for the rest of the course. Numerically integrating ODE's is something you should have mastered in your pre-requisite 3E4 course.
