www.gusucode.com > matlab 案例源码 matlab代码程序 > matlab/SolveNonstiffODEsExample.m
%% Solve Nonstiff ODEs % This page contains two examples of solving nonstiff ordinary differential % equations using |ode45|. MATLAB(R) has three solvers for nonstiff ODEs. % % * |ode45| % * |ode23| % * |ode113| % % For most nonstiff problems, |ode45| performs best. However, |ode23| is % recommended for problems that permit a slightly cruder error tolerance or % in the presence of moderate stiffness. Likewise, |ode113| can be more % efficient than |ode45| for problems with stringent error tolerances. % % If the nonstiff solvers take a long time to solve the problem or % consistently fail the integration, then the problem might be stiff. See % <docid:matlab_math.bu22m86> for more information. %% Example: Nonstiff van der Pol Equation % The van der Pol equation is a second order ODE % % $$y''_1 - \mu \left( 1 - y_1^2\right) y'_1+y_1=0,$$ % % where $\mu > 0$ is a scalar parameter. Rewrite this equation as a system % of first-order ODEs by making the substitution $y'_1 = y_2$. The % resulting system of first-order ODEs is % % $$ % \begin{array}{cl} % y'_1 &= y_2\\ % y'_2 &= \mu (1-y_1^2) y_2 - y_1.\end{array} % $$ % %% % The system of ODEs must be coded into a function file that the ODE solver % can use. The general functional signature of an ODE function is % % dydt = odefun(t,y) % % That is, the function must accept both |t| and |y| as inputs, even if it % does not use |t| for any computations. % % The function file |vdp1.m| codes the van der Pol equation using $\mu = % 1$. The variables $y_1$ and $y_2$ are represented by |y(1)| and |y(2)|, % and the two-element column vector |dydt| contains the expressions for % $y'_1$ and $y'_2$. % % <include>vdp1.m</include> % %% % Solve the ODE using the |ode45| function on the time interval |[0 20]| % with initial values |[2 0]|. The output is a column vector of time points % |t| and a solution array |y|. Each row in |y| corresponds to a time % returned in the corresponding row of |t|. The first column of |y| % corresponds to $y_1$, and the second column to $y_2$. [t,y] = ode45(@vdp1,[0 20],[2; 0]); %% % Plot the solutions for $y_1$ and $y_2$ against |t|. plot(t,y(:,1),'-o',t,y(:,2),'-o') title('Solution of van der Pol Equation (\mu = 1) using ODE45'); xlabel('Time t'); ylabel('Solution y'); legend('y_1','y_2') %% % The |vdpode| function solves the same problem, but it accepts a % user-specified value for $\mu$. The van der Pol equations become stiff as % $\mu$ increases. For example, with the value $\mu = 1000$ you need to use % a stiff solver such as |ode15s| to solve the system. %% Example: Nonstiff Euler Equations % The Euler equations for a rigid body without external forces are a % standard test problem for ODE solvers intended for nonstiff problems. % % The equations are % % $$ \begin{array}{cl} y'_1 &= y_2y_3 \\ y'_2 &= -y_1y_3 \\ y'_3 &= % -0.51y_1y_2. \end{array}$$ % %% % The function file |rigidode| defines and solves this first-order system % of equations over the time interval |[0 12]|, using the vector of initial % conditions |[0; 1; 1]| corresponding to the initial values of $y_1$, % $y_2$, and $y_3$. The local function |f(t,y)| encodes the system of % equations. % % |rigidode| calls |ode45| with no output arguments, so the solver uses the % default output function |odeplot| to automatically plot the solution % points after each step. % % <include>rigidode.m</include> % %% % Solve the nonstiff Euler equations by calling the |rigidode| function. rigidode title('Solution of Rigid Body w/o External Forces using ODE45') legend('y_1','y_2','y_3','Location','Best')