National Instruments NI MATRIXx Xmath Specification Sheet Download Page 68

Page: 68 / 127

Chapter 3

Multiplicative Error Reduction

Xmath Model Reduction Module

3-22

ni.com

The values of

(

) along the

-axis are the same as the values of

around a circle with diameter defined by [

–

–1

0] on the positive

real axis (refer to Figure 3-2). Also, the values of

along the

-axis

are the same as the values of

(

) around a circle with diameter defined by

[–

–1

0, –

0].

We can implement an arbitrary bilinear transform using the

subsys( )

function, which substitutes a given transfer function for the s- or z-domain
operator, as previously shown.

To implement

use:

gtildesys=subsys(gsys,makep([-b,1]/makep([1,-a])

To implement

use:

gsys=subsys(gtildesys,makep([b,1]/makep([1,a])

Note

The systems substituted in the previous calls to subsys invert the function

specification because these functions use backward polynomial rotation.

Any zero (or rank reduction) on the

-axis of

(

) becomes a zero (or rank

reduction) in

[

] > 0 of

, and if

(

) has a zero (or rank reduction)

at infinity, this is shifted to a zero (or rank reduction) of

at the point

–1

, again in

[

] > 0. If all poles of

(

) are inside the circle of diameter

[–

–1

0], all poles of

will be in

[

] < 0, and if

(

) has no

zero (or rank reduction) on this circle,

will have no zero (or rank

reduction) on the

-axis, including

∞

(

) is nonsingular for almost all values of

, it will be nonsingular or

have no zero or rank reduction on the circle of diameter [–

–1

0, –

for almost all choices of

. If

and

are chosen small enough,

(

) will

have all its poles inside this circle and no zero or rank reduction on it, while

then will have all poles in

[

] < 0 and no zero or rank reduction on

the

-axis, including

∞

The steps of the algorithm, when

(

) has a zero on the

-axis or at

∞

are as follows:

For small

with 0 <

–1

, form

as shown for

gtildesys

Reduce

, this being possible because

is stable and

has full rank on

, including

∞

Form

as shown for

gsys

( )

s a

–

-------------------

⎝

⎠

⎛

⎞

G s

( )

s a

---------------

⎝

⎠

⎛

⎞

( )

s a

–

-------------------

⎝

⎠

⎛

⎞

( )

s a

---------------

⎝

⎠

⎛

⎞

Summary of Contents for NI MATRIXx Xmath

Page 1: ...NI MATRIXx TM Xmath TM Model Reduction Module Xmath Model Reduction Module April 2007 370755C 01 ...

Page 2: ...3400 Lebanon 961 0 1 33 28 28 Malaysia 1800 887710 Mexico 01 800 010 0793 Netherlands 31 0 348 433 466 New Zealand 0800 553 322 Norway 47 0 66 90 76 60 Poland 48 22 3390150 Portugal 351 210 311 210 Russia 7 495 783 6851 Singapore 1800 226 5886 Slovenia 386 3 425 42 00 South Africa 27 0 11 805 8197 Spain 34 91 640 0085 Sweden 46 0 8 587 895 00 Switzerland 41 56 2005151 Taiwan 886 02 2377 2222 Thail...

Page 3: ...ration National Instruments respects the intellectual property of others and we ask our users to do the same NI software is protected by copyright and other intellectual property laws Where NI software may be used to reproduce software or other materials belonging to others you may use NI software only to reproduce materials that you may reproduce in accordance with the terms of any applicable lic...

Page 4: ...alic text denotes variables emphasis a cross reference or an introduction to a key concept Italic text also denotes text that is a placeholder for a word or value that you must supply monospace Text in this font denotes text or characters that you should enter from the keyboard sections of code programming examples and syntax examples This font is also used for the proper names of disk drives path...

Page 5: ...ankel Singular Values 1 8 Internally Balanced Realizations 1 10 Singular Perturbation 1 11 Spectral Factorization 1 13 Low Order Controller Design Through Order Reduction 1 15 Chapter 2 Additive Error Reduction Introduction 2 1 Truncation of Balanced Realizations 2 2 Reduction Through Balanced Realization Truncation 2 4 Singular Perturbation of Balanced Realization 2 5 Hankel Norm Approximation 2 ...

Page 6: ...uring Zero Error at DC 3 8 Hankel Singular Values of Phase Matrix of Gr 3 9 Further Error Bounds 3 9 Reduction of Minimum Phase Unstable G 3 9 Imaginary Axis Zeros Including Zeros at 3 10 Related Functions 3 14 mulhank 3 14 Restrictions 3 14 Algorithm 3 14 right and left 3 15 Consequences of Step 5 and Justification of Step 6 3 18 Error Bounds 3 20 Imaginary Axis Zeros Including Zeros at 3 21 Rela...

Page 7: ...ditional Background 4 20 Related Functions 4 21 Chapter 5 Utilities hankelsv 5 1 Related Functions 5 2 stable 5 2 Algorithm 5 2 Related Functions 5 4 compare 5 4 Chapter 6 Tutorial Plant and Full Order Controller 6 1 Controller Reduction 6 5 ophank 6 9 wtbalance 6 12 fracred 6 20 Appendix A Bibliography Appendix B Technical Support and Professional Services Index ...

Page 8: ...se they can include additional information that is useful for simulation plotting and signal labeling Document Organization This manual includes the following chapters Chapter 1 Introduction starts with an outline of the manual and some useful notes It also provides an overview of the Model Reduction Module describes the functions in this module and introduces nomenclature and concepts used throug...

Page 9: ...ectors are represented in lowercase G s is used to denote a transfer function of a system where s is the Laplace variable G q is used when both continuous and discrete systems are allowed H s is used to denote the frequency response over some range of frequencies of a system where s is the Laplace variable H q is used to indicate that the system can be continuous or discrete A single apostrophe fo...

Page 10: ...s all Model Reduction functions Each topic explains a function s inputs outputs and keywords in detail Refer to Chapter 2 MATRIXx Publications Online Help and Customer Support of the MATRIXx Getting Started Guide for complete instructions on using the help feature Overview The Xmath Model Reduction Module MRM provides a collection of tools for reducing the order of systems Many of the functions ar...

Page 11: ...with multiplicative errors Model reduction with frequency weighting of an additive error including controller reduction Utility functions Functions Figure 1 1 MRM Function balmoore redschur ophank balance mreduce bst mulhank wtbalance fracred hankelsv stable compare truncate Utility Functions Additive Error Model Reduction Multiplicative Model Reduction Frequency Weighted Model Reduction ...

Page 12: ...tem must be linear continuous time and stable with full rank along the jω axis including infinity compare Must be a state space system fracred A state space system must be linear and continuous hankelsv A system must be linear and stable mreduce A submatrix of a matrix must be nonsingular for continuous systems and variant for discrete systems mulhank A state space system must be linear continuous...

Page 13: ...ipted asterisk denotes matrix transposition and complex conjugation λmax A for a square matrix A denotes the maximum eigenvalue presuming there are no complex eigenvalues Reλi A and λi A for a square matrix A denote an arbitrary real part and an arbitrary magnitude of an eigenvalue of A for a transfer function X denotes An all pass transfer function W s is one where for all ω to each pole there co...

Page 14: ... grammians of the realization defined by A B C D Sometimes in the code WC is used for P and WO for Q They have a number of properties P 0 with P 0 if and only if A B is controllable Q 0 with Q 0 if and only if A C is observable and With vec P denoting the column vector formed by stacking column 1 of P on column 2 on column 3 and so on and denoting Kronecker product The controllability grammian can...

Page 15: ...he transfer function matrix has a finite impulse response I A A vec P vec BB lyapunov can be used to evaluate P and Q Hankel Singular Values If P Q are the controllability and observability grammians of a transfer function matrix in continuous or discrete time the Hankel Singular Values are the quantities λi 1 2 PQ Notice the following All eigenvalues of PQ are nonnegative and so are the Hankel si...

Page 16: ...el singular values in the various functions These procedures are summarized in Table 1 2 Table 1 2 Calculating Hankel Singular Values balance For a discussion of the balancing algorithm refer to the Internally Balanced Realizations section the Hankel singular values are given by diag R1 2 HSV balmoore For a discussion of the balancing algorithm refer to the Internally Balanced Realizations section...

Page 17: ...ned at the expense of difficulty of observation and conversely Balanced realizations are those when ease of control has been balanced against ease of observation Given an arbitrary realization there are a number of ways of finding a state variable coordinate transformation bringing it to balanced form A good survey of the available algorithms for balancing is in LHPW87 One of these is implemented ...

Page 18: ...ar value decomposition thus L cQLc VRU with actually V U Finally the coordinate basis transformation is given by T LcVR 1 4 resulting in T QT T 1P T 1 R1 2 Singular Perturbation A common procedure for approximating systems is the technique of Singular Perturbation The underlying heuristic reasoning is as follows Suppose there is a system with the property that the unforced responses of the system ...

Page 19: ...excited by the same u with initial condition given by x1 to and x2 to A 1 22A21x1 to A22 1B2u to then x1 and y computed from Equation 1 1 and from Equation 1 2 should be similar If Equation 1 1 and Equation 1 2 are excited with the same u have the same x1 to and Equation 1 1 has arbitrary x2 then x1 and y computed from Equation 1 1 and Equation 1 2 should be similar after a possible initial transi...

Page 20: ...y zero mean unit intensity white noise Then the output of S is a stationary process with a spectrum Φ s related to W s by 1 3 Evidently so that Φ jω is nonnegative hermitian for all ω when W jω is a scalar so is Φ jω with Φ jω W jω 2 In the matrix case Φ is singular for some ω only if W does not have full rank there and in the scalar case only if W has a zero there Spectral factorization as shown ...

Page 21: ... s 0 or in Re s 0 if Φ jω is nonsingular for all ω Moreover the particular W s previously defined is unique to within right multiplication by a constant orthogonal matrix In the scalar case this means that W s is determined to within a 1 multiplier Example 1 1 Example of Spectral Factorization Suppose Then Equation 1 3 is satisfied by which is stable and minimum phase Also Equation 1 3 is satisfie...

Page 22: ...speaking in any design procedure it is better to postpone approximation to a late step of the procedure if approximation is done early the subsequent steps of the design procedure may have unpredictable effects on the approximation errors Hence the scheme based on high order controller design followed by reduction is generally to be preferred Controller reduction should aim to preserve closed loop...

Page 23: ...odule 1 16 ni com multiplicative reduction as described in Chapter 4 Frequency Weighted Error Reduction is a sound approach Chapter 3 Multiplicative Error Reduction and Chapter 4 Frequency Weighted Error Reduction develop these arguments more fully ...

Page 24: ... will sit inside a closed loop or if one is reducing a controller that again is sitting in a closed loop focus on additive error model reduction may not be appropriate It is however extremely appropriate in considering reducing the transfer function of a filter One pertinent application comes specifically from digital filtering a great many design algorithms lead to a finite impulse response FIR f...

Page 25: ...Help truncate redschur One only can speak of internally balanced realizations for systems which are stable if the aim is to reduce a transfer function matrix G s which contains unstable poles one must additively decompose it into a stable part and unstable part reduce the stable part and then add the unstable part back in The function stable described in Chapter 5 Utilities can be used to decompos...

Page 26: ... flat it is zero at ω and may take its largest value at ω 0 so that there is in general no matching of DC gains of the original and reduced system The actual error may be considerably less than the error bound at all frequencies so that the error bound formula can be no more than an advance guide However the bound is tight when the dimension reduction is 1 and the reduction is of a continuous time...

Page 27: ... for a reduction to order n 1 when it is not just a bound but exact a maximum error of 2 Another situation where poor approximation can arise is when a highly oscillatory system is to be replaced by a system with a real pole Reduction Through Balanced Realization Truncation This section briefly describes functions that reduce balance and truncate to achieve reduction balmoore Computes an internall...

Page 28: ...iscrete time but otherwise it looks the same Figure 2 1 Different Approaches for Obtaining the Same Reduced Order Model Singular Perturbation of Balanced Realization Singular perturbation of a balanced realization is an attractive way to produce a reduced order model Suppose G s is defined by Full Order Realization Balanced Realization of Reduced Order Model in continuous time Nonbalanced Realizat...

Page 29: ...ference is that balanced truncation provides exact matching at ω but does not match at DC while singular perturbation is exactly the other way round Perfect matching at DC can be a substantial advantage especially if input signals are known to be band limited Singular perturbation can be achieved with mreduce Figure 2 1 shows the two alternative approaches For both continuous time and discrete tim...

Page 30: ...ithm available for obtaining Further the optimum which is minimizing does a reasonable job of minimizing because it can be shown that where n deg G with this bound subject to the proviso that G and are allowed to be nonzero and different at s The bound on is one half that applying for balanced truncation However It is actual error that is important in practice not bounds The Hankel norm approximat...

Page 31: ...ctive like that of balance is to generate an internally balanced state variable realization The implementations are not identical balmoore only can be applied on systems that have a stable A matrix and are controllable and observable that is minimal Checks which are rather time consuming are included The computation is intrinsically not well conditioned if Sys is nearly nonminimal The first part o...

Page 32: ... If it is larger than then the constant transfer function matrix D achieves the bound For continuous systems the actual approximation error depends on frequency but is always zero at ω In practice it is often greatest at ω 0 if the reduction of state dimension is 1 the error bound is exact with the maximum error occurring at DC The bound also is exact in the special case of a single input single o...

Page 33: ...t normally the approximation of a strictly proper system through mreduce will not be strictly proper in contrast to the situation with balmoore For discrete systems the D matrix is also normally changed so that for example a system which was strictly causal or guaranteed to contain a delay that is D 0 will be approximated by a system SysR without this property The presentation of the Hankel singul...

Page 34: ...e Alternatively it may be used after an initial application of balance or redschur If Sys was calculated from redschur and VA VD were posed as arguments then SysR is calculated as in redschur refer to the redschur section truncate should be contrasted with mreduce which achieves a reduction through a singular perturbation calculation If Sys is balanced the same error bound formulas apply though no...

Page 35: ...t that there is no balancing offers numerical advantages especially if Sys is nearly nonminimal Sys should be stable and this is checked by the algorithm In contrast to balmoore minimality of Sys that is controllability and observability is not required If the Hankel singular values of Sys are ordered as then those of SysR in the continuous time case are A restriction of the algorithm is that is r...

Page 36: ... submatrices are obtained as follows and then a singular value decomposition is found From these quantities the transformation matrices used for calculating SysR are defined and the reduced order system is An error bound is available In the continuous time case it is as follows Let G jω and GR jω be the transfer function matrices of Sys and SysR respectively For the continuous case V AWcWoVA Sasc ...

Page 37: ...hen the poles and zeros approximately alternate along the jω axis It is not normally tight in the discrete time case and for both continuous time and discrete time cases it is not tight if there are repeated singular values The presentation of the Hankel singular values may suggest a logical dimension for the reduced order system thus if it may be sensible to choose nsr k Related Functions ophank ...

Page 38: ... Sys nsr passes There seems to be no general rule to suggest which setting produces the more accurate approximation Gr Therefore if accuracy of approximation for a given order is critical both should be tried As noted previously if an approximation involving an unstable system is desired the default onepass 1 is specified Behaviors The following explanation deals first with the keyword onepass Sup...

Page 39: ...rt restricted to having degree ni 1 and with no restriction on the degree of the unstable part one can never obtain a lower bound on the approximation error than σni in the scalar or SISO G s case the Gr s which achieves the previous bound is unique while in the matrix or MIMO G s case the Gr s which achieves the previous bound may not be unique Glo84 The algorithm we use to find Gr s and Gu s how...

Page 40: ... balmoore or effectively by redschur it also is one half that obtained when applying mreduce to a balanced realization In general the D matrices of G and Gr are different that is G Gr in contrast to balmoore and redschur Similarly G 0 Gr 0 in general in contrast to the result when mreduce is applied to a balanced realization The price paid for obtaining a smaller error bound overall through Hankel...

Page 41: ...oincide with one of 0 n1 n2 an error message is obtained generally all the σi are different so the occurrence of error messages will be rare The next step of the algorithm is to calculate the sum G s Gr s Gu s following SCL90 A separate function ophred is called for this purpose The controllability and observability grammians P and Q are found in the usual way AP PA BB QA A Q C C and then a singul...

Page 42: ...property The particular choice in the algorithm of Glo84 and flagged in corollary 7 3 of Glo84 is equivalent to the previous construction in the sense of yielding the same though the actual formulas used here and in Glo84 for the construction procedure are quite different In a number of situations including the case of scalar SISO G s this is the only choice The next step of the algorithm is to ca...

Page 43: ...t is larger than then one chooses This ensures satisfaction of the error bound for G Gr given previously because Multipass Algorithm We now explain the multipass algorithm For simplicity in first explaining the idea suppose that the Hankel singular values at every stage or pass are distinct 1 Find a stable order ns 1 approximation Gn 1 s of G s with This can be achieved by the algorithm already gi...

Page 44: ...ns k No separate optimization of the D matrix of Gnsr is required The approximation error bound is the same as for the first algorithm The actual approximation error may be different Notice that this second algorithm does not calculate an unstable Gu s such that Discrete Time Systems No special algorithm is presented for discrete time systems Any stable discrete time transfer function matrix G z c...

Page 45: ...s domain equivalent with sysS subsys sysZ bilinsys There is no simple rule for choosing α the choice α 1 is probably as good as any The orders of Gr and Gu are the same as those of Hr and Hu respectively The error formulas are as follows Impulse Response Error If Gr is determined by the first single pass algorithm the impulse response error for t 0 between the impulse responses of G and Gr can be ...

Page 46: ... of approximation of the impulse response too Unstable System Approximation A transfer function G s with stable and unstable poles can be reduced by applying stable to separate G s into a stable and unstable part The former is reduced and then the unstable part can be added back on For additional information about the ophank function refer to the Xmath Help Related Functions stable redschur balmoo...

Page 47: ... The more specific reason arises in considering the problem of plant approximation with a high order possibly very high order plant being initially prescribed with no controller having been designed and with a requirement to provide a simpler model of the plant possibly to allow controller design Consider the arrangement of Figure 3 1 controller C s designed for G s j with G s I Δ G s Figure 3 1 C...

Page 48: ...incipal one is that if we are reducing the plant without knowledge of the controller we cannot calculate the measure because we do not know C jω Nevertheless one could presume that for a well designed system will be close to I over the operating bandwidth of the system and have smaller norm than 1 tending to zero as ω in fact outside the operating bandwidth of the system This suggests that in the ...

Page 49: ...m approximation Some of the similarities and differences are as follows When the errors are small the error bound formula for bst is about one half of that for bst With bst the actual multiplicative error as a function of frequency goes to zero as ω or after using an optional transformation given in the algorithm description to zero as ω 0 with mulhank the error tends to be more uniform as a funct...

Page 50: ...ms are accepted for discrete systems use makecontinuous before calling bst then discretize the result Sys bst makecontinuous SysD SysD discretize Sys Algorithm The modifications described in this section allow you to circumvent the previous restrictions The objective of the algorithm is to approximate a high order stable transfer function matrix G s by a lower order G r s with in the square G s ca...

Page 51: ...eros of G s in Re s 0 State variable realizations of W G and the stable part of W 1 s G s can be connected in a nice way and when the stable part of W 1 s G s has a balanced realization we say that the realization of G is stochastically balanced Truncating the balanced realization of the stable part of W 1 s G s induces a corresponding truncation in the realization of G s and the truncated realiza...

Page 52: ... the nonsingularity condition on G jω G jω will normally result in an error message that no stabilizing solution exists To obtain the best numerical results singriccati is invoked with the keyword method schur Although DW CW are not needed for the remainder of the algorithm they are simply determined in the square case by with minor modification in the nonsquare case The real point of the algorith...

Page 53: ...ep 4 is identical with that used in redschur except the matrices P Q which determine VA VD and so forth are the controllability and observability grammians of CW sI A 1B rather than of C sI A 1B the controllability grammian of G s and the observability grammian of W s The error formula WaS90 is 3 2 All νi obey νi 1 One can only eliminate νi where νi 1 Hence if nsr is chosen so that νnsr 1 1 the al...

Page 54: ...ever the error bound is not guaranteed to be tight except when nsr ns 1 Securing Zero Error at DC The error G 1 G Gr as a function of frequency is always zero at ω When the algorithm is being used to approximate a high order plant by a low order plant it may be preferable to secure zero error at ω 0 A method for doing this is discussed in GrA90 for our purposes 1 We need a bilinear transformation ...

Page 55: ...f the error measures or for plant reduction as opposed to or The BST algorithm ensures that in addition to Equation 3 2 there holds WaS90a which also means that for a scalar system and if the bound is small Reduction of Minimum Phase Unstable G For square minimum phase but not necessarily stable G it also is possible to use this algorithm with minor modification to try to minimize for Gr of a cert...

Page 56: ...stable and of full rank on the jω axis 2 Form Hr of the desired order to minimize approximately 3 Set Gr H 1 r Observe that The reduced order Gr will have the same poles in Re s 0 as G and be minimum phase Imaginary Axis Zeros Including Zeros at We shall now explain how to handle the reduction of G s which has a rank drop at s or on the jω axis The key is to use a bilinear transformation Saf87 Con...

Page 57: ...round a circle with diameter defined by b 1 j0 a j0 Figure 3 3 Bilinear Mapping from G s to Case 2 We can implement an arbitrary bilinear transform using the subsys function which substitutes a given transfer function for the s or z domain operator To implement use gtildesys subsys gsys makep b 1 makep 1 a To implement use gsys subsys gtildesys makep b 1 makep 1 a Note The systems substituted in t...

Page 58: ...ecause is stable and has full rank on s jω including ω 3 Form as shown for gsys The error will be overbounded by the error and Gr will contain the same zeros in Re s 0 as G If there is no zero or rank reduction of G s at the origin one can take a 0 and b 1 bandwidth over which a good approximation of G s is needed and at the very least b 1 sufficiently large that the poles of G s lie in the circle...

Page 59: ... along the jω axis corresponds to multiplicative approximation of G s around a circle in the right half plane touching the jω axis at the origin For those points on the jω axis near the circle there will be good multiplicative approximation of G jω If it is desired to have a good approximation of G s over an interval jΩ jΩ then a choice such as ε 1 5 Ω or 10 Ω needs to be made Reduction then proce...

Page 60: ... following restrictions The user must ensure that the input system is stable and nonsingular at s infinity The algorithm may be problematic if the input system has a zero on the jω axis Only continuous systems are accepted for discrete systems use makecontinuous before calling mulhank then discretize the result Sys mulhank makecontinuous SysD SysD discretize Sys Algorithm The objective of the algo...

Page 61: ...f the stable part of W 1 s G s and multiplicative reduction of G s The reduction procedure then can be repeated on the new phase function of the just found approximation to obtain a further reduction again in G s right and left A description of the algorithm for the keyword right follows It is based on ideas of Glo86 in part developed in GrA86 and further developed in SaC88 The procedure is almost...

Page 62: ...ankel singular values νi of are computed by calling hankelsv The value of nsr is obtained if not prespecified either by prompting the user or by the error bound formula GrA89 Gre88 Glo86 3 3 with νi νi 1 being assumed If νk νk 1 νk r for some k that is νk has multiplicity greater than unity then νk appears once only in the previous error bound formula In other words the number of terms in the prod...

Page 63: ... is slightly different from that used in ophank Construct an SVD of with Σ1 of dimension ns r ns r and nonsingular Also obtain an orthogonal matrix T satisfying where and are the last r rows of and the state variable matrices appearing in a balanced realization of It is possible to calculate T without evaluating as it turns out refer to AnJ and the algorithm does this Now with there holds F s Cw s...

Page 64: ...gh the state variable formulas and Continue the reduction procedure starting with and repeating the process till Gr of the desired degree nsr is obtained For example in the second iteration is given by 3 4 Consequences of Step 5 and Justification of Step 6 A number of properties are true is of order ns r with 3 5 F p s F s vns 1 F s F s G W G s G s W s F s F s W s I vnsT I vnsT 1 W s F s F s G s G...

Page 65: ...e first as r Hankel singular values of F has the same zeros in Re s 0 as G s These properties mean that one is immediately positioned to repeat the reduction procedure on with almost all needed quantities being on hand W s G s W s W s G s G s P A F A FP B FB F BW P CG BG DG B FD V1C P DC F B WU1Σ1 B F I vnsT D Q A F A F Q C FC F CW DG 1 CG B W Q I vnsT I vnsT 1 C F D I vnsT 1 DC F B WU1 Σ1 B FD V1...

Page 66: ...hen The error bound Equation 3 3 is only exact when there is a single reduction step Normally this algorithm has a lower error bound than bst in particular if the νi are all distinct and the error bounds are approximately G G G 2 G 3 G 1 G G 3 G 1 G G G 1 G G 1 G G 2 G 1 G G 1 G 2G 2 1 G 2 G 3 G 1 G G 1 G G I 1 vns G 1 G 2 1 vns 1 G 2 1 G 3 1 vns 2 G 1 G G 3 vns 1 vns vns 1 1 vns 1 vns 2 1 vns 1 v...

Page 67: ...d whereas no such particular property of the error holds for mulhank Imaginary Axis Zeros Including Zeros at When G jω is singular or zero on the jω axis or at reduction can be handled in the same manner as explained for bst The key is to use a bilinear transformation Saf87 Consider the bilinear map defined by where 0 a b 1 and mapping G s into through 86 9 vi dB 20log10 G nsr G i nsr 1 ns 8 69 vi...

Page 68: ...0 of and if G s has a zero or rank reduction at infinity this is shifted to a zero or rank reduction of at the point b 1 again in Re s 0 If all poles of G s are inside the circle of diameter b 1 j0 a j0 all poles of will be in Re s 0 and if G s has no zero or rank reduction on this circle will have no zero or rank reduction on the jω axis including ω If G s is nonsingular for almost all values of ...

Page 69: ... thus a 0 the zeros of approach those of G s Thus for sufficiently small b one or more zeros of may be identified as lying on the imaginary axis The remedy is to increase a and or b above the desirable values The previous procedure for handling jω axis zeros or zeros at infinity will be deficient if the number of such zeros is the same as the order of G s for example if G s 1 d s for some stable d...

Page 70: ...1 makep 1 0 Notice that the number of zeros of G s in the circle of diameter 0 ε 1 j0 sets a lower bound on the degree of Gr s for such zeros become right half plane zeros of and must be preserved by bst Zeros at s are never in this circle so a procedure for reducing G s 1 d s is available There is one potential source of failure of the algorithm Because G s is stable certainly will be as its pole...

Page 71: ...atrices Their presence reflects a desire that the approximation process be more accurate at certain frequencies where V or W have large singular values than at others where they have small singular values For scalar G jω all the indices above are effectively the same with the effective weight just V jω W jω or W jω V jω When the system G is processing signals which do not have a flat spectrum and ...

Page 72: ...closed loop systems formed from the original and reduced order controller together with the plant This behavior is dependent in part on the plant and so one would expect that a procedure for approximating controllers ought in some way to reflect the plant This can be done several ways as described in the Controller Robustness Result section The following result is a trivial variant of one in Vid85...

Page 73: ... frequency than well away from it a fact very much in accord with classical control where one learns the importance of good loop shaping round this frequency The above measures EIS and EOS are advanced after a consideration of stability and the need for its preservation in approximating C by Cr If one takes the viewpoint that the important thing to preserve is the closed loop transfer function mat...

Page 74: ...s themselves deserve certain comments The two stability based measures EIS and EOS are derived from a sufficiency condition for stability rather than a necessity and sufficiency condition and so capture stability a little crudely For any constant nonsingular N the error measure EIS can be replaced by and the robustness result remains valid Use of an N may improve or worsen the quality of the appro...

Page 75: ... take a stable polynomial of the same degree as d and then represent the controller as a ratio of two stable transfer functions namely Now is the numerator and the denominator Write as Then we have the equivalence shown in Figure 4 1 Figure 4 1 Controller Representation Through Stable Fractions Evidently C s can be formed by completing the following steps 1 Construction of the one input two output...

Page 76: ... s has been determined by an LQG optimal design this spectrum turns out to be white that is independent of frequency so that no weight apart perhaps from scaling is needed A second possibility is to use a weight based on a stability robustness measure These points are now discussed in more detail To understand the construction of a natural fractional representation for C s suppose that and let KR ...

Page 77: ...elated values all stable transfer function matrices In particular For matrix C s the left and right matrix fraction descriptions are distinct entities It is the right MFD which corresponds to Figure 4 1 refer to Figure 4 2 Figure 4 2 C s Implemented to Display Right MFD Representation C s DL 1 s NL s NR s DR 1 s DL I KR sI A KEC 1 B NL KR sI A KEC 1 KE NR KR sI A BKR 1 KE DR I C sI A BKR 1 KE C s ...

Page 78: ... motivates the possibility of using no frequency dependent weighting in approximating G s but observe that after approximating the signal will no longer be white noise so that argument is questionable Simple appeal to duality motivates using no frequency dependent weighting for H s These are two of the options offered by fracred Two more fracred options depend on examining stability robustness the...

Page 79: ...ations for frequency weighted balanced truncation of G and subsequent construction of Cr s are exceptionally easy using this weight The second fracred option is the dual of this The error measure is where It is possible to argue heuristically the relevance of these error measures from a second point of view It turns out that KR C sI A BKR 1 KE C s P ˆ s P s I I P G 1 P CsI A KEC 1 B I C sI A KEC 1...

Page 80: ...efer to LiA90 More precisely if the LTR design is performed with input noise or process noise weighting tending to infinity reduction with fracred and type left stab which uses the error measure leads to effectively the same reduction as wtbalance with the type input stab If the LTR design is performed with state or output weighting tending to infinity in the index determining the state feedback l...

Page 81: ...ures through frequency weighted balanced truncation These are shown in Table 4 1 These error measures have certain interpretations as shown in Table 4 2 In case C s is not a compensator in a closed loop and the error measure is of interest you can work with type input spec and C V in lieu of C and V There is no restriction on the stability of C s or indeed of P s in the algorithm though if C s is ...

Page 82: ...n Cr is stabilizing for P The smaller the error measure is the greater the stability robustness output stab A similar stability robustness argument but based on breaking the loop at the controller input indicates that if C is stabilizing for P and the error measure is less that 1 then Cr is stabilizing for P The smaller the error measure is the greater the stability robustness match If T PC I PC 1...

Page 83: ...Cr s 6 Check closed loop stability with Cr s introduced in place of C s at least in case C s is a compensator More details of steps 3 and 4 will be given for the case when there is an input weight only The case when there is an output weight only is almost the same and the case when both weights are present is very similar refer to Enn84a for a treatment Let be a stable transfer function matrix to...

Page 84: ... colored noise that is the output of the shaping filter W s which is excited by white noise Small weighted Hankel singular values are a pointer to the possibility of eliminating states from Cs s without incurring a large error in No error bound formula is known however The actual reduction procedure is virtually the same as that of redschur except that Pcc is used Thus Schur decompositions of PccQ...

Page 85: ...se since it is stability of the closed loop system which is the key issue of concern except for type input spec but here there is only a single weight and so the theory guarantees preservation of stability Related Functions balance redschur stable fracred fracred SysCR HSV fracred Sys Kr Ke type nscr Qyy The fracred function uses fractional representations to calculate a reduction of a continuous ...

Page 86: ...d by with u the plant input and y the plant output The associated series compensator under unity negative feedback is and this may be written as a left or right MFD as follows 4 5 4 6 The reduction procedures right perf and left perf have similar rationales We shall describe right perf refer to AnM89 and LiA86 The first rationale involves observing that to reduce C s one might as well reduce its n...

Page 87: ... Figure 4 5 Internal Structure of Controller Recognize that the controller C s shown within the hazy rectangle in Figure 4 5 can be constructed by implementing and and then applying an interconnection rule connect the output of the second transfer function matrix back to the input at point X in Figure 4 5 Er s KR C sI A KE C s KR sI A 1 KE I C sI A 1 KE 1 KR sI A KEC 1 KE KE P s 1 s KR C A BKR X K...

Page 88: ... Compute the square roots of the eigenvalues of PQ Hankel singular values of the fractional representation of Equation 4 5 The maximum order permitted is the number of nonzero eigenvalues of PQ that are larger than ε 3 Introduce the order of the reduced order controller possibly by displaying the Hankel singular values HSVs to the user Broadly speaking one can throw away small HSVs but not large o...

Page 89: ... of C s will in general be quite different from those coming from the right MFD Equation 4 6 It may be possible to reduce much more with the left MFD than with the right MFD or vice versa before closed loop stability is lost As noted in the fracred input listing type left stab and right stab focus on a stability robustness measure in conjunction with Equation 4 5 and Equation 4 6 respectively Leav...

Page 90: ... controllability grammian for E s remains unaltered while the observability grammian is altered Hence Equation 4 5 at least with Qyy I and Equation 4 12 are the same while Equation 4 6 and Equation 4 13 are quite different The calculations leading to Equation 4 13 are set out in LAL90 The argument for type left perf is dual Another insight into Equation 4 14 is provided by relations set out in NJB...

Page 91: ...ariance is very large and right stab will probably be attractive if the output weighting in the performance index is very large LiA90 The reduced controllers will then actually be very similar to those obtained using wtbalance with the option input stab in the first case and output stab in the second case One example gives the HSVs summarized in Table 4 3 for an eighth order controller The most at...

Page 92: ...e and unstable parts that is given G s the function determines Gs s and Gu s the first with all poles in Re s 0 the second with all poles in Re s 0 such that The function is used within some of the other functions of the Model Reduction Module It should also be used when reduction of an unstable G s is contemplated The normal reduction functions for example balmoore or redschur require stability o...

Page 93: ...that Sys SysS SysU Continuous systems have unstable poles if real parts tol Discrete systems have unstable poles if magnitudes 1 tol The direct term D matrix is included in SysS If Sys has poles clustered near tol or 1 tol then SysS and SysU might be ill conditioned To avoid this problem choose tol to a value that is not close to the majority of poles Algorithm The algorithm begins by transforming...

Page 94: ...esent they are assigned to the unstable part of Sys In this circumstance you get the message The system has poles near or upon the jw axis for continuous systems and the following for discrete systems The system has poles near the unit circle Note If A has eigenvalues clustered near tol 1 tol in discrete time then X is likely to be ill conditioned and consequently SysS and SysU will also be ill co...

Page 95: ...table parts and then for example reduces the stable part with redschur the reduced stable part is added to the original unstable part to provide the desired system reduction Related Functions sylvester schur redschur balmoore compare respdiff compare Sys SysRed FTvec Fmin Fmax npts radians type The compare function provides a number of different graphical tests which can be used to compare two sta...

Page 96: ...e file specification appropriate to your operating system from the Xmath Commands area For example load XMATH demos mr_disc Plant and Full Order Controller The plant in question comprises four spinning disks connected by a flexible shaft A motor applies torque to the third disk and the output variable of interest is the angular displacement of the first disk The plant transfer function which is no...

Page 97: ...o guarantee stability in the presence of unstructured uncertainty More specifically the loop gain has to lie outside the shaded region shown in Figure 6 1 Figure 6 1 Loop Gain Constraints A diag 0 1 0 0 0 015 0 765 0 765 0 015 0 028 1 410 1 410 0 028 0 04 1 85 1 85 0 04 B 0 026 0 251 0 033 0 886 4 017 0 145 3 604 0 280 C 0 996 0 105 0 261 0 009 0 001 0 043 0 002 0 026 40 dB decade Frequency rad se...

Page 98: ...sc lqgcomp sys Kr Ke poles sysc ans a column vector 0 296674 0 292246 j 0 296674 0 292246 j 0 15095 0 765357 j 0 15095 0 765357 j 0 239151 1 415 j 0 239151 1 415 j 0 129808 1 84093 j 0 129808 1 84093 j The compensator itself is open loop stable A brief explanation of how Q and Wt are chosen is as follows First Q is chosen to ensure that the loop gain which would be relevant were the state measurab...

Page 99: ...s well as the constraints evidently the compensated plant meets the constraints You can enter the following commands to create a plot equivalent to Figure 6 2 sysol sys sysc svals svplot sys w radians svalsc svplot sysc w radians svalsol svplot sysol w radians plot svals x_log grid ylab line_width 2 hold plot svalsc keep plot svalsol keep f2 plot wc constr keep legend plant compensator compensated...

Page 100: ...ransfer function are 6 264 10 2 4 901 10 2 2 581 10 2 2 474 10 2 1 545 10 2 1 335 10 2 9 467 10 3 9 466 10 3 A reduction to order 2 is attempted The ending Hankel singular values that is σ3 σ4 σ8 have a sum that is not particularly small with respect to σ1 and σ2 this is an indication that problems may arise in the reduction syscr hsv redschur sysc 2 svalsRol svplot sys syscr w radians plot svalso...

Page 101: ...es in contrast to the full order controller Two are such as to cause violation of the specifications The closed loop gains differ by some 4 to 5 dB between the full order and reduced order controller in the vicinity of 0 1 radians per second The step response has overshoot of 50 as opposed to 40 and the ripple persists for longer We use the compare function refer to the compare section of Chapter ...

Page 102: ...Tutorial National Instruments Corporation 6 7 Xmath Model Reduction Module Generate Figure 6 4 compare syscl sysclr w radians type 5 f4 plot keep legend original reduced Figure 6 4 Closed Loop Gain with redschur ...

Page 103: ...apter 6 Tutorial Xmath Model Reduction Module 6 8 ni com Generate Figure 6 5 tvec 0 140 99 140 compare syscl sysclr tvec type 7 f5 plot keep legend original reduced Figure 6 5 Step Response with redschur ...

Page 104: ...hank is next used to reduce the controller with the results shown in Figures 6 6 6 7 and 6 8 Generate Figure 6 6 syscr sysu hsv ophank sysc 2 svalsrol svplot sys syscr w radians plot svalsol keep f6 plot wc constr keep grid title Open loop gain using ophank Figure 6 6 Open Loop Gain Using ophank ...

Page 105: ...del Reduction Module 6 10 ni com Generate Figure 6 7 syscl feedback sysol sysolr sys syscr sysclr feedback sysolr compare syscl sysclr w radians type 5 f7 plot keep legend original reduced Figure 6 7 Closed Loop Gain with ophank ...

Page 106: ... keep legend original reduced Figure 6 8 Step Response with ophank The open loop gain closed loop gain and step response are all inferior to those obtained with redschur This emphasizes the point that one cannot automatically assume that because the error bound formula for ophank is more attractive than that for redschur the error itself will be better for ophank ...

Page 107: ...oop gain Notice though that The error measure for wtbalance does not reflect the open loop gain it reflects the closed loop gain While the error in dB looks large as an absolute value it is not extremely so wtbalance works with additive not multiplicative error Hence it cannot be concluded that the algorithm is not working Use of the option match spec with wtbalance might be conjectured as a devic...

Page 108: ...odel Reduction Module The following function calls produce Figure 6 9 svalsrol svplot sys syscr w radians plot svalsol keep f9 plot wc constr keep grid legend reduced original constrained title Open Loop Gain Using wtbalance Figure 6 9 Open Loop Gain with wtbalance ...

Page 109: ... Reduction Module 6 14 ni com Generate Figure 6 10 syscl feedback sysol sysolr sys syscr sysclr feedback sysolr compare syscl sysclr w radians type 5 f10 plot keep legend original reduced Figure 6 10 Closed Loop Gain with wtbalance ...

Page 110: ...ith wtbalance Figures 6 9 6 10 and 6 11 are obtained for wtbalance with the option input spec Evidently there is little difference between this and the result with the option match One notices marginally better matching in the region of interest 0 1 to 5 rad per second at the expense of matching at other frequencies The weighted Hankel singular values again indicate that it is reasonable to seek a...

Page 111: ...er 6 Tutorial Xmath Model Reduction Module 6 16 ni com Generate Figure 6 12 vtf poly 0 1 10 poly 1 1 4 sysv check vtf ss convert svalsv svplot sysv w radians Figure 6 12 Frequency Response of the Weight V jω ...

Page 112: ...erate Figure 6 13 syscr sysclr hsv wtbalance sys sysc input spec 2 sysv svalsrol svplot sys syscr w radians plot svalsol keep f13 plot wc constr keep grid legend reduced original constrained title Open Loop Gain with wtbal input spec Figure 6 13 Open Loop Gain from wtbalance with input spec ...

Page 113: ... Module 6 18 ni com Generate Figure 6 14 syscl feedback sysol sysolr sys syscr sysclr feedback sysolr compare syscl sysclr w radians type 5 f14 plot keep legend original reduced Figure 6 14 System Singular Values of wtbalance with input spec ...

Page 114: ...nal Instruments Corporation 6 19 Xmath Model Reduction Module Generate Figure 6 15 tvec 0 140 99 140 compare syscl sysclr tvec type 7 f15 plot keep legend original reduced Figure 6 15 Step Response of wtbalance with input spec ...

Page 115: ...ft perf all produce instability Given the relative magnitudes of the Hankel singular values this is perhaps not surprising Figures 6 16 6 17 and 6 18 illustrate the results using right stab Generate Figure 6 16 svalsrol svplot sys syscr w radians plot svalsol keep f16 plot wc constr keep grid legend reduced original constrained title Open Loop Gain Using fracred Figure 6 16 Open Loop Gain Using fr...

Page 116: ...oration 6 21 Xmath Model Reduction Module Generate Figure 6 17 syscl feedback sysol sysolr sys syscr sysclr feedback sysolr compare syscl sysclr w radians type 5 f17 plot keep legend original reduced Figure 6 17 Closed Loop Response with fracred ...

Page 117: ... that from wtbalance with option match We can create a table to examine the values of the Hankel singular values based on different decompositions approaches set precision 3 Optional set format fixed we set a smaller precision here so we could fit the table in the manual syscr hsvrs fracred sys Kr Ke right stab 2 syscr hsvls fracred sys Kr Ke left stab 2 syscr hsvrp fracred sys Kr Ke right perf 2 ...

Page 118: ...ring hsvls right perf string hsvrp left perf string hsvlp hsvtable a rectangular matrix of strings right stab 3 308 0 728 0 112 0 078 0 024 0 018 0 011 0 010 left stab 1 403 1 331 1 133 1 092 0 965 0 549 0 526 0 313 right perf 0 034 0 016 0 013 0 010 0 004 0 004 0 000 0 000 left perf 4 907 4 874 3 846 3 781 1 225 1 175 0 505 0 041 ...

Page 119: ...nd Control Las Vegas November 1984 pp 127 132 Enn84a D F Enns Model reduction for control systems design PhD Thesis Dept of Aeronautics and Astronautics Stanford University CA USA 1984 GCP88 K Glover R F Curtain and J R Partington Realisation and approximation of linear infinite dimensional systems with error bounds SIAM J Controls and Optimization Vol 26 1988 pp 863 898 Glo84 K Glover All optimal...

Page 120: ...and R C Ward Computation of system balancing transformations and other applications of simultaneous diagonalizing algorithms IEEE Transactions on Automatic Control Vol AC 32 1987 pp 115 122 LiA86 Y Liu and BDO Anderson Controller reduction via stable factorization and balancing Int J Control Vol 44 1986 pp 507 531 LiA89 Y Liu and BDO Anderson Singular perturbation approximation of balanced systems...

Page 121: ...n continuous and discrete time model truncation Proceedings for the 29th CDC 1990 BBK88 S Boyd V Balakrishnan and P Kabamba A bisection method for computing the L norm of a transfer matrix and related problems Mathematical Controls Signals and Systems Vol 2 No 3 pp 207 219 1989 BeP79 A Berman and R J Plemmons Nonnegative Matrices in the Mathematical Sciences Computer Science and Applied Mathematic...

Page 122: ...lover and J C Doyle State space formulae for all stabilizing controllers that satisfy an L norm bound and relations to risk sensitivity Systems and Control Letters Vol 11 pp 167 172 1988 DGKF89 J C Doyle K Glover P K Khargonekar and B Francis State space solutions to standard H2 and L control problems IEEE Transactions on Automatic Control Vol AC 34 No 8 pp 831 847 August 1989 Gu80 N K Gupta Frequ...

Page 123: ...c Control Vol AC 26 February 1981 SA88 G Stein and M Athans The LQG LTR Procedure for Multivariable Control Design IEEE Transactions on Automatic Control Vol AC 32 No 2 February 1987 pp 105 114 Za81 G Zames Feedback and optimal sensitivity model reference transformations multiplicative semi norms and approximate inverses IEEE Transactions on Automatic Control Vol AC 26 pp 301 320 1981 KS72 H Kwake...

Page 124: ...ure every question receives an answer For information about other technical support options in your area visit ni com services or contact your local office at ni com contact Training and Certification Visit ni com training for self paced training eLearning virtual classrooms interactive CDs and Certification program information You also can register for instructor led hands on courses at locations...

Page 125: ...bst 1 5 1 14 3 3 C compare 1 5 5 4 controller reduction 4 2 with fractional representations 4 5 controller robustness 4 2 conventions used in the manual iv D diagnostic tools NI resources B 1 documentation conventions used in the manual iv NI resources B 1 drivers NI resources B 1 E equality bounds tight 1 7 error bound 2 7 for balanced stochastic truncation 3 8 for balanced truncation 2 2 for imp...

Page 126: ... 5 1 13 2 10 mulhank 1 5 1 14 3 14 multiplicative error 1 1 3 1 N National Instruments support and services B 1 NI support and services B 1 nomenclature 1 2 nomenclature for MRM 1 6 numerical conditioning 2 8 O ophank 1 5 2 14 discrete time systems 2 21 error formulas 2 22 impulse response error 2 22 multipass 2 20 onepass 2 18 unstable system approximation 2 23 P Padé approximation 1 5 perturbati...

Page 127: ...ility requirements 1 5 stable 1 5 5 2 sup 1 6 support technical B 1 T technical support B 1 tight equality bounds 1 7 training and certification NI resources B 1 transfer function allpass 1 6 troubleshooting NI resources B 1 truncate 1 5 2 4 2 11 U unstable zeros 2 3 W Web resources B 1 wtbalance 1 5 4 10 ...

Search results

Summary of Contents for NI MATRIXx Xmath

Reviews:

Brands by name

Popular brands

National Instruments NI MATRIXx Xmath, Specification Sheet

Search results

Summary of Contents for NI MATRIXx Xmath

Reviews:

Brands by name

Popular brands