药物计算分析导论-第二部分-11-12章Matlab.ppt

2003年 浙江大学中药科学与工程学系 Chap 11MoreAdvancedTimeIntegrationOfAnODE 2003年 浙江大学中药科学与工程学系 ThisLecture MoreaccurateschemesMorecomplicatedODEsVariabletimestepandembeddedmethodsusedtomakesureerrorsarewithinatolerance 2003年 浙江大学中药科学与工程学系 Adams BashforthSchemes IntheforwardEulerschemeweonlyusedthevalueoftherighthandsideattheprevioustimestep i e weonlyusedalinearapproximationtothetimederivative 2003年 浙江大学中药科学与工程学系 ABSchemes Idea wesetwhereIfinterpolatesfn fn 1 fn 2 fn 1 Nstagesi e 2003年 浙江大学中药科学与工程学系 f0 f1 f2 f3 Weinterpolatethefunctionthroughthefirst4points Thenweintegrateunderthecurvebetweent 3andt 4 2003年 浙江大学中药科学与工程学系 ABSchemes EssentiallyweuseinterpolationandaNewton Cotesquadratureformulatoformulate 2003年 浙江大学中药科学与工程学系 Runge KuttaSchemes SeevanLoanforderivationofRunge Kutta2andRunge Kutta4 Thefollowing simple schemeduetoJameson SchmidtandTurkel 1981 2003年 浙江大学中药科学与工程学系 Runge KuttaSchemes Beware itonlyworkswhenfisafunctionofyandnottheresistheorderofthescheme 2003年 浙江大学中药科学与工程学系 ErrorEstimate Matlabhasanumberoftimeintegratorsbuiltin ode23ode45andothers 2003年 浙江大学中药科学与工程学系 ode23 Forn 1 timestepsode23usestwoestimatesforyn 1 A2ndorderRKschemeanda3rdorderRKschemeareusedtobuildtwoguessesforyn 1 Ifthedifferencebetweenthesetwoestimatesarewithinatoleranceode23progressesontocalculatingyn 2Ifthedifferenceisgreaterthanthespecifiedtolerance ode23reducesthedtandtriesagain Itrepeatsuntilthedifferenceislowerthanthetolerance End 2003年 浙江大学中药科学与工程学系 T Y ode23 odefun tspan y0 withtspan t0tf integratesthesystemofdifferentialequationsfromtimet0totfwithinitialconditionsy0 Y0Avectorofinitialconditions OdefunAfunctionthatevaluatestheright handsideofthedifferentialequations 2003年 浙江大学中药科学与工程学系 PlanetsExampleUsingode23 Idea replaceourhomegrownEulerForwardschemewithMatlab sode23schemeintheplanets1 mscript 2003年 浙江大学中药科学与工程学系 Parametersfortheplanetsasbefore 2003年 浙江大学中药科学与工程学系 Initialvelocities GatherX Y VX VYintoonevector 2003年 浙江大学中药科学与工程学系 T CoordVel ode23 forcing TimeLimits CoordVel CalltoMatlabode23function functionwhichcalculatesf y 0FinalTime VectorholdsX Y VX VY WorksinMatlabv6 1 0 atleast 2003年 浙江大学中药科学与工程学系 Calltoode23 Findoutthetimeateachtimestep Extractfinalcoordinatesandvelocities Plotplanetorbits d size X returnsthesizesofeachdimensionofarrayXinavectordwithndims X elements 2003年 浙江大学中药科学与工程学系 Theroutinewhichcalculatestheforcingfunction 2003年 浙江大学中药科学与工程学系 TeamExercise Grabplanets2 mandforcing mRunthescriptUsetheTstepsvectortofindoutthetimestepforeachintegrationstage Plotagraphshowingthetimestep dt ateachtimestep Usehelptofindouthowtochangethetoleranceusedbyode23 hintyouwillneedtouseodeset Rerunthesimulationwithatoleranceof0 1 2003年 浙江大学中药科学与工程学系 2003年 浙江大学中药科学与工程学系 Application One DimensionalElectrostaticMotion 2003年 浙江大学中药科学与工程学系 ChargeRepulsion NowwewillconsiderthecaseofchargedparticleswiththesamesignchargeInsteadofattractingeachother thechargesrepeleachother 2003年 浙江大学中药科学与工程学系 ParticleAccelerations a1 a2 Eachparticlehasavectoraccelerationdirectlyawayfromtheotherparticle 2003年 浙江大学中药科学与工程学系 TeamProject Q1 Modifytheplanets2 mandforcing mtosimulatethefollowing ThereareNelectricallychargedparticlesconfinedtomoveinthex directiononlyDistributethechargesinitiallyatequispacedpointsin 1 1 Theequationsofmotionofthechargesare 2003年 浙江大学中药科学与工程学系 TeamProject Q1cont Includethefollowinginoneprojectwriteupperteam Ascatterplot withxasthehorizontalaxisandtastheverticalaxis showingthepathsofallthechargedparticlesusingode23Replaceode23withode45andrerunPlottheode23andode45 t x pathsonthesamegraph Plottheode23andode45timestepsizesonthesamegraphNamesofteammembers 2003年 浙江大学中药科学与工程学系 Chap 12Electrostatics CollidingDisksProject 2003年 浙江大学中药科学与工程学系 BigTeamProject Ok wehavebeenwarmingupwithsomesmallprojectsNowforsomethingalittlemoreinvolved 2003年 浙江大学中药科学与工程学系 ThePhysics Thedomainofinterestisatwo dimensionalperiodicboxofsidelengthLItispopulatedwithNcirclesofmassMeachandradiusR 2003年 浙江大学中药科学与工程学系 x y x 0 y 0 x L y L Domain 2003年 浙江大学中药科学与工程学系 x y x 0 y 0 x L y L Particles R 2003年 浙江大学中药科学与工程学系 DefiningASetofMotionRules Rule1 ThedisksareinitiallyrandomlydistributedThediskshaverandominitialvelocity withtheabsolutevalueofeachcomponentofvelocityboundedby1 Rule2 EffectivelyNewton sfirst allparticlesdonotaccelerateuntiltheycollide 2003年 浙江大学中药科学与工程学系 ParticlesCollide 2 1 Before 2 1 After 2003年 浙江大学中药科学与工程学系 Rule4 Momentum 动量 isconservedinacollision Rule5 Angularmomentum 动量矩 isconservedinthecollision Rule6 Totalkineticenergy 动能 isconservedThetotalkineticenergyoftheparticlesisthesamebeforeandafterthecollision 2003年 浙江大学中药科学与工程学系 ModelSimplification WearegoingtouseRules4 5 and6todeterminethechangeinvelocityoftwodiskswhentheycollide Wecansimplifythisproblembyconsideringacoordinateframewhichis co movingwiththecenterofoneofthedisks saydisk1 atitsoriginwithx axisalignedindirectionofthesegmentconnectingthetwodiskcenters 2003年 浙江大学中药科学与工程学系 InCo MovingFrame 2 1 Before 2003年 浙江大学中药科学与工程学系 MathVersionof4 inco movingframe Conservationoflinearmomentum Whichreducesto 2003年 浙江大学中药科学与工程学系 MathVersionof5 inco movingframe Conservationofangularmomentum Andsincethecenteroftwocannotbeattheoriginthisimplies 2003年 浙江大学中药科学与工程学系 MathVersionof6 inco movingframe Conservationofkineticenergy Reducesto 2003年 浙江大学中药科学与工程学系 Summaryof4 5 6 4 5 6 2003年 浙江大学中药科学与工程学系 Solutions Collision Nocollision 2003年 浙江大学中药科学与工程学系 InCo MovingFrame 2 1 Before 2 1 After 2003年 浙江大学中药科学与工程学系 RevertingToNon movingFrame Therulesessentiallysaythatiftwoequalmassdiskscollidethen theyexchangevelocityintheircenter to centervectortheyretaintheirownvelocityinthedirectionofthetangenttotheircenter to centervector 2003年 浙江大学中药科学与工程学系 GeneralFormula 2003年 浙江大学中药科学与工程学系 SummaryofCollisionFormula Giventherelativepositionoftwodisksbeforethecollisionandtheirvelocitieswecandeterminethepost collisionvelocities 2003年 浙江大学中药科学与工程学系 Changeinve