a77b88022f6398dc0bd8897a0c875a7db769f0e9
[scilab.git] / scilab / modules / cacsd / macros / ss2ss.sci
1 // Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
2 // Copyright (C) INRIA -
3 //
4 // Copyright (C) 2012 - 2016 - Scilab Enterprises
5 //
6 // This file is hereby licensed under the terms of the GNU GPL v2.0,
7 // pursuant to article 5.3.4 of the CeCILL v.2.1.
8 // This file was originally licensed under the terms of the CeCILL v2.1,
9 // and continues to be available under such terms.
10 // For more information, see the COPYING file which you should have received
11 // along with this program.
12
13 function [Sl1,right,left]=ss2ss(Sl,T,F,G,flag)
14     // State-space to state-space conversion
15     // Returns the linear system Sl1=[A1,B1,C1,D1]
16     // where A1=inv(T)*A*T,B1=inv(T)*B,C1=C*T,D1=D.
17     // Optional parameters F and G are state feedback
18     // and output injection respectively. For example,
19     // Sl1=ss2ss(Sl,T,F) returns Sl1=[A1,B1,C1,D1] with
20     // A1=inv(T)*(A+B*F)*T;B1=inv(T)*B;C1=(C+D*F)*T;D1=D;
21     // If F is given as input then right is a non singular
22     // linear system such that Sl1=Sl*right.
23     // Sl1*invsyslin(right) is a factorization of Sl.
24     // Idem for left: if F and G are given, Sl1=left*Sl*right.
25     // Example: Sl=ssrand(2,2,5); trzeros(Sl);
26     // Sl1=ss2ss(Sl,rand(5,5),rand(2,5),rand(5,2));
27     // trzeros(Sl1), trzeros(rand(2,2)*Sl1*rand(2,2))
28     // See also : projsl
29
30     [A,B,C,D]=abcd(Sl);
31     [LHS,RHS]=argn(0);
32     select RHS
33     case 2 then
34         Sl1=syslin(Sl(7),inv(T)*A*T,inv(T)*B,C*T,D);
35         right=eye(A);left=right;
36     case 3 then
37         A1=A+B*F;C1=C+D*F;
38         A1=inv(T)*A1*T;B1=inv(T)*B;C1=C1*T;D1=D
39         Sl1=syslin(Sl(7),A1,B1,C1,D1);
40         right=syslin(Sl(7),A+B*F,B,F,eye(F*B));
41         left=eye(size(C,1),size(C,1));
42     case 4 then
43         A1=A+B*F+G*C+G*D*F;C1=C+D*F;B1=B+G*D
44         A1=inv(T)*A1*T;B1=inv(T)*B1;C1=C1*T;D1=D
45         Sl1=syslin(Sl(7),A1,B1,C1,D1);
46         right=syslin(Sl(7),A+B*F,B,F,eye(F*B));
47         // Warning left is computed as [ At + G*Ct,G;Ct,I]
48         // where [At Bt; Ct Dt] is Sl1*right
49         At=A+B*F; Ct=C+D*F
50         left=syslin(Sl(7),At+G*Ct,G,Ct,eye(Ct*G));
51     case 5 then
52         if flag==1 then
53             // x in R^n , y in R^p, u in R^m
54             // output injection [ A + GC, (B+GD,-G)]
55             //                  [   C   , (D   , 0)]
56             // then feeback ( since output injection increase the
57             // size of the feedback the F matrix must be of size
58             // (m+p,n) --> F=[F1;F2] with F1 (m,n) and F2 (p,n)
59             //
60             // Sl1= [ A+GC +BF1+G*D*F1 -GF2, (B+GD,-G)]
61             //  [ C+D*F1               , (D   , 0)]
62             //
63             // We have then the following property
64             // Sl1 equiv  left*sysdiag(sys,eye(p,p))*right
65             //
66             //
67             n=size(A,"r");p=size(C,"r");
68             A1=A+G*C+[B+G*D,-G]*F;B1=[B+G*D,-G];C1=C+[D,zeros(p,p)]*F;
69             D1=[D,zeros(p,p)];
70             A1=inv(T)*A1*T;B1=inv(T)*B1;C1=C1*T;D1=D1
71             Sl1=syslin(Sl(7),A1,B1,C1,D1);
72             left=syslin(Sl(7),A+G*C,[G,-G],C,[eye(p,p),zeros(p,p)]);
73             // Now we compute the right associated to left*Sl1
74             A1=A+G*C;B1=[B+G*D,-G];C1=C;D1=[D,zeros(p,p)];
75             right=syslin(Sl(7),A1+B1*F,B1,F,eye(F*B1));
76             return
77         end
78         if flag==2 then
79             // x in R^n , y in R^p, u in R^m
80             // feedback first F of size(m,n)
81             //                      [ A+BF,B]
82             //                  [ C+DF,D]
83             // then output injection
84             // Sl1= [ A+GC +BF+G*D*F, (B+GD,-G)]
85             //  [ C+D*F         , (D   , 0)]
86             // this is a generalization of the case 4
87             // We have then the following property
88             // Sl1 equiv left*sysdiag(sys*right,eye(p,p)))
89             //
90             p = size(C,"r");
91             A1=A+B*F+G*C+G*D*F;
92             B1=[B+G*D,-G];
93             C1=C+ D*F;
94             D1=[D,zeros(p,p)];
95             A1=inv(T)*A1*T;B1=inv(T)*B1;C1=C1*T;D1=D1
96             Sl1=syslin(Sl(7),A1,B1,C1,D1);
97             right=syslin(Sl(7),A+B*F,B,F,eye(F*B));
98             // Warning left is computed as [ At + G*Ct,(G,-G);
99             //                             [    Ct    ,(I, 0)]
100             // where [At Bt; Ct Dt] is Sl1*right
101             At=A+B*F; Ct=C+D*F
102             left=syslin(Sl(7),At+G*Ct,[G,-G],Ct,[eye(Ct*G),zeros(Ct*G)]);
103         end
104     end
105 endfunction