9  Gauss-based Method for separating the algebraic part from the differential part

Author
Affiliations

Enrico Bertolazzi

University of Trento, Department of Industrial Engineering

We are now prepared to explore the utility of the Gauss-like method Theorem 8.1 in determining the index of a Differential-Algebraic Equation (DAE) and transforming a generic linear DAE into an Ordinary Differential Equation (ODE).

9.1 Formulation of the DAE

Consider a generic linear DAE expressed in the following form:

\bm{E}\bm{x}^\prime + \bm{A}\bm{x} = \bm{f}(t). \tag{9.1}

The matrix \bm{E} and \bm{A}\in \Bbb{R}^{n \times n}. Define the rank of \bm{E} as:

\text{rank}~\bm{E} = r \leq n.

9.2 Change of Coordinates

Now, consider a pair of matrices \bm{L}^{(0)} \in \Bbb{R}^{n \times n} and \bm{R}^{(0)} \in \Bbb{R}^{n \times n} such that:

\bm{L}^{(0)} \bm{E} \bm{R}^{(0)} = \bm{E}^{(0)} = \begin{bmatrix} \bm{\Lambda}^{(0)} & 0 \\ 0 & 0 \end{bmatrix}.

As in the previous section, the superscript ^{(0)} indicates the iteration number, which corresponds to the Kronecker index of the DAE. At this first iteration, \bm{\Lambda}^{(0)} \in \Bbb{R}^{r \times r}.

If r = n, we are done since the DAE has index \nu = 0, and the matrix \bm{E} is full rank and invertible. If r < n, we introduce the following change of coordinates:

\bm{x}^{(0)} = \bm{R}^{(0)}\bm{x}, \qquad (\bm{x}^{(0)})' = \bm{R}^{(0)}\bm{x}',

and define the matrix \bm{A}^{(0)} and vector \bm{f}^{(0)} as follows:

\begin{aligned} \bm{A}^{(0)} &:= \bm{L}^{(0)}\bm{A}\bm{R}^{(0)}, \\ \bm{f}^{(0)} &:= \bm{L}^{(0)}\bm{f}. \end{aligned}

so that

\bm{E}^{(0)}(\bm{x}^{(0)})' + \bm{A}^{(0)}\bm{x}^{(0)} = \bm{f}^{(0)}(t). \tag{9.2}

obtained by multiply the left-hand side of Equation 9.1 by \bm{L}^{(0)}.

Special Structure of the Matrix

Notice that by partitioning

(\bm{x}^\prime)_{1}^{(0)},\bm{x}_{1}^{(0)} \in \Bbb{R}^{r}, \quad (\bm{x}^\prime)_{2}^{(0)},\bm{x}_{2}^{(0)} \in \Bbb{R}^{n-r}, \quad \text{and} \quad \bm{x}^{(0)} = \begin{bmatrix} \bm{x}_{1}^{(0)}\\ \bm{x}_{2}^{(0)} \end{bmatrix}.

then Equation 9.2 has the following special structure:

\begin{bmatrix} \bm{\Lambda}^{(0)} & \bm{0} \\ \bm{0} & \bm{0} \end{bmatrix} \begin{bmatrix} (\bm{x}^\prime)_{1}^{(0)}(t)\\ (\bm{x}^\prime)_{2}^{(0)}(t) \end{bmatrix} + \begin{bmatrix} \bm{A}_{11}^{(0)} & \bm{A}_{12}^{(0)} \\ \bm{A}_{21}^{(0)} & \bm{A}_{22}^{(0)} \end{bmatrix} \begin{bmatrix} \bm{x}_{1}^{(0)}\\ \bm{x}_{2}^{(0)} \end{bmatrix} = \begin{bmatrix} \bm{f}_{1}^{(0)}(t)\\ \bm{f}_{2}^{(0)}(t) \end{bmatrix}. \tag{9.3}

The equation Equation 9.3 can also be expressed as a pair of matrix equations:

\begin{aligned} \bm{\Lambda}^{(0)}(\bm{x}^\prime)_{1}^{(0)} + \bm{A}_{11}^{(0)}\bm{x}_{1}^{(0)} + \bm{A}_{12}^{(0)}\bm{x}_2^{(0)} &= \bm{f}_{1}^{(0)}(t), \qquad & (a) \\ \bm{A}_{21}^{(0)}\bm{x}_{1}^{(0)} + \bm{A}_{22}^{(0)}\bm{x}_{2}^{(0)} &= \bm{f}_{2}^{(0)}(t), \qquad & (b) \end{aligned} \tag{9.4}

Decoupling of Components

Notice that

Hence, in these coordinates, the differential and algebraic components of the DAE are decoupled.

Differentiation of Algebraic Constraints

Now we differentiate Equation 9.4(b) with respect to time:

\bm{A}_{21}^{(0)}(\bm{x}^\prime)_{1}^{(0)} + \bm{A}_{22}^{(0)}(\bm{x}^\prime)_{2}^{(0)} = (\bm{f}_2^\prime(t))^{(0)}. \tag{9.5}

Writing Equation 9.5 with Equation 9.4(a) in matrix form yields:

\begin{bmatrix} \bm{\Lambda}^{(0)} & \bm{0} \\ \bm{A}_{21}^{(0)} & \bm{A}_{22}^{(0)} \end{bmatrix} \begin{bmatrix} (\bm{x}^\prime)_{1}^{(0)}\\ (\bm{x}^\prime)_{2}^{(0)} \end{bmatrix} + \begin{bmatrix} \bm{A}_{11}^{(0)} & \bm{A}_{12}^{(0)} \\ \bm{0} & \bm{0} \end{bmatrix} \begin{bmatrix} \bm{x}_{1}^{(0)}\\ \bm{x}_{2}^{(0)} \end{bmatrix} = \begin{bmatrix} \bm{v}_{1}^{(0)}\\ (\bm{f}_2^\prime(t))^{(0)} \end{bmatrix}. \tag{9.6}

9.3 Rank Increase Observation

As observed, the differentiation with respect to time of the algebraic constraints in the matrix representation produces merely a shift to the left of the two blocks \bm{A}_{21}^{(0)} and \bm{A}_{22}^{(0)} while preserving the original structure.

Thus, in the transition from the DAE Equation 9.4 to the differentiated version Equation 9.6, we observe an increase in rank at each iteration while maintaining the structure of the equations. The final rank gives the Kronecker index of the original DAE.

Definition 9.1 The minimum integer such that \textrm{dim}~\bm{\Lambda}^{(\nu)} = n is the Kronecker index (or differential index) of the DAE.

In this way, we are able to compute the index of a DAE systematically through a sequence of transformations, eventually leading to a formulation where the components are decoupled, thus facilitating the analysis and solution of the DAE.

NoteExample

Consider the following set of linear differential equation,

\begin{aligned} x^\prime_{1}(t) + x^\prime_{2}(t) + x^\prime_{3}(t) &= 2x_{1}(t) + x_{3}(t) + t, \\[0.25em] x^\prime_{1}(t) + x^\prime_{2}(t) + x^\prime_{3}(t) &= x_{1}(t) + x_{2}(t) + t^{3}, \\[0.25em] 0 &= x_{1}(t) + x_{2}(t) + x_{3}(t) + \sin(t), \end{aligned}

that we can write in a more compact matrix form as follows,

\begin{bmatrix} 1 & 1 & 1 \\ 1 & 1 & 1 \\ 0 & 0 & 0 \end{bmatrix} \begin{bmatrix} x^\prime_{1}\\ x^\prime_{2}\\ x^\prime_{3} \end{bmatrix} - \begin{bmatrix} 2 & 0 & 1\\ 1 & 1 & 0\\ 1 & 1 & 1 \end{bmatrix} \begin{bmatrix} x_{1} \\ x_{2} \\ x_{3} \end{bmatrix} = \begin{bmatrix} t \\ t^3 \\ \sin(t) \end{bmatrix}.

Which corresponds to the usual form \bm{E}\bm{x}^\prime(t)+\bm{A}\bm{x}(t)=\bm{f}(t) obvious meaning for the symbols. Now consider to form the table

\begin{bmatrix} \bm{E} & \bm{A} & \bm{f}(t) \end{bmatrix}\in\Bbb{R}^{n\times 2 n+1}

in such a way that,

\begin{bmatrix} \bm{E} & \bm{A} & \bm{f}(t) \end{bmatrix} \begin{bmatrix} \bm{x}^\prime \\ \bm{x} \\ -1 \end{bmatrix} = \bm{0}

First it is necessary to check if the pair (\bm{E},\bm{A}) is a regular pencil,

\begin{vmatrix} \bm{A}-\lambda\bm{E} \end{vmatrix} = \begin{vmatrix} -2-\lambda & -\lambda & -1-\lambda \\ -1-\lambda & -1-\lambda & -\lambda \\ -1 & -1 & -1 \end{vmatrix} = -2

the determinant is different from zero, so the pair (\bm{E},\bm{A}) is a regular pencil.

Now it is easy to manipulate the following table

\begin{bmatrix} 1 & 1 & 1 & -2 & 0 & -1 & t\\ 1 & 1 & 1 & -1 & -1 & 0 & t^3\\ 0 & 0 & 0 & -1 & -1 & -1 & \sin(t) \end{bmatrix}

subtract first row to the second row

\begin{bmatrix} 1 & 1 & 1 & -2 & 0 & -1 & t\\ 0 & 0 & 0 & 1 & -1 & 1 & t^3-t\\ 0 & 0 & 0 & -1 & -1 & -1 & \sin(t) \end{bmatrix}

constraints:

\color{blue} \begin{aligned} x_1 - x_2 + x_3 - t^3 + t & = 0 \\ x_1 + x_2 + x_3 + \sin(t) & = 0 \end{aligned}

differentiate second and third row (reduce index by one)

\begin{bmatrix} 1 & 1 & 1 & 2 & 0 & 1 & t\\ -1 & 1 & -1 & 0 & 0 & 0 & 3t^2-1\\ -1 & -1 & -1 & 0 & 0 & 0 & \cos(t) \end{bmatrix}

add first row to the second and the third row

\begin{bmatrix} 1 & 1 & 1 & 2 & 0 & 1 & t\\ 0 & 2 & 0 & 2 & 0 & 1 & 3t^2+t-1\\ 0 & 0 & 0 & 2 & 0 & 1 & \cos(t)+t \end{bmatrix}

constraint:

\color{blue} 2x_1 + x_3 - \cos(t)-t = 0

derive third row (reduce index by one)

\begin{bmatrix} 1 & 1 & 1 & 2 & 0 & 1 & t\\ 0 & 2 & 0 & 2 & 0 & 1 & 3t^2+t-1\\ 2 & 0 & 1 & 0 & 0 & 0 & -\sin(t)+1 \end{bmatrix}

subtract 2 times first row to third

\begin{bmatrix} 1 & 1 & 1 & 2 & 0 & 1 & t\\ 0 & 2 & 0 & 2 & 0 & 1 & 3t^2+t-1\\ 0 & -2 & -1 & -4 & 0 & -2 & -\sin(t)+1-2t \end{bmatrix}

add second row to third

\begin{bmatrix} 1 & 1 & 1 & 2 & 0 & 1 & t\\ 0 & 2 & 0 & 2 & 0 & 1 & 3t^2+t-1\\ 0 & 0 & -1 & -2 & 0 & -1 & 3t^2-t-\sin(t) \end{bmatrix}

add third row to first and then multiply third row by (-1) and second row by (1/2)

\begin{bmatrix} 1 & 1 & 0 & 0 & 0 & 0 & 3t^2-\sin(t)\\ 0 & 1 & 0 & 1 & 0 & 1 & \frac{1}{2}(3t^2+t-1)\\ 0 & 0 & 1 & 2 & 0 & 1 & t+\sin(t)-3t^2 \end{bmatrix}

subtract second row to first

\begin{bmatrix} 1 & 0 & 0 & -1 & 0 & -1 & \frac{3}{2}t^2+\frac{1}{2}(1-t)-\sin(t)\\ 0 & 1 & 0 & 1 & 0 & 1 & \frac{1}{2}(3t^2+t-1)\\ 0 & 0 & 1 & 2 & 0 & 1 & 3t^2-t-\sin(t) \end{bmatrix}

then we obtain the ODE

\color{DarkGreen} \left\{ \begin{aligned} x_1^\prime & = x_1+x_3 + \frac{3}{2}t^2+\frac{1}{2}(1-t)-\sin(t)\\ x_2^\prime & = -x_1 - x_3 +\frac{1}{2}(3t^2+t-1) \\ x_3^\prime & =-2x_1 - x_3+3t^2-t-\sin(t) \end{aligned} \right.

with the hidden invariants (or constraints)

\color{blue} \begin{aligned} x_1 - x_2 + x_3 - t^3 + t & = 0 \\ x_1 + x_2 + x_3 + \sin(t) & = 0 \\ 2x_1 + x_3 - \cos(t)-t & = 0 \end{aligned}