www.gusucode.com > matlab 案例源码 matlab代码程序 > matlab/SolveStiffDAEExample.m
%% Solve Stiff Differential Algebraic Equation % This example shows how to use |ode23t| to solve a stiff differential % algebraic equation (DAE) problem that arises from an electrical circuit. % This problem is originally from p.377 of E. Hairer and G. Wanner, % _Solving Ordinary Differential Equations II: Stiff and % Differential-Algebraic Problems_, 2nd ed. (Berlin: Springer, 1996). % % The one-transistor amplifier problem coded by <matlab:edit('amp1dae.m') % amp1dae.m> can be rewritten in semi-explicit form, but this example % solves it in its original form $M u' = \phi(u)$. The problem includes a % constant, singular mass matrix that is not diagonal. The transistor % amplifier circuit contains six resistors, three capacitors, and a % transistor. % % <<../transistor_amplifier.png>> % % The initial voltage signal is $U_e(t)=0.4 \sin\left( 200 \pi t \right)$, % and the other parameters are constant. The goal is to solve for the % output voltage through node 5. % % Set the values of the initial and operating voltages. Ue = @(t) 0.4*sin(200*pi*t); Ub = 6; %% % Using Kirchoff's law to equalize the current through each node (1 through % 5), you obtain a system of five equations. The mass matrix of this system % has the form % % $$ M = \left( \begin{array}{ccccc} -c_1 & c_1 & 0 & 0 & 0\\ c_1 & -c_1 & % 0 & 0 & 0\\ 0 & 0 & -c_2 & 0 & 0\\ 0 & 0 & 0 & -c_3 & c_3\\ 0 & 0 & 0 & % c_3 & -c_3 \end{array} \right),$$ % % where $c_k = k \times 10^{-6}$ for $k=1,2,3$. % % Use the |odeset| function to pass in the mass matrix. Leave the % |'MassSingular'| option at its default value of |'maybe'| to test the % automatic detection of a DAE problem by the solver. c = 1e-6 * (1:3); M = zeros(5,5); M(1,1) = -c(1); M(1,2) = c(1); M(2,1) = c(1); M(2,2) = -c(1); M(3,3) = -c(2); M(4,4) = -c(3); M(4,5) = c(3); M(5,4) = c(3); M(5,5) = -c(3); opts = odeset('Mass',M); %% % The function |transampdae.m| contains the system of equations for this % example. Save this function in your current folder to run the example. % % <include>transampdae.m</include> % %% % Set the initial conditions. This example uses the consistent initial % conditions computed by Hairer and Wanner. u0(1) = 0; u0(2) = Ub/2; u0(3) = Ub/2; u0(4) = Ub; u0(5) = 0; %% % Solve the DAE system over the time interval |[0 0.05]| using |ode23t|. tspan = [0 0.05]; [t,u] = ode23t(@transampdae,tspan,u0,opts); %% % Plot the initial voltage $U_e(t)$ and output voltage $U_5(t)$. plot(t,Ue(t),'o',t,u(:,5),'.') axis([0 0.05 -3 2]); legend('Input Voltage U_e(t)','Output Voltage U_5(t)', 'Location', 'NorthWest'); title('One transistor amplifier DAE problem solved by ODE23T'); xlabel('t');