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