Main Content

odeMassMatrix

ODE mass matrix

Since R2023b

Description

An odeMassMatrix object represents the mass matrix for a system of ordinary differential equations.

ode objects can represent problems of the form M(t,y)y'=f(t,y), where M(t,y) is a mass matrix that can be full or sparse. The mass matrix encodes linear combinations of derivatives on the left side of the equation.

Create an ode object to represent the ODE problem, and specify the odeMassMatrix object as the value of the MassMatrix property to incorporate a mass matrix into the problem.

Creation

Description

M = odeMassMatrix creates an odeMassMatrix object with empty properties.

example

M = odeMassMatrix(Name=Value) specifies one or more property values using name-value arguments.

Properties

expand all

Mass matrix, specified as a matrix or handle to a function that evaluates the mass matrix.

ode objects can represent problems of the form M(t,y)y'=f(t,y), where M(t,y) is a mass matrix that can be full or sparse. The mass matrix encodes linear combinations of derivatives on the left side of the equation.

  • When the mass matrix is nonsingular, the equation simplifies to y'=M1f(t,y) and the ODE has a solution for any initial value. However, it is often more convenient and natural to express the model in terms of the mass matrix directly using M(t,y)y'=f(t,y), and avoiding the computation of the matrix inverse reduces the storage and execution time needed to solve the problem.

  • When the mass matrix is singular, the problem is a system of differential algebraic equations (DAEs). A DAE has a solution only when the initial values are consistent; that is, you must specify the initial slope y'0 using the InitialSlope property of the ode object such that the initial conditions are all consistent, M(t0,y0)y'0=f(t0,y0). For more information, see Solve Differential Algebraic Equations (DAEs).

You can specify the value of the MassMatrix property as:

  • A constant matrix with calculated values for M(t,y).

  • A handle to a function that computes the matrix elements and that accepts two input arguments, M = Mass(t,y). To give the function access to parameter values in the Parameters property, specify a third input argument in the function definition, M = Mass(t,y,p).

Example: M = odeMassMatrix(MassMatrix=@Mass) specifies the function Mass that evaluates the mass matrix.

Example: M = odeMassMatrix(MassMatrix=[1 0; -2 1]) specifies a constant mass matrix.

Data Types: single | double | function_handle
Complex Number Support: Yes

Singular mass matrix toggle, specified as "maybe", "yes", or "no". The default value of "maybe" causes the solver to test whether the problem is a differential algebraic equation (DAE) by testing whether the mass matrix is singular. Avoid this check by specifying "yes" if you know the system is a DAE with a singular mass matrix, or "no" if it is not.

Example: M = odeMassMatrix(MassMatrix=X,Singular="no")

State dependence of the mass matrix, specified as "weak", "strong", or "none".

  • For ODE problems of the form M(t)y'=f(t,y), set the StateDependence property to "none". This value ensures that the solver calls the mass matrix function with a single argument as in Mass(t) (or two arguments as in Mass(t,p) if you use the Parameters property of the ode object to pass parameters).

  • If the mass matrix depends on y, then set StateDependence to either "weak" (default) or "strong". In both cases, the solver calls the mass matrix function with two arguments as in Mass(t,y) (or three arguments as in Mass(t,y,p) if you use the Parameters property of the ode object to pass parameters). However, a value of "weak" results in implicit solvers using approximations when solving algebraic equations.

Example: M = odeMassMatrix(MassMatrix=@Mass,StateDependence="none") specifies that the mass matrix function Mass depends only on t.

Mass matrix sparsity pattern, specified as a sparse matrix. Use this property to specify the sparsity pattern of the matrix y[M(t,y)v]. The sparse matrix S has S(i,j) = 1 if, for any k, the (i,k) component of M(t,y) depends on component j of y. This property is similar to the MvPattern option of odeset.

Note

The SparsityPattern property is used by the ode15s, ode23t, and ode23tb solvers when StateDependence is "strong".

Example: M = odeMassMatrix(MassMatrix=@Mass,StateDependence="strong",SparsityPattern=S)

Data Types: double

Examples

collapse all

Consider this system of first-order equations.

y1=y1y3-y2y2=y1-10=y1+y2+y3

The left side of the equations contain time derivatives, yn=dyndt. However, because the derivative for y3 does not appear in the system, the equations define a system of differential algebraic equations. Rewriting the system in the form My=f(t,y) shows a constant, singular mass matrix on the left side.

[100010000][y1y2y3]˙=[y1y3-y2y1-1y1+y2+y3]

Solve the system of equations using the initial values [1 1 -2] by creating an ode object to represent the problem.

  • Specify the initial values in the InitialValue property.

  • Specify the system of equations as an anonymous function in the ODEFcn property.

  • Use an odeMassMatrix object to specify the constant, singular mass matrix in the MassMatrix property.

F = ode;
F.InitialValue = [1 1 -2];
F.ODEFcn = @(t,y) [y(1)*y(3)-y(2); 
                   y(1)-1; 
                   y(1)+y(2)+y(3)];
F.MassMatrix = odeMassMatrix(MassMatrix=[1 0 0; 0 1 0; 0 0 0],Singular="yes");

Display the ode object. The SelectedSolver property shows that the ode15s solver was automatically chosen for this problem.

F
F = 
  ode with properties:

   Problem definition
               ODEFcn: @(t,y)[y(1)*y(3)-y(2);y(1)-1;y(1)+y(2)+y(3)]
          InitialTime: 0
         InitialValue: [1 1 -2]
             Jacobian: []
           MassMatrix: [1x1 odeMassMatrix]

   Solver properties
    AbsoluteTolerance: 1.0000e-06
    RelativeTolerance: 1.0000e-03
               Solver: auto
       SelectedSolver: ode15s

  Show all properties


Solve the system of equations over the time interval [0 10] by using the solve method. Plot all three solution components.

S = solve(F,0,10);
plot(S.Time,S.Solution,"-o")
legend("y_1","y_2","y_3",Location="southeast")

Version History

Introduced in R2023b