# Second-Order Model Reduction Based on Gramians

Authors: Cong Teng
Journal: Journal of Control Science and Engineering
 (2012)

---

## Abstract

Some new and simple Gramian-based model order reduction algorithms are presented on second-order linear dynamical systems, namely, SVD methods. Compared to existing Gramian-based algorithms, that is, balanced truncation methods, they are competitive and more favorable for large-scale systems. Numerical examples show the validity of the algorithms. Error bounds on error systems are discussed. Some observations are given on structures of Gramians of second order linear systems.


---

## Body

## 1. Introduction

Models of linear dynamical systems are often given in second-order form(1)Mq̈(t)+Gq̇(t)+Kq(t)=B0u(t)y(t)=C0q(t)+D0q̇(t),where M,G,K∈ℝn×n, B0∈ℝn×m, C0, and D0∈ℝp×n are given matrices, u(t)∈ℝm is the given vector of inputs, y(t)∈ℝp is the unknown vector of outputs, q(t)∈ℝn is the unknown vector of internal variables, n is the dimension of the system, and M is assumed to be invertible.Models of mechanical systems in particular are usually of second-order form (1). For such a system, M=MT, G and K=KT are respectively the mass, damping, and stiffness matrices, u(t)∈ℝm is the input, B0u(t)=f(t)∈ℝn is the vector of external forces, and q(t)∈ℝn is the vector of internal variables.Many applications lead to very large models wheren is very large, while m,p≪n. Due to limitations on time and computer storage, it is often difficult or impossible to directly simulate or control these large-scale systems. Therefore, it is often desirable to reduce the system to a smaller-order model:(2)M̂q̈(t)+Ĝq̇(t)+K̂q(t)=B̂0u(t)ŷ(t)=Ĉ0q(t)+D̂0q̇(t),
such that  (1)
∥y-ŷ∥/∥u∥<tol, that is, outputs corresponding to the same inputs are close.(2)
The reduced system preserves important system properties such as stability and passivity.(3)
The procedure is computationally efficient.By lettingx=[qq̇], second-order system (1) is then transformed to a corresponding 2n-dimensional first-order system (3)Σ:{ẋ(t)=Ax(t)+Bu(t)y(t)=Cx(t)+Du(t),where(4)A=[0I-M-1K-M-1G]=[0I-KM-GM]B=[0M-1B0]=[0BM]C=[C0D0]D=0.Since algorithms and theory are well developed for first-order model reduction, it is natural to use these techniques to develop algorithms for second-order systems.Closely related to first-order system are the two Lyapunov equations(5)AP+PA*+BB*=0,A*Q+QA+C*C=0.Under assumption that first-order system is stable, it is well known that the reachability and observability Gramians 𝒫 and 𝒬 are the unique solutions to Lyapunov equation (5), where both 𝒫 and 𝒬 are symmetric positive definite. This Gramian 𝒫 has also the following variational interpretation proposed in [1]. Let J(u,t1,t2)=∫t1t2u*(t)u(t)dt=∫t1t2∥u(t)∥2dt. The optimal value to optimization problem (6)min⁡uJ(u,-∞,0)s.t.ẋ=Ax+Bu,x(0)=x0isx0*𝒫-1x0, that is, the minimal energy required to steer the state of the system from t=-∞ to state x0 at time t=0 is x0*𝒫-1x0. Similarly, the optimal value to its dual system is x0*𝒬-1x0.Partition𝒫 and 𝒬 into four equal blocks: (7)P=(P11P12P21P22),Q=(Q11Q12Q21Q22).In [2–4], it is shown that the optimal value to optimization problem(8)min⁡q̇0∈Rnmin⁡uJ(u,-∞,0)s.t.Mq̈(t)+Gq̇(t)+Kq(t)=Bu(t),q(0)=q0
is q0*𝒫11-1q0, that is, the minimal energy required to reach the given position q0 over all past inputs and initial velocities. And the optimal value to optimization problem(9)min⁡q0∈Rnmin⁡uJ(u,-∞,0)s.t.Mq̈(t)+Gq̇(t)+Kq(t)=Bu(t),q̇(0)=q̇0
is q̇0*𝒫22-1q̇0, that is, the minimal energy required to reach the given velocity q̇0 over all past inputs and initial positions. In [3, 4], 𝒫11 and 𝒫22 are defined to be position Gramian and velocity Gramian, respectively. By duality, 𝒬11 and 𝒬22 are position Gramian and velocity Gramian, respectively corresponding to 𝒬.It is well known that(10)P=∫0∞eAτBBTeATτdτ,Q=∫0∞eATτCTCeAτdτ.Supposex(t) is the solution of first-order system for impulse input u(t)=δ(t). Then(11)P=∫0∞x(t)x*(t)dt.
For second-order system and its corresponding first-order system, (12)P=∫0∞x(t)x*(t)dt=∫0∞(q(t)q̇(t))(q*(t)q̇*(t))dt=∫0∞(q(t)q*(t)q(t)q̇*(t)q̇(t)q*(t)q̇(t)q̇*(t))dt.Therefore,(13)P11=∫0∞q(t)q*(t)dt,P22=∫0∞q̇(t)q̇*(t)dt,where 𝒫11 and 𝒫22 are position and velocity Gramians, respectively. Similar results can be applied to 𝒬, 𝒬11, and 𝒬22 through dual systems.For Gramian based structure preserving mode order reduction of second-order systems, in [3, 5], some balanced truncation methods were derived. In this paper, we present two new algorithms based on SVD rather than balanced truncation. The algorithms are derived from second-order Gramians and system ℋ2-norm. Compared to existing techniques, they are simple and more suitable for large-scale settings. In Section 4, we apply these methods to three numerical examples, the computational results show the validity of our algorithms. Error bounds and structures of Gramians on second-order systems are also discussed.


## 2. Balanced Truncation Method

For first-order system(14)Σ∶{ẋ(t)=Ax(t)+Bu(t),y(t)=Cx(t).
balanced truncation method is to transform the system to another coordinate system such that both 𝒫 and 𝒬 are equal and diagonal: (15)P=Q=diag⁡(σ1,σ2,…,σN),whereσ1≥σ2≥⋯≥σN. It then truncates the model by keeping the states corresponding to k largest eigenvalues of the product 𝒫𝒬, that is, k largest Hankel singular values of the system.The main problem is to find balancing projection matrices. A numerically stable way to get the balancing truncated system is given in [6]: suppose the Cholesky factorization of 𝒫 and 𝒬 are is (16)P=UCUCT,Q=LOLOT,whereUC is upper triangular and LO is lower triangular matrices. Take SVD of UCTLO, (17)UCTLO=UΣVT.LetΣk be the first k by k principle submatrix of Σ, Uk, and Vk consist of the first k columns of U and V, respectively. Then the balanced projection matrices are (18)W=LOVkΣk-1/2,V=UCUkΣk-1/2.The reduced system is obtained in:(19)Â=WTAV,B̂=WTB,Ĉ=CV.Balanced truncation method has a global error bound [1]: (20)‖Σerr⁡‖H∞=‖Σ-Σred‖H∞≤2(σk+1+⋯+σn).For second-order system(21)Mq̈(t)+Gq̇(t)+Kq(t)=B0u(t)y(t)=C0q(t)+D0q̇(t),there are two balanced truncation methods. One was given in [5] which is to balance both 𝒫11 and 𝒬11 to form the projection matrices W1r and V1r and then get the reduced system by keeping r largest eigenvalues of 𝒫11𝒬11. It is equivalent to performing projections to its corresponding first-order system as following:(22)Â=WTAV=(W1rTW1rT)(0I-M-1K-M-1G)(V1rV1r)=(0I-W1r-1M-1KV1r-W1r-1M-1GV1r),B̂=WB=(W1rTW1rT)[0M-1B0]=[0W1rTM-1B0]Ĉ=CV=[C0D0](V1rV1r)=[C0V1rD0V1r].The reduction rules are then(23)M̂=Ir,K̂=W1rTM-1KV1r,Ĝ=W1rTM-1GV1r,B̂0=W1rTM-1B0,Ĉ0=C0V1r,D̂0=D0V1r.Similar results can be applied to 𝒫22 and 𝒬22. This gives the following algorithm.Algorithm 1 (Balanced truncation method based on𝒫11 and 𝒬11 (or 𝒫22 and 𝒬22, resp.) [5]).
(1) Compute𝒫11 and 𝒬11.(2) Compute the balanced truncation matricesW1r,V1r∈ℝn×r such that W1rTV1r=Ir and W1rT𝒫11W1r=V1rT𝒬11V1r=diag(σ1,σ2,…,σr), where σ12,σ22,…,σr2 are the r largest eigenvalues of   𝒫11𝒬11.(3) Perform projection to(M,G,K,B0,C0,D0) and get (M̂,Ĝ,K̂,B̂0,Ĉ0,D̂0):(24)M̂=Ir,K̂=W1rTM-1KV1r,Ĝ=W1rTM-1GV1r,B̂0=W1rTM-1B0,Ĉ0=C0V1r,D̂0=D0V1r.In [5], it also gave two other similar methods. One is to balance 𝒫11 and 𝒬22 to get the projection matrices, the other is to balance 𝒫22 and 𝒬11 to get the projection matrices.Another second-order balanced truncation method is given in [3] which is to balance both 𝒫11 and 𝒬11 to get the projection matrices W1r and V1r by keeping r largest eigenvalues of 𝒫11𝒬11, balance both 𝒫22 and 𝒬22 to get the projection matrices W2r and V2r, then use W=(W1rW2r) and V=(V1rV2r) as projection matrices for the corresponding first-order system: (25)Â=WTAV=(W1rTW2rT)(0I-M-1K-M-1G)(V1rV2r)=(0W1rTV2r-W2rTM-1KV1r-W2rTM-1GV2r),B̂=WTB=(W1rTW2rT)[0M-1B0]=[0W2rTM-1B0]Ĉ=CV=[C0D0](V1rV2r)=[C0V1rD0V2r].In order to let reduced system have the companion form of second-order system, it then takes the following transformation by lettingH-1=W1rTV2r,(26)Ã=(IH-1)Â(IH)=(IH-1)(0W1rTV2r-W2rTM-1KV1r-W2rTM-1GV2r)(IH)=(0I-H-1W2rTM-1KV1r-H-1W2rTM-1GV2rH),B̃=(IH-1)B̂=(IH-1)[0W2rTM-1B0]=[0H-1W2rTM-1B0],C̃=Ĉ(IH)=[C0V1rD0V2r](IH)=[C0V1rD0V2rH].The corresponding algorithm is as follows:Algorithm 2 (Balanced truncation method based on(𝒫11𝒫22) and (𝒬11𝒬22) [3]).
(1) Compute𝒫11, 𝒫22, 𝒬11, and 𝒬22.(2) Compute the balanced truncation matrices for𝒫11 and 𝒬11 to get W1r,V1r∈ℝn×r, and balanced transformation matrices for 𝒫22 and 𝒬22 to get W2r,V2r∈ℝn×r.(3) Perform projection to get reduced system(M̂,Ĝ,K̂,B̂0,Ĉ0,D̂0):(27)M̂=Ir,K̂=H-1(W2rTM-1KV1r),Ĝ=H-1(W2rTM-1GV2r)H,B̂0=H-1(W2rTM-1B0),Ĉ0=C0V1r,D̂0=(D0V2r)H,
where H-1=W1rTV2r.


## 3. Some New Algorithms Based on SVD

In this section, we propose some new algorithms based on SVD method directly.First,ℋ2 norm is an important quantity for bounding linear systems.Lemma 3 (see [7]).
Given a first-order system(28)Σ=(ABC),suppose 𝒫 and 𝒬 are reachability and observability Gramians. Then, ℋ2 norm of the system can be expressed as
(29)‖Σ‖H22=trace{CPCT}=trace{BTQB}.ℋ2 norm for second-order system (1) with D0=0 can then be stated as in the following proposition.Proposition 4.
Given a second-order system (1) with D0=0, suppose
(30)Σ=(ABC)is the transformed first-order system, where
(31)A=[0I-M-1K-M-1G],B=[0M-1B0],C=[C00].Assume reachability and observability Gramians 𝒫 and 𝒬 are divided into four equal blocks, and 𝒫=[𝒫11𝒫12𝒫12T𝒫22], 𝒬=[𝒬11𝒬12𝒬12T𝒬22]. Then the following holds:
(32)‖Σ‖H22=trace{B0TM-TQ22M-1B0}=trace{C0P11C0T}.Proof.
From Lemma3,
(33)‖Σ‖H22=trace{BTQB}=trace{[0(M-1B0)T][Q11Q12Q12TQ22][0M-1B0]}=trace{B0TM-TQ22M-1B0},(34)‖Σ‖H22=trace{[C00][P11P12P12TP22][C0T0]}=trace{C0P11C0T}.From (8), we know the positions which are difficult to reach are spanned by the eigenvectors corresponding to small eigenvalues of 𝒫11. We would like to eliminate those positions in reduced systems. Proposition 4 shows that 𝒫11 can independently decide the system ℋ2 norm. Therefore, we may get the reduced system by keeping the positions corresponding to r largest eigenvalues of 𝒫11. This motivates the following procedure for model order reduction. Suppose SVD of 𝒫11 is 𝒫11=W1SV1T, where S=diag(σ1,σ2,…,σn), and σ1≥σ2≥⋯≥σn>0. So W1,V1 are orthogonal matrices and W1=V1 since 𝒫11 is symmetric positive definite. Let the transformation matrices for the corresponding first-order system be(35)W=(W1W1),V=W,and the projection matrices be(36)Wk=(W1rW1r),Vk=Wk,whereW1r consists of the first r columns of W1. For simplicity, we denote Wk and Vk by W and V, respectively. Then use W and V as projection matrices for the corresponding first-order system, that is, perform projections as in (22) and (23). Same results can be applied to 𝒬22. From the dual system, symmetrically 𝒬11 and 𝒫22 are also crucial in weighting the system ℋ2 norms independently. This results in the following algorithms.Algorithm 5 (Second-order model reduction—SVD on𝒫11 (or, 𝒫22, 𝒬11, and 𝒬22, resp.)).
(1) Compute𝒫11 (or 𝒫22, 𝒬11, and 𝒬22, resp.), and denote it as P for simplicity.(2) Take SVD onP and get P=W1S1W1T. Form the orthogonal projection matrices W1r which consists of the first r columns of W1.(3) Perform projection to get the reduced system(M̂,Ĝ,K̂,B̂0,Ĉ0,D̂0):(37)M̂=Ir,K̂=W1rTM-1KW1r,Ĝ=W1rTM-1GW1r,B̂0=W1rTM-1B0,Ĉ0=C0W1r,D̂0=D0W1r.Besides eliminating the positions which are difficult to reach in reduced system, it is also desirable to delete the velocities that are difficult to reach. From (9), these velocities are spanned by the eigenvectors of 𝒫22 corresponding to small eigenvalues. This motivates the following procedure for model order reduction. Suppose SVD of 𝒫11 is 𝒫11=W1S1V1T, where S1=diag(σ1,σ2,…,σn), σ1≥σ2≥⋯≥σn>0. And SVD of 𝒫22 is 𝒫22=W2S2V2T, where S2=diag(δ1,δ2,…,δn), δ1≥δ2≥⋯≥δn>0. So W1, V1, W2, and V2 are all orthogonal matrices, and W1=V1, W2=V2. Let the transformation matrices for the corresponding first-order system be(38)W=(W1W2),V=W,and the projection matrices be(39)Wk=(W1rW2r),Vk=Wk,whereW1r and W2r consist of the first r columns of W1 and W2, respectively. By using similar idea to Algorithm 2, in order to let the reduced system Â have the companion form of second-order model, we adopt a matrix H such that H-1=W1rTV2r, and then perform the projections as in (26) and (27). Similar results can be applied to 𝒬11 and 𝒬22. This results in the following algorithm.Algorithm 6 (Second-order model reduction—SVD on𝒫11 and 𝒫22 (or 𝒬11 and 𝒬22, resp.)).
(1) Compute𝒫11 and 𝒫22.(2) Take SVD on𝒫11 and 𝒫22 to get 𝒫11=W1S1W1T and 𝒫22=W2S2W2T,  respectively.(3) Perform projection to get reduced system(M̂,Ĝ,K̂,B̂0,Ĉ0,D̂0):(40)M̂=Ir,K̂=H-1(W2rTM-1KW1r),Ĝ=H-1(W2rTM-1GW2r)H,B̂0=H-1(W2rTM-1B0),Ĉ0=C0W1r,D̂0=(D0W2r)H,
where H-1=W1rTW2r.Now, consider Algorithm5 and SISO second-order systems. Note that by taking SVD,  𝒫11=W1S1W1T. Denote S11=diag(σ1,σ2,…,σr) and S12=diag(σr+1,σr+2,…,σn). So S1=diag(S11,S12). Suppose z∶=W1TC0T, that is,
(41)(z1z2⋮zn)=W1T(c1c2⋮cn).Then from Proposition 4, (42)‖Σ‖H22=C0P11C0T=C0W1SW1TC0T=zTSz=σ1z12+σ2z22+⋯+σnzn2.In reduced system,Ĉ0=C0W1r, 𝒫̂11=S1. Therefore, (43)‖Σred‖H22=Ĉ0P̂11Ĉ0T=(C0W1r)S1(W1rTC0T)=σ1z12+σ2z22+⋯+σrzr2.In many cases, the eigenvalues of𝒫11 decay rapidly. Therefore, the reduced system produced by Algorithm 5 takes the main part of the system ℋ2 norm. Actually we also have, (44)‖Σ‖H22-‖Σred‖H22=σr+1zr+12+σr+2zr+22+⋯+σnzn2≤σr+1‖z‖22=σr+1‖C0‖22.This is not theℋ2-norm of error system ∥Σ-Σred∥ℋ2, but it gives some information to support Algorithm 5. This result for SISO second-order systems can easily extends to MIMO systems. Similar results can be applied to methods in Algorithm 6.Now, consider complexity of the algorithms. In our SVD method on𝒫11 in Algorithm 5, after computing 𝒫11 which takes O(n3) time, it then uses one-time O(n3) operation which is SVD on 𝒫11 to get the projection matrices. In our Algorithm 6, after computing 𝒫11 and 𝒫22 each of which takes O(n3) time, it then uses twice O(n3) operations which are SVD on 𝒫11 and SVD on 𝒫22 to get the projection matrices. In balanced truncation method in Algorithm 1, after computing 𝒫11 and 𝒬11, it also has one-time SVD. Besides this, it uses two Cholesky factorizations 𝒫11=UCUCT and 𝒬11=LOLOT which all take time O(n3) in order to get the projection matrices. Note that Algorithm 2 has more flops than Algorithm 1. Therefore, our SVD algorithms are simpler in time cost than existing balanced truncation methods. There are some other techniques in getting the balanced truncation systems, for example, balance-free truncation [8]. But it is not hard to see that they are also more complex than our SVD methods. This makes our SVD methods more suitable to solve large-scale problems.


## 4. Computational Results

In this section, we apply the model reduction methods to three numerical examples: the building model, the clamped beam model, and international space station model, see [9] for detailed descriptions. In Table 1, four model reduction methods are compared: SVD method on 𝒫11 in Algorithm 5 (“SVD1”), SVD method on 𝒫11 and 𝒫22 in Algorithm 6 (“SVD12”), balanced truncation on 𝒫11 and 𝒬11 in Algorithm 1 (“BT1”), and balanced truncation on 𝒫11 and 𝒬11, 𝒫22 and 𝒬22 in Algorithm 5 (“BT12”). The comparison is based on the relative error of ℋ2 norm, where (45)‖ɛ‖2‖H‖2=∶‖Σ-Σ̂‖H2‖Σ‖H2,Table 1
Errors of reduced models produced by four different methods.
Modelnrmp∥ℋ∥2∥ɛSVD1∥2/∥ℋ∥2∥ɛSVD12∥2/∥ℋ∥2∥ɛBT1∥2/∥ℋ∥2∥ɛBT12∥2/∥ℋ∥2Building model244114.5e−034.1e−015.1e−013.2e−013.2e−01Clamped beam17429113.3e+023.3e−021.3e−011.3e−035.7e−04International state space13527331.0e−046.7e−041.0e−043.2e−063.2e−06n is the size of M, r is the size of M̂ in reduced model, and m is the size of input vector u(t), p is the size of output vector y(t). Figures 1 and 2 show the amplitude of frequency response of original and reduced systems for clamped beam model for r=20 and r=29, respectively.Figure 1
Frequency response on clamped beam model with sizen=174 and the reduced model of size r=20 by four different methods.Figure 2
Frequency response on clamped beam model with sizen=174 and the reduced model of size r=29 by four different methods.


## 5. Error Bounds and Discussions

Error bounds for first-order model reduction do not apply for second-order systems. Consider second-order system (1) with D0=0∶(46)Σ∶{Mq̈(t)+Gq̇(t)+Kq(t)=B0u(t)y(t)=C0q(t).Suppose the reduced systemΣ̂ is obtained by keeping r largest eigenvalues of the orthogonal eigenspace of 𝒫11. Sorensen and Antoulas in [10] provided a priori error bound showing that the ℋ2 norm of error system is bounded by a constant times the summation of neglected singular values of 𝒫11.It considers structured systems in [10]. Let Q(s) and P(s) be polynomial matrices in s: (47)Q(s)=∑j=1rQjsj,Qj∈Rn×n,P(s)=∑j=1r-1Pjsj,sPj∈Rn×m,whereQ is invertible, and Q-1P is a strictly proper rational matrix. Denote by Q(d/dt), P(d/dt) the differential operators (48)Q(ddt)=∑j=1rQjdjdtj,P(ddt)=∑j=1r-1Pjdjdtj.The structured systems are defined by the following equations:(49)Σ∶{Q(ddt)x(t)=P(ddt)u(t)y(t)=Cx(t),where C∈ℝp×n. The Gramian is defined as the Gramian of x(t) when the input is an impulse: (50)P∶=∫0∞x(t)x(t)*dt.Let(51)P=VSV*be the eigensystem of𝒫, where V is orthogonal, the diagonal elements of S are in decreasing order: (52)S=(σ1σ2⋱σn),whereσ1≥σ2≥⋯≥σn. Partition V=[V1,V2] where V1 consists of the first r columns of V. The reduced model is derived from(53)Q̂j=V1*QjV1,P̂j=V1*Pj,Ĉ=CV1,and the reduced system is then(54)Σ̂∶{Q̂(ddt)x̂(t)=P̂(ddt)u(t)ŷ(t)=Ĉx̂(t).
Partition accordingly, (55)[C1,C2]=CV,S=(S1S2).Theorem 7 (see[10]).
The reduced modelΣ̂ derived from the dominant eigenspace of the Gramian 𝒫 for Σ as described above satisfies
(56)‖Σ-Σ̂‖H22≤trace{C2S2C2*}+Ktrace{S2},where 𝒦 is a modest constant depending on Σ and Σ̂, and the diagonal elements of S2 are the neglected smallest eigenvalues of 𝒫.Specially for second-order system (1) with D0=0:(57)Σ∶{Mq̈(t)+Gq̇(t)+Kq(t)=B0u(t)y(t)=C0q(t),it is easy to see that it can be described by (49) with Q(s)=Ms2+Gs+K and P(s)=B0. So there is the following corollary.Corollary 8.
For second-order system (57), suppose 𝒫11=(S1S2) is diagonal, where S1=diag(σ1,σ2,…,σr), S2=diag(σr+1,σr+2,…,σn), and σ1≥σ2≥⋯≥σn. Suppose the reduced system Σ̂ is obtained by keeping r largest eigenvalues of 𝒫11 as described in (53) and (54). Then
(58)‖Σerr‖H22=‖Σ-Σ̂‖H22≤c0trace{S2}=c0(σr+1+σr+2+…+σn),where c0 is a modest constant depending on Σ and Σ̂.Now, we explain the corollary and what the constantc0 is. Note Q(s)=Ms2+Gs+K∈ℂn×n and Q(s)=[Q11(s)Q12(s)Q21(s)Q22(s)] with Q11(s)∈ℂr×r. Let W(s)=Q11(s)-1Q12(s). Partition accordingly C0=[C01,C02]. Then, (59)c0=sup⁡ω∈R‖(C01W(iω))*(C01W(iω)-2C02)‖2+‖C02‖22.If the reduced system has no poles on the imaginary axis,supω∈ℝ∥W(iω)∥2 is finite, and so supω∈ℝ∥(C01W(iω))*(C01W(iω)-2C02)∥2 is finite. Therefore, c0 is a constant depending on Σ and Σ̂. There is a similar result for 𝒬11 by duality.Consider Algorithm5, SVD method on 𝒫11. Note that (60)P11=∫0∞q(t)q*(t)dt,and𝒫11=VSV*, [V1,V2]=V. So the reduction rules in Algorithm 5 are exactly the same as (53). Similar result applies for 𝒬11. Therefore, the error systems produced by our Algorithm 5, SVD method on 𝒫11 or 𝒬11, all have good global error bounds. So far, there is no priori error bound for error systems based on 𝒫22 or 𝒬22.Even though balanced truncation methods in Algorithms1 and 2, and SVD method on both 𝒫11 and 𝒫22 (or 𝒬11 and 𝒬22) in Algorithm 6, all have good performance in those three numerical examples, so far no priori error bounds are provided for these methods.For second-order system (1) with D0=0, the following theorem provides some observations on structures of gramians 𝒫 and 𝒬 in corresponding first-order system.Theorem 9.
Given second-order system (1) with D0=0, then in any coordinate system, there is no state space transformation which gives a block diagonal 𝒬.Proof.
From the observability Lyapunov equation(61)ATQ+QA+CTC=0,
one can get
(62)(0-KMTI-GMT)(Q11Q12Q12TQ22)+(Q11Q12Q12TQ22)(0I-KM-GM)+[C0T0][C00]=0.Equating each block gives
(63)-KMTQ12T-Q12KM+C0TC0=0,-KMTQ22+Q11-Q12GM=0,Q12-GMTQ22+Q12T-Q22GM=0.
From the first equation in (63), one can get 𝒬12 is not zero, otherwise C0=0, and then from (61) the solution for this Lyapunov equation is 𝒬=∫0∞eATτCTCeAτdτ=0. This contradicts with 𝒬 being symmetric positive definite. So 𝒬 is not block diagonal.Algorithm2 [3], that is, balanced truncation method on (𝒫11,𝒫22) and (𝒬11,𝒬22), is to balance both 𝒫11 and 𝒬11 to get the projection matrices W1r and V1r by keeping r largest eigenvalues of 𝒫11𝒬11, balance both 𝒫22 and 𝒬22 to get the projection matrices W2r and V2r, then use W=(W1rW2r) and V=(V1rV2r) as projection matrices for the corresponding first-order system. This is equivalent to assuming both 𝒫 and 𝒬 are block diagonal in corresponding first-order system, that is, (64)P=(P11P22),Q=(Q11Q22),in order to keepr largest eigenvalues of 𝒫11𝒬11 and 𝒫22𝒬22 in reduced corresponding first-order system. Similar results also apply for methods in Algorithm 6, SVD method on 𝒬11 and 𝒬22. From Theorem 9, 𝒬 can not be block diagonal, therefore the reduced system obtained from Algorithm 2 or Algorithm 6 may lose some information. But Algorithm 1 [5] (balanced truncation based on (𝒫11,𝒬11) (or (𝒫22,𝒬22) resp.)) and Algorithm 5 (SVD on 𝒫11 (or 𝒫22, 𝒬11, 𝒬22 resp.)), avoid this drawback in assuming 𝒬 being block diagonal, but they only work for either position Gramians or velocity Gramians, while Algorithms 2 and 6 have the advantage in taking into account both position Gramians and velocity Gramians.


## 6. Conclusions

Two Gramian based second-order model reduction methods are proposed in this paper based on SVD on diagonal blocks of Gramians. Even though they are the well known SVD methods, to our knowledge, they are the simplest ways in Gramian based second-order model reduction. Two existing balanced truncation methods [3, 5] were presented in 2006 and 2008, respectively. When compared to these two existing techniques, our SVD methods are competitive in numerical examples and more suitable for large-scale setting problems. By using the result given in [10], there is global error bound for one method, SVD on 𝒫11 or 𝒬11. Except this, so far no priori error bounds are provided for other methods on second-order model reduction, even through they work good in numerical examples. Interesting future work would be to provide global error bounds for these methods.



---

*Source: 302498-2012-04-11.xml*