{VERSION 2 3 "IBM INTEL NT" "2.3" } {USTYLETAB {CSTYLE "Maple Input" -1 0 "Helvetica" 1 9 255 0 0 1 0 1 0 0 1 0 0 0 0 }{CSTYLE "2D Math" -1 2 "Times" 1 10 0 0 0 0 0 0 2 0 0 0 0 0 0 }{CSTYLE "2D Comment" 2 18 "" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 } {CSTYLE "2D Output" 2 20 "" 0 1 0 0 255 1 0 0 0 0 0 0 0 0 0 }{CSTYLE " " -1 256 "Courier" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 }{CSTYLE "" -1 257 "C ourier" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 }{CSTYLE "" -1 258 "Courier" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 }{CSTYLE "" -1 259 "Courier" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 }{CSTYLE "" -1 260 "Courier" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 }{PSTYLE "Normal" -1 0 1 {CSTYLE "" -1 -1 "Times" 1 12 0 0 0 0 2 2 2 0 0 0 0 0 0 }0 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }{PSTYLE "Heading 1" 0 3 1 {CSTYLE "" -1 -1 "" 1 18 0 0 0 0 0 1 0 0 0 0 0 0 0 }1 0 0 0 6 6 0 0 0 0 0 0 -1 0 }{PSTYLE "Heading 2" 3 4 1 {CSTYLE "" -1 -1 "" 1 14 0 0 0 0 0 0 0 0 0 0 0 0 0 }0 0 0 -1 4 4 0 0 0 0 0 0 -1 0 }{PSTYLE " Heading 3" 4 5 1 {CSTYLE "" -1 -1 "" 1 12 0 0 0 0 1 0 0 0 0 0 0 0 0 } 0 0 0 -1 0 0 0 0 0 0 0 0 -1 0 }{PSTYLE "Maple Output" 0 11 1 {CSTYLE " " -1 -1 "" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 }3 3 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }{PSTYLE "" 11 12 1 {CSTYLE "" -1 -1 "" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 }1 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }{PSTYLE "Maple Plot" 0 13 1 {CSTYLE "" -1 -1 "" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 }3 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }{PSTYLE "Title" 0 18 1 {CSTYLE "" -1 -1 "" 1 18 0 0 0 0 0 1 1 0 0 0 0 0 0 }3 0 0 -1 12 12 0 0 0 0 0 0 19 0 }{PSTYLE "R3 Font 0 " -1 256 1 {CSTYLE "" -1 -1 "Helvetica" 0 12 0 0 0 0 2 1 2 0 0 0 0 0 0 }0 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }{PSTYLE "R3 Font 2" -1 257 1 {CSTYLE "" -1 -1 "Courier" 0 10 0 0 0 0 2 2 2 0 0 0 0 0 0 }0 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }{PSTYLE "R3 Font 3" -1 258 1 {CSTYLE "" -1 -1 "Courier" 0 14 0 0 0 0 2 2 2 0 0 0 0 0 0 }0 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }{PSTYLE "R3 Font 4" -1 259 1 {CSTYLE "" -1 -1 "Courier" 0 14 0 0 0 0 2 2 2 0 0 0 0 0 0 }0 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }{PSTYLE "R3 \+ Font 5" -1 260 1 {CSTYLE "" -1 -1 "Courier" 0 14 0 0 0 0 2 2 2 0 0 0 0 0 0 }0 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }{PSTYLE "R3 Font 6" -1 261 1 {CSTYLE "" -1 -1 "Courier" 0 14 0 0 0 0 2 2 2 0 0 0 0 0 0 }0 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }{PSTYLE "R3 Font 7" -1 262 1 {CSTYLE "" -1 -1 "Courier" 0 14 0 0 0 0 2 2 2 0 0 0 0 0 0 }0 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }{PSTYLE "R3 Font 8" -1 263 1 {CSTYLE "" -1 -1 "Courier" 0 14 0 0 0 0 2 2 2 0 0 0 0 0 0 }0 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }{PSTYLE "R3 \+ Font 0" -1 264 1 {CSTYLE "" -1 -1 "" 1 10 255 0 0 1 2 1 2 0 0 0 0 0 0 }0 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }{PSTYLE "R3 Font 2" -1 265 1 {CSTYLE "" -1 -1 "Courier" 1 10 0 0 0 1 2 2 2 0 0 0 0 0 0 }0 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }} {SECT 0 {EXCHG {PARA 18 "" 0 "" {TEXT -1 65 "Reversible, Exothermic, G as Phase reaction in a Catalytic Reactor" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 0 "" }}}{EXCHG {PARA 5 "" 0 "" {TEXT -1 17 "Model Formulation " }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 0 "" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 79 "The component material balances for a tubular reactor may be written as follows" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 8 "re start:" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 35 "CMB := Diff(X,W)= -r[A]/F[A,0]: CMB;" }}{PARA 11 "" 1 "" {XPPMATH 20 "6#/-%%DiffG6$%\"XG %\"WG,$*&&%\"rG6#%\"AG\"\"\"&%\"FG6$F.\"\"!!\"\"F4" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 29 "The reaction rate is given by" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 45 "rateeqn:=r[A]=-k*(C[A]^2-C[C]/K[C]): rateeq n;" }}{PARA 11 "" 1 "" {XPPMATH 20 "6#/&%\"rG6#%\"AG,$*&%\"kG\"\"\",&* $&%\"CGF&\"\"#F+*&&F/6#F/F+&%\"KGF3!\"\"F6F+F6" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 50 "with the reaction rate coefficient calculated from" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 48 "keqn:=k=k[ref]*exp(E[A]/R *(1/T[ref]-1/T)): keqn;" }}{PARA 11 "" 1 "" {XPPMATH 20 "6#/%\"kG*&&F$ 6#%$refG\"\"\"-%$expG6#*(&%\"EG6#%\"AGF)%\"RG!\"\",&*$&%\"TGF'F3F)*$F7 F3F3F)F)" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 48 "The reaction equilibr ium coefficient is given by" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 60 "Keqn:=K[C]=K[C,ref]*exp(Delta(H[R])/R*(1/T[ref]-1/T)): Keqn;" }} {PARA 11 "" 1 "" {XPPMATH 20 "6#/&%\"KG6#%\"CG*&&F%6$F'%$refG\"\"\"-%$ expG6#*(-%&DeltaG6#&%\"HG6#%\"RGF,F7!\"\",&*$&%\"TG6#F+F8F,*$F " 0 "" {MPLTEXT 1 0 55 "CAeqn:=C[A]=C[A,0]*(1-X)/(1+epsilon*X)*y*T[0]/T: CAeqn;" }} {PARA 11 "" 1 "" {XPPMATH 20 "6#/&%\"CG6#%\"AG*.&F%6$F'\"\"!\"\"\",&F, F,%\"XG!\"\"F,,&F,F,*&%(epsilonGF,F.F,F,F/%\"yGF,&%\"TG6#F+F,F5F/" }}} {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 60 "CCeqn:=C[C]=-epsilon*C[A,0]* X/(1+epsilon*X)*y*T[0]/T: CCeqn;" }}{PARA 11 "" 1 "" {XPPMATH 20 "6#/& %\"CG6#F%,$*0%(epsilonG\"\"\"&F%6$%\"AG\"\"!F*%\"XGF*,&F*F**&F)F*F/F*F *!\"\"%\"yGF*&%\"TG6#F.F*F5F2F2" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 21 "The energy balance is" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 69 "EB := Diff(T,W)=(U[a]*(T[a]-T)+r[A]*Delta(H[R]))/(F[A,0]*C[P,A]): \+ EB;" }}{PARA 11 "" 1 "" {XPPMATH 20 "6#/-%%DiffG6$%\"TG%\"WG*(,&*&&%\" UG6#%\"aG\"\"\",&&F'F.F0F'!\"\"F0F0*&&%\"rG6#%\"AGF0-%&DeltaG6#&%\"HG6 #%\"RGF0F0F0&%\"FG6$F8\"\"!F3&%\"CG6$%\"PGF8F3" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 29 "The pressure drop is given by" }}}{EXCHG {PARA 0 "> \+ " 0 "" {MPLTEXT 1 0 56 "dpeqn:=Diff(y,W)=-alpha*(1+epsilon*X)/2/y*T/T[ 0]: dpeqn;" }}{PARA 11 "" 1 "" {XPPMATH 20 "6#/-%%DiffG6$%\"yG%\"WG,$* ,%&alphaG\"\"\",&F,F,*&%(epsilonGF,%\"XGF,F,F,F'!\"\"%\"TGF,&F26#\"\"! F1#F1\"\"#" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 6 "where " }{XPPEDIT 18 0 "y=P/P[0]" "/%\"yG*&%\"PG\"\"\"&F%6#\"\"!!\"\"" }{TEXT -1 23 " is the pressure ratio." }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 31 "The param eters in the model are" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 150 " params := \{C[P,A]=40,C[P,C]=80,Delta(H[R])=-40000,E[A]=41800,k[ref]=0 .5,K[C,ref]=25000,R=8.314,U[a]=0.8,T[a]=500,alpha=0.015,T[ref]=450,eps ilon=-0.5\};" }}{PARA 12 "" 1 "" {XPPMATH 20 "6#>%'paramsG<./&%\"UG6#% \"aG$\"\")!\"\"/&%\"TGF)\"$+&/%&alphaG$\"#:!\"$/&F06#%$refG\"$]%/%(eps ilonG$!\"&F-/&%\"CG6$%\"PG%\"AG\"#S/-%&DeltaG6#&%\"HG6#%\"RG!&++%/&FB6 $FDFB\"#!)/&%\"kGF9$\"\"&F-/&%\"EG6#FE\"&+=%/&%\"KG6$FBF:\"&+]#/FN$\"% 9$)F6" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 0 "" }}}{EXCHG {PARA 5 "" 0 "" {TEXT -1 28 "Method 1 - As an ODE problem" }}}{EXCHG {PARA 0 "" 0 " " {TEXT -1 0 "" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 211 "We are going t o use the BESIRK code to solve this particular problem. BESIRK is an i mplementation of a semi-implicit Runge-Kutta method (see Schwalbe et a l, 1996) and is much faster than the methods built in to " }{TEXT 257 14 "dsolve/numeric" }{TEXT -1 62 " that we used in another example. We read the code into Maple." }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 38 "read `c:/maple/numerics/integ/besirk`:" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 175 "BESIRK needs three arguments. The first is a list of the differential (and algebraic equations), the second is a list of initi al values and the third is the integration range." }}}{EXCHG {PARA 0 " " 0 "" {TEXT -1 0 "" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 42 "Values at \+ the start of the reactor include" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 87 "startvalues := \{F[A,0]=5.0,P[0]=10,C[A,0]=0.271,T[0] =450,y[0]=1,C[C,0]=0\}: startvalues;" }}{PARA 11 "" 1 "" {XPPMATH 20 " 6#<(/&%\"FG6$%\"AG\"\"!$\"#]!\"\"/&%\"TG6#F)\"$]%/&%\"yGF0\"\"\"/&%\"C G6$F8F)F)/&%\"PGF0\"#5/&F8F'$\"$r#!\"$" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 54 "The initial condition is given as a list of equations." } }}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 37 "Start := [C[A]=C[A,0],C[C] =0]: Start;" }}{PARA 11 "" 1 "" {XPPMATH 20 "6#7$/&%\"CG6#%\"AG&F&6$F( \"\"!/&F&6#F&F+" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 31 "DEqns := [CMB,EB,dpeqn]: DEqns;" }}{PARA 11 "" 1 "" {XPPMATH 20 "6#7%/-%%DiffG 6$%\"XG%\"WG,$*&&%\"rG6#%\"AG\"\"\"&%\"FG6$F/\"\"!!\"\"F5/-F&6$%\"TGF) *(,&*&&%\"UG6#%\"aGF0,&&F9F?F0F9F5F0F0*&F,F0-%&DeltaG6#&%\"HG6#%\"RGF0 F0F0F1F5&%\"CG6$%\"PGF/F5/-F&6$%\"yGF),$*,%&alphaGF0,&F0F0*&%(epsilonG F0F(F0F0F0FRF5F9F0&F96#F4F5#F5\"\"#" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 59 "The algebraic equations are collected in a list as follows." }} }{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 46 "AEqns:=[rateeqn,keqn,Keqn,C Aeqn,CCeqn]: AEqns;" }}{PARA 12 "" 1 "" {XPPMATH 20 "6#7'/&%\"rG6#%\"A G,$*&%\"kG\"\"\",&*$&%\"CGF'\"\"#F,*&&F06#F0F,&%\"KGF4!\"\"F7F,F7/F+*& &F+6#%$refGF,-%$expG6#*(&%\"EGF'F,%\"RGF7,&*$&%\"TGF;F7F,*$FGF7F7F,F,/ F5*&&F66$F0F6#*(-%&DeltaG6#&%\"HG6#FCF,FCF7FDF,F,/F/*.&F06$F(\" \"!F,,&F,F,%\"XGF7F,,&F,F,*&%(epsilonGF,FfnF,F,F7%\"yGF,&FG6#FZF,FGF7/ F3,$*0FinF,FXF,FfnF,FgnF7FjnF,F[oF,FGF7F7" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 8 "Maple's " }{TEXT 258 4 "subs" }{TEXT -1 180 " command is q uite limited when it comes to making substitutions involving indexed a nd unindexed variables with the same root name. To deal with this limi tation we have written the " }{TEXT 259 4 "Subs" }{TEXT -1 56 " comman d that is read into the worksheet as part of the " }{TEXT 256 10 "util s.mpl " }{TEXT -1 127 "package. We have chosen to use this utulity pac kage here largely because the workaround in standard Maple is even mor e awkward." }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 32 "read `c:/mapl e/utils/utils.mpl`:" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 65 "Subs (params,Start,startvalues,AEqns,right):\nAIC:=Subs(\",\",right);" }} {PARA 12 "" 1 "" {XPPMATH 20 "6#>%$AICG7'/&%\"rG6#%\"AG,$-%$expG6#,&$ \"+2(es6\"!\")\"\"\"*$%\"TG!\"\"$!+\"=kw-&!\"'$!'0sO!\"(/%\"kG,$F,$\" \"&F6/&%\"KG6#%\"CG,$-F-6#,&$!+)3Z\"p5F2F3F4$\"+'*=;6[F9\"&+]#/&FFF),$ **,&F3F3%\"XGF6F3,&F3F3FU$!\"&F6F6%\"yGF3F5F6$\"']>7!\"$/&FFFE,$**FUF3 FVF6FYF3F5F6$\"'](4'!\"%" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 142 "The \+ next step is to eliminate the algebraic variables and to provide the p arameters with their numerical values in the differential equations." }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 32 "DEqns2:=subs(AEqns,AEqns, DEqns):" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 48 "DEqns3:=subs(par ams,startvalues,DEqns2): DEqns3;" }}{PARA 12 "" 1 "" {XPPMATH 20 "6#7% /-%%DiffG6$%\"XG%\"WG,$*&-%$expG6#,&$\"+2(es6\"!\")\"\"\"*$%\"TG!\"\"$ !+\"=kw-&!\"'F3,&**,&F3F3F(F6\"\"#,&F3F3F($!\"&F6!\"#%\"yGF=F5FA$\"+]- =([\"F@*,F(F3F>F6FBF3F5F6-F-6#,&$!+)3Z\"p5F2F3F4$\"+'*=;6[F9F6$!++++RC !#7F3$\"+++++5!#5/-F&6$F5F),($\"+++++?!\"*F3F5$!+++++SFOF+$FQ!\"(/-F&6 $FBF),$*(F>F3FBF6F5F3$!+nmmm;!#9" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 52 "The initial values of the differential variables are" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 30 "IC := [W=0,X=0,T=450,y=1]: IC;" }} {PARA 11 "" 1 "" {XPPMATH 20 "6#7&/%\"WG\"\"!/%\"XGF&/%\"TG\"$]%/%\"yG \"\"\"" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 46 "We integrate the system of ODEs using BESIRK. " }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 34 " result := BESIRK(DEqns3,IC,0..20);" }}{PARA 12 "" 1 "" {XPPMATH 20 "6# >%'resultG-%&ARRAYG6$7#;\"\"!\"\"\"7$/6#F*7&F*F*$\"$]%F*$F+F*/6#F+7&$ \"#?F*$\"1n)=K!Q(*\\s!#;$\"1$y(=@qj\\6!#7$\"1m!Q\\Bd!owF:" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 81 "The time taken to compute the result can \+ be determined using the global variable " }{TEXT 260 10 "BESIRKtime" } {TEXT -1 37 " which is computed while BESIRK runs." }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 11 "BESIRKtime;" }}{PARA 11 "" 1 "" {XPPMATH 20 "6#$\"%]L!\"$" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 323 "Note that BE SIRK integrated all the way to the end of the reactor in a single step . We will not be able to obtain much useful information about what hap pened along the length of the tube using this result so we shall repea t the integration using an optional argument to BESIRK that allows us \+ to limit the maximum step size." }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 43 "result := BESIRK(DEqns3,IC,0..20,hmax=0.2):" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 38 "The time taken for this computation is" }}} {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 11 "BESIRKtime;" }}{PARA 11 "" 1 "" {XPPMATH 20 "6#$\"%7'*!\"$" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 493 "We do not print the result in this case since to do so requires c onsiderably more space. To plot the results we use the makelist functi on in the BESIRKpackage that is designed to create lists of points to \+ plot from BESIRK integrations. The arguments are the name of the outpu t table and the columns corresponding to the x and y axis variables. H ere we plot the conversion (the second column) as a function of reacto r length (the first column). Other variables can be displayed in simil ar ways." }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 40 "plot(makelist(r esult,1,2),labels=[W,X]);" }}{PARA 13 "" 1 "" {INLPLOT "6$-%'CURVESG6$ 7bq7$\"\"!F(7$$\"1+++++++?!#;$\"1VifqZT)[\"!#=7$$\"1+++++++SF,$\"1!))Q 2X;s,$F/7$$\"1,++++++gF,$\"1TN+'z*G)e%F/7$$\"1+++++++!)F,$\"1^ETMfk.iF /7$$\"\"\"F($\"1(RlORLa'yF/7$$\"1+++++++7!#:$\"1Cs(RJ]fd*F/7$$\"1+++++ ++9FG$\"1mF\")*zlP8\"!#<7$$\"1+++++++;FG$\"1T#zsf>`J\"FO7$$\"1+++++++= FG$\"1)G]h()RD]\"FO7$$\"\"#F($\"1%)[.g9t&p\"FO7$$\"1+++++++AFG$\"1oaPO EA&*=FO7$$\"1+++++++CFG$\"1j.)4]n85#FO7$$\"1+++++++EFG$\"1cKe!f[XJ#FO7 $$\"1+++++++GFG$\"1'HA`%)z^`#FO7$$\"\"$F($\"1ii=]*4Pw#FO7$$\"1,++++++K FG$\"1Pg%G6E1+$FO7$$\"1,++++++MFG$\"1d9LX$ekC$FO7$$\"1,++++++OFG$\"1Y \"H!HVy,NFO7$$\"1,++++++QFG$\"1l])zBv'FO7$$\"1-++++++eFG$\"1RGwJu _rrFO7$$\"1.++++++gFG$\"1S&Gn2!*Rh(FO7$$\"1.++++++iFG$\"1Nu%[dh>3)FO7$ $\"1.++++++kFG$\"1&3tCKAzd)FO7$$\"1.++++++mFG$\"1;x2$R'o/\"*FO7$$\"1.+ +++++oFG$\"1j:?#odam*FO7$$\"1/++++++qFG$\"1$[#Q/(*QE5F,7$$\"1/++++++sF G$\"1Yu#HF?/4\"F,7$$\"1/++++++uFG$\"1kaY)o?\"f6F,7$$\"1/++++++wFG$\"1p ZQL\"\\IB\"F,7$$\"1/++++++yFG$\"1ZDJJE&GJ\"F,7$$\"1/++++++!)FG$\"1Eu%3 h$G*R\"F,7$$\"1.++++++#)FG$\"1(HYSa>K\\\"F,7$$\"1-++++++%)FG$\"1%e.,&f o&f\"F,7$$\"1,++++++')FG$\"1&)39(z#)yq\"F,7$$\"1,++++++))FG$\"176=tP@J =F,7$$\"\"*F($\"1-3]!z;t'>F,7$$\"1*************>*FG$\"1,6Fd[l$\"1R q#f,\"y\"*[F,7$$\"1************R6Fd[l$\"1wg4ze,x_F,7$$\"1************f 6Fd[l$\"1jp,fYhYcF,7$$\"1************z6Fd[l$\"1Eu?^#z])fF,7$$\"1****** *******>\"Fd[l$\"1?$))>(Q\\!G'F,7$$\"1************>7Fd[l$\"1?#p_X*oElF ,7$$\"1************R7Fd[l$\"19\"R,>3Ns'F,7$$\"1************f7Fd[l$\"1u BI_kSvoF,7$$\"1************z7Fd[l$\"1Ql=V@O*)pF,7$$\"1*************H\" Fd[l$\"1)H$R&[/I2(F,7$$\"1************>8Fd[l$\"1-YDy\")RLrF,7$$\"1)*** ********R8Fd[l$\"10&pCovk<(F,7$$\"1)***********f8Fd[l$\"1e1KEGsF,7$$\"1)************R\"Fd[l$\"1kw kQ,7VsF,7$$\"1)***********>9Fd[l$\"12dV1iP`sF,7$$\"1)***********R9Fd[l $\"1rXOarPgsF,7$$\"1)***********f9Fd[l$\"1uT#)=!y]E(F,7$$\"1)********* **z9Fd[l$\"1`;:Fd[l$\"1]nAG>=rsF,7$$\"1)***********R:Fd[l$\"1hrI+YqrsF, 7$$\"1)***********f:Fd[l$\"1\"*Gn3%4=F(F,7$$\"1)***********z:Fd[l$\"1g %y=n7;F(F,7$$\"1)************f\"Fd[l$\"1^N&[U'>rsF,7$$\"1(***********> ;Fd[l$\"1qOI1'=1F(F,7$$\"1(***********R;Fd[l$\"1W:z&>?*psF,7$$\"1(**** *******f;Fd[l$\"1i\\2<-8psF,7$$\"1(***********z;Fd[l$\"1QZc\\#p#osF,7$ $\"1(************p\"Fd[l$\"1<.M'*=NnsF,7$$\"1(***********>isF,7$$\"1(***********>=Fd[l$\"1#e\\n:v5E(F,7$$\"1 (***********R=Fd[l$\"17I\\fH$*fsF,7$$\"1(***********f=Fd[l$\"1&>o/Mn(e sF,7$$\"1(***********z=Fd[l$\"1Saw0)yvD(F,7$$\"1'*************=Fd[l$\" 1sBKWwOcsF,7$$\"1'***********>>Fd[l$\"1C%=g)R8bsF,7$$\"1'***********R> Fd[l$\"16$)\\Vy(QD(F,7$$\"1'***********f>Fd[l$\"1f_LU\"*f_sF,7$$\"1'** *********z>Fd[l$\"1)f>Gu(H^sF,7$$\"1'*************>Fd[l$\"1M\">dXt*\\s F,7$$\"1'4++++++#Fd[l$\"1G\">dXt*\\sF,-%'COLOURG6&%$RGBG$F^[l!\"\"F(F( -%+AXESLABELSG6$%\"WG%\"XG" 2 293 180 180 2 0 1 0 2 9 0 4 2 1.000000 45.000000 45.000000 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 1 0 0 0 -30718 1 0 0 0 0 0 0 }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 0 "" }}} {EXCHG {PARA 5 "" 0 "" {TEXT -1 27 "Method 2 - As a DAE problem" }}} {EXCHG {PARA 0 "" 0 "" {TEXT -1 0 "" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 778 "In the approach demonstrated above the original system of diff erential and algebraic equations was converted to a system of differen tial equations only. This was easy to do in this case because the vari ables on the right hand sides of the ODEs could be eliminated using th e algebraic equations. The advantage of doing so is that the numerical solution can be computed quite rapidly. A minor disadvantage is that \+ we no longer have values for the other variables to hand and we need t o compute them using the results of the integration. An alternative ap proach that provides values for all of the variables at every step of \+ the integration is to solve the system of model equations as the DAE p roblem it really is. BESIRK is quite capable of solving DAE problems a s we show below." }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 32 "DAEqns \+ := [op(DEqns),op(AEqns)];" }}{PARA 12 "" 1 "" {XPPMATH 20 "6#>%'DAEqns G7*/-%%DiffG6$%\"XG%\"WG,$*&&%\"rG6#%\"AG\"\"\"&%\"FG6$F1\"\"!!\"\"F7/ -F(6$%\"TGF+*(,&*&&%\"UG6#%\"aGF2,&&F;FAF2F;F7F2F2*&F.F2-%&DeltaG6#&% \"HG6#%\"RGF2F2F2F3F7&%\"CG6$%\"PGF1F7/-F(6$%\"yGF+,$*,%&alphaGF2,&F2F 2*&%(epsilonGF2F*F2F2F2FTF7F;F2&F;6#F6F7#F7\"\"#/F.,$*&%\"kGF2,&*$&FNF 0FhnF2*&&FN6#FNF2&%\"KGFboF7F7F2F7/F\\o*&&F\\o6#%$refGF2-%$expG6#*(&% \"EGF0F2FLF7,&*$&F;FhoF7F2*$F;F7F7F2F2/Fco*&&Fdo6$FNFioF2-F[p6#*(FFF2F LF7F`pF2F2/F_o*.&FNF5F2,&F2F2F*F7F2FXF7FTF2FenF2F;F7/Fao,$*0FZF2F]qF2F *F2FXF7FTF2FenF2F;F7F7" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 56 "The par ameters are replaced with their numerical values." }}}{EXCHG {PARA 0 " > " 0 "" {MPLTEXT 1 0 41 "DAEqns2:=Subs(params,startvalues,DAEqns);" } }{PARA 12 "" 1 "" {XPPMATH 20 "6#>%(DAEqns2G7*/-%%DiffG6$%\"XG%\"WG,$& %\"rG6#%\"AG$!+++++?!#5/-F(6$%\"TGF+,($\"+++++?!\"*\"\"\"F7$!+++++S!#7 F-$F2!\"(/-F(6$%\"yGF+,$*(,&F7!\"$/ FY,$**F*F " 0 "" {MPLTEXT 1 0 53 "temp:=[op(IC),op(AIC)]:\nDAEIC:=Subs(temp,temp,right) ;" }}{PARA 11 "" 1 "" {XPPMATH 20 "6#>%&DAEICG7+/%\"WG\"\"!/%\"XGF(/% \"TG\"$]%/%\"yG\"\"\"/&%\"rG6#%\"AG$!'0sO!\"(/%\"kG$\"\"&!\"\"/&%\"KG6 #%\"CG\"&+]#/&FBF4$\"++++5F!#5/&FBFAF(" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 56 "We carry out essentially the same integration as before. " }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 39 "DAEresult:=BESIRK(DAEqn s2,DAEIC,0..20):" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 102 "The time nee ded for the integration is much longer than before (but this is a more difficult problem)." }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 11 "BES IRKtime;" }}{PARA 11 "" 1 "" {XPPMATH 20 "6#$\"''\\;\"!\"$" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 131 "The concentrations are variable numbers \+ 8 and 9 in each line of the output table; we can plot the concentratio n profiles as follows" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 58 "pl ot(\{seq(makelist(DAEresult,1,k),k=[8,9])\},labels=[W,C]);" }}{PARA 13 "" 1 "" {INLPLOT "6%-%'CURVESG6$7?7$\"\"!$\"1++++++5F!#;7$$\"1+++++ +DJF+$\"1833O\\5'o#F+7$$\"1++++++]iF+$\"1))Q2a%=$ \"1P$=G$zDoDF+7$$\"1+++]PM?AF>$\"1k)zLM=$GDF+7$$\"1+++D\"GOs#F>$\"1\"e 2i$fc#[#F+7$$\"1++]P4CxKF>$\"1Cm^S%f(HCF+7$$\"1++voC\\vPF>$\"1tISN&3(z BF+7$$\"1,+++SutUF>$\"1Xq\"oH=pK#F+7$$\"1,]P%o?=#[F>$\"1,!pZ/'>lAF+7$$ \"1,DJ+(*3:`F>$\"1)*zAUo$e?#F+7$$\"1,+D;(e$3eF>$\"1f>c#QQA9#F+7$$\"1^7 yLY&4N'F>$\"1fO7:jYm?F+7$$\"1w$f&f4HRoF>$\"1$o%RlM-#*>F+7$$\"1,vL&GFwK (F>$\"1c1$Qa,/\">F+7$$\"1RMpVszkyF>$\"18z?.EW5=F+7$$\"1q\\)y>%ob%)F>$ \"1E#zel)*\\o\"F+7$$\"1[tlcC[()*)F>$\"1mQEUc[a:F+7$$\"1E(Har!G>&*F>$\" 1o`cRL9.9F+7$$\"1Mz+!)eU55!#9$\"1I^`xRq37F+7$$\"1#R>$pMxu5F^r$\"1<8TX* 4qk*!#<7$$\"1'*>cI6e4C\"F^r$\"1+&4\" GyfJZFfr7$$\"1#eA!Q1!fM\"F^r$\"1%>,\"y[_$*RFfr7$$\"18bk(>[hf\"F^r$\"1, J^Y]oGPFfr7$$\"1P2B80Ur=F^r$\"1&)*y(f/;!e$Ffr7$$\"#?F($\"15d/$fH)3NFfr -%'COLOURG6&%$RGBG$\"#5!\"\"F(F(-F$6$7?7$F(F(7$F-$\"1jw8Aq`aJ!#>7$F2$ \"1B)Q:FnCS'Fcu7$F7$\"1IWLVs/45!#=7$F<$\"1*)HnUztH9Fju7$FB$\"1IPRHDC7> Fju7$FG$\"1uIr:tGpCFju7$FL$\"1sHyH>bnI,[1aFju7$Fjn$\"1szz-b,^jFju7$F_o$\"1l$e[*HuzsF ju7$Fdo$\"1'4eoL'R(H)Fju7$Fio$\"1UGDI11T&*Fju7$F^p$\"1[q8SWoz5Ffr7$Fcp $\"1X!=%*)RO@7Ffr7$Fhp$\"1z$yc!oz+9Ffr7$F]q$\"1V'=dfpdj\"Ffr7$Fbq$\"17 DcVS$G*=Ffr7$Fgq$\"1ZE#\\ly!4AFfr7$F\\r$\"1l>(=w]%\\EFfr7$Fbr$\"1*)3`' )39sKFfr7$Fhr$\"1&e00(R-tSFfr7$F]s$\"1\"y6*QZ&H([Ffr7$Fbs$\"1\"3l6wW05 &Ffr7$Fgs$\"1Iy66;)z'\\Ffr7$F\\t$\"1SL\"pJ7&QZFfr7$Fat$\"1m0$*)49_i%Ff r-Fft6&FhtF(FitF(-%+AXESLABELSG6$%\"WG%\"CG" 2 304 202 202 2 0 1 0 2 9 0 4 2 1.000000 45.000000 45.000000 10030 10061 10056 10074 0 0 0 20030 0 12020 0 0 0 0 0 0 0 1 1 0 0 0 0 16 0 0 0 0 0 0 }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 89 "The higher curve on the left corresponds \+ to the concentration of A, the lower curve to C." }}}{EXCHG {PARA 0 " " 0 "" {TEXT -1 0 "" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 71 "The temper ature is variable number 3 in the output list; its profile is" }}} {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 60 "plot(\{seq(makelist(DAEresul t,1,k),k=[3])\},labels=[`W`,`T`]);" }}{PARA 13 "" 1 "" {INLPLOT "6$-%' CURVESG6$7?7$\"\"!$\"$]%F(7$$\"1++++++DJ!#;$\"1rI;aG/CX!#87$$\"1++++++ ]iF.$\"1-$oJ$p1\\XF17$$\"1+++++](o*F.$\"1/-Jwh\"yd%F17$$\"1++++](oM\"! #:$\"1*[?u?\\5h%F17$$\"1++++D\"Gw\"F?$\"1tm_4uu\\YF17$$\"1+++]PM?AF?$ \"1d2vAo?&p%F17$$\"1+++D\"GOs#F?$\"1$3`?(4=\\ZF17$$\"1++]P4CxKF?$\"1F) 3FV6T\"[F17$$\"1++voC\\vPF?$\"15\")ea[Qy[F17$$\"1,+++SutUF?$\"1h$ow=*> \\\\F17$$\"1,]P%o?=#[F?$\"1>%yl9Nh.&F17$$\"1,DJ+(*3:`F?$\"1D\"4Vt!=C^F 17$$\"1,+D;(e$3eF?$\"1jvC0\"=OA&F17$$\"1^7yLY&4N'F?$\"17D(y?-%\\`F17$$ \"1w$f&f4HRoF?$\"1BY\\%fI7[&F17$$\"1,vL&GFwK(F?$\"1ma\"G+Rej&F17$$\"1R MpVszkyF?$\"1>%f!y>mSeF17$$\"1q\\)y>%ob%)F?$\"1cR_r\\:ChF17$$\"1[tlcC[ ()*)F?$\"1@')zVn#RX'F17$$\"1E(Har!G>&*F?$\"1$*=Zy_s')oF17$$\"1Mz+!)eU5 5!#9$\"1e>i(=%oMvF17$$\"1#R>$pMxu5F_r$\"1]([$=4+@&)F17$$\"1'*>cI6e4C\"F_r$\"1N[*)\\3x<6!#77$$\"1#eA!Q1!fM\" F_r$\"1+;#G)eag6Fas7$$\"18bk(>[hf\"F_r$\"1=5L6wOi6Fas7$$\"1P2B80Ur=F_r $\"1n;>%*3!Q:\"Fas7$$\"#?F($\"1bncbrj\\6Fas-%'COLOURG6&%$RGBG$\"#5!\" \"F(F(-%+AXESLABELSG6$%\"WG%\"TG" 2 313 199 199 2 0 1 0 2 9 0 4 2 1.000000 45.000000 45.000000 10030 10061 10056 10074 0 0 0 20030 0 12020 0 0 0 0 0 0 0 1 1 0 0 0 -16400 0 0 0 0 0 0 0 }}}{EXCHG {PARA 0 " " 0 "" {TEXT -1 28 "The reaction rate profile is" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 47 "plot(makelist(DAEresult,1,5),labels=[`W`,`r`] );" }}{PARA 13 "" 1 "" {INLPLOT "6$-%'CURVESG6$7?7$\"\"!$!1+++++0sO!#< 7$$\"1++++++DJ!#;$!1,E@85FGQF+7$$\"1++++++]iF/$!1hHWmn/'*RF+7$$\"1++++ +](o*F/$!1CpK;GV&>%F+7$$\"1++++](oM\"!#:$!1*3*HM=#\\V%F+7$$\"1++++D\"G w\"F?$!1$)p+$4>is%F+7$$\"1+++]PM?AF?$!1*ya6?!)e3&F+7$$\"1+++D\"GOs#F?$ !1?3nNv.QbF+7$$\"1++]P4CxKF?$!1cr'f(49>hF+7$$\"1++voC\\vPF?$!1EbuoeYNn F+7$$\"1,+++SutUF?$!1@R&)4X:juF+7$$\"1,]P%o?=#[F?$!1p;#)\\s1G%)F+7$$\" 1,DJ+(*3:`F?$!1&*y.!H7w[*F+7$$\"1,+D;(e$3eF?$!1'ya3v\"fy5F/7$$\"1^7yLY &4N'F?$!1P&=\\^:&e7F/7$$\"1w$f&f4HRoF?$!1%)oe&3_gY\"F/7$$\"1,vL&GFwK(F ?$!1)HU#Rn4M%ob%)F?$!1q,bc dxYFF/7$$\"1[tlcC[()*)F?$!1ft&)f]maNF/7$$\"1E(Har!G>&*F?$!1!)fq!*>j?ZF /7$$\"1Mz+!)eU55!#9$!1Jcu;w^RlF/7$$\"1#R>$pMxu5F_r$!1e-#y(46p))F/7$$\" 1'*>cI6e4C\"F_r$!1-XL1>VzUF/7$$\"1#eA !Q1!fM\"F_r$!1\"eY%f_)z8)F+7$$\"18bk(>[hf\"F_r$\"1e%*p#pVM=\"!#=7$$\"1 P2B80Ur=F_r$\"1-?M_sdvHF[t7$$\"#?F($\"1w&4-Q+'RLF[t-%'COLOURG6&%$RGBG$ \"#5!\"\"F(F(-%+AXESLABELSG6$%\"WG%\"rG" 2 336 213 213 2 0 1 0 2 9 0 4 2 1.000000 45.000000 45.000000 10030 10061 10056 10074 0 0 0 20030 0 12020 0 0 0 0 0 0 0 1 1 0 0 0 -126 10 0 0 0 0 0 0 }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}}{MARK "0 0 0" 49 }{VIEWOPTS 1 1 0 1 1 1803 }