matlab中ode45函数编写.docx
《matlab中ode45函数编写.docx》由会员分享，可在线阅读，更多相关《matlab中ode45函数编写.docx（23页珍藏版）》请在冰豆网上搜索。
matlab中ode45函数编写
functionvarargout=ode45（ode,tspan,y0,options,varargin）
%ODE45Solvenon-stiffdifferentialequations,mediumordermethod.
%[TOUT,YOUT]=ODE45（ODEFUN,TSPAN,Y0）withTSPAN=[T0TFINAL]integrates
%thesystemofdifferentialequationsy'=f（t,y）fromtimeT0toTFINAL
%withinitialconditionsY0.ODEFUNisafunctionhandle.ForascalarT
%andavectorY,ODEFUN（T,Y）mustreturnacolumnvectorcorresponding
%tof（t,y）.EachrowinthesolutionarrayYOUTcorrespondstoatime
%returnedinthecolumnvectorTOUT.Toobtainsolutionsatspecific
%timesT0,T1,...,TFINAL（allincreasingoralldecreasing）,useTSPAN=
%[T0T1...TFINAL].
%
%[TOUT,YOUT]=ODE45（ODEFUN,TSPAN,Y0,OPTIONS）solvesasabovewithdefault
%integrationpropertiesreplacedbyvaluesinOPTIONS,anargumentcreated
%withtheODESETfunction.SeeODESETfordetails.Commonlyusedoptions
%arescalarrelativeerrortolerance'RelTol'（1e-3bydefault）andvector
%ofabsoluteerrortolerances'AbsTol'（allcomponents1e-6bydefault）.
%Ifcertaincomponentsofthesolutionmustbenon-negative,use
%ODESETtosetthe'NonNegative'propertytotheindicesofthese
%components.
%
%ODE45cansolveproblemsM（t,y）*y'=f（t,y）withmassmatrixMthatis
%nonsingular.UseODESETtosetthe'Mass'propertytoafunctionhandle
%MASSifMASS（T,Y）returnsthevalueofthemassmatrix.Ifthemassmatrix
%isconstant,thematrixcanbeusedasthevalueofthe'Mass'option.If
%themassmatrixdoesnotdependonthestatevariableYandthefunction
%MASSistobecalledwithoneinputargumentT,set'MStateDependence'to
%'none'.ODE15SandODE23Tcansolveproblemswithsingularmassmatrices.
%
%[TOUT,YOUT,TE,YE,IE]=ODE45（ODEFUN,TSPAN,Y0,OPTIONS）withthe'Events'
%propertyinOPTIONSsettoafunctionhandleEVENTS,solvesasabove
%whilealsofindingwherefunctionsof（T,Y）,calledeventfunctions,
%arezero.Foreachfunctionyouspecifywhethertheintegrationis
%toterminateatazeroandwhetherthedirectionofthezerocrossing
%matters.ThesearethethreecolumnvectorsreturnedbyEVENTS:
%[VALUE,ISTERMINAL,DIRECTION]=EVENTS（T,Y）.FortheI-theventfunction:
%VALUE（I）isthevalueofthefunction,ISTERMINAL（I）=1iftheintegration
%istoterminateatazeroofthiseventfunctionand0otherwise.
%DIRECTION（I）=0ifallzerosaretobecomputed（thedefault）,+1ifonly
%zeroswheretheeventfunctionisincreasing,and-1ifonlyzeroswhere
%theeventfunctionisdecreasing.OutputTEisacolumnvectoroftimes
%atwhicheventsoccur.RowsofYEarethecorrespondingsolutions,and
%indicesinvectorIEspecifywhicheventoccurred.
%
%SOL=ODE45（ODEFUN,[T0TFINAL],Y0...）returnsastructurethatcanbe
%usedwithDEVALtoevaluatethesolutionoritsfirstderivativeat
%anypointbetweenT0andTFINAL.ThestepschosenbyODE45arereturned
%inarowvectorSOL.x.ForeachI,thecolumnSOL.y（:
I）contains
%thesolutionatSOL.x（I）.Ifeventsweredetected,SOL.xeisarowvector
%ofpointsatwhicheventsoccurred.ColumnsofSOL.yearethecorresponding
%solutions,andindicesinvectorSOL.iespecifywhicheventoccurred.
%
%Example
%[t,y]=ode45（@vdp1,[020],[20]）;
%plot（t,y（:
1））;
%solvesthesystemy'=vdp1（t,y）,usingthedefaultrelativeerror
%tolerance1e-3andthedefaultabsolutetoleranceof1e-6foreach
%component,andplotsthefirstcomponentofthesolution.
%
%ClasssupportforinputsTSPAN,Y0,andtheresultofODEFUN（T,Y）:
%float:
double,single
%
%Seealso
%otherODEsolvers:
ODE23,ODE113,ODE15S,ODE23S,ODE23T,ODE23TB
%implicitODEs:
ODE15I
%optionshandling:
ODESET,ODEGET
%outputfunctions:
ODEPLOT,ODEPHAS2,ODEPHAS3,ODEPRINT
%evaluatingsolution:
DEVAL
%ODEexamples:
RIGIDODE,BALLODE,ORBITODE
%functionhandles:
FUNCTION_HANDLE
%
%NOTE:
%TheinterpretationofthefirstinputargumentoftheODEsolversand
%somepropertiesavailablethroughODESEThavechangedinMATLAB6.0.
%Althoughwestillsupportthev5syntax,anynewfunctionalityis
%availableonlywiththenewsyntax.Toseethev5help,typein
%thecommandline
%moreon,typeode45,moreoff
%NOTE:
%Thisportiondescribesthev5syntaxofODE45.
%
%[T,Y]=ODE45（'F',TSPAN,Y0）withTSPAN=[T0TFINAL]integratesthe
%systemofdifferentialequationsy'=F（t,y）fromtimeT0toTFINALwith
%initialconditionsY0.'F'isastringcontainingthenameofanODE
%file.FunctionF（T,Y）mustreturnacolumnvector.Eachrowin
%solutionarrayYcorrespondstoatimereturnedincolumnvectorT.To
%obtainsolutionsatspecifictimesT0,T1,...,TFINAL（allincreasing
%oralldecreasing）,useTSPAN=[T0T1...TFINAL].
%
%[T,Y]=ODE45（'F',TSPAN,Y0,OPTIONS）solvesasabovewithdefault
%integrationparametersreplacedbyvaluesinOPTIONS,anargument
%createdwiththeODESETfunction.SeeODESETfordetails.Commonly
%usedoptionsarescalarrelativeerrortolerance'RelTol'（1e-3by
%default）andvectorofabsoluteerrortolerances'AbsTol'（all
%components1e-6bydefault）.
%
%[T,Y]=ODE45（'F',TSPAN,Y0,OPTIONS,P1,P2,...）passestheadditional
%parametersP1,P2,...totheODEfileasF（T,Y,FLAG,P1,P2,...）（see
%ODEFILE）.UseOPTIONS=[]asaplaceholderifnooptionsareset.
%
%ItispossibletospecifyTSPAN,Y0andOPTIONSintheODEfile（see
%ODEFILE）.IfTSPANorY0isempty,thenODE45callstheODEfile
%[TSPAN,Y0,OPTIONS]=F（[],[],'init'）toobtainanyvaluesnotsupplied
%intheODE45argumentlist.Emptyargumentsattheendofthecalllist
%maybeomitted,e.g.ODE45（'F'）.
%
%ODE45cansolveproblemsM（t,y）*y'=F（t,y）withamassmatrixMthatis
%nonsingular.UseODESETtosetMassto'M','M（t）',or'M（t,y）'ifthe
%ODEfileiscodedsothatF（T,Y,'mass'）returnsaconstant,
%time-dependent,ortime-andstate-dependentmassmatrix,respectively.
%ThedefaultvalueofMassis'none'.ODE15SandODE23Tcansolveproblems
%withsingularmassmatrices.
%
%[T,Y,TE,YE,IE]=ODE45（'F',TSPAN,Y0,OPTIONS）withtheEventspropertyin
%OPTIONSsetto'on',solvesasabovewhilealsolocatingzerocrossings
%ofaneventfunctiondefinedintheODEfile.TheODEfilemustbe
%codedsothatF（T,Y,'events'）returnsappropriateinformation.See
%ODEFILEfordetails.OutputTEisacolumnvectoroftimesatwhich
%eventsoccur,rowsofYEarethecorrespondingsolutions,andindicesin
%vectorIEspecifywhicheventoccurred.
%
%SeealsoODEFILE
%ODE45isanimplementationoftheexplicitRunge-Kutta（4,5）pairof
%DormandandPrincecalledvariouslyRK5（4）7FM,DOPRI5,DP（4,5）andDP54.
%Itusesa"free"interpolantoforder4communicatedprivatelyby
%DormandandPrince.Localextrapolationisdone.
%DetailsaretobefoundinTheMATLABODESuite,L.F.Shampineand
%M.W.Reichelt,SIAMJournalonScientificComputing,18-1,1997.
%MarkW.ReicheltandLawrenceF.Shampine,6-14-94
%Copyright1984-2009TheMathWorks,Inc.
%$Revision:
5.74.4.10$$Date:
2009/04/2103:
24:
15$
solver_name='ode45';
%Checkinputs
ifnargin<4
options=[];
ifnargin<3
y0=[];
ifnargin<2
tspan=[];
ifnargin<1
error（'MATLAB:
ode45:
NotEnoughInputs',...
'Notenoughinputarguments.SeeODE45.'）;
end
end
end
end
%Stats
nsteps=0;
nfailed=0;
nfevals=0;
%Output
FcnHandlesUsed=isa（ode,'function_handle'）;
output_sol=（FcnHandlesUsed&&（nargout==1））;%sol=odeXX（...）
output_ty=（~output_sol&&（nargout>0））;%[t,y,...]=odeXX（...）
%Theremightbenooutputrequested...
sol=[];f3d=[];
ifoutput_sol
sol.solver=solver_name;
sol.extdata.odefun=ode;
sol.extdata.options=options;
sol.extdata.varargin=varargin;
end
%Handlesolverarguments
[neq,tspan,ntspan,next,t0,tfinal,tdir,y0,f0,odeArgs,odeFcn,...
options,threshold,rtol,normcontrol,normy,hmax,htry,htspan,dataType]=...
odearguments（FcnHandlesUsed,solver_name,ode,tspan,y0,options,varargin）;
nfevals=nfevals+1;
%Handletheoutput
ifnargout>0
outputFcn=odeget（options,'OutputFcn',[],'fast'）;
else
outputFcn=odeget（options,'OutputFcn',@odeplot,'fast'）;
end
outputArgs={};
ifisempty（outputFcn）
haveOutputFcn=false;
else
haveOutputFcn=true;
outputs=odeget（options,'OutputSel',1:
neq,'fast'）;
ifisa（outputFcn,'function_handle'）
%WithMATLAB6syntaxpassadditionalinputargumentstooutputFcn.
outputArgs=varargin;
end
end
refine=max（1,odeget（options,'Refine',4,'fast'））;
ifntspan>2
outputAt='RequestedPoints';%outputonlyattspanpoints
elseifrefine<=1
outputAt='SolverSteps';%computedpoints,norefinement
else
outputAt='RefinedSteps';%computedpoints,withrefinement
S=（1:
refine-1）/refine;
end
printstats=strcmp（odeget（options,'Stats','off','fast'）,'on'）;
%Handletheeventfunction
[haveEventFcn,eventFcn,eventArgs,valt,teout,yeout,ieout]=...
odeevents（FcnHandlesUsed,odeFcn,t0,y0,options,varargin）;
%Handlethemassmatrix
[Mtype,M,Mfun]=odemass（FcnHandlesUsed,odeFcn,t0,y0,options,varargin）;
ifMtype>0%non-trivialmassmatrix
Msingular=odeget（options,'MassSingular','no','fast'）;
ifstrcmp（Msingular,'maybe'）
warning（'MATLAB:
ode45:
MassSingularAssumedNo',['ODE45assumes'...
'MassSingularis''no''.SeeODE15SorODE23T.']）;
elseifstrcmp（Msingular,'yes'）
error（'MATLAB:
ode45:
MassSingularYes',...
['MassSingularcannotbe''yes''forthissolver.SeeODE15S'...
'orODE23T.']）;
end
%IncorporatethemassmatrixintoodeFcnandodeArgs.
[odeFcn,odeArgs]=odemassexplicit（FcnHandlesUsed,Mtype,odeFcn,odeArgs,Mfun,M）;
f0
升级会员 