gusucode.com > demos工具箱matlab源码程序 > demos/fem2ode.m
function fem2ode(N) %FEM2ODE Stiff problem with a constant mass matrix, M*y' = f(t,y). % The parameter N controls the discretization, and the resulting system % consists of N equations. For example, to solve a 20x20 system, use % fem2ode(20) % By default, N is 9. % % F(T,Y) returns the derivatives vector for a finite element % discretization of a partial differential equation. By default, the % solvers of the ODE Suite solve systems of the form y' = f(t,y). To solve % a system My' = f(t,y), use ODESET to set the property 'Mass' to a constant % matrix M. % % See also ODE23S, ODE15S, ODE23T, ODE23TB, ODESET, FUNCTION_HANDLE. % Mark W. Reichelt and Lawrence F. Shampine, 11-11-94. % Copyright 1984-2014 The MathWorks, Inc. % Problem parameter, shared with the nested function. if nargin < 1 N = 9; end tspan = [0; pi]; y0 = sin((pi/(N+1))*(1:N)'); % Constant mass matrix e = repmat(pi/(6*(N+1)),N,1); % h = pi/(N+1); e = (h/6)+zeros(N,1); M = spdiags([e 4*e e], -1:1, N, N); % mass matrix options = odeset('Mass',M); [t,y] = ode23s(@f,tspan,y0,options); figure; surf(1:N,t,y); zlim([0 1]); view(142.5,30); title(['Finite element problem with constant mass matrix, ' ... 'solved by ODE23S']); xlabel('space'); ylabel('time'); zlabel('solution'); % ----------------------------------------------------------------------- % Nested function -- N is provided by the outer function. % function dydt = f(t,y) % Derivative function. N is provided by the outer function. e = repmat(exp(t)*(N+1)/pi,N,1); % h = pi/(N+1); e = (exp(t)/h)+zeros(N,1); d = repmat(-2*exp(t)*(N+1)/pi,N,1); R = spdiags([e d e], -1:1, N, N); dydt = R * y; end % ----------------------------------------------------------------------- end % fem2ode