Choose a Jacobian Method for an Implicit Solver

For implicit solvers, Simulink® must compute the solver Jacobian, which is a submatrix of the Jacobian matrix associated with the continuous representation of a Simulink model. In general, this continuous representation is of the form:

`$\begin{array}{l}\stackrel{˙}{x}=f\left(x,t,u\right)\\ y=g\left(x,t,u\right).\end{array}$`

The Jacobian, J, formed from this system of equations is:

`$J=\left(\begin{array}{cc}\frac{\partial f}{\partial x}& \frac{\partial f}{\partial u}\\ \frac{\partial g}{\partial x}& \frac{\partial g}{\partial u}\end{array}\right)=\left(\begin{array}{cc}A& B\\ C& D\end{array}\right).$`

In turn, the solver Jacobian is the submatrix, ${J}_{x}$.

`${J}_{x}=A=\frac{\partial f}{\partial x}.$`

Sparsity of Jacobian

For many physical systems, the solver Jacobian Jx is sparse, meaning that many of the elements of Jx are zero.

Consider the following system of equations:

`$\begin{array}{l}{\stackrel{˙}{x}}_{1}={f}_{1}\left({x}_{1},{x}_{3}\right)\\ {\stackrel{˙}{x}}_{2}={f}_{2}\left({x}_{2}\right)\\ {\stackrel{˙}{x}}_{3}={f}_{3}\left({x}_{2}\right).\end{array}$`

From this system, you can derive a sparsity pattern that reflects the structure of the equations. The pattern, a Boolean matrix, has a 1 for each${x}_{i}$ that appears explicitly on the right-hand side of an equation. Therefore, you obtain:

`${J}_{x,pattern}=\left(\begin{array}{ccc}1& 0& 1\\ 0& 1& 0\\ 0& 1& 0\end{array}\right)$`

The sparse perturbation and the sparse analytical methods may be able to take advantage of this sparsity pattern to reduce the number of computations necessary and improve performance.

Solver Jacobian Methods

When you choose an implicit solver from the Solver pane of the configuration parameters dialog box, a parameter called ```Solver Jacobian method``` and a drop-down menu appear. This menu has five options for computing the solver Jacobian.

Note

If you set Automatic solver parameter selection to `error` in the Solver Diagnostics pane, and you choose a different solver than that suggested by Simulink, you may receive an error.

Limitations

The solver Jacobian methods have these limitations associated with them.

• If you select an analytical Jacobian method, but one or more blocks in the model do not have an analytical Jacobian, then Simulink applies a perturbation method.

• If you select sparse perturbation and your model contains data store blocks, Simulink applies the full perturbation method.

Heuristic 'auto' Method

The default setting for the `Solver Jacobian method` is `auto`. Selecting this choice causes Simulink to determine which of the remaining four methods best suits your model. This flowchart depicts the algorithm.

Sparse methods are beneficial for models that have a large number of states. If 50 or more states exist in your model, `auto` chooses a sparse method. Unlike other implicit solvers, `ode23s` is a sparse method because it generates a new Jacobian at every time step. A sparse analytical or a sparse perturbation method is, therefore, advantageous. Selecting `auto` also ensures that the analytical methods are used only if every block in your model can generate an analytical Jacobian.

Full and Sparse Perturbation Methods

The full perturbation method solves the full set of perturbation equations and uses LAPACK for linear algebraic operations. This method is costly from a computational standpoint, but remains a recommended method for establishing baseline results.

The sparse perturbation method attempts to improve the run-time performance by taking mathematical advantage of the sparse Jacobian pattern. Returning to the sample system of three equations and three states,

`$\begin{array}{l}{\stackrel{˙}{x}}_{1}={f}_{1}\left({x}_{1},{x}_{3}\right)\\ {\stackrel{˙}{x}}_{2}={f}_{2}\left({x}_{2}\right)\\ {\stackrel{˙}{x}}_{3}={f}_{3}\left({x}_{2}\right).\end{array}$`

The solver Jacobian is:

 `$\begin{array}{c}{J}_{x}=\left(\begin{array}{ccc}\frac{\partial {f}_{1}}{\partial {x}_{1}}& \frac{\partial {f}_{1}}{\partial {x}_{2}}& \frac{\partial {f}_{1}}{\partial {x}_{3}}\\ \frac{\partial {f}_{2}}{\partial {x}_{1}}& \frac{\partial {f}_{2}}{\partial {x}_{2}}& \frac{\partial {f}_{2}}{\partial {x}_{3}}\\ \frac{\partial {f}_{3}}{\partial {x}_{1}}& \frac{\partial {f}_{3}}{\partial {x}_{2}}& \frac{\partial {f}_{3}}{\partial {x}_{3}}\end{array}\right)\\ =\left(\begin{array}{ccc}\frac{{f}_{1}\left({x}_{1}+\Delta {x}_{1},{x}_{2},{x}_{3}\right)-{f}_{1}}{\Delta {x}_{1}}& \frac{{f}_{1}\left({x}_{1},{x}_{2}+\Delta {x}_{2},{x}_{3}\right)-{f}_{1}}{\Delta {x}_{2}}& \frac{{f}_{1}\left({x}_{1},{x}_{2},{x}_{3}+\Delta {x}_{3}\right)-{f}_{1}}{\Delta {x}_{3}}\\ \frac{{f}_{2}\left({x}_{1}+\Delta {x}_{1},{x}_{2},{x}_{3}\right)-{f}_{2}}{\Delta {x}_{1}}& \frac{{f}_{2}\left({x}_{1},{x}_{2}+\Delta {x}_{2},{x}_{3}\right)-{f}_{2}}{\Delta {x}_{2}}& \frac{{f}_{2}\left({x}_{1},{x}_{2},{x}_{3}+\Delta {x}_{3}\right)-{f}_{2}}{\Delta {x}_{3}}\\ \frac{{f}_{3}\left({x}_{1}+\Delta {x}_{1},{x}_{2},{x}_{3}\right)-{f}_{3}}{\Delta {x}_{1}}& \frac{{f}_{3}\left({x}_{1},{x}_{2}+\Delta {x}_{2},{x}_{3}\right)-{f}_{3}}{\Delta {x}_{2}}& \frac{{f}_{3}\left({x}_{1},{x}_{2},{x}_{3}+\Delta {x}_{3}\right)-{f}_{3}}{\Delta {x}_{3}}\end{array}\right)\end{array}$`

It is, therefore, necessary to perturb each of the three states three times and to evaluate the derivative function three times. For a system with n states, this method perturbs the states n times.

By applying the sparsity pattern and perturbing states x1 and x 2 together, this matrix reduces to:

 `${J}_{x}=\left(\begin{array}{ccc}\frac{{f}_{1}\left({x}_{1}+\Delta {x}_{1},{x}_{2}+\Delta {x}_{2},{x}_{3}\right)-{f}_{1}}{\Delta {x}_{1}}& 0& \frac{{f}_{1}\left({x}_{1},{x}_{2},{x}_{3}+\Delta {x}_{3}\right)-{f}_{1}}{\Delta {x}_{3}}\\ 0& \frac{{f}_{2}\left({x}_{1}+\Delta {x}_{1},{x}_{2}+\Delta {x}_{2},{x}_{3}\right)-{f}_{2}}{\Delta {x}_{2}}& 0\\ 0& \frac{{f}_{3}\left({x}_{1}+\Delta {x}_{1},{x}_{2}+\Delta {x}_{2},{x}_{3}\right)-{f}_{3}}{\Delta {x}_{2}}& 0\end{array}\right)$`

The solver can now solve columns 1 and 2 in one sweep. While the sparse perturbation method saves significant computation, it also adds overhead to compilation. It might even slow down the simulation if the system does not have a large number of continuous states. A tipping point exists for which you obtain increased performance by applying this method. In general, systems having a large number of continuous states are usually sparse and benefit from the sparse method.

The sparse perturbation method, like the sparse analytical method, uses UMFPACK to perform linear algebraic operations. Also, the sparse perturbation method supports both RSim and rapid accelerator mode.

Full and Sparse Analytical Methods

The full and sparse analytical methods attempt to improve performance by calculating the Jacobian using analytical equations rather than the perturbation equations. The sparse analytical method, also uses the sparsity information to accelerate the linear algebraic operations required to solve the ordinary differential equations.

For details on how to access and interpret the sparsity pattern in MATLAB®, see Exploring the Solver Jacobian Structure of a Model.

Code Generation Support

While the sparse perturbation method supports RSim, the sparse analytical method does not. Consequently, regardless of which sparse method you select, any generated code uses the sparse perturbation method. This limitation applies to rapid accelerator mode as well.