Unformatted text preview: ChBE 2120 Fall 2010 Homework 6 ‐ Answer Part1 (1) function fVal = hw6CstrFo_solutions(T) F = 0.1; V = 1; Ca0 = 1; Ea = 50000; deltaH = ‐5420000; R = 8.314; k0 = 1; rho = 1000; cp = 4184; T0 = 298; Ca = F*Ca0/(F + k0*exp(‐Ea/R/T)*V); deltaH*k0*exp(‐Ea/R/T)*V*Ca; fVal = F*rho*cp*(T0 ‐ T) + deltaH*k0*exp(‐Ea/R/T)*V*Ca; (2), (3) function root = secant_solutions(rootFun,x_0,x_neg1) tol = 1e‐5; guessCurr = x_0; guessPrev = x_neg1; guessCurrFval = rootFun(guessCurr); guessPrevFval = rootFun(guessPrev); while(abs(guessCurrFval ‐ guessPrevFval) > tol) newGuess = guessCurr ‐ (guessPrev ‐ guessCurr)/ ... (guessPrevFval ‐ guessCurrFval)*guessCurrFval; newGuessFval = rootFun(newGuess); guessPrev = guessCurr; guessPrevFval = guessCurrFval; guessCurr = newGuess; guessCurrFval = newGuessFval; end root = guessCurr; (4) % Solutions % ChBE 2120 % Homework 6 % Part 1, question 3 % A simple script to execute our secant method clear all xRange = [200:.1:1000]; for i=1:length(xRange) yRange(i) = hw6CstrFo_solutions(xRange(i)); end plot(xRange,yRange) title('Preliminary analysis') xlabel('Tempearture') ylabel('Value of the function, = 0 for correct temperature') % So we see here that our function definitely crosses zero sometime around % 300 K, so we'll arbitrarily choose 320 and 325 as our guesses x_0 = 325; x_neg1 = 320; secantResult = secant_solutions(@hw6CstrFo_solutions,x_0,x_neg1) Part2 (1) function fVal = hw6CstrSo_solutions(x) F = 1; V = 1; Ca0 = 10; Cb0 = 5; EaR = 273; deltaH = 1; k0 = 1; rho = 1; cp = 1; T0 = 273; fVal = zeros(4,1); Ca = x(1); Cb = x(2); Cc = x(3); T = x(4); fVal(1) = F*Ca0 ‐ F*Ca ‐ 2*k0*exp(‐EaR/T)*V*Ca^2*Cb; fVal(2) = F*Cb0 ‐ F*Cb ‐ k0*exp(‐EaR/T)*V*Ca^2*Cb; fVal(3) = ‐F*Cc + k0*exp(‐EaR/T)*V*Ca^2*Cb; fVal(4) = F*rho*cp*(T0 ‐ T) + deltaH*k0*exp(‐EaR/T)*V*Ca^2*Cb; function fVal = hw6CstrJac_solutions(x) F = 1; V = 1; Ca0 = 10; Cb0 = 5; EaR = 273; deltaH = 1; k0 = 1; rho = 1; cp = 1; T0 = 273; fVal = zeros(4,4); Ca = x(1); Cb = x(2); Cc = x(3); T = x(4); fVal(1,1) = ‐F ‐ 2*k0*exp(‐EaR/T)*V*2*Ca*Cb; fVal(1,2) = ‐2*k0*exp(‐EaR/T)*V*Ca^2; fVal(1,3) = 0; fVal(1,4) = ‐2*k0*exp(‐EaR/T)*V*Ca^2*Cb*(EaR/T^2); fVal(2,1) = ‐2*k0*exp(‐EaR/T)*V*Ca*Cb; fVal(2,2) = ‐F ‐ k0*exp(‐EaR/T)*V*Ca^2; fVal(2,3) = 0; fVal(2,4) = ‐k0*exp(‐EaR/T)*V*Ca^2*Cb*(EaR/T^2); fVal(3,1) = 2*k0*exp(‐EaR/T)*V*Ca*Cb; fVal(3,2) = k0*exp(‐EaR/T)*V*Ca^2; fVal(3,3) = ‐F; fVal(3,4) = k0*exp(‐EaR/T)*V*Ca^2*Cb*(EaR/T^2); fVal(4,1) = 2*k0*deltaH*exp(‐EaR/T)*V*Ca*Cb; fVal(4,2) = k0*deltaH*exp(‐EaR/T)*V*Ca^2; fVal(4,3) = 0; fVal(4,4) = ‐F*rho*cp + k0*deltaH*exp(‐EaR/T)*V*Ca^2*Cb*(EaR/T^2); (2) function root = newtonRaphson_solutions(rootFun,jacFun,initGuess) tol = 1e‐10; guess = initGuess; guessFval = rootFun(guess); while(max(abs(guessFval) > tol)) jacVal = jacFun(guess); newGuess = jacVal\(jacVal*guess ‐ guessFval); newGuessFval = rootFun(newGuess); guess = newGuess; guessFval = newGuessFval; end root = guess; (3) % Solutions % ChBE 2120 % Homework 6 % Part 2, question 3 % A simple script to execute our newton‐raphson method clear all initialGuess = [10; 5; 0; 500]; newtonRaphsonResult = newtonRaphson_solutions(@hw6CstrTo_solutions,... @hw6CstrJac_solutions,initialGuess) Part3 We are trying to get the outlet temperature of the reactor to be 300K. We are able to vary the inlet feed temperature. Varying the inlet feed temperature has a direct impact on the outlet temperature. So, we could say that T = f(T0). But we want T to equal 300. So, in order to make this a problem of root‐ finding (finding zeros), we should instead say that f(T0) = T ‐ 300. Then, the inlet temperature that makes that outlet be 300 will return a value of zero for our function. It is now a root‐solving problem. We can use a secant approach to solve this problem. We start with two values of T0 that are fairly close to each other. We then perform the secant method in order to find the temperature that causes the outlet temperature to be closest to 300. Of course, this leaves the question: how do we find T for these secant steps? We are given a T0. So, given a T0, we would call our Newton‐Raphson function to find the steady‐state temperature and concentration for that T0. The steady‐state temperature is then T, which is what we need for the function defined above. So, we are nesting the Newton‐Raphson method inside an execution of the bisection method ‐‐‐ for every time that we bisect the interval, we are doing a Newton‐Raphson solution of the problem. ALTERNATIVE: you could also take your 4 equations (the energy balance and the balance on components A, B, and C), and add a fifth equation: T ‐ 300. That way, you are requiring that another equation be made to equal zero, which requires that T = 300. Then you would just execute Newton‐Raphson on the five‐equation system, and you'd get a result back that would give you concentrations, the correct T0, and a T that equals 300. ANOTHER ALTERNATIVE: you could instead change your equations such that 300 is the established outlet temperature, and T0 and the concentrations are the only variables. This would mean one of two things: you could solve the four‐equation system for energy and component balances using a multi‐ dimensional N‐R approach, or you could notice that equation 4 would allow T0 to be defined only in terms of cA and cB, and equation 2 can define cB only in terms of cA, and so you can substitute those into equation 1 and solve it to find the Ca that makes equation 1 equal 0 when T is 300. This could be done using secant or Newton‐Raphson, and then you identify the T0 required based on equation 4. ...
View
Full Document
 Spring '07
 Gallivan
 function fVal

Click to edit the document details