Mercurial > octave-nkf

comparison scripts/control/hinf/hinfsyn_ric.m @ 7125:f084ba47812b

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
[project @ 2007-11-08 02:29:23 by jwe]
author jwe
date Thu, 08 Nov 2007 02:29:24 +0000
parents a1dbe9d80eee
children
comparison
equal deleted inserted replaced
7124:d07cb867891b 7125:f084ba47812b
45 ## @end table 45 ## @end table
46 ## @end deftypefn 46 ## @end deftypefn
47 47
48 function [Xinf, x_ha_err] = hinfsyn_ric (A, BB, C1, d1dot, R, ptol) 48 function [Xinf, x_ha_err] = hinfsyn_ric (A, BB, C1, d1dot, R, ptol)
49 49
50 if (nargin != 6)
51 print_usage ();
52 endif
53
50 x_ha_err = 0; # assume success 54 x_ha_err = 0; # assume success
51 Xinf = []; # default return value 55 Xinf = []; # default return value
52 n = issquare(A); 56 n = issquare (A);
53 nw = issquare(R); 57 nw = issquare (R);
54 if(rank(R) != nw) x_ha_err = 6; 58 if (rank (R) != nw)
59 x_ha_err = 6;
55 else # build hamiltonian Ha for X_inf 60 else # build hamiltonian Ha for X_inf
56 xx = ([BB; -C1'*d1dot]/R) * [d1dot'*C1, BB']; 61 xx = ([BB; -C1'*d1dot]/R) * [d1dot'*C1, BB'];
57 Ha = [A, 0*A; -C1'*C1, -A'] - xx; 62 Ha = [A, 0*A; -C1'*C1, -A'] - xx;
58 x_ha_err = 0; 63 x_ha_err = 0;
59 [d, Ha] = balance(Ha); 64 [d, Ha] = balance (Ha);
60 [u, s] = schur(Ha, "A"); 65 [u, s] = schur (Ha, "A");
61 rev = real(eig(s)); 66 rev = real (eig (s));
62 67
63 if (any(abs(rev) <= ptol)) # eigenvalues near the imaginary axis 68 if (any (abs (rev) <= ptol)) # eigenvalues near the imaginary axis
64 x_ha_err = 1; 69 x_ha_err = 1;
65 elseif (sum(rev > 0) != sum(rev < 0)) 70 elseif (sum (rev > 0) != sum (rev < 0))
66 ## unequal number of positive and negative eigenvalues 71 ## unequal number of positive and negative eigenvalues
67 x_ha_err = 2; 72 x_ha_err = 2;
68 else 73 else
69 ## compute positive Riccati equation solution 74 ## compute positive Riccati equation solution
70 u = d * u; 75 u = d * u;
71 Xinf = u(n+1:2*n,1:n) / u(1:n,1:n); 76 Xinf = u(n+1:2*n,1:n) / u(1:n,1:n);
72 if (!all(all(finite(Xinf)))) 77 if (! all (all (finite (Xinf))))
73 x_ha_err = 3; 78 x_ha_err = 3;
74 elseif (norm(Xinf-Xinf') >= 10*ptol) 79 elseif (norm (Xinf-Xinf') >= 10*ptol)
75 ## solution not symmetric 80 ## solution not symmetric
76 x_ha_err = 4; 81 x_ha_err = 4;
77 else 82 else
78 ## positive semidefinite? 83 ## positive semidefinite?
79 ## force symmetry (faster, avoids some convergence problems) 84 ## force symmetry (faster, avoids some convergence problems)
80 Xinf = (Xinf + Xinf')/2; 85 Xinf = (Xinf + Xinf')/2;
81 rev = eig(Xinf); 86 rev = eig (Xinf);
82 if (any(rev <= -ptol)) 87 if (any (rev <= -ptol))
83 x_ha_err = 5; 88 x_ha_err = 5;
84 endif 89 endif
85 endif 90 endif
86 endif 91 endif
87 endif 92 endif