Data-Driven Identification of Rayleigh-Damped Second-Order Systems

Data-Driv en Identification of Ra y leigh-Damped Second-Or der Sy stems Igor P ontes Duff, Pa w an Go yal, and Peter B enner Ab stract In this p a p er , we presen t a data-d r iven approach to identify s econd - order sys tems, ha ving internal Ra yleigh damping. This means that th e dampin g matr ix is given as a linear combin a tion of the mass and stiffness matr ices. These systems typically appear when per forming v ario us engineer ing s tudies, e.g., vib rational and structur al analy sis. In an exper im ental setup, the frequen cy response of a s ys tem can be measured via v ario us approaches, f or instance, by measurin g the vib rations using an accelerometer. As a consequence, giv en frequency samples, the iden tifi- cation of the underl ying system r e lies on rational approximation. T o that aim, we propo se an identification of the cor r esponding second-or der sy stem, e xtendin g the Loewner frame w ork f or this class of sy stems. The efficien cy of the proposed method is demonstrated b y means of v ar ious numer ical b enchmarks. 1 Introdu ction In this paper , we discuss a data-d r iven identification fra m ew ork f or a class of second- order ( SO ) sys tems of the f or m: Σ SO : =  M Ü x ( t ) + D Û x ( t ) + Kx ( t ) = Bu ( t ) , y ( t ) = Cx ( t ) , (1) Igor P ontes Duff, Pa w an Go yal are with the Max Plan ck Ins titute f or Dynamics of Comple x T ec hnical Sys tems, Sandtors traße 1, 39106 Magdeburg, German y; (e-mail: \protect\T1\textbraceleftpontes, goyalp\protect\T1\textbraceright@mpi- magdeburg.mpg.de ). P eter Benner i s wit h the M ax Planck Institute f or Dy namics of Complex T echnical Sys t ems, Sandtors tr . 1, 39106 Magdebu rg, German y , and also with t he T echn ische U niv ersität Che m- nitz, F aculty of Mathematics, Reich enhainer Straße 41, 09126 Chemnitz, Germany (email: benner@mp i-magdebur g.mpg.de). 1 2 Igor Pontes Duff, Pa wan Go y al, and P eter Benn er where x ( t ) ∈ R n is th e state vector , u ( t ) ∈ R m are th e inputs, y ( t ) ∈ R p are the out- puts or measurements, and M , D , K ∈ R n × n are, respectivel y , the mass matr ix, the damping matr ix and the stiffness ma tr ix , B ∈ R n × m and C ∈ R p × n . F or simp licity , we address the problem f or single-inpu t s ingle-ou tput (SISO) systems, i.e., m = p = 1. The multiple-inpu t multip le-output (MIMO) g eneralization is straightforward and can be done by follo win g th e lines of the MIMO e xten sion of t he classical Loewner framew o rk [16] based on tangential in terp olation. Su ch systems arise in man y engi- neer ing ap plications, inclu ding vibratio n analy sis [17], st r uctural d ynamics [9] and electr ic circuits. W e denote the S O sys tems (1) b y Σ SO = ( M , D , K , B , C ) . Moreo v er , we assum e a zero inhom ogeneous co ndition, i.e., x ( 0 ) = Û x ( 0 ) = 0. Hence, by me ans of the Laplac e transform , the in put-ou tput behavior of th e system Σ SO is associated with the transfer fu nction as f ollow s: H SO ( s ) = C  s 2 M + s D + K  − 1 B . (2) Furth e r m ore, thro ughou t the paper, we assum e the prop or tional Ray leigh damping hypothesis, i.e., the dampin g matrix D is giv en by a linear combinatio n of the mass and stiffness matrices: D = α M + β K , (3) f or α , β ≥ 0. This hypothesis is of ten consider ed in sev eral engin e e rin g applica- tion, where the d amper is numer ically constructed in order to a v oid non-d ampened oscillations, see [17] for mo re details. In the past twenty years, model order reduction of SO sys tems has been inv estigated e xtensiv ely ; see f or instance [18, 7, 21] for balancin g -type methods, and [6, 4, 26, 3] f or moment matching and H 2 -optimality based metho ds. R ecen tl y , the authors in [22] pro vided an extensiv e comp a riso n among comm on method s f or SO model order reduction applied to a larg e - scale mechanic a l ar tificial fishtail model. In all of th e abov e- mentioned w orks, the authors supp ose that the y hav e ac c e ss to the matr ices, defining th e or iginal sy stems and the reduced-order sy stems are constructed via Petro v -Galerkin projections. Thus, the m ain goal is to fin d projec tio n matrices V , W ∈ R n × r , leading to the SO red u ced-ord er sys tem ˆ H SO ( s ) = ˆ C  s 2 ˆ M + s ˆ D + ˆ K  − 1 ˆ B , (4 ) with ˆ M = W T MV , ˆ D = W T D V , ˆ K = W T KV , ˆ B = W T B and ˆ C = CV . Ho we ver , it is n ot necessary that the realizatio n is given or is f easible to ob tain ; thus, w e suppose that the or iginal system realization may not b e av ailab le. Instead, we assume to h av e ac cess only to f requency domain data, e.g ., ar ising from e xper iments or numer ical sim u lations. More precisel y , w e are inter e sted in sol v ing the follo wing problem . Data-Driv en I dentification of Ra y leigh-Damped Second-Order Systems 3 Problem 1 ( SO data-driven identification) Given in ter polation data {( σ i , ω i ) | σ i ∈ C and ω i ∈ C , i = 1 , . . . , ρ } , (5) constr u ct a S O realization Σ SO = ( M , D , K , B , C ) of app ropr iate dimensions, satisfying the prop o r tional Ray leigh damping h ypothesis, i. e., D = α M + β K , whose transf er function H SO ( s ) : = C ( s 2 M + s D + K ) − 1 B satisfies the interpo- lation condition s, i.e., H SO ( σ i ) = ω i , i = 1 , . . . ρ. (6) Problem 1 c o r responds to an identification problem wh ich aims at d eter mining a SO re alization th at not on l y inter p olates at given measuremen ts, but also satisfies the Rayleigh damping hypothesis. A similar pro blem f or time-d elay systems was studied in [2 0] and [23]. F ur ther more, we would lik e to mention that a data-dr iv en approa ch for structured n on-para m etr ic systems has bee n studied in [24]. Ho w ev er, the con str uction of th e s tru ctured reduced -order system is not a straightforward task. The pur po se of this pap er is thus to e xtend the app lica tio n do m ain of the Lo e wner framew o rk es tablished in [15, 16] to S O systems. With this aim, a ne w SO Loe wner framew o rk is dev elo ped, yield ing a Rayl eigh damped SO system of the f or m (2) that inter polates at giv en frequency measurements. The rest of the p a per is organized as f o llo w s. Section 2 reca lls some prelimin ar y results on the rational interpo lation Loewner fr amew ork p roposed in [16]. Section 3 presents an e xtension of these results to the class of Ra yleigh damp ed SO s ys tems. The section is div ided into tw o par ts. Th e first on e assumes th e know led ge of the Rayleigh damping par ameters, α and β , an d d e rives the Lo ewner matr ices f or SO sy stems. The seco n d part presents a heur istic pro cedure, o r iginally prop osed in [23] in the conte xt of time-delay sy stems, enablin g us to estimate the parame ters α and β . F inally , Section 4 illustrates the pro p osed fra m ew ork b y numer ical examples and Section 5 conclu d es the paper . 2 Classical Loe wner frame w ork In this section, we b r iefly recall th e Lo ewner fr amew ork [ 1 6]. A first-order ( F O ) sys tem Σ FO = ( E , A , B , C ) is a dynamica l sys tem of the f or m: Σ FO : =  E Û x ( t ) = Ax ( t ) + Bu ( t ) , x ( 0 ) = 0 , y ( t ) = Cx ( t ) , (7) 4 Igor Pontes Duff, Pa wan Go y al, and P eter Benn er with E , A ∈ R n × n , B ∈ R n × m and C ∈ R p × n , and the leading dimension n is th e order of the sys tem. For clarity of exposition, we focus f or n o w on the single-input single-outp ut (SISO) case, i.e., wh e n m = p = 1. The sy stem (7) is associated with the transf er function giv en b y H FO ( s ) = C ( s E − A ) − 1 B . (8) There e xist se v eral M OR techniques f or fir st-order systems such as e xp licit m o- ment matching [27, 25], imp licit mo m ent matching u sing Kr ylo v sub spaces [10, 13], Syl vester equations b ased m ethod [12], extensions for MI M O sys tems [11]. W e refer the read e r to the boo k s [2, 5] for more details. Ho we v er , our goal lies in the iden- tification of linear systems using only the frequency data. H e n ce, th e identification problem , in its SISO form , is s tated as f ollo w s. Problem 2 (First-order da ta-driven model reduction) Giv en interp olation data {( σ i , ω i ) | σ i ∈ C and ω i ∈ C , i = 1 , . . . , ρ } (9) constr u ct a minimal-or der realization Σ = ( E , A , B , C ) of appropr iate dim en- sions, whose transf er function H FO ( s ) = C ( s E − A ) − 1 B satisfies the inter p o la- tion condition s H FO ( σ i ) = ω i , i = 1 , . . . ρ. (10) A wide ran ge o f methods has b e en dev elop ed to sol v e Proble m 2, e .g., v ector fitting [ 14], the AAA algor ithm [1 9] and the Loewner framew o rk [16]. In this pap er , we focus o n the latter approach and, in what f ollow s, w e recall some of the results contained therein. Firstl y , w e assume that the nu mber of inter polatio n data is ev en, i.e., ρ = 2 ℓ , and as a result, the data can b e par titioned in two d isjoint sets as follo w s: r ight inter polation set P r : { ( λ i , w i ) | λ i ∈ C and w i , i = 1 , . . . , ℓ } , and (11a) left inter polatio n set P l : { ( µ j , v j ) | µ j ∈ C and v j ∈ C , j = 1 , . . . , ℓ } . (11b) Using this par tition, we associate the f ollowing Loewner matrices. Definition 1 (Loewner matrices [16]) Giv en the right P r and left P l inter po- lation sets, w e associate them with the Loewner matr ix L and shifted Loe wner matr ix L σ giv en b y L = ©   « v 1 − w 1 µ 1 − λ 1 · · · v 1 − w ℓ µ 1 − λ ℓ . . . . . . . . . v ℓ − w 1 µ ℓ − λ 1 · · · v ℓ − w ℓ µ ℓ − λ ℓ ª ® ® ® ¬ , L σ = ©    « µ 1 v 1 − λ 1 w 1 µ 1 − λ 1 · · · µ 1 v 1 − λ ℓ w ℓ µ 1 − λ ℓ . . . . . . . . . µ ℓ v ℓ − λ 1 w 1 µ ℓ − λ 1 · · · µ ℓ v ℓ − λ ℓ w ℓ µ ℓ − λ ℓ ª ® ® ® ¬ . (12) Data-Driv en I dentification of Ra y leigh-Damped Second-Order Systems 5 Remark 1 T he L oewner matrix L was intr o duced in [1]. As shown therein, its useful- ness der ives from the fact that its rank is equal to the or der of th e minimal realization H FO satisfying the inter polation cond itions in (10). Hence, it re v eals the comple xity of the reduced -order model s olving Problem 2 . Ne x t, let us introduc e the f ollo wing matr ices associated with the inter po lation problem as f ollo ws :  Λ = diag ( λ 1 , . . . , λ ℓ ) ∈ C ℓ × ℓ ˆ H ( Λ ) =  w 1 . . . w ℓ  T ∈ C ℓ × 1 and  M = diag ( µ 1 , . . . , µ ℓ ) ∈ C ℓ × ℓ ˆ H ( M ) =  v 1 . . . v ℓ  T ∈ C ℓ × 1 (13) Also, let 1 ∈ R ℓ × 1 be th e column vector with all entr ies equal to o ne. Hence, the Loewner matr ices satisfy the follo win g Sylv ester equations M L − L Λ = ˆ H ( M ) 1 T − 1 ˆ H ( Λ ) T , and (14a) M L σ − L σ Λ = M ˆ H ( M ) 1 T − 1 ˆ H ( Λ ) Λ . (14b) An elegant solution f o r Problem 2 based on the Loewner pair ( L , L σ ) w as prop osed in [16]. This is summar ized in the f o llo wing theorem. Theorem 1 (Loewner f ramew ork [ 16]) Le t L and L σ be the Lo ewner matri- ces associated w ith the partition in (13) . If ( L σ , L ) is a regular pencil with no µ i or λ j being an eig envalue, then the matrices ˆ E = − L , ˆ A = − L σ , ˆ B = V , ˆ C = W , pro vides a r e alization ˆ Σ FO = ( ˆ E , ˆ A , ˆ B , ˆ C ) for a minimal or der interpo lant of Problem 2, i.e. , the transf er function ˆ H FO ( s ) = W ( s L σ − L ) − 1 V satisfies the interpolation conditions in (10) . Theorem 1 allow s to ob tain a FO system ˆ H = ( ˆ E , ˆ A , ˆ B , ˆ C ) whose transfer function inter polates rig ht and left data as stated in Prob le m 2. H ow ev er, when mor e data than necessar y are pro vided, th en th e h ypothe sis o f Theorem 1 may no t be satisfied. Hence, a singular -value d ecompo sition (S VD) based procedu re h as been proposed in [16] to find an FO sys tem inter polating the frequency data. Ne x t, recall that a SO sys tem Σ SO = ( M , D , K , B , C ) can be wr itten as a first-order realization, f or instance, as follo w s: H SO _ FO ( s ) = C ( s E − A ) − 1 B , where 6 Igor Pontes Duff, Pa wan Go y al, and P eter Benn er E =  I 0 0 M  , A =  0 I − K − D  , B =  0 B  and C =  C 0  . As a consequence, the classical Loewner fr amew ork presented in Section 2 can be employ ed to find a first-orde r realization . How ev er, the intr insic SO str ucture will no t be preserved in the id entified mod el. But the classical Loewner f ramew ork yields an inf or mation about th e or der of a SO realization fitting the data, wh ich is outlined in the f ollo wing remark. Remark 2 ( Order of SO model) Let u s sup pose that the fr equency data in Problem 2 an d let L be a Loewner ma tr ix giv en in (12) co nstructed with this data. Then, the order of the SO sys tem fitting the data equals 1 2 rank ( L ) . In the follo win g section, we discuss an e xten sion o f the Loewner framew o r k f or the class of Ra yleigh damped SO systems. 3 Second-or der Loe wner F ramew ork This section contains ou r main contribu tion, which presents an e x tension of the Loewner frame w ork to the class of SO Ra yleigh da m ped sys tems (1). Here, we also assume tha t th e num ber of inter polatio n d ata is ev en, i.e., ρ = 2 ℓ , an d the data is p a rtitio n ed into tw o disjoint sets as in (1 1a) and ( 11b). Moreov er , the data is organized into the matr ices Λ , ˆ H ( Λ ) , M , ˆ H ( M ) as in (13). This section is d ivided into two par ts. In the first one, we assume t o ha v e a prior i kn ow led ge of th e Ra yleigh damping p arameters α and β and we der ive the equiv alents of the Loewner ma tr ices (12) and the Theorem 1 to the class of SO Ra yleigh damped sy s tems. The second par t is d edicated to p r oposing a heu r istic procedu re to es timate the parameters α and β using the frequency d ata a vailable. 3.1 Seco n d-order Loe w n er matrices In what follo w s, we assume that Problem 1 has a minimal order r solution H ⋆ SO , g iv en by H ⋆ SO ( s ) = C ⋆  s 2 M ⋆ + s D ⋆ + K ⋆  − 1 B ⋆ , (15) with D ⋆ = α M ⋆ + β K ⋆ . Her e, we also assume th a t the coefficients α and β from the Ray leigh-Dam ped hypothesis are kno wn. Then , later in this section, we will sho w how to constr uct a realizatio n equivalent to H ⋆ SO ( s ) that only depend s o n the frequency data. T o that aim, let us first recall a result from [4] enabling projection- based structu red preservin g model reduction. Data-Driv en I dentification of Ra y leigh-Damped Second-Order Systems 7 Theorem 2 (Structure preserving SO model reduction [4]) Con sider the SO tr ans- fer function H SO ( s ) as given in (2) . F or given in terpolation po ints λ i and µ i , i ∈ { 1 , . . . , ℓ } , le t the pr ojection matrices V and W be as follow s: V = h  λ 2 1 M + λ 1 D + K  − 1 B , . . . ,  λ 2 ℓ M + λ ℓ D + K  − 1 B i (16a) W = h  µ 2 1 M + µ 1 D + K  − T C T , . . . ,  µ 2 ℓ M + µ ℓ D + K  − T C T i (16b) Hence, the r ed uced-order model ˆ H SO ( s ) con structed b y P etro v-Galerkin pr ojection as in (4) satisfies the interpolation conditions H SO ( λ i ) = ˆ H SO ( λ i ) and H SO ( µ i ) = ˆ H SO ( µ i ) , for i = 1 , . . . , ℓ . The ab o v e theorem allow s us to co nstruct a SO r educed-o rder m odel by inter polatio n . Let us ap ply this theor em to the SO system H ⋆ SO ( s ) (1 5). For this, we will constr u ct th e matr ix V using the inter polation po ints in Λ , a n d the matr ix W u sing the inter p olation points in M . As a co nsequence, V and W are, respectivel y , the solutions of the f ollo wing matrix equations M ⋆ V Λ 2 + D ⋆ V Λ + K ⋆ V = B ⋆ 1 T , and (17a) M 2 W T M ⋆ + M W T D ⋆ + W T M = 1 C ⋆ , (17b) Multiplying the equatio ns on the lef t (17a) and (17b) on the left by W T and V T , respectivel y , on e obtains W T M ⋆ V Λ 2 + W T D ⋆ V Λ + W T K ⋆ V = W T B ⋆ 1 T , M 2 W T M ⋆ V + WD ⋆ V M + WK ⋆ V = C ⋆ V 1 T . If we set ˆ M = W T M ⋆ V , ˆ D = W T D ⋆ V , ˆ K = W T K ⋆ V , (18a) ˆ B = W T B ⋆ = ˆ H ( M ) , and ˆ C = C ⋆ V = ˆ H ( Λ ) T , (18b) then the SO system ˆ H SO = ( ˆ M , ˆ D , ˆ K , ˆ B , ˆ C ) is the reduced-or der m odel ob tained by Theorem 2, satisfyin g the inter polatio n conditio ns from Problem 1 . Hence, w e can rewrite the abo v e eq uations as follo ws: ˆ M Λ 2 + ˆ D Λ + ˆ K = ˆ H ( M ) 1 T , M 2 ˆ M + M ˆ D + ˆ K = 1 ˆ H ( Λ ) T . Moreov er, if we apply the Ray legh-Damped hypothesis, i.e., ˆ D = α ˆ M + β ˆ K , w e obtain 8 Igor Pontes Duff, Pa wan Go y al, and P eter Benn er ˆ M  Λ 2 + α Λ  + ˆ K ( β Λ + I ) = ˆ H ( M ) 1 T , (20a)  M 2 + α M  ˆ M + ( β M + I ) ˆ K = 1 ˆ H ( Λ ) T . (20b) Notice that the abo v e equation s can be so l v ed for ˆ M and ˆ K . Ho w e v er , in order to hav e an analytic expression for the matrices of the reduced - order system in a similar wa y as f or the Loe wner matrices (12), w e ne e d to introduce the f o llo wing c hange of variables: L SO : = −( I + β M ) ˆ M ( I + β Λ ) , L SO σ : = ( I + β M ) ˆ K ( I + β Λ ) , (21a) B SO : = ( I + β M ) ˆ H ( M ) , and C SO : = ˆ H ( Λ ) T ( I + β Λ ) . (21b) Notice t hat the two r ealizations ˆ Σ SO = ( ˆ M , α ˆ M + β ˆ K , ˆ K , ˆ B , ˆ C ) and ˆ Σ Loe w SO = (− L SO , − α L SO + β L SO σ , L SO σ , B SO , C SO ) are equivalent, i.e., they represent the s ame tran sf er function. Hence, the rea liza tio n ˆ H Loe w SO also satisfies th e inter polation conditions fro m Problem 1. Additionally , b y a simple comp utation, w e obtain that the m atr ice s L SO and L SO σ satisfy the follo wing equations L SO F ( Λ ) + L SO σ = −D ( M ) ˆ H ( M ) 1 T , F ( M ) L SO + L SO σ = − 1 ˆ H ( Λ ) T D ( Λ ) , where , f or a given matr ix Ω , F ( Ω ) : = ( I + β Ω ) − 1 ( Ω 2 + α Ω ) an d D ( Ω ) : = ( I + β Ω ) . As a consequen c e , L SO F ( Λ ) − F ( M ) L SO = 1 ˆ H ( Λ ) T D ( Λ ) − D ( M ) ˆ H ( M ) 1 T , (22a) L SO σ F ( Λ ) − F ( M ) L SO σ = 1 ˆ H ( Λ ) T N ( Λ ) − N ( M ) ˆ H ( M ) 1 T , (22b) where, f or a giv en matr ix Ω , N ( Ω ) : = ( Ω 2 + α Ω ) . Notice that the Sy lv ester equa- tions (22) are equivalent to (14a) f or the case of SO sy stems. Hence, using those equation, one can der ive a nalytic e xpression s o f L SO and L SO σ . Definition 2 ( SO L oewner matrices) Let us suppose α and β are known and let d ( s ) : = 1 + s β , n ( s ) : = s 2 + α s , an d f ( s ) : = n ( s ) d ( s ) , be scalar fun c tio ns. Then the SO Loewner matr ices, na m ely , the SO Loewner matr ix L SO and the shifted Loewner m atr ix L SO σ are given by Data-Driv en I dentification of Ra y leigh-Damped Second-Order Systems 9 L SO = ©    « d ( µ 1 ) v 1 − d ( λ 1 ) w 1 f ( µ 1 )− f ( λ 1 ) · · · d ( µ 1 ) v 1 − d ( λ ℓ ) w ℓ f ( µ 1 )− f ( λ ℓ ) . . . . . . . . . d ( µ ℓ ) v ℓ − d ( λ 1 ) w 1 f ( µ ℓ )− f ( λ 1 ) · · · d ( µ ℓ ) v ℓ − d ( λ ℓ ) w ℓ f ( µ ℓ )− f ( λ ℓ ) ª ® ® ® ¬ , (23) L SO σ = ©    « n ( µ 1 ) v 1 − n ( λ 1 ) w 1 f ( µ 1 )− f ( λ 1 ) · · · n ( µ 1 ) v 1 − n ( λ ℓ ) w ℓ f ( µ 1 )− f ( λ ℓ ) . . . . . . . . . n ( µ ℓ ) v ℓ − n ( λ 1 ) w 1 f ( µ ℓ )− f ( λ 1 ) · · · n ( µ ℓ ) v ℓ − n ( λ ℓ ) w ℓ f ( µ ℓ )− f ( λ ℓ ) ª ® ® ® ¬ . (24) Moreov er, by const r uction L SO = −( I + β M ) ˆ M ( I + β Λ ) = −( I + β M ) W T M ⋆ V ( I + β Λ ) . Thus, the f ollowing remark holds. Remark 3 I f w e hav e sufficient inter polation data, then rank ( V ) = rank ( W ) = r . As a consequence, the rank of the SO L oewner matr ix L SO giv es us the ord er of th e Ra yleigh damped SO m inimal realization inter polating the points, since rank  L SO  = rank  W T M ⋆ V  = rank  M ⋆  = order of the minimal SO inter polant . W e are now able to state the an a lo gue result to Theorem 1 f or Ra yleigh-d a m ped SO sys tems. Theorem 3 ( SO dat a-driven identificatio n) Assume t hat µ i , λ j for a ll i , j = 1 , . . . , ℓ . Additionally, supp ose that ( s 2 + α s ) L SO + ( β s + 1 ) L SO σ is i nvertible for all s = { λ 1 , . . . , λ ℓ } ∪ { µ 1 , . . . , µ ℓ } . Then ˆ M = − L SO , ˆ K = L SO σ , ˆ B = ( I + β Λ ) − 1 V SO ˆ C = W SO ( I + β M ) , and ˆ K = α ˆ M + β ˆ K satisfy the interpola tion conditions fr om Pr o b lem 1. W e no w consider the case where more data than necessar y are pro vided, w h ich is realistic for applications. In this case, the ass umption s of the abov e theor em are not satisfied; thu s, one needs to project onto the colum n span and the row span of a linear combin ation of the tw o Loewner m atr ices. More precisely , let the f ollo wing assumption be satisfied: rank   L SO L SO σ   = rank   L SO L SO σ   = r (25) Then, we co nsider the compact S VDs 10 Igor Pontes Duff, Pa wan Go y al, and P eter Benn er  L SO L SO σ  = Y ρ Σ l ˜ V T and  L SO L SO σ  = ˜ W Σ r X T ρ . (26) Using the projection matr ices V ρ and W ρ , w e are able to remo ve the r edundan cy in the data by means of the f ollowing re sult. Theorem 4 The SO r ealization ˆ Σ SO = ( ˆ M , ˆ D , ˆ K , ˆ B , ˆ C ) of a minimal interpo la nt of Problem 1 is giv en as: ˆ M = − Y T ρ L SO X ρ , ˆ K = − Y T ρ L SO σ X ρ , ˆ D = α ˆ M + β ˆ K , (27a) ˆ B = Y T ρ ˆ B L o e w , and ˆ C = ˆ C L o e w X ρ . (27b) Dependin g on whet her r in (25) is the e xact or ap pr oximate rank , we o btain either an interpolan t or a n appro ximate interpola nt o f t he data, r espectively . 3.2 Optimizing par ameters In the pre vious section, we ha v e shown how to co n struct a SO realization for given transfer f u nction measuremen ts and a pr ior i kno w ledge of th e parameter s α and β from the Ra yleigh-damped h ypothe sis. How e v er , there a re sev eral cases, where e xact values of α and β are not known bu t we rather can hav e a hint of the range for the parameters, i.e., α ∈ R α and β ∈ R β . Theref ore, as done for delay sy stems in [23], we also pr o pose a heur istic o p timization ap proach to obta in the param eters α and β for SO sy stems, satisfying the Ra yleigh-d a mped h ypothesis. For this pur po se, we split the data train in g D training and tes t set D test , e.g., in the ratio 80:20. Hen ce, w e ideally aim at sol v ing the optimization as f ollow s: min α ∈ R α , β ∈ R β J ( α , β ) (28) where J ( α , β ) : = Õ ( σ k , v k ) ∈ D test   ˆ H SO ( σ k α , β ) − v k   2 + Õ ( µ k , w k ) ∈ D test   ˆ H SO ( µ k α , β ) − w k   2 , where ˆ H SO is constr ucted using o nly the trainin g data. How ev er, the optimization problem (28) is no n-conv e x, and sol ving it is a challenging task. Therefore, we seek to so lv e a relaxed problem . For this pu rp ose, in the paper, w e make a 2-D gr id for the parameters α an d β in given inter vals. Th en, we seek to deter mine the para m eters on the gr id where the functio n J ( α, β ) is minimized. Nonetheless, solving the o p timization pr oblem (28) needs futu re inv estigation and so we leav e it as a possible future research problem. Data-Driv en I dentification of Ra y leigh-Damped Second-Order Systems 11 4 Numerical re sults In this section, we illus trate the efficiency of the prop o sed meth ods via se v eral numer ical ex amples, arising in various ap plications. All the simulation s are done on a CPU 2.6 GHz Intel ® Core™i5, 8 GB 1 6 00 MHz DDR3, MA TL AB ® 9.1.0. 44165 5 (R2016b ). 4.1 Demo e xample At first, we discu ss an ar tificial example to illus trate the pro posed me th od. Let us consider a SO system of orde r n = 2, Σ SO = ( M , D , K , B , C ) whose matr ices are given by : M =  1 0 0 1  , K =  1 0 0 2  , D = α M + β M , and B T = C =  2 3  , with α = 0 . 0 1 a nd β = 0 . 02. W e collect 20 samples ( σ j , ˆ H SO ( σ j )) , f or σ j ∈ ι [ 10 − 1 , 10 1 ] logarithmically spaced. Then , w e constr uct the F O and SO Loe wner matr ices in (12) and (23), receptivel y . In Figure 1 , w e plot the deca y o f the singular v alues of the L and L SO matr ices. It can be obser ved that rank ( L ) = 4 and r ank  L SO  = 2, as expected. In d eed, the demo sy stem has a minimal SO realization of order 2 and a min imal FO realization of order 4. By app l ying th e SVD p rocedur e, we co nstruct t w o r educed-o rder m odels of order 2, one f o r F O an d the other for SO . W e com pare the transfer function s of the o r ig inal and reduced-ord e r sy stems, and the r e sults ar e plotted in Figure 2 . The figure sho ws that the er ror between th e or iginal a n d SO reduced-o rder sys tem is of the lev el of machine precision, which mean s that the SO ap proach has r eco v ered an equivalent re a liza tion of the or iginal model. Additionally , the FO reduced sys tem of order 2 was not able to m imic the same behavior of the or ig inal sys tem, sh o wing that a larger ord er is required in this case. 4.2 Building e xample Let us no w consider the building model from the SLI COT librar y [8]. It descr ibes the displacement of a multi-storey buildin g , for example, d ur ing an ear thquake. It is a FO sy stem of order r = 48, wh o se dynam ics comes fro m a mechan ical system. The Ra yleigh damping co efficients here are α ≈ 0 . 4 947 and β ≈ 0 . 0 011. F or this e xamp le , w e collect 200 samples H ( i ω ) , with ω ∈ [ 1 0 0 , 10 2 ] . Then, w e build the F O and S O Lo ewner matr ices in (12) and (23), r eceptivel y . Additionally , using the heuristic p rocedur e in Subsection 3 .2, w e constructed th e reduced mod el 12 Igor Pontes Duff, Pa wan Go y al, and P eter Benn er Original sy stem FO Loewner SO Loewner 2 4 6 8 10 10 − 17 10 − 9 10 − 1 Re lative singular v alues Fig. 1 Demo e xample: Decay of the sing ular values f or the FO and SO Loe wner matrices. 10 − 1 10 0 10 1 10 0 10 2 Freq(s ) Magnitude 10 − 1 10 0 10 1 10 − 20 10 − 8 10 4 Freq(s ) Fig. 2 Demo ex ample: The figure on the left sho ws the Bode plot of the original sys tem and the FO and SO reduce d-order models. The figure on t he r ight show s the Bode plot of the er ror between t he original and reduced-ord er sys tems. assuming w e do not know a p r io ri the parameter s α and β . After this procedur e , we obtain α ∗ = 0 . 4 95 a nd β ∗ = 0 . 00 1, which are fairl y close to th e or ig inal pa rameters. In Figure 3, we plot the decay o f the singular values of the FO Lo e wner matr ix, th e SO Loewner matr ix fo r the or iginal parameters α and β , an d the SO Loewner matrix f or the estimated param eters α and β . Th e d ecay of the singular value s for the SO Loewner matr ix with origin al parameters is f aster than f or the FO Lo ewner matr ix. How ev er , f or th e SO Lo e wner m atrix with es timated parameters, th e decay of singular v alues starts fast an d then becomes slo w er . This show s that if the parameters α and β ar e not w ell identified, a higher reduced-or d er will be ne e ded to interp olate th e data. By applying the S VD pr ocedure, w e constr uct three red uced-or d er mod els of order 16. W e com pare the transf er fu nctions of the or ig inal an d red uced-or d er systems, and the results are plo tted in Figure 4 . This figu re show s th at f or the SO Loewner appr o ach (or iginal parameter s or with estimated par ameters) outper form the classical Loewner framew o rk. Data-Driv en I dentification of Ra y leigh-Damped Second-Order Systems 13 Original s y stem FO Loewner SO Loewner SO Loewner opt. par . 0 20 40 60 80 10 − 15 10 − 8 10 − 1 Re lative singular v alues Fig. 3 Build exa mple: Deca y of t he singular values f or t he FO Loewner matrix and fo r SO Loewn er matrices. 10 0 10 1 10 2 10 − 5 10 − 3 Freq(s ) Magnitude 10 0 10 1 10 2 10 − 7 10 − 4 10 − 1 Freq(s ) Fig. 4 Build ex ample: The figure on the left show s the Bode plot of the origi nal sys tem and the FO and SO reduce d-order models. The figure on t he r ight show s the Bode plot of the er ror between t he original and reduced-ord er sys tems. 4.3 Artificial Fishtail As the last ex ample, we consider the ar tificial fishtail mo del presented in [2 2]. This model comes fro m a finite-elem ent discretization of the contin u ous mechan ics model of an ar tificial fishtail. After discretization, the finite- dimensional sys tem has a S O realization of order 779 , 232. F o r this model, th e Ra yleigh damping is chosen with parameters α = 1 . 0 · 10 − 4 , β = 2 · 10 − 4 . It is a MIMO sy stem, bu t f o r the nu mer ical application, here we co nsider only th e first tran sf er function, i.e., from u 1 to y 1 . F or this e xamp le , w e collect 200 samples H ( i ω ) , with ω ∈ [ 1 0 1 , 10 4 ] . Then, w e build FO and SO Loe wner m a trice s in (12) and (23), receptiv ely . A dditionally , we also compute the reduced mod e l using the heur istic p rocedur e in Subsection 3 . 2, for which w e obtain the es timated parameters α ∗ ≈ 1 . 19 · 10 − 4 and β ∗ ≈ 2 · 10 − 4 . In Figure 3, we plot th e de cay of the singular v alues of the F O Loewner matr ix, the S O Lo ewner matr ix f or the or iginal par ameters α and β , and the SO Lo ewner matr ix f o r the estimated parame te r s α ∗ and β ∗ . By applying the SVD procedure, w e constr u ct thr ee reduced-o rder mode ls of order 8. W e co mpare the transf er functio ns of the or igin a l and re duced-o rder sys tems, and the results ar e p lotted in Figure 4. This 14 Igor Pontes Duff, Pa wan Go y al, and P eter Benn er Original s y stem FO Loewner SO Loewner SO Loewner opt. par . 0 20 40 60 80 10 − 10 10 − 5 10 0 Re lative singular v alues Fig. 5 Fishtail ex ample: Deca y of the singular values f or the FO Loe wner matr ix and fo r SO Loe wner matrices. 10 1 10 2 10 3 10 − 1 10 0 10 1 Freq(s ) Magnitude 10 1 10 2 10 3 10 − 5 10 − 4 10 − 3 10 − 2 Freq(s ) Fig. 6 F ishtail e xample: Th e figure on the lef t sho ws the Bode plot of the o riginal sy stem a nd the FO and SO r educed- order models. The fig ure on the right sho w s the Bo de plot of the er ror betw een the original and reduced- order sy stems . figure show s that the SO Loewner approach with or ig inal p arameters and S O Loewner with estimated parameters outper f or m the classical Loewner frame w ork. 5 Conclusions In this paper, we ha v e s tu died the problem of the identification of Ra y leigh-d a m ped second-o rder sys tems fr om frequency data. T o that aim , w e prop o se modified SO L o ewner ma tr ices which a re t he key tools to constru ct a realization inter po- lating the given data. Additionally , in the case o f red undant data, an SVD-based scheme is presented to construc t reduced-ord er m o dels. Moreov er , a heu r istic opti- mization p roblem is sketched to es timate the damping parameters. Finally , w e hav e illustrated the efficiency of the p roposed approach in some nu mer ical e xamples, an d we co mpared the results with the classical Loewner fr amew ork . Data-Driv en I dentification of Ra y leigh-Damped Second-Order Systems 15 A c kno wledg ement This w or k was supp or ted by Deutsche F orsc hungsgemeinsc h af t (DFG) , Collabor a - tiv e R esearch Center CR C 96 "T h er mo-energetic Design of Machine T oo ls". Re f eren ces 1. A. Anto ulas and B. Anderson , On the scalar r ational inte rpolation problem , IMA Jou rnal of Mathematical C ontrol and In f ormation, 3 (1986), pp. 61 –88. 2. A. C. A ntoulas , Appr o ximation of Larg e-Scale Dynamical Sys tems , SIAM Publications, Philadelphia, P A, 2005. 3. C. Be attie and P . Benner , H 2 -optimality conditions f or s tructure d dynamical syst ems , Preprint MPIMD/14-18, Max Planck Ins titute Magdeburg, Oct. 2014. A va ilable from http://www.mpi- magdeburg.mpg.de/preprints/ . 4. C. A. Beattie and S. Gugercin , Interpolat ory projection methods for structur e-pr eser ving model redu ction , Sy stems Co ntrol Lett. , 58 (2009), pp. 225–2 32. 5. P . Benner, A. C ohen, M. O hlber ger, and K. Willcox , Model Reduc tion and Appr oxi- mation: Theo ry and Algorithms , Computational Science & Engineering, SIAM P ublications, Philadelphia, P A, 2017. 6. Y . Cha hla oui, K. A. Galliv an, A. V an dendorpe, an d P . V an Dooren , Model r eduction of second-or der sys tems , in Dimension Reduction of Lar ge-Scale Sy stems, Spring er , 2005, pp. 149–172. 7. Y . Cha hla oui, D. Lemonnier, A. V andendorpe, a nd P . V an Dooren , Second-o rd er bal- anced truncation , Linea r Alg ebra Appl., 415 (2006 ), pp. 373–384 . 8. Y . Chahlaoui and P . V an Dooren , A collection of benc hmark examples for model reduction of linear ti m e inv ariant dy namical sys tems , SLICOT W orking N ote 20 02–2, U nive rsity of Manche ster , Fe b. 2002. A vailable from ww w.slicot.org . 9. R. R. Crai g and A. J. Kurdila , F undamentals of s t ructur al dynamics , John Wiley & Sons, 2006. 10. K. Galliv an, E. Grimme, and P . V an Dooren , A r ati onal Lanczos algorit hm for model r eduction , Numer . Algorithms, 12 (1996), pp. 33– 63. 11. K. Galliv an, A. V andendorpe, and P . V an Dooren , Model reduction of MIMO sys tems via tang enti al inte rpolation , S IAM J. Matrix Anal. Appl., 26 (2004), pp. 328 –349. 12. , Sylv es ter eq uations and projection- based model r eduction , J. Comput. Appl. Math., 162 (2004), p p. 213–229. 13. E. J. Grimme , Krylo v projection methods f or mod el r eduction , PhD thesis, U niv . of Illinois at Ur bana-Cha mpaign, US A, 1997. 14. B. Gusta vsen and A. Sem l yen , Rational appr o ximation of fr equ ency do main respo nses by v ector fitting , IEEE T rans. Po w er Delivery , 14 (1999), pp. 1052–1061. 15. C. Ionita , Lagr ang e rational interpolation and its applications to model r eduction and syst em identification , PhD thes is, Ri ce U niv ersity , 2013. 16. A. J . Ma y o and A. C. Antoulas , A fr amew or k f or the solution of the gene r alized r eali zation pr oblem , Li near Algebr a App l., 425 (2007 ), pp. 634–662. Special Issue in hon or of P . A. Fuhrmann , Edited by A. C. Antoulas, U . Helmke, J. Rose nthal, V . Vinniko v , an d E. Zerz . 17. L. Meir o vi tch , Principles and t ec hnique s of vibr ati ons , v ol. 1, Prentice Hall Ne w Jerse y , 1997. 18. D. G. Meyer and S. Sriniv as an , Balancing and model re duction for second -or der form li near sys tems . , IE EE T rans. Autom. Control, 41 (1996), pp. 1632–1644. 19. Y . Nakats ukasa, O. Sète, and L. N. Trefethen , The AAA algorithm for r ational appr o xi- mation , 40 (2018), pp. A1494–A1522. 20. I. Pontes Duff , C. Poussot - V assal, and C. Seren , Realization independent sing l e time-delay dynamical model interpolation and H 2 -optimal appr oximation , in 54th IEEE Conf erence on Decision and Co ntrol (CDC), I EEE, 20 15, pp. 466 2–4667. 16 Igor Pontes Duff, Pa wan Go y al, and P eter Benn er 21. T . Reis and T . Stykel. , Balanced t runcation model r eduction of second-or der sys tems , M ath. Comput. Model. Dyn. Sy st., 14 (2008) , pp. 391–406. 22. J. Saak, D. Siebel ts, and S. W . Werner , A comparison of second-or der model order reduction meth ods for an ar t ificial fishtail , a t-A utomatisierungstechn ik, 67 ( 2019), pp. 64 8–667. 23. P . Schulze and B. Unger , Data-driv en interpolation of dynamical systems with delay , Systems & Control Le tters, 97 ( 2016), pp. 125 –131. 24. P . Sc hulze, B. Unger, C. Bea ttie, and S. Gugercin , Data-driv en structure d r ealization , Linear Algebra Appl., 537 ( 2018), p p. 250–286. 25. D. C. Villemag ne and R. E. Skel ton , Model r eduction using a projection fo rmulation , Internat. J . C ontrol, 46 ( 1987), pp. 214 1–2169. 26. S. A. Wy att , Issues in interpolat ory model reduc tion: I nexact so lv es, second-or der syst ems and DAEs , PhD thesis, Virginia T ec h, 201 2. 27. A. Y ousuff , D. W agie, and R . S kel ton , Linear sys tem appr oximation via cov ariance equiv- alent realizations , J. Math. An al. Appl., 106 (1985), pp. 91–115.