1 CHAPTER 28 28 1 The solution with the 2nd order RK Heun without corrector can be laid out as t k1 c 0 10 20 30 40 50 10 25 34 375 40 23438 43 89648 46 1853 k2 c 2 1 25 0 78125 0 488281 0 305176 0 190735 30 37 5 42 1875 45 11719 46 94824 48 09265 1 0 625 0 390625 0 244141 0 152588 0 095367 1 5 0 9375 0 585938 0 366211 0 228882 0 143051 For the 4th order RK the solution is t c 0 10 20 30 40 50 10 25 72917 35 27317 41 06419 44 57801 46 71009 k1 2 1 213542 0 736342 0 446791 0 2711 0 164495 k2 cmid 20 31 79688 38 95487 43 29814 45 93351 47 53257 1 5 0 910156 0 552256 0 335093 0 203325 0 123371 cmid 17 5 30 27995 38 03445 42 73965 45 59463 47 32695 k3 1 625 0 986003 0 598278 0 363017 0 220268 0 133652 cend 26 25 35 58919 41 25594 44 69436 46 78069 48 04662 k4 1 1875 0 72054 0 437203 0 265282 0 160965 0 097669 A plot of both solutions along with the analytical result is displayed below 50 Analytical RK 4 RK 2 25 0 0 20 40 28 2 The mass balance equations can be written as dc1 dt dc2 dt dc3 dt dc4 dt dc5 dt 0 16c1 0 06c3 1 0 2c1 0 2c2 0 05c2 0 25c3 4 0 0875c3 0 125c4 0 0375c5 0 04c1 0 02c2 0 06c5 Selected solution results Euler s method are displayed below along with a plot of the results t 0 1 2 c1 0 0000 1 0000 2 0800 c2 0 0000 0 0000 0 2000 c3 0 0000 4 0000 7 0000 c4 0 0000 0 0000 0 3500 c5 0 0000 0 0000 0 0400 dc1 dt dc2 dt 1 0000 0 0000 1 0800 0 2000 1 0872 0 3760 dc3 dt 4 0000 3 0000 2 2600 dc4 dt dc5 dt 0 0000 0 0000 0 3500 0 0400 0 5703 0 0848 PROPRIETARY MATERIAL The McGraw Hill Companies Inc All rights reserved No part of this Manual may be displayed reproduced or distributed in any form or by any means without the prior written permission of the publisher or used beyond the limited distribution to teachers and educators permitted by McGraw Hill for their individual course preparation If you are a student using this Manual you are using it without permission 1 572917 0 9544 0 579102 0 351382 0 213208 0 129369 2 3 4 5 76 77 78 79 80 3 1672 4 2160 5 1999 0 5760 1 0942 1 7186 9 2600 10 9738 12 2851 0 9203 1 6201 2 3874 0 1248 0 2555 0 4307 1 0488 0 9839 0 9051 0 5182 0 6244 0 6963 1 7138 1 3113 1 0147 0 6999 0 7673 0 7927 0 1307 0 1752 0 2165 13 2416 13 2418 13 2419 13 2421 13 2422 13 2395 13 2399 13 2403 13 2406 13 2409 18 6473 18 6475 18 6476 18 6477 18 6478 16 8515 16 8621 16 8720 16 8813 16 8901 12 9430 12 9609 12 9777 12 9935 13 0084 0 0002 0 0002 0 0001 0 0001 0 0001 0 0004 0 0004 0 0003 0 0003 0 0003 0 0001 0 0001 0 0001 0 0001 0 0001 0 0106 0 0099 0 0093 0 0088 0 0082 0 0179 0 0168 0 0158 0 0149 0 0140 20 15 10 5 0 0 10 20 30 c1 40 c2 50 c3 c4 60 70 80 c5 Finally MATLAB can be used to determine the eigenvalues and eigenvectors a 16 06 0 0 0 2 2 0 0 0 0 05 25 0 0 0 0 0875 125 0375 04 02 0 0 06 a 0 1600 0 0600 0 0 0 0 2000 0 2000 0 0 0 0 0 0500 0 2500 0 0 0 0 0 0875 0 1250 0 0375 0 0400 0 0200 0 0 0 0600 v d eig a v 0 0 0 0 0 0 8192 1 0000 0 5735 0 0 d 0 1250 0 0 0 2500 0 0 0 0 0 0 0 0 0 0 4997 0 8662 0 1039 0 1582 0 0436 0 4956 0 8466 0 2604 0 5701 0 6892 0 3635 0 0043 0 0 0 0600 0 0 0 0 0 0 0686 0 0 0 0 0 0 2914 28 3 Substituting the parameters into the differential equation gives dc 14 583333 0 0833333c 0 15c 2 dt The mid point method can be applied with the result PROPRIETARY MATERIAL The McGraw Hill Companies Inc All rights reserved No part of this Manual may be displayed reproduced or distributed in any form or by any means without the prior written permission of the publisher or used beyond the limited distribution to teachers and educators permitted by McGraw Hill for their individual course preparation If you are a student using this Manual you are using it without permission 3 The results are approaching a steady state value of 9 586267 Challenge question The steady state form i e dc dt 0 of the equation is 0 14 58333 0 08333c 0 15c2 which can be solved for 9 586267 and 10 1418 Thus there is a negative root If we put in the initial y value as 10 1418 or higher precision the solution will stay at the negative root for awhile until roundoff errors may lead to the solution going unstable If we use an initial value less than the negative root the solution will also go unstable However if we pick a value that is slightly higher as per machine precision it will gravitate towards the positive root For example if we use 10 14 15 10 5 0 5 0 10 5 10 15 20 15 This kind of behavior will also occur for higher initial conditions For example using an initial condition of 16 gives 20 15 10 5 0 0 5 10 15 20 However if we start to use even higher initial conditions the solution will again go unstable For example if we use an initial condition of 17 the result is PROPRIETARY MATERIAL The McGraw Hill Companies Inc All rights reserved No part of this Manual may be displayed reproduced or distributed in any form or by any means without the prior written permission of the publisher or used beyond the limited distribution to teachers and educators permitted by McGraw Hill for their individual course preparation If you are a student using this Manual you are using it without permission 4 50 40 30 20 10 0 0 5 10 15 20 Therefore it looks like initial conditions roughly in the range from 10 14 up to about 16 5 will yield stable solutions that converge on the steady state solution of 9 586 Note that …
View Full Document
Unlocking...