search

SLICOT / SLICOT-Reference

SLICOT - A Fortran subroutines library for systems and control
BSD 3-Clause "New" or "Revised" License
48 stars 23 forks source link

Incorrect results MB03WD with increasing period P #3

Open iwoodsawyer opened 3 years ago

iwoodsawyer commented 3 years ago

The transformation matrices from the MB03WD periodic Schur decomposition are not computed correctly when the period P increases. The slicot example TMB03WD can be used to reproduce the problem, see attached data and result files. When P becomes 18 or larger, the transformation matrices getting incorrectly computed. With P=18 the calculated result NORM (Z'AZ - Aout) = 1.899 becomes very large. While with P=17 the periodic Schur decomposition is still computed correctly, because result NORM (Z'AZ - Aout) = 2.01839D-14 is still very small.

The problem seem to be with periodic Schur decomposition MB03WD and not with the periodic Hessenberg decomposition MB03VD/MB03VY, because using the slicot example TMB03VD with P=18 data shows that the periodic Hessenberg is computed correctly as the result NORM (Q'AQ - Aout) = 8.11581D-15 is very small.

Note, that the MB03WD function did not return any convergence error indicating the computation went wrong.

examples.zip

andreasvarga commented 2 years ago

For P = 17, the results seems to be correct. The computed eigenvalues are just the eigenvalues of the matrix

1.5 -.7 3.5 -.7 
1.  0.  2.  3. 
1.5 -.7 2.5 -.3 
1.  0.  2.  1. 
  -0.457125956138543 + 0.0im
 -0.3021830717156699 + 0.0im
  2.8796545139271057 - 1.3574054272339942im
  2.8796545139271057 + 1.3574054272339942im

raised to the power 17.

However, I was not able to reproduce your results for P = 18. The computed eigenvalues should be the above ones raised to the power 18, i.e.,

   7.597065506151983e-7 - 0.0im
 4.4143053112332784e-10 - 0.0im
   -8.433316071305162e7 - 1.1250751633904827e9im
   -8.433316071305162e7 + 1.1250751633904827e9im

and these are the eigenvalues computed with MB03WD

   -8.433316071306074e7 + 1.1250751633904915e9im
   -8.433316071306074e7 - 1.1250751633904915e9im
 4.4143053112340214e-10 + 0.0im
    7.59706550615205e-7 + 0.0im

which agree quite well with the above values.

I checked also using MB03BD (an alternative, to be prefered, forMB03WD), which provides essentially the same results. Also, the check of transformations is satisfactory in all cases.

andreasvarga commented 2 years ago

I also encountered failures for period 17 and above. I used matrices of the form

    rand*[1.5 -.7 3.5 -.7 
         1.  0.  2.  3. 
         1.5 -.7 2.5 -.3 
         1.  0.  2.  1.]