Data-Driven Estimation of Vinnicombe metric

Quantifying model mismatch in a control-relevant manner is fundamental in robust control. A well-known metric for this purpose is the $ν$-gap, or Vinnicombe metric, which measures the discrepancy between a nominal model and the real system from a clo…

Authors: Margarita A. Guerrero, Henrik S, berg

Data-Driven Estimation of Vinnicombe metric
Data-Driv en Estimation of V innicombe metric Mar ga rita A. Guerrero, Henrik Sa ndberg, and Cristian R. Rojas Abstract — Quantifying model mismatch in a control-rele vant manner is fundamental in robust contr ol. A well-k nown metric fo r this purpose is the ν -gap, or V innicombe metric, which measures the discrepancy between a n ominal model and the re al system from a closed-loop viewpoint. Howe ver , its computation typically requires explicit kn owledge of the tru e system. In th is letter , we p ropose an ident ification-free, data-driven method to estimate the ν -gap between d iscrete-time SISO systems directly from input-outp ut experiments. The method is complemented by a data-dri ven winding-number test, based on W elch-type a veraging, to ver ify a re qu ired topological condition fo r the computation of the metric. Numerical simulations on hea vy- duty gas-turbine models and a textbook example show that the pr oposed estimate closely matches MA TLAB © gapmetric , while correctly detecting cases in which the admissibility conditions fail. Index T e rms — Data-driven modeling, System id en tification, Robust control. I . I N T RO D U C T I O N A c o re aim o f auto matic control is to design ro bust controller s that operate reliably despite mo deling errors. In this context, robustness means that perfo rmance doe s not deteriora te significantly when the real system d iffers from th e mod el used for the controller synth esis. Several frameworks hav e been de veloped to address this problem , including small-ga in argum e n ts [1] , integral q u adratic con- straints (IQCs) [2] , and H ∞ control [ 3]. Essentially , these approa c h es pr ovide stability o r perf ormance guarantees for all plants who se un certainty satisfies a p rescribed de scr iption, such as norm bou nds, s truc tu ral assumptions, or multiplier- based conditions. A control-relevant way to quan tify mo del mismatc h is provided by the V innicomb e m etric [4 , 5], also k n own as th e ν -gap metric, defined fo r linear time-inv ariant systems (L TI). This metric is particularly appea ling because it measur e s discrepancies between plants fro m a closed-loop viewpoint and indu ces the weakest to pology in which clo sed-loop stability is a robust pro p erty . Thu s, if the ν -gap between a no m inal mo del and the true plant is suf ficiently small, then a c ontroller that robustly stabilizes th e nomin a l model will also stabilize the true p lant, un der suitable admissibility condition s. Moreover , while top o logically equivalent to the This work was partly supported by the Swedish Research Council under contrac t number 2023-05170, by the W all enberg AI, Autonomous Systems and Software Program (W ASP), and by the Swedish Civi l Defence and Resilie nce Agency (Project MAD-V AMCHS). T he authors are with the Di vision o f Decision a nd Control Systems, KTH Royal Institut e of T ech- nology , 100 44 Stockholm, Sweden (e-mails: mags3@kth.se, hsan@kth.se, crro@kth.se). Code av ailable at https://githu b.c om/mags-ono/nu-gap . gap metric [6 , 7] , the ν -gap ad m its a direct freq uency- response interpretation. In practice , however , th e u sefulness of b o th ro bust-control guaran tee s and ν -g a p-based reason ing d e p ends on how accu- rately th e p lan t uncer tainty or model mismatch can b e ch arac- terized. Wh en th e tru e plant is unknown and prior k nowledge is ina ccurate, con structing a mean ingful uncer tain ty descrip- tion may be difficult, conservati ve, or infeasible. This p r actical difficulty has motivated research o n esti- mating contro l-relev ant robustness quantities directly fro m data, in settings where on ly input-o utput measurements are av ailable and no explicit plan t k nowledge is assumed [ 8, 9]. One promin ent line of work estimates the H ∞ norm of the mo deling erro r directly f rom experimen ts b y mean s of iterativ e proce dures, such as power iterations [10, 11] and Thomp so n sam p ling [12]. More recen tly , a power-iteration- based app roach for estimating the structured singular value has been p r oposed in [13]. T o the best o f our k nowledge, th e re is no data-dr i ven method for estimating the V innico mbe metric. Mo ti vated by this, we pr opose an id entification-fr ee, iterative ap proach that estimates the ν -gap betwee n two plants directly from tim e - domain input-o utput data. T he proposed method builds on the core idea of power itera tions to compu te a n umerical estimate of the d istan ce between two d ynamical system s, while an iter ati ve winding- number calculation is used to verify the admissibility c o nditions und er which the ν -gap is well defined . Our main contr ibutions are as fo llows : • W e propo se a novel iden tification-free, data-d r i ven method fo r estimating the V innicomb e me tr ic directly from ti me- domain input-output experiments. • W e introd uce an identification- free p rocedur e for verify - ing th e adm issibility cond itions of the ν - gap throug h an iterativ e wind ing-nu m ber computa tio n based on W elch - inspired frequency-do main e stima tes. • W e demo nstrate the effecti veness of the prop osed method thro ugh numerical simulations, showing that the resulting estimates closely agree with MA TLAB ® gapmetric . • W e illustrate the practical applicability of th e p roposed framework to Gener a l Electric heavy-duty gas turbin es, showing h ow the V innicombe me tric can be interpreted and implemented in an industrial setting. The rem a inder of this p aper is organ ized as follows: In Section II, w e state ou r prob lem setup. Section III revie ws the V innicomb e m etric and th e power iteration m ethod. In Section IV, we outline our proposed a p proach , and demon- strate its efficacy in Section V. Th e paper is conclud ed in Σ G 0 Σ C v 1 u y v 2 − + − + Fig. 1: Standard feedback configurati on. Section VI. I I . P RO B L E M S TA T E M E N T Consider a stable un ity-feedba c k intercon nection of a nominal plan t G 0 with a controller C , deno ted by [ G 0 , C ] . Suppose that the nominal plant differs fro m an unkn own tru e system G . Our objective is to gua rantee that the co ntroller C will perfo rm as inten ded o n the true process regard less of the un m odeled d y namics. Throug hout the pap er , ( G 0 , G, C ) are d iscrete-time, linear time-in variant, sing le-input single- output systems. The aforem entioned p roblem can be h andled with a w e ll- known distance metric, namely th e ν -gap, also called the V innicombe metric [ 4, 5]. The V in n icombe metric, denoted by δ ν , is an appropriate measure for assessing how two systems differ in their c lo sed-loop beh avior , satisfying 0 6 δ ν ( G 1 , G 2 ) 6 1 (see Section III- A). This quantity admits a clo sed-loop interpretation: it mea- sures, in a n ormalized worst-case sense, the discrepancy between the complemen tary sensitivity beha vio r s induc e d by un ity feedbac k ar ound G an d G 0 . In particular, taking a sufficiently sm a ll value of δ ν implies that any co ntroller that stabilizes G 0 with adequ ate margin also stabilizes G (under additional con ditions summ arized in Section III-A). When the real plant is unknown, computing δ ν ( G, G 0 ) becomes unten able, as it requ ires knowledge of the true process. A classical workar o und is to postulate that G lies in a neighbo rhood o f the nomin al model, i.e., B τ ( G 0 ) : = { ˜ G : δ ν ( ˜ G, G 0 ) ≤ τ } [5 , Proposition 1.2] . Howev er, speci- fying such a set typically relies on prior m odeling kn owledge about unmo deled dynam ics an d may be c onservati ve o r inaccurate. Motiv ated b y th ese limitation s, we propose in Section IV an iden tificatio n-free, data-driv en procedure to compu te the V innicombe metric f or the single-in put sing le-outpu t (SISO) case using a power-iteration p rinciple [1 4], req uiring only input-o u tput experiments on the b lack-box system G and an iterative r edesign of the excitation sig n al. While we primarily assume that G 0 is an av ailable nomin al model, th e same methodolo gy can be used as a data-driven similarity test whe n both system s are unknown, e.g ., to compare th e opening /closing dynamics of two valves fro m closed-loo p measuremen ts. Before presen tin g the prop o sed app roach, we introd uce the technical preliminaries req uired th r ougho ut the manuscrip t. I I I . P R E L I M I N A R I E S A. The V inn icombe Metric Consider the standard negati ve-feedb ack inter connection of a nomin al d iscrete-time L TI multi-inpu t m u lti-output (MIMO) plant G 0 ( z ) and a compen sator C ( z ) as shown in Fig. 1. Let v 1 and v 2 denote additive d istu r bances at the controller input a n d o utput, respectiv ely , and le t ( y , u ) be the p lant o utput an d inp ut. The closed-lo op transfer matrix from [ v 1 v 2 ] T to [ y u ] T is T ( G 0 , C ) :=  G 0 I   I + C G 0  − 1  − C I  , (1) where I = 1 is the SISO case. W e defin e the stability mar gin of th e feedback intercon- nection [ G 0 , C ] as b G 0 ,C := ( k T ( G 0 , C ) k − 1 ∞ , if [ G 0 , C ] is st able , 0 , otherwise. (2) The n umber b G 0 ,C can b e regarded as a m easure o f the perfor mance of the feedback system co mprising G 0 and C , with larger values of (2) cor respondin g to better per f ormance. Note that (2) lies in the rang e [0 , 1] [4]. Definition 1: Let G 0 ( z ) an d G ( z ) b e discrete-time single- input single-o u tput ( SISO) transfer fu n ctions. The V i nn i- combe metric ( o r ν -g ap) [5] between G 0 and G is defin ed as δ ν ( G 0 , G ) =           G 0 − G p (1 + | G 0 | 2 )(1 + | G | 2 )      ∞ , if ( G 0 , G ) ∈ C , 1 , otherwise. (3) where ( G 0 , G ) ∈ C if the fu nctions f 1 ( z ) : = 1 + G 0 ( z ) G 0 ( z − 1 ) a nd f 2 ( z ) : = 1 + G ( z ) G 0 ( z − 1 ) have the same num ber of zer os in E : = { z ∈ C : | z | > 1 } . The con dition ( G 0 , G ) ∈ C is an index (win ding-nu mber) requirem ent that en sures a meanin g ful no tion of co ntinuity: if G is perturb ed co ntinuou sly away from G 0 , it r ules out situations in which an intermed iate plant ˜ G would necessarily yield δ ν ( ˜ G, G 0 ) = 1 e ven though δ ν ( G, G 0 ) is small [1 5 , Section 12.1]. In the SISO case, this con dition can be ch ecked via a Nyquist/windin g -numb er test, wh ere the winding number is den oted b y wno( · ) ; in pa r ticular , the Nyquist plot of G ( e j ω ) G 0 ( e − j ω ) sho uld n ot encir cle the critical po int − 1 for ω ∈ ( − π , π ] . Addition ally , in this same SISO setting, the pointwise quantity underly in g the ν -gap admits a geome tric interp r etation as the cho rdal distance betwee n th e pr o jected Nyquist diagram s of two plants onto the Riemann sphere [5 , Sec. 3]. Theref ore, δ ν may be viewed as the m aximum normalize d separation over frequen cy , which explains its r e le vance as a closed-loo p similarity metric. Pr oposition 1 ( [5, Pr op. 1.2]) : Given a n ominal plant G 0 and a pertu rbed plant G , th e interco nnection [ G, C ] is stable for all com p ensators C satisfying b G 0 ,C > β if and o nly if δ ν ( G 0 , G ) 6 β . The ν -gap metric δ ν ( G 0 , G ) q uantifies the distance be- tween two plan ts from a closed- loop p erspective. I n p artic- ular , δ ν ( G 0 , G ) ∈ [0 , 1] , with δ ν ( G 0 , G ) = 0 if an d only if G 0 = G . Hence, when G represen ts the (unkn own) true plant, a small value of δ ν ( G 0 , G ) indicates that the mismatch between G 0 and G is n egligible in terms o f closed- loop stability . B. P o wer I teration Metho d In or der to comp ute ( 3), one need s to co mpute the H ∞ norm of a tr ansfer fu nction. T o this end , con sid e r a linear discrete-time stable causal a n d proper sy stem of the form y k = G ( q ) u k + e k where { u k } is an input signal, G is a transfe r fun ction, { y k } is th e outp ut o f the system, an d { e k } is G a ussian, not necessarily white noise. Next, recall that the H ∞ norm is an ind uced norm, i.e., k G k ∞ = k G k i 2 ; see [ ? , App. A.5.7] for a pro of. Assume that, at iter ation t , the signals have a finite len gth N ∈ N , i.e., u t := [ u t 0 , . . . , u t N − 1 ] T , y t := [ y t 0 , . . . , y t N − 1 ] T ∈ R N × 1 , hence k G k i 2 ,N = sup u t ∈ R N × 1 \{ 0 } k y t k 2 k u t k 2 , (4) where k G k i 2 ,N → k G k ∞ for N → ∞ [1 1, Th. 3 ]. L et us assume tha t, at iteration t , th e system starts fr om zer o initial condition s. If there is no n oise affecting the outpu t of th e system, the map ping f r om u t to y t can be written as y t = G N u t , where G N =      g 0 0 0 · · · 0 g 1 g 0 0 · · · 0 . . . . . . . . . . . . . . . g N − 1 g N − 2 · · · g 1 g 0      is a triang u lar T oeplitz matrix. Th us, ( 4) is equiv alent to k G k 2 i 2 ,N ≈ sup u t 6 =0 ( u t ) T G T N G N u t ( u t ) T u t = λ max  G T N G N  . (5) The right hand side, kn own as the Rayleigh quotient of G T G , corresp o nds to the largest eigenv alue of G T G , and it can b e co mputed using th e power iteration me th od. Pr o c edur e 1 (P ower Iteration Metho d): Consider a square matrix A ∈ R n × n . The power iteration method perfor ms the following steps to estimate the eigenv ector of A correspo nding to its largest eigenv alue: 1. Choose x 0 ∈ R n random ly . 2. Set i ← 1 . 3. Compu te x i ← 1 k Ax i − 1 k 2 Ax i − 1 . 4. Set i ← i + 1 and g o to Step 3. For normalize d iter ates, a natural estimate of the magni- tude of the d ominant eige nvalue is k Ax i − 1 k 2 1 . In the SI SO frequen cy-domain setting, the H ∞ norm co- incides with the maxim um m a gnitude of the frequency response over all freq uencies, i.e., with the peak of the Bod e magnitud e plot. Thus, the sam e power-iteration principle can be interp r eted in the fr equency domain as an iterative search for the peak gain, where the maximizin g frequ ency is implicitly identified by th e iteratio ns. 1 If A is symmetric, the dominant eigen val ue may also be estimated through the Rayleigh quotient x T i − 1 Ax i − 1 . The two estimates may differ numerica lly under noise or finite-it eration ef fects. Our p r oposed proced ure retain s this templa te to ap p rox- imate the ν -gap in ( 3) b ased solely on experimen ts o n the unknown system G . W e assume e k = 0 throug hout the analysis, wh ile noise is included in th e simulations of Section V to evaluate robustness. I V . P RO P O S E D A P P RO AC H In th is section , we p ropose a data-driven proced u re to estimate the V innicom b e ν -gap metric from time-domain experiments. Assuming that a nom inal m odel G 0 is av ailable and that o nly input- output d ata can be co llected, we ad apt the m odel-based constructio n in [4] to operate directly on frequen cy-domain qua n tities estimated fro m time-d omain data. A. ν -gap estimatio n W e f o cus on the pr actically r e levant setting wher e G denotes the (u nknown) plan t and G 0 is a nominal mo del. For a given inpu t u , let U ( e j ω ) be its discrete-time Fourier transform (DFT) , and let Y ( e j ω ) a nd Y 0 ( e j ω ) d enote the correspo n ding outp ut spectra ob tained from experiments on G an d G 0 , r espectiv ely , so that Y ( e j ω ) = G ( e j ω ) U ( e j ω ) an d Y 0 ( e j ω ) = G 0 ( e j ω ) U ( e j ω ) . When ( G, G 0 ) ∈ C , the ν -g ap admits the equ i valent freq uency-dom ain expression δ ν ( G, G 0 ) = sup ω ∈ ( − π ,π ] | U ( e j ω ) | 2 | Y ( e j ω ) − Y 0 ( e j ω ) | p | U ( e j ω ) | 2 + | Y ( e j ω ) | 2 p | U ( e j ω ) | 2 + | Y 0 ( e j ω ) | 2 . As in H ∞ norm estimation , the ν -gap is the suprem um of th e modulu s of a freq uency function . In the no iseless case, this motiv ates a p ower -iteration- type inpu t r edesign that p r ogres- si vely concentrates excitation at f r equencies th a t dom inate the a bove ratio. Specifically , given spe c tr a ( U n , Y n , Y n 0 ) at iteration n , we u pdate the next input in the freq uency do main accordin g to ˜ U n +1  e j ω  =   U n  e j ω    2  Y n  e j ω  − Y n 0  e j ω  q | U n ( e j ω ) | 2 + | Y n ( e j ω ) | 2 q | U n ( e j ω ) | 2 + | Y n 0 ( e j ω ) | 2 . (6) The r ole of (6) is the n analog ous to th at o f the matrix - vector p roduct in the classical p ower iteration: each upd ate amplifies the freq uency com ponents th at con tribute m ost to the norm alized discrepancy betwe e n G and G 0 . The updated spectrum is map ped back to the time domain v ia an inverse DFT , no rmalized, and used in the n ext exp eriment. Thus, the input prog ressi vely con centrates n ear the frequen cy a t wh ich the ν -gap is attained . B. wno estimation The update rule in (6) suggests a straigh tforward im- plementation ; howe ver , the ind ex cond ition ( G 0 , G ) ∈ C cannot be overlooked. Th is condition r ules out patho logical situations where, despite a small closed-lo op mismatch be- tween G 0 and G , any con tin uous p ath fro m G 0 to G would necessarily cross an interm ediate transfer function e G for which the ν -gap saturates to one. Accordin gly , we verify that ( G 0 , G ) ∈ C in a data - driven fashion by co mputing winding number s of the index func tio ns f 1 and f 2 , i.e., by ch ecking that the Nyquist cur ve o f G ( e j ω ) G 0 ( e − j ω ) does not encircle the critical poin t − 1 f or ω ∈ ( − π , π ] . Using th e same infor mation from the Fourier transform of the input and outp u ts of (6), and inspire d by W elch ’ s method [16] we com pute f 2 as fo llows : f 2  e j ω k  ≈ 1 +  P N acc n =1 Y n [ k ] U n [ k ] ∗   P N acc n =1 Y 0 ,n [ k ] U n [ k ] ∗  ∗  max  P N acc n =1 | U n [ k ] | 2 , ε 0  2 , (7) where n is the curr ent iteration, ω k := 2 π k / N , k ∈ { 0 , 1 , . . . , N − 1 } is the DFT grid on the unit circle, N acc is th e nu mber of iterations u sed to appr oximate f 2 , and ε 0 is a small number to av oid d i vision by 0 . An analog o us expression ca n be written fo r f 1 by r eplacing ( Y n , Y 0 ,n ) with ( Y 0 ,n , Y 0 ,n ) in (7). Inspired by W elch’ s a veraged pe- riodog ram estimato r [17, Ch. 6, Eq . (6 .70)], we estimate the index fu n ctions by averaging acro ss p ower-iteration exper- iments. Unlike classical W elch’ s m ethod, which par titions a single r e c ord into m ultiple sub-b atches an d reco mputes the DFT for each b a tch, our appro ach reuses th e DFTs that ar e alrea d y c omputed at each power -iteratio n step . W e therefor e average over the first N acc ≪ M iter ations: as the power m e thod prog resses, the redesigned inpu t b ecomes increasingly narrowband ar ound the peak frequen cy , wh ich reduces spectral coverag e an d makes later iterations less informa tive fo r index verification [ 1 1]. Fin ally , for i ∈ { 1 , 2 } , let f i [ k ] := f i ( e j ω k ) deno te the samp les o f f i on the DFT grid. W e comp ute the wind ing num bers wno( f 2 ) and wno( f 1 ) by summin g discrete phase incremen ts of f i around the unit cir cle, where we set f i [ N ] := f i [0] to close the sampled path. Spe c ifica lly , we define Θ i = N − 1 X k =0 arg  f i [ k + 1] f i [ k ]  , wno( f i ) = round  Θ i 2 π  . ( 8) The d a ta - driven win d ing-nu mber comp utation is summarized in Algor ithm 2. These co nsiderations lea d to the pseu do-cod e for the data- driven comp utation of the V inn icombe metr ic in Algor ithm 1. The pro cedure term inates early if wno( f 1 ) 6 = wno( f 2 ) within the first N acc iterations (ind icating ( G 0 , G ) / ∈ C ); otherwise, it ru ns fo r th e fu ll M iter ations. V . E X P E R I M E N T S This section presents simulation studies to assess Algo- rithms 1 – 2. The study in cludes: (i) two heavy-duty gas- turbine case studies, wh e re ν - gap estimates are b e nch- marked again st MA TLAB’ s gapmetr ic ; and (ii) a syn- thetic/textbook examp le for wh ich the winding nu m ber con- dition is analy tica lly verifiab le. A. Case Stud y: Heavy- Duty Gas T urbine Heavy-duty gas tur bines (HDGT s) are large industrial gas turbines used mainly fo r power gen eration. They consist of an axial-flow com pressor, combustor chamb e r s where fuel is burned, and a tur bine section that expand s th e ho t Algorithm 1 Data-d riv en estimation o f th e ν -gap Requir e: A random input vector u 0 ∈ R N 1: fo r n = 0 , 1 , 2 , . . . , M − 1 do 2: Apply u n to G and G 0 , and colle ct outputs y n and y n 0 , respec- ti vely . 3: Compute the DFT of u n , y n and y n 0 , and calcul ate (6), for each frequenc y ω . 4: if n < N acc , update accumulat ors via S y u + = Y n [ k ] U n [ k ] ∗ , S y 0 u + = Y 0 ,n [ k ] U n [ k ] ∗ , S u + = | U n [ k ] | 2 end if 5: if n = N acc then 6: Run Algorithm 2 7: if inC is false, T erminate and set δ ν ← 1 end if 8: end if 9: Compute ˜ u n +1 as the IDFT of ˜ U n +1 ( e j ω ) 10: Normalize ˜ u n +1 by its 2 -norm to obtain u n +1 11: end for 12: Estimate δ ν as the 2 -norm of ˜ u M Algorithm 2 Data-d riv en ind ex condition ( G 1 , G 2 ) ∈ C Requir e: Accumula tors ( S y u , S y 0 u , S u ) , N acc batche s; tolerance ε 0 > 0 ; threshold tol f > 0 . 1: fo r n = 0 , 1 , . . . , N acc , Construct f 1 and f 2 via (7), end for 2: Compute minimum distance s to the origin: m 1 = min k | f 1 [ k ] | and m 2 = m in k | f 2 [ k ] | . 3: Compute winding numbers via (8), with f i [ N ] = f i [0] , i ∈ { 1 , 2 } 4: Set inC ←  m 1 > t ol f  ∧  m 2 > tol f  ∧  wno( f 1 ) = wno( f 2 )  gases to p roduce shaft work, typica lly driving an e lectrical generato r . A key con trol variable is the Ga s Control V alve (GCV), wh ich regulates the fuel flow into th e com bustor and therefor e dir ectly affects the generated p ower . In th is work, we adop t the well-known Rowen represen tation for heavy- duty gas tu rbines [18] an d co nsider a reduced inp ut-outpu t model fro m th e ( m easured) GCV p osition u GCV to gen erated power P , ob tained by cascadin g: ( i) a fu el dy namics G f , cl , (ii) a combustor delay , ( iii) a co mpressor-discharge dynam ics G cd , an d (iv) the f uel-to-torq ue (p ower) map F 2 . Specifically , the f uel-path dynamics W f ( s ) are mod eled as W f ( s ) = G cd ( s ) e − sT C R G f , cl ( s ) u GCV ( s ) , with G f , cl ( s ) = 1 T f s + 1 + K F , G cd ( s ) = 1 T cd s + 1 , where T f is the fuel-system time constant, K F the fuel- system feedbac k gain, T C R the co mbustion reactio n time delay , and T cd the com p ressor-discharge time constan t. The (affine) torqu e/power ma p is given by P = F 2 ( W f , N ) = A + B W f + C (1 − S ) , (9) where S is th e per-unit rotor speed an d { A, B , C } are Rowen torqu e -block parameter s. For grid-con n ected op eration, sp eed d eviations are small and we take S ≈ 1 over th e consider ed ope rating windows. Linearizing (9) ar ound an ope r ating point and using incre- mental signals (den o ted by ∆( · ) ) yields the L TI incr emental plant G ( s ) := ∆ P ( s ) ∆ u GCV ( s ) = B e − sT C R ( T f s + 1 + K F )( T cd s + 1) . (10 ) A.1. Nominal Mod e l In [1 9], the par ameters of Rowen’ s model fo r HDGT are estimated f or a 172- M W simple cycle, single shaf t p ower unit by using operation al d ata a t a nomin a l p oint. Following the setup in Sec tio n II, we tre a t the p arameter values in [19, T able V I I] as the closest av ailable d escription of the true plant G . Next, we assume th at a nomin a l tran sfer f unction G 0 from the (measured ) GCV po sition u GCV to generated power P has be e n identified by means of a multi-ob jec ti ve differential ev olution op timization algor ithm as described in [2 0]. Th e parameters in [20, T able II], are estimated based on real data acquired from a ga s turbine of 162 M W nominal power . Assume tha t we have derived a contr oller C that stabilizes G 0 on a hig h-fidelity in d ustrial simulator . T he objective is to assess wh e ther d eploying C on the true oper ating power station is justified from a closed-loo p ro bustness viewpoint, by estimating δ ν ( G 0 , G ) fro m data. T o emu late realistic o p erating constra in ts, we co nsider that the plant is under dispatch instructions 2 , so that only small variations in generated power ar e permitted du r ing testing. Un d er these constraints, we apply a bo unded excitation signal u 0 ∈ R N (small-signal variations arou nd the op erating point) to both G and G 0 , co llect the cor respondin g input-o utput data, a n d run A lg orithm 1 to obtain a data-driven estimate of the ν -gap and the associated in dex-conditio n check. As sh own in Fig. 2, the Monte Carlo (MC) average of the pr oposed estimate ˆ δ ν , computed over 100 run s, beco m es nearly settled after roug hly 60 iterations an d appr o aches MA TLAB’ s refere n ce value δ ν = 0 . 21 7 2 , in dicating conver - gence of the p r oposed p ower -iteration sch eme. Th e e stimate is c omputed using N acc = 1 0 , N = 9000 , T s = 0 . 05 , an d measuremen t-noise variance σ 2 y = 0 . 01 o n b oth plants. T o f urther interpret this r esult, we inspect the dom inant frequen cy in th e last iteration and estimate the n ominal local gain as p 0 ( ω p eak ) ≈ Y 0 ( ω peak ) U ( ω peak ) , wh ere Y 0 and U are the DFTs of the n ominal output an d input, respectively . In this c ase, the d ominant p eak oc curs at ω p eak ≈ 0 . 13 6 rad/samp le, with | p 0 ( ω p eak ) | ≈ 0 . 54 . Acco rding to V innicomb e’ s Riemann- sphere interpr e tation [5, Section 3], a ν -g ap ar ound 0 . 2 correspo n ds to a relativ ely small discrepancy near un ity g ain, namely , at freq u encies wh ere | p 0 ( ω ) | ≈ 1 . In the presen t case, however , the maximizing freq uency lies below unity gain, since | p 0 ( ω p eak ) | ≈ 0 . 54 . Still, the value ˆ δ ν = 0 . 217 2 suggests that the no minal and real plants are reason ably clo se from a clo sed-loop viewpoint. In par ticular , if the d eployed controller C satisfies b G 0 ,C > 0 . 2172 , then Pro p osition 1 guaran tee s that it also stabilizes the true plan t. A.2. Unkno wn Models Rowen’ s model provides simplified d ynamic r e presenta- tions fo r a br oad rang e of General E lectric (GE) heavy-duty , single-shaft g as turbin es and accommo d ates bo th liqu id - and gas-fuel systems. In pa r ticular , GE ’ s Fram e 6F an d 2 In power systems, “dispatch” refers to grid-operator instruct ions on po wer setpoints and allow able load varia tions. 0 20 40 60 80 100 0.16 0.17 0.18 0.19 0.2 0.21 0.22 0.23 0.24 δ ν = 0 . 2172 Iteration n ˆ δ ν ˆ δ ν Estimate Mean of 100 MC δ ν MA TLAB Fig. 2: Estimate of the ν -gap ˆ δ ν (black, dashed), avera ge of 100 Monte Carlo (MC) simulations (blue, solid), and actual δ ν (red, dashed). Frame 9 F are repr esentativ e mid- and hig h -power units (appro ximately 88 MW and 29 4 MW , respectively) for wh ich Rowen repo rts “simplified” para meter sets and associated fuel/dyn amic characteristics; see [18, Notes 4-7 ]. W ith th e rapid g r owth of weather-dependen t solar and wind gen er- ation, gas-fired power plants are increasingly required to provide renewable-balancing services, making opera tional flexibility a critical req uirement. G E ’ s Dry L ow NO x (DLN) combustor systems [2 1], including recent F-class up g rades such a s the DLN2.6+ platform [ 2 2], suppor t such flexible operation while meeting tightening emissions regulations. Maintaining th is performan c e, however , typically requires periodic tuning b y a techn ical advisor (T A), often throu gh on-site cam paigns la stin g 1–3 days de p ending on the scope. For Frame 6F/9F units, tun ing the gas contro l valves (GCVs) is central to DLN perfo rmance, since these valves regulate the fuel flow a nd, to g ether with other actuator s, d e termine th e combustor op erating c o ndition, em issions, and loa d response. W e consider a n atural-gas power station equipped with both a Frame 9 F and a Frame 6F un it (no t necessarily within the sam e train). Sup pose the 9F unit h a s b een successfully tuned to comply with NO x constraints, an d th e op e rator aim s to ach ie ve comp arable emissions co mpliance on the 6F un it. In practice, however , limited T A availability and d ispatch constraints may p rev ent executing a f ull emissions test cam - paign for the DL N2.6+ pro cedure. As a result, the p lant seeks an alternative, low-impact p rocedur e to assess wh ether the existing GCV tunin g/controller settings d ev elope d f or th e 9F are com patible with th e 6F , while remaining within allow able small power variations a p proved by d ispatch. Let G 1 ( z ) a nd G 2 ( z ) den ote the (unkn own) discre te- time transfer functions fr o m the GCV inp ut u GCV to gen erated power P for the 9F and 6 F units, re spectiv ely . By applyin g Algorithm 1 to short, bo unded tests on ea c h unit, the main - tenance team can o b tain a d ata-driven estimate of the ν -gap between G 1 and G 2 and thus q uantify the degree of closed - loop mismatch. This enab les an assessment of whethe r a GCV contro ller that stabilizes G 1 is likely to d eli ver co m- parable performa n ce on G 2 , i.e. , whether the p lant-to-plan t mismatch is suf ficiently small for safe deployment. Using the p arameter sets in [18, Notes 4- 7], we pr oceed to estima te δ ν ( G 1 , G 2 ) . 0 20 40 60 80 100 120 140 0.28 0.29 0.3 0.31 0.32 0.33 0.34 0.35 δ ν = 0 . 3323 Iteration n ˆ δ ν ˆ δ ν Estimate Mean of 100 MC δ ν MA TLAB Fig. 3: Estimate of the ν -gap ˆ δ ν (black, dashed), aver age of 100 MC simulatio ns (blue, solid), and actual δ ν (red, dashed). Using the same par ameters as b efore ( N acc = 10 , N = 9000 , T s = 0 . 0 5 , and σ 2 y = 0 . 0 1 ), the MC mean of the propo sed estimate over 1 00 runs stabilizes after abo ut 100 iterations a nd conv erges in Fig. 3 to MA TLAB’ s value δ ν = 0 . 3323 . The dom inant peak is attained at ω p eak ≈ 0 . 0126 rad/sample, with | p 0 ( ω p eak ) | ≈ 0 . 644 , which places the maximizing fr equency closer to nom inal unity g ain. From V innico mbe’ s Riemann-sph ere interpr e ta tio n, this makes the geometric intuition beh ind the ν -gap m o re d irectly relev ant here; n ev ertheless, the larger value ˆ δ ν = 0 . 3 323 still ind icates a no n -negligible degree of closed- loop similarity . Fr om a practical viewpoint, this sugg e sts th at a controller tuned f o r the GE 9F unit may also be app licable to the GE 6F un it. A formal guar antee, however , would require estimating the correspo n ding stability margin b G 9 F ,C , f or example throu g h a power-iteration-based estimate of k T ( G 9 F , C ) k ∞ . In th at case, Pro position 1 would app ly whenever b G 9 F ,C > 0 . 3 323 . B. Case Stud y: P oles in RHP W e ev aluate the ν - gap me tr ic on a discrete-time SISO case study with samp ling time T s = 1 s. Consider the stab le plant G 1 ( z ) = 1 − z − 1 1 − 0 . 8 z − 1 , G 2 ( z ) = 1 . 8 z − 1 G 1 ( z ) . (11) Plant G 2 differs from G 1 by a gain of 1 . 8 and an a dditional one-samp le delay z − 1 , which introdu ces extra phase lag. W e can obtain explicit expr essions for the index functions on the unit circle z = e j ω . Sinc e the plants ar e SISO with real coefficients, G 1 ( e − j ω ) = G 1 ( e j ω ) , and thu s f 1 ( e j ω ) = 1 + G 1 ( e j ω ) G 1 ( e − j ω ) = 1 +   G 1 ( e j ω )   2 , f 2 ( e j ω ) = 1 + G 2 ( e j ω ) G 1 ( e − j ω ) = 1 + G 2 ( e j ω ) G 1 ( e j ω ) . By (11), we ca n simp lify f 2 ( e j ω ) as f 2 ( e j ω ) = 1 + 1 . 8 e − j ω   G 1 ( e j ω )   2 . Note th a t f 1 ( e j ω ) is real and strictly positive fo r all ω , hence wno( f 1 ) = 0 . In contr a st, the factor e − j ω induced by the one-samp le delay r otates th e complex term in f 2 ( e j ω ) as ω varies, so the cu rve f 2 ( e j ω ) may encircle the origin, yieldin g wno( f 2 ) 6 = 0 . 0 2 4 6 8 10 12 14 0.7 0.75 0.8 0.85 0.9 0.95 1 δ ν = 0 . 9308 Iteration n ˆ δ ν ˆ δ ν Estimate Mean of 100 MC δ ν MA TLAB Fig. 4: Estimate of the ν -gap ˆ δ ν (black, dashed), aver age of 100 MC simulatio ns (blue, solid), and actual δ ν (red, dashed). Using N acc = 10 , N = 10 00 , T s = 1 , and σ 2 y = 0 . 01 , Fig. 4 shows that the p ower-iteration estimate initially approa c h es MA TLAB’ s v alue δ ν = 0 . 9 308 . Ho wever , a t n = N acc , Alg orithm 2 detects wno( f 1 ) 6 = wno( f 2 ) , so the pair fails the set C condition. Hence, ou r p r ocedure ap plies a hard stop an d d oes n ot r efine th e estimate fu rther . V I . C O N C L U S I O N W e ha ve presented an ide ntification-fre e , data- driven method to estimate the V inn icombe ν -gap between two discrete-time SISO systems directly f r om time-do main input- output exp eriments. The propo sed procedu r e combines it- erative ν -gap e stima tio n with a data- driven admissibility check, pr oviding a control-r elev ant m easure of model mis- match for closed -loop analysis. Num erical stud ies, in cluding an indu strially motivated case, show clo se a greement with MA TLAB’ s g apmetric and corre c t detection o f ind ex- condition violations. Future work will addr ess extension to MIMO systems. R E F E R E N C E S [1] C. A. Desoer and M. V idyasagar , F eedback Systems: Input-Output Pr operties . SIAM, 2009. [2] A. Megretski and A. Rantzer , “System analysis via integra l quadrat ic constrai nts, ” IEEE T ransaction s on Automatic Contr ol , vol. 42, pp. 819–830, 1997. [3] K. Zhou, J. C. Doyle, and K. Glov er , R obust and Optimal Contr ol . Prentic e-Hall, 1996. [4] G. V innicombe , “Structure d uncertainty and the graph topology , ” in Pr oceedings of the 30th IE EE CDC , pp. 541–542, 1991. [5] G. V innicombe, “On the frequency response interpret ation of an inde xed L 2 -gap metric, ” in American Contr ol Conferen ce , pp. 1133– 1137, 1992. [6] A. El-Sakkary , “The gap m etric : Robustness of stabilizat ion of feed- back systems, ” IEEE T ransactions on Automati c Contro l , vol. 30, pp. 240–247, 1985. [7] M. Cantoni and H. Pfifer , “Gap metric computation for time-v arying linea r systems on finite horizons, ” IF AC -P apersOnLin e , vol. 50, no. 1, pp. 14513–14518, 2017. [8] A. Koch, Determining input-output pr operties of linear time-i nv ariant systems fro m data . Logos V erlag, 2022. [9] T . Koening s, M. Kruege r, H . L uo, and S. X. Ding, “ A data- driv en computat ion method for the gap metric and the optimal stabil ity margin , ” IEE E T ransactio ns on Automatic Contr ol , vol. 63, no. 3, pp. 805–810, 2017. [10] B. W ahlber g, M. B. Syber g, and H. Hjalmarsson, “Non-para metric methods for L 2 -gain estimation using itera tiv e experiment s, ” Auto- matica , vol. 46, no. 8, pp. 1376–1381, 2010. [11] C. R. Rojas, T . Oomen, H. Hjalmarsson, and B. W ahlber g, “ Analyzi ng iterat ions in identific ation with applicati on to nonparamet ric H ∞ - norm estimation, ” Automatica , vol. 48, no. 11, pp. 2776–2790, 2012. [12] M. I. M ¨ uller and C. R. Rojas, “Gain estimation of linear dynamical systems using Thompson sampling, ” in Proc eedings of the 22nd Inter- national Confer ence on Artificial Intellig ence and Statistic s (A IST ATS) , pp. 1535–1543, 2019. [13] M. A. Guerrero, B. L akshminara yanan, and C. R. Rojas, “Dat a- dri ven estimation of s tructu red singula r value s, ” IEEE Contr ol Systems Letter s , vol. 9, pp. 1976–1981, 2025. [14] G. H. Golub and C. F . V an L oan, Matrix Computations . John Hopkins Uni versity Press, 2013. [15] K. J. ˚ Astr ¨ om and R. Murray , F eedback Systems: An Intr oduction for Scient ists and Engineers . Princeton Univ ersity Press, 2021. [16] P . W elch, “The use of fast Fourier transform for the estimat ion of po wer spectra: A method based on time av eraging over short, modified periodog rams, ” IEEE T ransactions on Audio and Electr oacoustics , vol. 15, no. 2, pp. 70–73, 1967. [17] L. Ljung, System Identi fication: Theory for the User , 2nd Edition . Prentic e Hall, 1999. [18] W . I. Ro wen, “Simplifie d mathematica l representa tions of heavy -duty gas turbines, ” Journal of E nginee ring for P ower , vol. 105, pp. 865– 869, 1983. [19] M. R. B. T av akoli, B. V ahidi, and W . Gawlik, “ An educational guide to extrac t the parameters of heavy duty gas turbines model in dynamic studies based on operational data, ” IE EE T ransactions on P ower Systems , vol. 24, no. 3, pp. 1366–1374, 2009. [20] A. Khormali, I. Y ousefi, H. Y ahyaei, and S. M. Aliyari, “Identificat ion of an industrial gas turbine based on rowe n’ s model and using multi- object iv e optimizat ion method, ” in 3rd R SI ICR OM , pp. 482–487, 2015. [21] L. B. Dav is and S. Black, “Dry lo w NO x combust ion systems for GE hea vy-duty gas turbines, ” ASME, 1996. [22] K. V enkatara man, S. E. Lewi s, J. Natarajan, S. R. Thomas, and J. V . Citeno , “F-cla ss DLN technolo gy adv ancements: DLN2. 6+, ” in T urbo Expo: P ower for Land, Sea, and Air , vol. 54631, pp. 587–594, 2011.

Original Paper

Loading high-quality paper...

Comments & Academic Discussion

Loading comments...

Leave a Comment