* Bug #4731 fixed - lqr() failed when time domain of input was a number. 37/11137/3
Charlotte HECQUET [Thu, 28 Mar 2013 13:05:40 +0000 (14:05 +0100)]
Change-Id: I31f0cc34b56ca171f567c96fe23aaae5a6a58ef6

scilab/CHANGES_5.4.X
scilab/CHANGES_5.5.X
scilab/modules/cacsd/macros/lqr.sci
scilab/modules/cacsd/tests/nonreg_tests/bug_4731.dia.ref [new file with mode: 0644]
scilab/modules/cacsd/tests/nonreg_tests/bug_4731.tst [new file with mode: 0644]

index db96392..729162b 100644 (file)
@@ -1,8 +1,6 @@
                     Changes between version 5.4.0 and 5.4.1
                     =======================================
 
-Syntax changes
-===============
 
 Improvements
 =============
index 50ad75a..5dfe367 100644 (file)
@@ -61,6 +61,8 @@ Bug fixes
 
 * paramfplot2d: When input argument theta was column vector, an error occurred.
 
+* Bug #4731 fixed - lqr() failed when time domain of input was a number.
+
 * Bug #5539 fixed - sylv() help page was wrong in the discrete-time case.
 
 * Bug #6693 fixed - modulo did not accept polynomial inputs. Help page was not updated.
index fd9f255..8dc5ac5 100644 (file)
 // http://www.cecill.info/licences/Licence_CeCILL_V2-en.txt
 
 function [K,X]=lqr(P12)
-  //lqr gain for full-state LQ problem
-  //(discrete or continuous)
-  //          discrete                        continuous
-  //      |I   0   0|   | A    0   B  |      |I   0   0|   | A    0    B  |
-  //     z|0   A'  0| - |-C'C  I   -S'|    s |0   I   0| - |-C'C -A'  -S' |
-  //      |0   B'  0|   | S    0   D'D|      |0   0   0|   | S   -B'   D'D|
-  if typeof(P12)<>'state-space' then
-     error(msprintf(gettext("%s: Wrong type for input argument #%d: Linear state space expected.\n"),"lqr",1))
-  end
-
-  [A,B2,C1,D12]=P12(2:5);
-  Q=C1'*C1;R=D12'*D12;S=D12'*C1;
-  [n,nu]=size(B2);
-  [ny,n]=size(C1);
-  select P12(7)
-  case [] then
-    error(msprintf(gettext("%s: Wrong value for input argument #%d: Time domain must be ''c'' or ''d''.\n"),"lqr",1))
-  case 'c' then
-    Z=0*A;I=eye(A);O=zeros(n,nu);
-    bigE=[I,Z,O; ...
-          Z,I,O; ...
-          zeros(nu,2*n+nu)];
-
-    bigA=[A ,Z,   B2; ...
-         -Q ,-A',-S'; ...
-          S ,B2', R];
-    Ri=inv(R);
-    Left=[I,  Z,  -B2*Ri;
-          Z,  I,  S'*Ri;
-          zeros(nu,2*n), Ri];
-    LA=Left*bigA;LE=Left*bigE;N=1:2*n;
-    //[wsmall,ks1]=schur(LA(N,N),LE(N,N),'c');
-    [wsmall,ks1]=schur(LA(N,N),'c');
-    if ks1<>n then 
-      error(msprintf(gettext("%s: Stable subspace is too small.\n"),"lqr"));
+    //lqr gain for full-state LQ problem
+    //(discrete or continuous)
+    //          discrete                        continuous
+    //      |I   0   0|   | A    0   B  |      |I   0   0|   | A    0    B  |
+    //     z|0   A'  0| - |-C'C  I   -S'|    s |0   I   0| - |-C'C -A'  -S' |
+    //      |0   B'  0|   | S    0   D'D|      |0   0   0|   | S   -B'   D'D|
+    if typeof(P12)<>'state-space' then
+        error(msprintf(gettext("%s: Wrong type for input argument #%d: Linear state space expected.\n"),"lqr",1))
     end
-    X12=wsmall(1:n,1:n);phi12=wsmall(n+1:$,1:n);X=phi12/X12;
-    if rcond(X12)< 1.d-5 then 
-      warning(msprintf(gettext("%s: Bad conditionning.\n"),"lqr"));
-    end
-    K=-Ri*(B2'*X+S);
-    return;
-    //////////////////////////////
-    // Other implementation ... //
-    //////////////////////////////
-    //[Q,Z,Qd,Zd,numbeps,numbeta]=kroneck(bigE,bigA);
-    //[w,ks]=schur(bigA,bigE,'c');
-    //if ks<>n then error('lqr: stable subspace too small!');end
-    //ws=w(:,1:n);
-    //X12=ws(1:n,:);
-    //phi12=ws(n+1:2*n,:);
-    //u12=ws(2*n+1:2*n+nu,:);
-    //if rcond(X12)< 1.d-5 then warning('lqr: bad conditionning!');end
-    //K=u12/X12;
-    //X=phi12/X12;
-    //return
-  case 'd' then
-    I=eye(A);Z=0*I;
+
+    [A,B2,C1,D12]=P12(2:5);
     Q=C1'*C1;R=D12'*D12;S=D12'*C1;
-    bigE=[I,Z,0*B2; ...
-          Z,A',0*B2; ...
-          0*B2',-B2',0*B2'*B2];
+    [n,nu]=size(B2);
+    [ny,n]=size(C1);
+    if P12(7) == [] then
+        error(msprintf(gettext("%s: Wrong value for input argument #%d: Time domain must be ''c'' or ''d''.\n"),"lqr",1))
+    elseif P12(7) == 'c'
+        Z=0*A;I=eye(A);O=zeros(n,nu);
+        bigE=[I,Z,O; ...
+        Z,I,O; ...
+        zeros(nu,2*n+nu)];
 
-    bigA=[A,Z, B2; ...
-         -Q ,I, -S'; ...
-          S, 0*B2', R];
+        bigA=[A ,Z,   B2; ...
+        -Q ,-A',-S'; ...
+        S ,B2', R];
+        Ri=inv(R);
+        Left=[I,  Z,  -B2*Ri;
+        Z,  I,  S'*Ri;
+        zeros(nu,2*n), Ri];
+        LA=Left*bigA;LE=Left*bigE;N=1:2*n;
+        //[wsmall,ks1]=schur(LA(N,N),LE(N,N),'c');
+        [wsmall,ks1]=schur(LA(N,N),'c');
+        if ks1<>n then 
+            error(msprintf(gettext("%s: Stable subspace is too small.\n"),"lqr"));
+        end
+        X12=wsmall(1:n,1:n);phi12=wsmall(n+1:$,1:n);X=phi12/X12;
+        if rcond(X12)< 1.d-5 then 
+            warning(msprintf(gettext("%s: Bad conditionning.\n"),"lqr"));
+        end
+        K=-Ri*(B2'*X+S);
+        return;
+        //////////////////////////////
+        // Other implementation ... //
+        //////////////////////////////
+        //[Q,Z,Qd,Zd,numbeps,numbeta]=kroneck(bigE,bigA);
+        //[w,ks]=schur(bigA,bigE,'c');
+        //if ks<>n then error('lqr: stable subspace too small!');end
+        //ws=w(:,1:n);
+        //X12=ws(1:n,:);
+        //phi12=ws(n+1:2*n,:);
+        //u12=ws(2*n+1:2*n+nu,:);
+        //if rcond(X12)< 1.d-5 then warning('lqr: bad conditionning!');end
+        //K=u12/X12;
+        //X=phi12/X12;
+        //return
+    elseif P12(7) == 'd' | type(P12(7))==1 
+        I=eye(A);Z=0*I;
+        Q=C1'*C1;R=D12'*D12;S=D12'*C1;
+        bigE=[I,Z,0*B2; ...
+        Z,A',0*B2; ...
+        0*B2',-B2',0*B2'*B2];
 
-    Ri=inv(R);
+        bigA=[A,Z, B2; ...
+        -Q ,I, -S'; ...
+        S, 0*B2', R];
 
-    Left=[I,  Z,  -B2*Ri; ...
-          Z,  I,  S'*Ri; ...
-          zeros(nu,2*n), Ri];
-    LA=Left*bigA;LE=Left*bigE;N=1:2*n;
-    [wsmall,ks1]=schur(LA(N,N),LE(N,N),'d');
-    if ks1<>n then 
-      error(msprintf(gettext("%s: Stable subspace is too small.\n"),"lqr"));
-    end
-    X12=wsmall(1:n,1:n);phi12=wsmall(n+1:$,1:n);X=phi12/X12;
-    if rcond(X12)< 1.d-5 then 
-      warning(msprintf(gettext("%s: Bad conditionning.\n"),"lqr"));
+        Ri=inv(R);
+
+        Left=[I,  Z,  -B2*Ri; ...
+        Z,  I,  S'*Ri; ...
+        zeros(nu,2*n), Ri];
+        LA=Left*bigA;LE=Left*bigE;N=1:2*n;
+        [wsmall,ks1]=schur(LA(N,N),LE(N,N),'d');
+        if ks1<>n then 
+            error(msprintf(gettext("%s: Stable subspace is too small.\n"),"lqr"));
+        end
+        X12=wsmall(1:n,1:n);phi12=wsmall(n+1:$,1:n);X=phi12/X12;
+        if rcond(X12)< 1.d-5 then 
+            warning(msprintf(gettext("%s: Bad conditionning.\n"),"lqr"));
+        end
+        K=-pinv(B2'*X*B2+R)*(B2'*X*A+S);
+        return
+        ////////////////////
+        // Other form ... //
+        ////////////////////
+        //[w,ks]=schur(bigA,bigE,'d');
+        //if ks<>n then error('lqr: stable subspace too small!');end
+        //ws=w(:,1:n);
+        //X12=ws(1:n,:);
+        //phi12=ws(n+1:2*n,:);
+        //u12=ws(2*n+1:2*n+nu,:);
+        //if rcond(X12)< 1.d-5 then warning('lqr: bad conditionning!');end
+        //K=u12/X12;
+        //X=phi12/X12;
+        //return
     end
-    K=-pinv(B2'*X*B2+R)*(B2'*X*A+S);
-    return
-    ////////////////////
-    // Other form ... //
-    ////////////////////
-    //[w,ks]=schur(bigA,bigE,'d');
-    //if ks<>n then error('lqr: stable subspace too small!');end
-    //ws=w(:,1:n);
-    //X12=ws(1:n,:);
-    //phi12=ws(n+1:2*n,:);
-    //u12=ws(2*n+1:2*n+nu,:);
-    //if rcond(X12)< 1.d-5 then warning('lqr: bad conditionning!');end
-    //K=u12/X12;
-    //X=phi12/X12;
-    //return
-  end
+
 endfunction
diff --git a/scilab/modules/cacsd/tests/nonreg_tests/bug_4731.dia.ref b/scilab/modules/cacsd/tests/nonreg_tests/bug_4731.dia.ref
new file mode 100644 (file)
index 0000000..d32e879
--- /dev/null
@@ -0,0 +1,25 @@
+// =============================================================================
+// Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
+// Copyright (C) 2013 - Scilab Enterprises - Charlotte HECQUET
+//
+//  This file is distributed under the same license as the Scilab package.
+// =============================================================================
+//
+// <-- Non-regression test for bug 4731 -->
+//
+// <-- Bugzilla URL -->
+// http://bugzilla.scilab.org/show_bug.cgi?id=4731
+//
+// <-- Short Description -->
+// lqr fails if dt argument of syslin is a number
+A=[0.5 0.3;0.3 0.9];B=[0.9; 0.3];
+Q=diag([2,5]);R=2;
+Big=sysdiag(Q,R);
+[w,wp]=fullrf(Big);C1=wp(:,1:2);D12=wp(:,3:$);
+P=syslin(1,A,B,C1,D12);
+[K,X]=lqr(P);
+S=spec(A+B*K);
+R=norm(A'*X*A-(A'*X*B)*pinv(B'*X*B+R)*(B'*X*A)+Q-X,1);
+assert_checkalmostequal(real(S),[0.24977769272823835; 0.47346570058644721]);
+assert_checkequal(imag(S),[0;0]);
+assert_checktrue(R<100*%eps);
diff --git a/scilab/modules/cacsd/tests/nonreg_tests/bug_4731.tst b/scilab/modules/cacsd/tests/nonreg_tests/bug_4731.tst
new file mode 100644 (file)
index 0000000..09cbbf8
--- /dev/null
@@ -0,0 +1,26 @@
+// =============================================================================
+// Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
+// Copyright (C) 2013 - Scilab Enterprises - Charlotte HECQUET
+//
+//  This file is distributed under the same license as the Scilab package.
+// =============================================================================
+//
+// <-- Non-regression test for bug 4731 -->
+//
+// <-- Bugzilla URL -->
+// http://bugzilla.scilab.org/show_bug.cgi?id=4731
+//
+// <-- Short Description -->
+// lqr fails if dt argument of syslin is a number
+
+A=[0.5 0.3;0.3 0.9];B=[0.9; 0.3];
+Q=diag([2,5]);R=2;
+Big=sysdiag(Q,R);
+[w,wp]=fullrf(Big);C1=wp(:,1:2);D12=wp(:,3:$);
+P=syslin(1,A,B,C1,D12);
+[K,X]=lqr(P);
+S=spec(A+B*K);
+R=norm(A'*X*A-(A'*X*B)*pinv(B'*X*B+R)*(B'*X*A)+Q-X,1);
+assert_checkalmostequal(real(S),[0.24977769272823835; 0.47346570058644721]);
+assert_checkequal(imag(S),[0;0]);
+assert_checktrue(R<100*%eps);