Transcription of cc - CaltechAUTHORS
1 Ph~:1s: 1984,byJohnWiley&Sons, ,JohnWiley& :Davis, '8'072483-21590 PrintedintheUnitedStatesofAmerica109 8 7 6 5 4 3 2 mathematicalprojectistypicallynotinteres tedinsophisticatedtheoreticaltreatments, butratherinthesolutionofa , ,a "black-box" , ,althoughthesubjectscoveredinthisbookare thesameasinothertexts,thetreatmentisdiff erentinthatit computationalpropertyofinterestoraretheu nderlyingmethodsofa basicknowledgeofmathematics, ,whenpossible, :Class1 :Problemsthatshouldbesolvedwithsoftwareo fthetypedescribedinthetext(designatedbya nasteriskaftertheproblemnumber).Thelevel ofthisbookisintroductory, textfora havesuccessfullyusedthismaterialfora ,Chapter1 to3, gratefullyacknowledgethefollowingindivid ualswhohaveeitherdirectlyorindirectlycon tributedtothisbook:KennethDenison,JulioD iaz,PeterMer-cure,KathleenRichter,PeterR ony,LayneWatson, alsothanktheDepartmentofChemicalEngineer ingatVirginiaPolytechnicInstituteandStat eUniversityforitssupport,andI ,andmostofall, dedicatethisbooktomywife,whouncomplainin glygaveupa.
2 DiscreteVariableMethodsIntroductionBackg roundInitial-ValueMethodsShootingMethods MultipleShootingSuperpositionfiniteDiffe renceMethodsLinearSecond-OrderEquationsF luxBoundaryConditionsIntegrationMethodNo nlinearSecond-OrderEquationsFirst-OrderS ystemsHigher-OrderMethodsMathematicalSof twareProblemsReferencesBibliographyBound ary-ValueProblemsforOrdinaryDifferential Equations:finiteElementMethodsIntroducti onBackgroundPiecewisePolynomialFunctions TheGalerkinMethodNonlinearEquationsInhom ogeneousDirichletandFluxBoundaryConditio nsMathematicalSoftwareCollocationMathema ticalSoftwareProblemsReferencesBibliogra phyContents53535354546365676871757983858 7919395979797991041091101111121191231251 26 ContentsChapter4 Chapter5 ParabolicPartialDifferentialEquationsinO neSpaceVariableIntroductionClassificatio nofPartialDifferentialEquationsMethodofL inesfiniteDifferencesLow-OrderTimeApprox imationsTheThetaMethodBoundaryandInitial ConditionsNonlinearEquationsInhomogeneou sMediaHigh-OrderTimeApproximationsFinite ElementsGalerkinCollocationMathematicalS oftwareProblemsReferencesBibliographyPar tialDifferentialEquationsinTwoSpaceVaria bles127127127128130130133135140142147154 154158162167172174177 Introduction177 EllipticPDEs-FiniteDifferences177 Background177 Laplace'sEquationina Square178 DirichletProblem178 NeumannProblem179 RobinProblem180 VariableCoefficientsandNonlinearProblems 184 NonuniformGrids185 IrregularBoundaries190 DirichletCondition190 NormalDerivativeConditions191
3 EllipticPDEs-finiteElements192 Background192 Collocation194 Gakr~n200xliAppendicesParabolicPDEsinTwo SpaceVariablesMethodofLinesAlternatingDi rectionImplicitMethodsMathematicalSoftwa reParabolicsEllipticsProblemsReferencesB ibliographyContents211211212214215219222 224227A:ComputerArithmeticandErrorContro l229 ComputerNumberSystem229 NormalizedFloatingPointNumberSystem230 Round-OffErrors230B:Newton'sMethod235C:G aussianElimination237 DenseMatrix237 TridiagonalMatrix241D:B-Splines243 , ,diffusion-reaction,mass-heattransfer, computationalpropertyofinterestoraretheu nderlyingmethodsofa (IVP)orboundary-valueproblems(BVP).Examp lesofthetwotypesare:IVP:y"=-yx( )yeO)=2,y'(O)1( )BVP:y"=-yx( )yeO)=2,y(l)=1( ) [Eqs.(LIb)and( )] ,theconditionsaregivenatthesamevalueofx, whereasinthecaseoftheBVP, , turnsoutthatthenumericalmethodsforeachty peofproblem,IVPorBVP, ,leavingBVPstoChapters2 (m)=f(x,Y,y',y".
4 ,y(m-1 )( )withinitialconditionsy(xo)=Yoy'(XO)=yby (m-1)(XO)=Y6m-1)wherefisa knownfunctionandYo,yb,..,Y6m-1) ( ) ,wedefinea newsetofdependentvariablesY1(X),Yz(x),.. ,Ym(x)byYl=YYz=Y'Y3=y"Ym=y(m-1)( )andtransform( )intoY;'=Yzy~=Y3=f1(X,Yl,Yz,,Ym)=fz(x,Y1 'Yz,,Ym)( )y:r,=f(x,YvYz,..,Ym)=fm(x,YvYz,..,Ym)wi thY1(XO)=YoyzCxo)=ybExplicitMethodsInvec tornotation( )becomesy'(x)=rex,y)y(xo)=Yowhere3( )[heX)]y(x)=Y2~X),Ym(X)[fleX,y)]rex,y)=f 2(X;y),fm(x,y)[Yo]Y_ybo ~m'-l)Itiseasytoseethat( )canrepresenteitheranmth-orderdifferenti alequation,a systemofequationsofmixedorderbutwithtota lorderofm, ,subroutinesforsolvingIVPsas-sumethatthe problemisintheform( ).Inordertosimplifytheanalysis,webeginby examininga singlefirst-orderIVP,afterwhichweextendt hediscussiontoincludesystemsoftheform( ).
5 Considertheinitial-valueproblemy'=f(x,y) ,Y(Xo)=Yo( )Weassumethataflayiscontinuousonthestrip Xo~x~XN'thusguaranteeingthat( )possessesa uniquesolution[1].Ify(x)istheexactsoluti onto( ),itsgraphisa curveinthexy-planepassingthroughthepoint (xo,Yo).Adiscretenumericalsolutionof( )isdefinedtobea setofpoints[(Xi'u;)]~o,whereUo=Yoandeach point(Xi'u;)isanapproximationtothecorres pondingpoint(Xi'Y(Xi)) setofpoints, numericalsolution[(Xi'Ui)]~O'EXPLICITMET HODSW eagainconsider( )asthemodeldifferentialequationandbeginb ydividingtheinterval[xo,XN]intoNequallys pacedsubintervalssuchthatXi=Xo+ih,i=0,1, 2,..,N( )Theparameterhiscalledthestep-sizeanddoe snotnecessarilyhavetobeuniformovertheint erval.(Variablestep-sizesareconsideredla ter.)4 Initial-ValueProblemsforOrdinaryDifferen tialEquationsIfy(x)istheexactsolutionof( ),thenbyexpandingy(x)aboutthepointXiusin gTaylor'stheoremwithremainderweobtain:Y( Xi+1)=y(xi)+(Xi+1-xi)y'(Xi)+(Xi+12~xyy"( ~J,Thesubstitutionof( )into( )gives( )( )Thesimplestnumericalmethodisobtainedbyt runcating( ) (xJ,Ui+1=Ui+hf(xi,uJ,Uo=Yoi=0,1.))))
6 ,N-1,( )Thismethodis ,wefindthattheapplicationof( )tocomputeUi+1createsanerrorinthevalueof Ui+ ,ei+ ,z(x),byZ'(X)=f(x,z),Z(Xi)=Ui( )Anexpressionforthelocaltruncationerror, ei+1=Z(Xi+1)Ui+1'canbeobtainedbycomparin gtheformulaforUi+1withtheTaylor' (xi+h)=z(xi)+hf(Xi'z(xi +~~z"(~JorZ(Xi+h)=Ui+hf(Xi'uJ+~~Z"(~i)'i tfollowsthat( )ei+1=~~Z"(~i)=0(h2)( )ThenotationO() denotestermsoforder(), ,f(h)=O(hL)ifIf(h)1~Ahlash~0,whereAandIa reconstants[1].Theglobalerrorisdefinedas ( )andisthusthedifferencebetweenthetruesol utionandthenumericalsolutionatX=Xi+ +1andc;i+ +1andc;i+ methodispth-orderaccurateifei+l=0(hP+1)5 ( )andfrom( )and( )theEulermethodis explicitsincethefunctionfis evaluatedwithknowninformation( ,attheleft-handsideofthesubinterval). ( ).Topartiallyanswerthisquestion,weconsid erExample1, [2] largeexcessofhydrogen,thereactionispseud o-first-orderattemperaturesbelow200 Cwiththerategivenbymole/(gofcatalyst s)whereRg=gasconstant, (mole'K)- Q -Ea=2700cal/molePH2=hydrogenpartialpress ure(torr)ko= (gcat s torr)Ko= (mole'K)T=absolutetemperature(K)CB=conce ntrationofbenzene(mole/cm3).))))
7 PriceandButt[3]studiedthisreactionina , typicalrun,PH2=685torrPB=densityoftherea ctorbed, , C6 SLOPE=f(xO'YO) (x)(X3'U3)IIISLOPE=f(x2,u2)IIIIIID efineC~=feedconcentrationofbenzene(mole/ cm3)z=axialreactorcoordinate(cm)Lreactor lengthydimensionlessconcentrationofbenze ne(CB/C~)x=dimensionlessaxialcoordinate( z/L).Theone-dimensionalsteady-statemater ialbalanceforthereactorthatexpressesthef actthatthechangeintheaxialconvectionofbe nzeneisequaltotheamountconvertedbyreacti oniswithC~atxoSinceeisconstant,:=-PBePH2 koKoTexp[(-~g-;'Ea)]yLetExplicitMethodsU singthedataprovided,wehave<!> ,thematerialbalancedy= ( )Nowwesolvethematerialbalanceequationusi ngtheEulermethod[Eq.( )]:whereUi+1= ;,i=0,1,2,..,N-1h= problemwheretheanalyticalsolutiondecreas eswithincreasingvaluesoftheindependentva riable,a numericalmethodisunstableif theglobalerrorgrowswithincreasingvalueso ftheindependentvariable(fora rigorousdefinitionofstability,see[4]).
8 Therefore,forthisproblemtheEulermethodis unstablewhenN= , , (forthisproblem), >20, practicalstandpoint,the"effective"reacti onzonewouldbeap-proximately0~x~ ,thena "short" ,weseethata largenumberofstepsarerequiredtoachievea "good" , reactant,thusbeingconvertedtoproductsast hefluidprogressestowardthereactoroutlet( x=1), ,a longerreactorwouldallowforgreaterconvers ion, ,smalleryvaluesatx=1.( ) ::= ,y=1atx=0 AnalyticalxSolutiontN=10N=20N=100N= (-1) (-2) (-1) (-1)- (-3) (-1) (-1) (-1) (-4) (-2) (-1) (-2)- (-5) (-2)0044837(-2) (-2) (-6) (-3) (-2) (-3) (-7) (-3) (- 3) (-3) (-8) (-4) (-3) (-4) (-9) (-4) (-4) (-4) (-10) (-5) (-4) (-5)- (-12) (-5) (-5) (-5) (-13)0045600(-6) (-5) (-6) (-14) (-6) (-6) (-6) (-15)oo40006(-7) (-6) (-7)- (-16) (-7) (-7) (-7) (-17) (-8) (-7) (-7)- (-18) (-8) (-7) (-8) (-19) (-9) (-8) (-8) (-20) (-10) (-8) (-9) (-21) (-10) (-9)t( -3) ,STABILITYInExample1 it wasseenthatforsomechoicesofthestep-size, theapproximatesolutionwasunstable, ,wewillexaminethequestionofstabilityusin gthetestequationdy=AydxyeO)=YowhereAisa ( )givesUt+1=Ut+AhutorUt+l=(1+hA)Ut=(1+hA) 2Ut_1=Theanalyticalsolutionof( )isy(xt+1)
9 =yoeAXi+l=yoe(i+l)hA( )( )( )Comparing( )with( )showsthattheapplicationofEuler'smethodt o( )isequivalenttousingtheexpression(1+hA) ,Y=1atx=0dxAnalyticalxSolutionN=100N=100 0N= machinenumber(seeAppendixA),theneo= ( ),withUoreplacedbyYo-eo,Ui+1(1+h1l.)i+1( yo-eo)andtheglobalerror(5i+1is(5i+1=Y(Xi +1)-Ui+1=yoe(i+1)hA-(1+hA)i+1(Yo-eo)or(5 i+1=[e(i+1)Ah-(1+hA)i+1]Yo+(1+hA)i+1eo( )Hence, ,thereisanerrorthatresultsfromtheEulerme thodapproximation(1+hA) thepropagationeffectoftheinitialerror, ,if11+hAl>1,thiscomponentwillgrowand,nom atterwhatthemagnitudeofeois,it willbecomethedominanttermin(5i+l'Therefo re,tokeepthepropagationeffectsofprevious errorsboundedwhenusingtheEulermethod,wer equire11+hAl<s;1( )Theregionofabsolutestabilityisdefinedby thesetofh(realnonnegative)andAvaluesforw hicha perturbationina singlevalueUiwillproducea changeinsubsequentvaluesthatdoesnotincre asefromsteptostep[4].))))
10 Thus,onecanseefrom( )thatthestabilityregionfor( )correspondstoa unitdiskinthecomplexhA-planecenteredat(- 1,0).IfAisreal,then- 2<S;hA<S;0( )Noticethatif thepropagationeffectisthedominanttermin( ),theglobalerrorwilloscillateinsignif-2<S;hA<S; 'YDifferentialEquationsReferringtoExample1,findthemaximumallowablestep-sizeforstabilityandfornonoscillatorybehaviorforthematerialbalance equationsofthe"long"and"short" :'AL= :'As= :0;;,:h'A;;,:- 2 Fornonoscillatoryerror:0;;,:h'A>-1(real)(real)UnstableStable,erroroscill ationsStable, < ~h~ < < ~h~ < , "reasonably"accuratesolutionisa > ,h= ,whileforN=20,h= ,whenN=20,theglobalerrorshouldoscillatei f thepropagationerroristhedominantterminEq .( ). ( ):6i+l=[e(i+l)Ah-(1+h'A)i+l]yo+(1+h'A)i+ leo=(A)yo+(B)eoForN=10,h= 'Ah= ,o12(A) (B) andeoissmall,theglobalerrorisdominatedby term(A)andnotthepropagationterm, ,term(B).
