Movatterモバイル変換


[0]ホーム

URL:


Skip to content

Navigation Menu

Sign in
Appearance settings

Search code, repositories, users, issues, pull requests...

Provide feedback

We read every piece of feedback, and take your input very seriously.

Saved searches

Use saved searches to filter your results more quickly

Sign up
Appearance settings

Commitb10a823

Browse files
committed
Balred() now will return a vector of reduced models if orders is passed as a vector.
1 parentf9584cd commitb10a823

File tree

1 file changed

+76
-60
lines changed

1 file changed

+76
-60
lines changed

‎control/modelsimp.py

Lines changed: 76 additions & 60 deletions
Original file line numberDiff line numberDiff line change
@@ -233,7 +233,7 @@ def balred(sys, orders, method='truncate'):
233233
234234
Examples
235235
--------
236-
>>> rsys = balred(sys,order, method='truncate')
236+
>>> rsys = balred(sys,orders, method='truncate')
237237
238238
"""
239239
ifmethod=='matchdc':
@@ -257,73 +257,89 @@ def balred(sys, orders, method='truncate'):
257257
# else:
258258
dico='C'
259259

260+
rsys= []#empty list for reduced systems
261+
262+
#check if orders is a list or a scalar
263+
try:
264+
order=iter(orders)
265+
exceptTypeError:#if orders is a scalar
266+
orders= [orders]
267+
260268
#first get original system order
261269
nn=sys.A.shape[0]#no. of states
262270
mm=sys.B.shape[1]#no. of inputs
263271
rr=sys.C.shape[0]#no. of outputs
264272
#first do the schur decomposition
265273
T,V,l=schur(sys.A,sort='lhp')#l will contain the number of eigenvalues in the open left half plane, i.e. no. of stable eigenvalues
266274

267-
rorder=orders- (nn-l)
268-
ifrorder<=0:
269-
raiseValueError("System has %i unstable states which is more than ORDER(%i)"% (nn-l,orders))
270-
271-
ifl>0:#handles the stable/unstable decomposition if unstable eigenvalues are found
272-
#Author: M. Clement (mdclemen@eng.ucsd.edu) 2016
273-
print("Unstable eigenvalues found, performing stable/unstable decomposition")
274-
As=np.asmatrix(T)
275-
Bs=V.T*sys.B
276-
Cs=sys.C*V
277-
#from ref 1 eq(1) As = [A_ Ac], Bs = [B_], and Cs = [C_ C+]; _ denotes stable subsystem
278-
# [0 A+] [B+]
279-
A_=As[0:l,0:l]
280-
Ac=As[0:l,l::]
281-
Ap=As[l::,l::]
282-
283-
B_=Bs[0:l,:]
284-
Bp=Bs[l::,:]
285-
286-
C_=Cs[:,0:l]
287-
Cp=Cs[:,l::]
288-
#do some more tricky math IAW ref 1 eq(3)
289-
B_tilde=np.bmat([[B_,Ac]])
290-
D_tilde=np.bmat([[np.zeros((rr,mm)),Cp]])
291-
292-
subSys=StateSpace(A_,B_tilde,C_,D_tilde)
293-
294-
job='B'# balanced (B) or not (N)
295-
equil='N'# scale (S) or not (N)
296-
n=np.size(subSys.A,0)
297-
m=np.size(subSys.B,1)
298-
p=np.size(subSys.C,0)
299-
Nr,Ar,Br,Cr,hsv=ab09ad(dico,job,equil,n,m,p,subSys.A,subSys.B,subSys.C,nr=rorder,tol=0.0)
300-
301-
rsubSys=StateSpace(Ar,Br,Cr,np.zeros((p,m)))
302-
303-
A_r=rsubSys.A
304-
#IAW ref 1 eq(4) B^{tilde}_r = [B_r, Acr]
305-
B_r=rsubSys.B[:,0:mm]
306-
Acr=rsubSys.B[:,mm:mm+(nn-l)]
307-
C_r=rsubSys.C
308-
309-
#now put the unstable subsystem back in
310-
Ar=np.bmat([[A_r,Acr], [np.zeros((nn-l,rorder)),Ap]])
311-
Br=np.bmat([[B_r], [Bp]])
312-
Cr=np.bmat([[C_r,Cp]])
313-
314-
rsys=StateSpace(Ar,Br,Cr,sys.D)
315-
275+
foriinorders:
276+
rorder=i- (nn-l)
277+
ifrorder<=0:
278+
raiseValueError("System has %i unstable states which is more than ORDER(%i)"% (nn-l,i))
279+
280+
foriinorders:
281+
ifl>0:#handles the stable/unstable decomposition if unstable eigenvalues are found
282+
#Author: M. Clement (mdclemen@eng.ucsd.edu) 2016
283+
print("Unstable eigenvalues found, performing stable/unstable decomposition")
284+
rorder=i- (nn-l)
285+
As=np.asmatrix(T)
286+
Bs=V.T*sys.B
287+
Cs=sys.C*V
288+
#from ref 1 eq(1) As = [A_ Ac], Bs = [B_], and Cs = [C_ C+]; _ denotes stable subsystem
289+
# [0 A+] [B+]
290+
A_=As[0:l,0:l]
291+
Ac=As[0:l,l::]
292+
Ap=As[l::,l::]
293+
294+
B_=Bs[0:l,:]
295+
Bp=Bs[l::,:]
296+
297+
C_=Cs[:,0:l]
298+
Cp=Cs[:,l::]
299+
#do some more tricky math IAW ref 1 eq(3)
300+
B_tilde=np.bmat([[B_,Ac]])
301+
D_tilde=np.bmat([[np.zeros((rr,mm)),Cp]])
302+
303+
subSys=StateSpace(A_,B_tilde,C_,D_tilde)
304+
305+
job='B'# balanced (B) or not (N)
306+
equil='N'# scale (S) or not (N)
307+
n=np.size(subSys.A,0)
308+
m=np.size(subSys.B,1)
309+
p=np.size(subSys.C,0)
310+
Nr,Ar,Br,Cr,hsv=ab09ad(dico,job,equil,n,m,p,subSys.A,subSys.B,subSys.C,nr=rorder,tol=0.0)
311+
312+
rsubSys=StateSpace(Ar,Br,Cr,np.zeros((p,m)))
313+
314+
A_r=rsubSys.A
315+
#IAW ref 1 eq(4) B^{tilde}_r = [B_r, Acr]
316+
B_r=rsubSys.B[:,0:mm]
317+
Acr=rsubSys.B[:,mm:mm+(nn-l)]
318+
C_r=rsubSys.C
319+
320+
#now put the unstable subsystem back in
321+
Ar=np.bmat([[A_r,Acr], [np.zeros((nn-l,rorder)),Ap]])
322+
Br=np.bmat([[B_r], [Bp]])
323+
Cr=np.bmat([[C_r,Cp]])
324+
325+
rsys.append(StateSpace(Ar,Br,Cr,sys.D))
326+
327+
else:
328+
job='B'# balanced (B) or not (N)
329+
equil='N'# scale (S) or not (N)
330+
n=np.size(sys.A,0)
331+
m=np.size(sys.B,1)
332+
p=np.size(sys.C,0)
333+
Nr,Ar,Br,Cr,hsv=ab09ad(dico,job,equil,n,m,p,sys.A,sys.B,sys.C,nr=i,tol=0.0)
334+
335+
rsys.append(StateSpace(Ar,Br,Cr,sys.D))
336+
337+
#if orders was a scalar, just return the single reduced model, not a list
338+
iflen(orders)==1:
339+
returnrsys[0]
340+
#if orders was a list/vector, return a list/vector of systems
316341
else:
317-
job='B'# balanced (B) or not (N)
318-
equil='N'# scale (S) or not (N)
319-
n=np.size(sys.A,0)
320-
m=np.size(sys.B,1)
321-
p=np.size(sys.C,0)
322-
Nr,Ar,Br,Cr,hsv=ab09ad(dico,job,equil,n,m,p,sys.A,sys.B,sys.C,nr=orders,tol=0.0)
323-
324-
rsys=StateSpace(Ar,Br,Cr,sys.D)
325-
326-
returnrsys
342+
returnrsys
327343

328344
defminreal(sys,tol=None,verbose=True):
329345
'''

0 commit comments

Comments
 (0)

[8]ページ先頭

©2009-2025 Movatter.jp