# Structured backward errors for eigenvalues of linear port-Hamiltonian descriptor systems

When computing the eigenstructure of matrix pencils associated with the passivity analysis of perturbed port-Hamiltonian descriptor system using a structured generalized eigenvalue method, one should make sure that the computed spectrum satisfies the symmetries that corresponds to this structure and the underlying physical system. We perform a backward error analysis and show that for matrix pencils associated with port-Hamiltonian descriptor systems and a given computed eigenstructure with the correct symmetry structure there always exists a nearby port-Hamiltonian descriptor system with exactly that eigenstructure. We also derive bounds for how near this system is and show that the stability radius of the system plays a role in that bound.

## Authors

• 15 publications
• 10 publications
01/17/2022

### Control of port-Hamiltonian differential-algebraic systems and applications

The modeling framework of port-Hamiltonian descriptor systems and their ...
10/06/2021

### A structure preserving shift-invert infinite Arnoldi algorithm for a class of delay eigenvalue problems with Hamiltonian symmetry

In this work we consider a class of non-linear eigenvalue problems that ...
11/09/2021

### Structure-Preserving Linear Quadratic Gaussian Balanced Truncation for Port-Hamiltonian Descriptor Systems

We present a new balancing-based structure-preserving model reduction te...
03/13/2020

### Finding the closest normal structured matrix

Given a structured matrix A we study the problem of finding the closest ...
07/05/2019

### Fast optical absorption spectra calculations for periodic solid state systems

We present a method to construct an efficient approximation to the bare ...
11/12/2018

### Reciprocal and Positive Real Balanced Truncations for Model Order Reduction of Descriptor Systems

Model order reduction algorithms for large-scale descriptor systems are ...
12/09/2019

### Structured backward errors in linearizations

A standard approach to compute the roots of a univariate polynomial is t...
##### This week in AI

Get the week's most popular data science and artificial intelligence research sent straight to your inbox every Saturday.

## 1 Introduction

We study the perturbation analysis of the eigenstructure (finite and infinite eigenvalues, left and right eigenvectors) of matrix pencils associated the passivity analysis of linear time-invariant descriptor systems (generalized state-space systems) of the form

 E˙x(t)=Ax(t)+Bu(t), x(0)=0,y(t)=Cx(t)+Du(t), (1)

where , , and

are vector-valued functions denoting, respectively, the

input, state, and output of the system. Denoting real and complex -vectors ( matrices) by , (, ), respectively, the coefficient matrices satisfy , , , and . Note that we require that input and output dimensions are both equal to ; and that is a square regular pencil , i. e.  does not vanish identically for all .

We will particularly focus on systems that are positive real or passive and their port-Hamiltonian realizations (see next section). Our work is motivated by two applications, the first is the perturbation analysis arising from computational methods to compute the eigenstructure [26, 41, 54, 55] and the second arises from the need to obtain small perturbations that bring the system back to this structure when it has been destroyed in the process of discretization, model reduction, or other computational techniques, [1, 12, 13, 20, 21, 27, 28, 31, 46]. In both applications the eigenstructure of an originally passive system is perturbed due to perturbations in the process. And then one either wants to determine a nearby passive system with the perturbed eigenstructure (if this exists) or one wants to perturb the eigenstructure so that it is that of a nearby passive system [25]. A similar problem arises in stability analysis and the computation of stability radii and smallest pertubations that make a system stable [23, 24, 35]. While most of these mentioned previous works are for standard passive systems, here we deal with descriptor systems, as they arise from the linearization around stationary solutions of general systems of differenential-algebraic equations [11, 17, 33, 41].

Throughout this article we will use the following notation. The Hermitian (or conjugate) transpose (transpose) of a vector or matrix is denoted by (

) and the identity matrix is denoted by

or

if the dimension is clear. We denote the set of Hermitian and skew-Hermitian matrices in

, respectively, by and . Positive definiteness (semi-definiteness) of is denoted by (). The set of all positive definite (positive semidefinite) matrices in is denoted (). With of a Hermitian matrix we denote the triple of integers of numbers of positive, negative and zero eigenvalues of . The real and imaginary parts of a complex matrix are written as and , respectively, and is the imaginary unit. The 2-norm of a matrix will be denoted by and the Frobenius norm by . The Frobenius norm of a list of matrices is defined as .

The eigenstructure of matrix pencils is characterized by the Kronecker canonical form.

###### Theorem 1.1.

Let . Then there exist nonsingular matrices and such that

 S(λE−A)T=diag(Lϵ1,…,Lϵp,L⊤η1,…,L⊤ηq,Jλ1ρ1,…,Jλrρr,Nσ1,…,Nσs), (2)

where the block entries have the following properties:

1. Every entry is a bidiagonal block of size , , of the form

 λ⎡⎢ ⎢⎣10⋱⋱10⎤⎥ ⎥⎦−⎡⎢ ⎢⎣01⋱⋱01⎤⎥ ⎥⎦.
2. Every entry is a bidiagonal block of size , , of the form

 λ⎡⎢ ⎢ ⎢ ⎢ ⎢ ⎢⎣10⋱⋱10⎤⎥ ⎥ ⎥ ⎥ ⎥ ⎥⎦−⎡⎢ ⎢ ⎢ ⎢ ⎢ ⎢⎣01⋱⋱01⎤⎥ ⎥ ⎥ ⎥ ⎥ ⎥⎦.
3. Every entry is a Jordan block of size , , , of the form

 λ⎡⎢ ⎢ ⎢ ⎢ ⎢ ⎢⎣1⋱⋱1⎤⎥ ⎥ ⎥ ⎥ ⎥ ⎥⎦−⎡⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢⎣λj1⋱⋱⋱1λj⎤⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥⎦.
4. Every entry is a nilpotent block of size , , of the form

 λ⎡⎢ ⎢ ⎢ ⎢ ⎢ ⎢⎣01⋱⋱⋱10⎤⎥ ⎥ ⎥ ⎥ ⎥ ⎥⎦−⎡⎢ ⎢ ⎢ ⎢ ⎢ ⎢⎣1⋱⋱1⎤⎥ ⎥ ⎥ ⎥ ⎥ ⎥⎦.

The Kronecker canonical form is unique up to permutation of the blocks.

A value is called a (finite) generalized eigenvalue of if , and is said to be an eigenvalue of if zero is an eigenvalue of . The blocks as in (iii) are associated with the finite eigenvalues of , and the blocks as in (iv) correspond to the eigenvalue . The size of the largest block is called the index of the pencil , where, by convention, if is invertible. The matrix pencil is called regular if and for some , otherwise it is called singular.

### 1.1 Positive-realness, passivity, and port-Hamiltonian systems

By applying the Laplace transform to (1) and eliminating the state, we obtain the transfer function

 T(s):=D+C(sE−A)−1B, (3)

mapping the Laplace transform of to the Laplace transform of . On the imaginary axis , describes the frequency response of the system. We have the following extensions of the concepts of positive realness and passivity to descriptor systems, see e. g., [12, 20].

###### Definition 1.1.
1. A transfer function as in (3) is positive real if it is i) analytic in the open right half complex plane (including ), and ii) for all in the closed right half complex plane. Moreover, is strictly positive real if for all in the closed right half complex plane.

2. A system of the form (1) is passive if there exists a state-dependent storage function, , such that for any the dissipation inequality

 H(x(t1))−H(x(t0))≤∫t1t0R(y(t)Hu(t))dt (4)

holds. If for all , inequality (4) is strict then the system is called strictly passive.

It is well-known, see e. g., [6, 20], that a system with regular pencil that is controllable ( for all ), and observable ( for all ) is passive if and only if it is positive real and stable (all finite eigenvalues of are in the closed left half complex plane, and those on the imaginary axis including are semisimple), and it is strictly passive if and only if it is strictly positive real and asymptotically stable (all finite eigenvalues of are in the open left half complex plane, and the infinite eigenvalues are semisimple).

In recent years, the special class of port-Hamiltonian (pH) realizations of passive systems has received a lot attention. PH systems are a tool for energy-based modeling, see

[50]; with the energy storage function , the dissipation inequality (4) holds and so pH systems are always passive. The (robust) representation of passive systems as pH systems has been analyzed in [6, 7], and in the extension to pH descriptor systems in [8, 42, 49, 51].

###### Definition 1.2.

A linear time-invariant port-Hamiltonian (pH) descriptor system has the generalized state-space form

 E˙x=(J−R)x+(G−P)u,y=(G+P)Hx+(S−N)u, (5)

where the coefficient matrices satisfy the symmetry conditions

 V:=[JG−GHN]=−VH, W:=[RPPHS]=WH≥0, E=EH≥0. (6)

The correspondence with the generalized state-space realization (1) is given via , , , and .

The pH representation seems to be a very robust representation [36, 37], it allows easy ways for structure preserving model reduction [4, 29, 48] and it greatly simplifies optimization methods for computing stability and passivity radii [23, 24, 25, 46].

### 1.2 Eigenstructure computation

To analyze whether a system is passive, one can compute the eigenstructure of the para-Hermitian matrix function (even matrix pencil) , where , which is given by

 S(s):=s⎡⎢⎣0E0−EH00000⎤⎥⎦−⎡⎢⎣0ABAH0CHBHCDH+D⎤⎥⎦. (7)

An often more advantageous representation of this pencil (in the context of pH systems) is obtained by applying a congruence transformation. Consider the unitary matrix

 X:=1√2[InInIn−In],

and , then one can form the specially structured even pencil

 ^S(s):=^XHS(s)^X=s⎡⎢⎣0E0−EH00000⎤⎥⎦−⎡⎢⎣−R−JG−JHR−PGH−PHS⎤⎥⎦. (8)

The system is passive if the pencil (7) (or equivalently the pencil (8)) is regular, has no purely imaginary eigenvalues and is of index at most one, see [20], so computing the eigenvalues and the structure at allows to check passivity.

In view of this fact it is important to understand the perturbation theory and the backward error analysis for the pencils (7) (or equivalently the pencil (8)). In this respect, the advantage of the form (8) is that perturbations can be mapped back directly to the data matrices , while in (7) this holds for the data matrices . In both cases, for the backward error analysis, we should also make sure that an arbitrary perturbation of the pencil can be mapped back in a structured sense to perturbations in the data matrices, meaning (i) that the zero blocks should not be perturbed, (ii) that the perturbed matrices , should remain Hermitian, positive semidefinite, and skew-Hermitian, and (iii) that the repeated block entries should have repeated perturbations as well.

There exist simple and well-conditioned transformations to go back and forth between the two representations (1) and (5), since

 [−RGSJ−P−N]=1√2X[ABDAHCHDH],

and

 [ABDAHCHDH]=√2XH[−RGSJ−P−N].

Thus for a perturbation analysis we can use either of the two system matrices. We will consider the perturbations of the system pencil (7) and then show how to extend them to the pencil (8).

### 1.3 Backward error analysis

Let us assume that we have determined (via a computational method) an approximate eigenstructure of the pencil . A backward error analysis yields that this eigenstructure corresponds to the exact eigenstructure of a perturbed pencil

 (S+ΔS)(s):=s(E+ΔE)−(A+ΔA),

where and is the perturbation level. If the eigenstructure is determined by a backward stable algorithm, then is a small multiple of the machine precision (round-off unit), but in other approximations it may be much larger, e. g., when the perturbation arises from model reduction or other approximations.) But even if the relative perturbation is small, it is likely to destroy the structure present in the original pencil .

In view of this, in this paper, we study the following questions.

1. Does the perturbed (computed) eigenstructure correspond to that of a pencil with the same block and symmetry structure, i. e. that of a port-Hamiltonian descriptor system.

2. If the answer to the first question is positive, then what is the nearest port-Hamiltonian descriptor system that has exactly this eigenstructure?

3. If the answer to the first question is negative, then what is the nearest port-Hamiltonian descriptor system.

Related questions have already been studied in [3, 6, 23, 25, 35, 44, 45, 46] in the context of finding best pH representations of stable and passive systems and the computation of stability and passivity radii of linear time-invariant dynamical systems. However, all these papers mainly deal with the classical port-Hamiltonian systems, i. e. the case ; here we study pH descriptor systems, which have extra properties that need to be incorporated [39, 40].

A lower bound for the backward errors that one can expect is the stability radius of the generalized eigenvalue problem , since pH systems are guaranteed to be stable. The stability radius of a pencil is defined as the smallest perturbation that causes to be on the border of the stability region [22]. In the descriptor case this happens when an eigenvalue reaches the imaginary axis, when the system has an index or when the pencil becomes singular [19].

In general, to characterize the smallest perturbation that makes a pencil singular is an open problem for unstructured descriptor systems [15, 30, 38] and requires very complex optimization methods even in special cases. However, for pH descriptor systems it has recently been shown in [40] that these distances are easily characterized. Actually the distance to singularity is given by the smallest perturbation that generates a common nullspace of , while actually the distance to instability and the structured distance to the nearest problem with an index are the same and are characterized by the smallest perturbation that generates a common nullspace of and under structure preservering perturbations.

The classical stability radius is given by

 ρ(E,A)=inf∥(ΔE,ΔA)∥F{Λ(E+ΔE,A+ΔA)∩ıR≠∅}=infωσn(A−ıωE)/√1+ω2, (9)

and the minimizing perturbation can be constructed from the

-th singular value triple

at the minimizing frequency ,

 ΔE:=ıωσnunvHn/(1+ω2),ΔA:=σnunvHn/(1+ω2). (10)

For large scale system with pH structure, recently a computational method to compute the stability radius has been derived in [2].

The paper is organized as follows. In the next Section 2, we construct a congruence transformation that restores the special structure of the pencil and we compute upper bounds for its departure from the identity. In Section 3 we illustrate the results of Section 2 with a number of numerical experiments. In Section 4 we end with a few concluding remarks.

## 2 Computing structured perturbation matrices that realize backward errors

In this section we address the first question whether an eigenstructure associated with a system of the form (1) corresponds to that of a pencil associated with a pH descriptor system. Assume that is a regular pencil and that i) for all , i. e. the system is controllable, and ii) for all , i. e. the system is observable, see [14, 41] for a detailed discussion. Before we characterize structured backward errors we need the following lemma.

###### Lemma 2.1.

Consider a controllable and observable descriptor system of the form (1) associated with a strictly passive pH descriptor system of the form (5) and with positive definite. If is a regular pencil, then the finite generalized eigenvalues of are symmetric with respect to the imaginary axis and there are exactly semisimple infinite generalized eigenvalues. Moreover,

 Inertia(ıE)={n,n,m},Inertia(A−ıωE)={n+m,n,0},forallω∈R.
###### Proof.

Since and since is regular, the pencil has exactly infinite eigenvalues. Since by the assumption of strict passivity, it follows that and hence the index is one and the finite eigenvalues of are the eigenvalues of the Hamiltonian matrix

 H:=[E−1A00−AHE−H]−[E−1B−CH](DH+D)−1[CBHE−H]

obtained by forming the Schur complement of with respect to the block . It is well-known, see [41, 47] that Hamiltonian matrices have a spectrum that is symmetric with respect to the imaginary axis. The inertia of the Hermitian matrix is clearly , since is invertible. Since we have assumed controllability and observability, it is also well-known [41, 55] that has no purely imaginary eigenvalues. ∎

A similar result as Lemma 2.1 can also be obtained for the case that and/or are only semidefinite. In this case one has to separate the differential and the algebraic equations and one has to make the stronger assumption that the pencil is of index one. This can be achieved via structured staircase forms, see e.g. [14] for general descriptor system and [5] for pH descriptor systems. In the following we treat the case discussed in Lemma 2.1, i. e. we assume that and are positive definite so that we are not on the boundary of the set of passive systems.

When we perturb the pencil , it is clear that we cannot allow for arbitrary perturbations. The symmetry of the finite spectrum follows from the fact that is Hermitian end is skew Hermitian. We will therefore require that the perturbation preserves this, and hence that the backward errors and are also Hermitian and skew-Hermitian, respectively, i. e. that stays an even pencil when perturbed.

If we start to perturb a matrix, then its inertia remains constant in an open neighborhood of the matrix only if it has no zero eigenvalues. Otherwise the inertia will change for arbitrarily small perturbations, unless we impose constraints on the type of perturbations that are allowed. Therefore we will need to impose that our perturbation preserves the rank of the matrix .

When computing the eigenstructure of even pencils such as (7) or (8), then there exist algorithms that guarantee these properties, see [9, 10] and the references therein. We will employ the even implicitly restarted Arnoldi method of [43], in which stays Hermitian, stays skew-Hermitian, and the null-space of is preserved. When our perturbation results from an eigenvalue algorithm, we can therefore assume that the perturbation of the pencil satisfies

 ΔE=⎡⎢ ⎢⎣ΔE11ΔE120ΔE21ΔE220000⎤⎥ ⎥⎦=−ΔHE,ΔA=⎡⎢ ⎢⎣ΔA11ΔA12ΔA13ΔA21ΔA22ΔA23ΔA31ΔA32ΔA33⎤⎥ ⎥⎦=ΔHA. (11)

If the perturbation arises from an approximation of the model (such as discretization or model reduction), then this approximation process needs to be done in such a way, that constraints (such as Kirchhoff’s conditions in networks, or position constraints as in mechanical systems) that result from the physical properties of the system are not destroyed, see [BeaGM19]. If this is done correctly then again the structure (LABEL:strut1) is typically preserved.

### 2.1 Bounds on the structured backward errors

If we use a backward stable algorithm structure preserving algorithm from [43] to compute the eigenstructure, then where is a small multiple of the machine precision (unit round-off), and and have the structure indicated in (11). To see whether the computed eigenstructure is that associated with a pH descriptor system and to compute the backward error, we need to find a transformation that preserves the computed eigenstructure, preserves the structure indicated in (11), annihilates the diagonal blocks and , and also restores the property .

To preserve the computed eigenstructure and the Hermitian character of , we will perform a congruence transformation; and in order to preserve the structure of indicated in (11), we will constrain it to be block lower triangular, i. e. we look for a transformation

 Z:=⎡⎢⎣Z11Z120Z21Z220Z31Z32Z33⎤⎥⎦

such that

and

with . We also require that is as close as possible to the identity matrix, such that remain as small as possible. This suggests that we choose and and look for a submatrix of

 [Z11Z12Z21Z22]:=I2n+Y=I2n+[Y11Y12Y21Y22]

near the identity matrix, and satisfying the matrix equations

 (I+YH)[ΔA11A+ΔA12AH+ΔA21ΔA22](I+Y) = [0A+ΔAAH+ΔHA0], (12) (I+YH)[ΔE11E+ΔE12−EH+ΔE21ΔE22](I+Y) = [0E+ΔE−EH−ΔHE0]. (13)

Removing common terms on both sides and using the notation and , we can rewrite these equations as

 YH[ΔA11AΔAHΔΔA22]+[ΔA11AΔAHΔΔA22]Y+[ΔA11ΔA12ΔA21ΔA22]=[0ΔAΔHA0]−YH[ΔA11AΔAHΔΔA22]Y,

and

 YH[ΔE11EΔ−EHΔΔE22]+[ΔE11EΔ−EHΔΔE22]Y+[ΔE11ΔE12ΔE21ΔE22]=[0ΔE−ΔHE0]−YH[[]cccΔE11EΔ−EHΔΔE22]Y,

in which we need to zero out the diagonal blocks. Considering these equations, it seems reasonable to choose and then solve the remaining quadratic equations

 AΔY21+YH21AHΔ = −ΔA11−YH21ΔA22Y21, (14) EΔY21−YH21EHΔ = −ΔE11−YH21ΔE22Y21, (15) AHΔY12+YH12AΔ = −ΔA22−YH12ΔA11Y12, (16) −EHΔY12+YH12EΔ = −ΔE22−YH12ΔE11Y12 (17)

for the unknowns and . If we decompose and in their Hermitian and skew Hermitian parts, , and , with , , , and , then, using the function which stacks the columns of a matrix in a large vector, we have

 vec(Y12) = vec(W12)+vec(V12),vec(YH12)=vec(W12)−vec(V12), vec(Y21) = vec(W21)+vec(V21),vec(YH21)=vec(W21)−vec(V21).

We can then rewrite the equations (14)–(17) using Kronecker products as

 [In⊗EΔ−¯¯¯¯EΔ⊗InIn⊗AΔ¯¯¯¯AΔ⊗In][vec(W21)+vec(V21)vec(W21)−vec(V21)]=[−vec(ΔE11)−vec(ΔA11)]+O(∥Y∥2), (18) (19)

If we ignore the quadratic terms on the right hand side, then we obtain linear systems that are solvable when the pencils and have no common eigenvalues, see e. g., [34], which is the case when these pencils come from a sufficiently small perturbation of a system which is strictly passive. We have the following result.

###### Lemma 2.2.

Consider the linear systems (18)–(19) with the quadratic terms set to , set

 K1(E,A):=[In⊗E−¯¯¯¯E⊗InIn⊗A¯¯¯¯A⊗In],K2(E,A):=[−In⊗EHET⊗InIn⊗AHAT⊗In], (20)

and let

 ^δ:=max(∥ΔA12∥2,∥ΔE12∥2)<12min{σ2n2(K1(E,A)),σ2n2(K2(E,A))},

where denotes the th singular value of the matrix . Then the solution satisfies the bound

 √2∥(Y21,Y12)∥F≤∥(ΔE11,ΔA11,ΔE22,ΔA22)∥Fmin{σ2n2(K1(E,A)),σ2n2(K2(E,A))}−2^δ. (21)
###### Proof.

Define for , then it follows from standard perturbation theory, see e. g., [32], that

 σ2n2(K1(EΔ,AΔ))≥σ2n2(K1(E,A))−2^δ,σ2n2(K2(EΔ,AΔ))≥σ2n2(K2(E,A))−2^δ.

The bound (21) then follows from the solutions of the linear systems (18)–(19), which can be written as

 √2X[vec(W21)vec(V21)]=[vec(Y21)vec(YH21)]=−^K−11[vec(ΔE11)vec(ΔA11)],√2X[vec(W12)vec(V12)]=[vec(Y12)vec(YH12)]=−^K−12[vec(ΔE22)vec(ΔA22)], (22)

and the fact that for any matrix . ∎

An estimate of the smallest singular values

and is obtained from considering the triple in (10) which yields

 [ıω(uTn⊗uHn)−(uTn⊗uHn)][In⊗E−¯¯¯¯E⊗InIn⊗A¯¯¯¯A⊗In] = σn[(uTn⊗vHn)(vTn⊗uHn)], [ıω(vTn⊗vHn)−(vTn⊗vHn)][−In⊗EHET⊗InIn⊗AHAT⊗In] = σn[(vTn⊗uHn)(uTn⊗vHn)].

It follows from these identities that the smallest singular value of the unperturbed block Kronecker products must be smaller or equal to . Defining the smallest structured singular value of a structured matrix as the smallest structured perturbation that makes it singular, one can expect that this is a very good estimate, since the smallest structured singular value equals the stability radius . The quality of this estimate is illustrated via numerical examples in Section 3.

### 2.2 An iteration solution procedure

The solution of the quadratic equations (14)–(15) in and (16)–(17) in , can be obtained using the iterative schemes

 AΔ[Y21]i+1+[YH21]i+1AHΔ=−ΔA11−[YH21]iΔA22[Y21]i,EΔ[Y21]i+1−[YH21]i+1EHΔ=−ΔE11−[YH21]iΔE22[Y21]i,

and

 AHΔ[Y12]i+1+[YH12]i+1AΔ=−ΔA22−[YH12]iΔA11[Y12]i,−EHΔ[Y12]i+1+[YH12]i+1EΔ=−ΔE22−[YH12]iΔE11[Y12]i.

Using an analysis similar to that of [52], we can show that these iterations converge to a solution of the quadratic equations (14), (15), (16), and (17), see [52, Theorem 2.11, p. 242] and [18]. We obtain the following main result.

###### Theorem 2.1.

Consider the system of matrix equations (14), (15), (16), (17). Let

 δ := min{σ2n2(K1(EΔ,AΔ)),σ2n2(K2(EΔ,AΔ))}−2max{∥ΔA12∥2,∥ΔE12∥2}, θ := ∥(ΔE11,ΔE22,ΔA11,ΔA22)∥F, ω := √2∥∥(ΔA11,ΔA22,ΔE11,ΔE22)∥∥F.

If and , then there exists a solution of these equations satisfying

 ∥(Y12,Y21)∥F≤2θ/δ. (23)
###### Proof.

Lemma 2.2 and the assumption guarantee that the linear system of matrix equations (22) is solvable. If we write its solution in terms of the representation with the matrices for and with for , then we obtain the bound

 ∥([W21]1,[V21]1,[W12]1,[V12]1)∥F≤∥(ΔE11,ΔE22,ΔA11,ΔA22)∥Fδ=θδ=:ρ0,

using Lemma 2.2. The iterative schemes can then be written as

 [vec([W21]i+1)vec([V21]i+1)] = [vec([W21]0)vec([V21]0)] (24) +(^K1√2X)−1[vec([W21−V21]iΔE