From: <par...@us...> - 2011-07-31 22:23:51
|
Revision: 8427 http://octave.svn.sourceforge.net/octave/?rev=8427&view=rev Author: paramaniac Date: 2011-07-31 22:23:44 +0000 (Sun, 31 Jul 2011) Log Message: ----------- control: add discrete-time support to ncfsyn (case D!=0) Modified Paths: -------------- trunk/octave-forge/main/control/devel/ncfsyn/makefile_ncfsyn.m trunk/octave-forge/main/control/devel/ncfsyn/ncfsyn.m trunk/octave-forge/main/control/devel/ncfsyn/slsb10kd.cc Added Paths: ----------- trunk/octave-forge/main/control/devel/ncfsyn/slsb10zd.cc Modified: trunk/octave-forge/main/control/devel/ncfsyn/makefile_ncfsyn.m =================================================================== --- trunk/octave-forge/main/control/devel/ncfsyn/makefile_ncfsyn.m 2011-07-31 21:57:09 UTC (rev 8426) +++ trunk/octave-forge/main/control/devel/ncfsyn/makefile_ncfsyn.m 2011-07-31 22:23:44 UTC (rev 8427) @@ -13,6 +13,7 @@ SB02OV.f SB02MV.f SB02OU.f SB02MR.f mkoctfile "-Wl,-framework" "-Wl,vecLib" \ + slsb10zd.cc \ SB10ZD.f MA02AD.f SB02OD.f select.f MB01RX.f \ MB02VD.f SB02OY.f SB02OW.f SB02OV.f SB02OU.f \ SB02MR.f MA02GD.f SB02MV.f \ No newline at end of file Modified: trunk/octave-forge/main/control/devel/ncfsyn/ncfsyn.m =================================================================== --- trunk/octave-forge/main/control/devel/ncfsyn/ncfsyn.m 2011-07-31 21:57:09 UTC (rev 8426) +++ trunk/octave-forge/main/control/devel/ncfsyn/ncfsyn.m 2011-07-31 22:23:44 UTC (rev 8427) @@ -52,7 +52,7 @@ if (isct (Gs)) # continuous-time [ak, bk, ck, dk, rcond] = slsb10id (a, b, c, d, factor); elseif (any (d(:))) # discrete-time, d != 0 - [ak, bk, ck, dk, rcond] = slsb10zd (a, b, c, d, factor); + [ak, bk, ck, dk, rcond] = slsb10zd (a, b, c, d, factor, 0.0); else # discrete-time, d == 0 [ak, bk, ck, dk, rcond] = slsb10kd (a, b, c, factor); endif @@ -269,7 +269,7 @@ %! %! FACTOR = 1.1; %! -%! G = ss (A, B, C, [], 1); +%! G = ss (A, B, C, [], 1); # value of sampling time doesn't matter %! K = ncfsyn (G, [], [], FACTOR); %! [AK, BK, CK, DK] = ssdata (K); %! @@ -299,3 +299,118 @@ %!assert (BK, BKe, 1e-4); %!assert (CK, CKe, 1e-4); %!assert (DK, DKe, 1e-4); + + +## discrete-time case D!=0, direct access to sb10zd +%!shared AK, BK, CK, DK, RCOND, AKe, BKe, CKe, DKe, RCONDe +%! A = [ 0.2 0.0 3.0 0.0 -0.3 -0.1 +%! -3.0 0.2 -0.4 -0.3 0.0 0.0 +%! -0.1 0.1 -1.0 0.0 0.0 -3.0 +%! 1.0 0.0 0.0 -1.0 -1.0 0.0 +%! 0.0 0.3 0.6 2.0 0.1 -0.4 +%! 0.2 -4.0 0.0 0.0 0.2 -2.0 ]; +%! +%! B = [ -1.0 -2.0 +%! 1.0 3.0 +%! -3.0 -4.0 +%! 1.0 -2.0 +%! 0.0 1.0 +%! 1.0 5.0 ]; +%! +%! C = [ 1.0 -1.0 2.0 -2.0 0.0 -3.0 +%! -3.0 0.0 1.0 -1.0 1.0 -1.0 +%! 2.0 4.0 -3.0 0.0 5.0 1.0 ]; +%! +%! D = [ 10.0 -6.0 +%! -7.0 8.0 +%! 2.0 -4.0 ]; +%! +%! FACTOR = 1.1; +%! +%! [AK, BK, CK, DK, RCOND] = slsb10zd (A, B, C, D, FACTOR, 0.0); +%! +%! AKe = [ 1.0128 0.5101 -0.1546 1.1300 3.3759 0.4911 +%! -2.1257 -1.4517 -0.4486 0.3493 -1.5506 -1.4296 +%! -1.0930 -0.6026 -0.1344 0.2253 -1.5625 -0.6762 +%! 0.3207 0.1698 0.2376 -1.1781 -0.8705 0.2896 +%! 0.5017 0.9006 0.0668 2.3613 0.2049 0.3703 +%! 1.0787 0.6703 0.2783 -0.7213 0.4918 0.7435 ]; +%! +%! BKe = [ 0.4132 0.3112 -0.8077 +%! 0.2140 0.4253 0.1811 +%! -0.0710 0.0807 0.3558 +%! -0.0121 -0.2019 0.0249 +%! 0.1047 0.1399 -0.0457 +%! -0.2542 -0.3472 0.0523 ]; +%! +%! CKe = [ -0.0372 -0.0456 -0.0040 0.0962 -0.2059 -0.0571 +%! 0.1999 0.2994 0.1335 -0.0251 -0.3108 0.2048 ]; +%! +%! DKe = [ 0.0629 -0.0022 0.0363 +%! -0.0228 0.0195 0.0600 ]; +%! +%! RCONDe = [ 0.27949D-03 0.66679D-03 0.45677D-01 0.23433D-07 0.68495D-01 0.76854D-01 ].'; +%! +%!assert (AK, AKe, 1e-4); +%!assert (BK, BKe, 1e-4); +%!assert (CK, CKe, 1e-4); +%!assert (DK, DKe, 1e-4); +%!assert (RCOND, RCONDe, 1e-4); + + +## discrete-time case D!=0 +%!shared AK, BK, CK, DK, RCOND, AKe, BKe, CKe, DKe, RCONDe +%! A = [ 0.2 0.0 3.0 0.0 -0.3 -0.1 +%! -3.0 0.2 -0.4 -0.3 0.0 0.0 +%! -0.1 0.1 -1.0 0.0 0.0 -3.0 +%! 1.0 0.0 0.0 -1.0 -1.0 0.0 +%! 0.0 0.3 0.6 2.0 0.1 -0.4 +%! 0.2 -4.0 0.0 0.0 0.2 -2.0 ]; +%! +%! B = [ -1.0 -2.0 +%! 1.0 3.0 +%! -3.0 -4.0 +%! 1.0 -2.0 +%! 0.0 1.0 +%! 1.0 5.0 ]; +%! +%! C = [ 1.0 -1.0 2.0 -2.0 0.0 -3.0 +%! -3.0 0.0 1.0 -1.0 1.0 -1.0 +%! 2.0 4.0 -3.0 0.0 5.0 1.0 ]; +%! +%! D = [ 10.0 -6.0 +%! -7.0 8.0 +%! 2.0 -4.0 ]; +%! +%! FACTOR = 1.1; +%! +%! G = ss (A, B, C, D, 1); # value of sampling time doesn't matter +%! K = ncfsyn (G, [], [], FACTOR); +%! [AK, BK, CK, DK] = ssdata (K); +%! +%! AKe = [ 1.0128 0.5101 -0.1546 1.1300 3.3759 0.4911 +%! -2.1257 -1.4517 -0.4486 0.3493 -1.5506 -1.4296 +%! -1.0930 -0.6026 -0.1344 0.2253 -1.5625 -0.6762 +%! 0.3207 0.1698 0.2376 -1.1781 -0.8705 0.2896 +%! 0.5017 0.9006 0.0668 2.3613 0.2049 0.3703 +%! 1.0787 0.6703 0.2783 -0.7213 0.4918 0.7435 ]; +%! +%! BKe = [ 0.4132 0.3112 -0.8077 +%! 0.2140 0.4253 0.1811 +%! -0.0710 0.0807 0.3558 +%! -0.0121 -0.2019 0.0249 +%! 0.1047 0.1399 -0.0457 +%! -0.2542 -0.3472 0.0523 ]; +%! +%! CKe = [ -0.0372 -0.0456 -0.0040 0.0962 -0.2059 -0.0571 +%! 0.1999 0.2994 0.1335 -0.0251 -0.3108 0.2048 ]; +%! +%! DKe = [ 0.0629 -0.0022 0.0363 +%! -0.0228 0.0195 0.0600 ]; +%! +%! RCONDe = [ 0.27949D-03 0.66679D-03 0.45677D-01 0.23433D-07 0.68495D-01 0.76854D-01 ].'; +%! +%!assert (AK, AKe, 1e-4); +%!assert (BK, BKe, 1e-4); +%!assert (CK, CKe, 1e-4); +%!assert (DK, DKe, 1e-4); Modified: trunk/octave-forge/main/control/devel/ncfsyn/slsb10kd.cc =================================================================== --- trunk/octave-forge/main/control/devel/ncfsyn/slsb10kd.cc 2011-07-31 21:57:09 UTC (rev 8426) +++ trunk/octave-forge/main/control/devel/ncfsyn/slsb10kd.cc 2011-07-31 22:23:44 UTC (rev 8427) @@ -86,7 +86,6 @@ int lddk = max (1, m); // arguments out - int nk; Matrix ak (ldak, n); Matrix bk (ldbk, np); Matrix ck (ldck, n); Added: trunk/octave-forge/main/control/devel/ncfsyn/slsb10zd.cc =================================================================== --- trunk/octave-forge/main/control/devel/ncfsyn/slsb10zd.cc (rev 0) +++ trunk/octave-forge/main/control/devel/ncfsyn/slsb10zd.cc 2011-07-31 22:23:44 UTC (rev 8427) @@ -0,0 +1,154 @@ +/* + +Copyright (C) 2011 Lukas F. Reichlin + +This file is part of LTI Syncope. + +LTI Syncope is free software: you can redistribute it and/or modify +it under the terms of the GNU General Public License as published by +the Free Software Foundation, either version 3 of the License, or +(at your option) any later version. + +LTI Syncope is distributed in the hope that it will be useful, +but WITHOUT ANY WARRANTY; without even the implied warranty of +MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +GNU General Public License for more details. + +You should have received a copy of the GNU General Public License +along with LTI Syncope. If not, see <http://www.gnu.org/licenses/>. + +Positive feedback controller for a discrete-time system (D != 0). +Uses SLICOT SB10ZD by courtesy of NICONET e.V. +<http://www.slicot.org> + +Author: Lukas Reichlin <luk...@gm...> +Created: August 2011 +Version: 0.1 + +*/ + +#include <octave/oct.h> +#include <f77-fcn.h> +#include "common.cc" + +extern "C" +{ + int F77_FUNC (sb10zd, SB10ZD) + (int& N, int& M, int& NP, + double* A, int& LDA, + double* B, int& LDB, + double* C, int& LDC, + double* D, int& LDD, + double& FACTOR, + double* AK, int& LDAK, + double* BK, int& LDBK, + double* CK, int& LDCK, + double* DK, int& LDDK, + double* RCOND, + double& TOL, + int* IWORK, + double* DWORK, int& LDWORK, + bool* BWORK, + int& INFO); +} + +DEFUN_DLD (slsb10zd, args, nargout, + "-*- texinfo -*-\n\ +Slicot SB10ZD Release 5.0\n\ +No argument checking.\n\ +For internal use only.") +{ + int nargin = args.length (); + octave_value_list retval; + + if (nargin != 6) + { + print_usage (); + } + else + { + // arguments in + Matrix a = args(0).matrix_value (); + Matrix b = args(1).matrix_value (); + Matrix c = args(2).matrix_value (); + Matrix d = args(3).matrix_value (); + + double factor = args(4).double_value (); + double tol = args(5).double_value (); + + int n = a.rows (); // n: number of states + int m = b.columns (); // m: number of inputs + int np = c.rows (); // np: number of outputs + + int lda = max (1, n); + int ldb = max (1, n); + int ldc = max (1, np); + int ldd = max (1, np); + + int ldak = max (1, n); + int ldbk = max (1, n); + int ldck = max (1, m); + int lddk = max (1, m); + + // arguments out + Matrix ak (ldak, n); + Matrix bk (ldbk, np); + Matrix ck (ldck, n); + Matrix dk (lddk, np); + ColumnVector rcond (6); + + // workspace + int liwork = 2 * max (n, m+np); + int ldwork = 16*n*n + 5*m*m + 7*np*np + 6*m*n + 7*m*np + + 7*n*np + 6*n + 2*(m + np) + + max (14*n+23, 16*n, 2*m-1, 2*np-1); + + OCTAVE_LOCAL_BUFFER (int, iwork, liwork); + OCTAVE_LOCAL_BUFFER (double, dwork, ldwork); + OCTAVE_LOCAL_BUFFER (bool, bwork, 2*n); + + // error indicator + int info; + + + // SLICOT routine SB10ZD + F77_XFCN (sb10zd, SB10ZD, + (n, m, np, + a.fortran_vec (), lda, + b.fortran_vec (), ldb, + c.fortran_vec (), ldc, + d.fortran_vec (), ldd, + factor, + ak.fortran_vec (), ldak, + bk.fortran_vec (), ldbk, + ck.fortran_vec (), ldck, + dk.fortran_vec (), lddk, + rcond.fortran_vec (), + tol, + iwork, + dwork, ldwork, + bwork, + info)); + + if (f77_exception_encountered) + error ("hinfsyn: slsb10zd: exception in SLICOT subroutine SB10ZD"); + + if (info != 0) + error ("hinfsyn: slsb10zd: SB10ZD returned info = %d", info); + + // resizing + ak.resize (n, n); + bk.resize (n, np); + ck.resize (m, n); + dk.resize (m, np); + + // return values + retval(0) = ak; + retval(1) = bk; + retval(2) = ck; + retval(3) = dk; + retval(4) = rcond; + } + + return retval; +} This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |