Create the following second-order differential algebraic equation.

syms y(t);
eqn = diff(y(t),t,2) == (1-y(t)^2)*diff(y(t),t) - y(t);

Use reduceDifferentialOrder to rewrite that equation as a system of two first-order differential equations. Here, vars is a vector of state variables of the system. The new variable Dy(t) represents the first derivative of y(t) with respect to t.

[eqs,vars] = reduceDifferentialOrder(eqn,y(t))

eqs =
 diff(Dyt(t), t) + y(t) + Dyt(t)*(y(t)^2 - 1)
                       Dyt(t) - diff(y(t), t)
 
vars =
   y(t)
 Dyt(t)

Set initial conditions for y(t) and its derivative Dy(t) to 2 and 0 respectively.

initConditions = [2 0];

Find the mass matrix M of the system and the right sides of the equations F.

[M,F] = massMatrixForm(eqs,vars)

M =
[  0, 1]
[ -1, 0]
 
F =
 - y(t) - Dyt(t)*(y(t)^2 - 1)
                      -Dyt(t)

M and F refer to the form $$M\left(t,x\left(t\right)\right)\dot{x}\left(t\right)=F\left(t,x\left(t\right)\right).$$. To simplify further computations, rewrite the system in the form $$\dot{x}\left(t\right)=f\left(t,x\left(t\right)\right)$$.

f = M\F

f =
                          Dyt(t)
 - Dyt(t)*y(t)^2 - y(t) + Dyt(t)

Convert f to a MATLAB function handle by using odeFunction. The resulting function handle is input to the MATLAB ODE solver ode15s.

odefun = odeFunction(f,vars);
ode15s(odefun, [0 10], initConditions)

Create the system of differential algebraic equations. Here, the symbolic functions x1(t) and x2(t) represent the state variables of the system. The system also contains constant symbolic parameters a, b, and the parameter function r(t). These parameters do not represent state variables. Specify the equations and state variables as two symbolic vectors: equations as a vector of symbolic equations, and variables as a vector of symbolic function calls.

syms x1(t) x2(t) a b r(t)
eqs = [diff(x1(t),t) == a*x1(t) + b*x2(t)^2,...
       x1(t)^2 + x2(t)^2 == r(t)^2];
vars = [x1(t) x2(t)];

Find the mass matrix M and vector of the right side F for this system. M and F refer to the form $$M\left(t,x\left(t\right)\right)\dot{x}\left(t\right)=F\left(t,x\left(t\right)\right).$$.

[M,F] = massMatrixForm(eqs,vars)

M =
[ 1, 0]
[ 0, 0]
 
F =
        b*x2(t)^2 + a*x1(t)
 r(t)^2 - x1(t)^2 - x2(t)^2

Use odeFunction to generate MATLAB function handles from M and F. The function handle F contains symbolic parameters.

M = odeFunction(M,vars)
F = odeFunction(F,vars,a,b,r(t))

M = 
  function_handle with value:
    @(t,in2)reshape([1.0,0.0,0.0,0.0],[2,2])

F = 
  function_handle with value:
    @(t,in2,param1,param2,param3)[param1.*in2(1,:)+...
    param2.*in2(2,:).^2;param3.^2-in2(1,:).^2-in2(2,:).^2]

Specify the parameter values.

a = -0.6;
b = -0.1;
r = @(t) cos(t)/(1+t^2);

Create the reduced function handle F.

F = @(t,Y) F(t,Y,a,b,r(t));

Specify consistent initial conditions for the DAE system.

t0 = 0;
y0 = [-r(t0)*sin(0.1); r(t0)*cos(0.1)];
yp0 = [a*y0(1) + b*y0(2)^2; 1.234];

Create an option set that contains the mass matrix M of the system and vector yp0 of initial conditions for the derivatives.

opt = odeset('mass',M,'InitialSlope',yp0);

Now, use ode15s to solve the system of equations.

ode15s(F, [t0, 1], y0, opt)

Define the system of differential equations. Find the mass matrix M and the right side F.

syms x(t) y(t)
eqs = [diff(x(t),t)+2*diff(y(t),t) == 0.1*y(t), ...
       x(t)-y(t) == cos(t)-0.2*t*sin(x(t))];
vars = [x(t) y(t)];
[M,F] = massMatrixForm(eqs,vars);

Write the MATLAB code for M and F to the files myfileM and myfileF. odeFunction overwrites existing files. Include the comment Version: 1.1 in the files You can open and edit the output files.

M = odeFunction(M,vars,'File','myfileM','Comments','Version: 1.1');

function expr = myfileM(t,in2)
%MYFILEM
%    EXPR = MYFILEM(T,IN2)

%    This function was generated by the Symbolic Math Toolbox version 7.3.
%    01-Jan-2017 00:00:00

%Version: 1.1
expr = reshape([1.0,0.0,2.0,0.0],[2, 2]);

F = odeFunction(F,vars,'File','myfileF','Comments','Version: 1.1');

function expr = myfileF(t,in2)
%MYFILEF
%    EXPR = MYFILEF(T,IN2)

%    This function was generated by the Symbolic Math Toolbox version 7.3.
%    01-Jan-2017 00:00:00

%Version: 1.1
x = in2(1,:);
y = in2(2,:);
expr = [y.*(1.0./1.0e1);-x+y+cos(t)-t.*sin(x).*(1.0./5.0)];

Specify consistent initial values for x(t) and y(t) and their first derivatives.

xy0 = [2; 1];    % x(t) and y(t)
xyp0 = [0; 0.05*xy0(2)];    % derivatives of x(t) and y(t)

Create an option set that contains the mass matrix M, initial conditions xyp0, and numerical tolerances for the numerical search.

opt = odeset('mass', M, 'RelTol', 10^(-6),...
               'AbsTol', 10^(-6), 'InitialSlope', xyp0);

Solve the system of equations by using ode15s.

ode15s(F, [0 7], xy0, opt)

Create the system of differential algebraic equations. Here, the symbolic functions x1(t) and x2(t) represent the state variables of the system. Specify the equations and state variables as two symbolic vectors: equations as a vector of symbolic equations, and variables as a vector of symbolic function calls.

syms x1(t) x2(t)

a = -0.6;
b = -0.1;
r = @(t) cos(t)/(1 + t^2);

eqs = [diff(x1(t),t) == a*x1(t) + b*x2(t)^2,...
       x1(t)^2 + x2(t)^2 == r(t)^2];
vars = [x1(t) x2(t)];

Find the mass matrix M and vector of the right side F for this system. M and F refer to the form $$M\left(t,x\left(t\right)\right)\dot{x}\left(t\right)=F\left(t,x\left(t\right)\right).$$.

[M,F] = massMatrixForm(eqs,vars)

M =
[ 1, 0]
[ 0, 0]
 
F =
               - (3*x1(t))/5 - x2(t)^2/10
 cos(t)^2/(t^2 + 1)^2 - x1(t)^2 - x2(t)^2

Generate MATLAB function handles from M and F. Because most of the elements of the mass matrix M are zeros, use the Sparse argument when converting M.

M = odeFunction(M,vars,'Sparse',true)
F = odeFunction(F,vars)

M = 
  function_handle with value:
    @(t,in2)sparse([1],[1],[1.0],2,2)

F = 
  function_handle with value:
    @(t,in2)[in2(1,:).*(-3.0./5.0)-in2(2,:).^2./1.0e+1;...
    cos(t).^2.*1.0./(t.^2+1.0).^2-in2(1,:).^2-in2(2,:).^2]

Specify consistent initial conditions for the DAE system.

t0 = 0;
y0 = [-r(t0)*sin(0.1); r(t0)*cos(0.1)];
yp0= [a*y0(1) + b*y0(2)^2; 1.234];

Create an option set that contains the mass matrix M of the system and vector yp0 of initial conditions for the derivatives.

opt = odeset('mass',M,'InitialSlope', yp0);

Solve the system of equations using ode15s.

ode15s(F, [t0, 1], y0, opt)