LAPACK  3.4.1
LAPACK: Linear Algebra PACKage
sdrgsx.f
Go to the documentation of this file.
00001 *> \brief \b SDRGSX
00002 *
00003 *  =========== DOCUMENTATION ===========
00004 *
00005 * Online html documentation available at 
00006 *            http://www.netlib.org/lapack/explore-html/ 
00007 *
00008 *  Definition:
00009 *  ===========
00010 *
00011 *       SUBROUTINE SDRGSX( NSIZE, NCMAX, THRESH, NIN, NOUT, A, LDA, B,
00012 *                          AI, BI, Z, Q, ALPHAR, ALPHAI, BETA, C, LDC, S,
00013 *                          WORK, LWORK, IWORK, LIWORK, BWORK, INFO )
00014 * 
00015 *       .. Scalar Arguments ..
00016 *       INTEGER            INFO, LDA, LDC, LIWORK, LWORK, NCMAX, NIN,
00017 *      $                   NOUT, NSIZE
00018 *       REAL               THRESH
00019 *       ..
00020 *       .. Array Arguments ..
00021 *       LOGICAL            BWORK( * )
00022 *       INTEGER            IWORK( * )
00023 *       REAL               A( LDA, * ), AI( LDA, * ), ALPHAI( * ),
00024 *      $                   ALPHAR( * ), B( LDA, * ), BETA( * ),
00025 *      $                   BI( LDA, * ), C( LDC, * ), Q( LDA, * ), S( * ),
00026 *      $                   WORK( * ), Z( LDA, * )
00027 *       ..
00028 *  
00029 *
00030 *> \par Purpose:
00031 *  =============
00032 *>
00033 *> \verbatim
00034 *>
00035 *> SDRGSX checks the nonsymmetric generalized eigenvalue (Schur form)
00036 *> problem expert driver SGGESX.
00037 *>
00038 *> SGGESX factors A and B as Q S Z' and Q T Z', where ' means
00039 *> transpose, T is upper triangular, S is in generalized Schur form
00040 *> (block upper triangular, with 1x1 and 2x2 blocks on the diagonal,
00041 *> the 2x2 blocks corresponding to complex conjugate pairs of
00042 *> generalized eigenvalues), and Q and Z are orthogonal.  It also
00043 *> computes the generalized eigenvalues (alpha(1),beta(1)), ...,
00044 *> (alpha(n),beta(n)). Thus, w(j) = alpha(j)/beta(j) is a root of the
00045 *> characteristic equation
00046 *>
00047 *>     det( A - w(j) B ) = 0
00048 *>
00049 *> Optionally it also reorders the eigenvalues so that a selected
00050 *> cluster of eigenvalues appears in the leading diagonal block of the
00051 *> Schur forms; computes a reciprocal condition number for the average
00052 *> of the selected eigenvalues; and computes a reciprocal condition
00053 *> number for the right and left deflating subspaces corresponding to
00054 *> the selected eigenvalues.
00055 *>
00056 *> When SDRGSX is called with NSIZE > 0, five (5) types of built-in
00057 *> matrix pairs are used to test the routine SGGESX.
00058 *>
00059 *> When SDRGSX is called with NSIZE = 0, it reads in test matrix data
00060 *> to test SGGESX.
00061 *>
00062 *> For each matrix pair, the following tests will be performed and
00063 *> compared with the threshhold THRESH except for the tests (7) and (9):
00064 *>
00065 *> (1)   | A - Q S Z' | / ( |A| n ulp )
00066 *>
00067 *> (2)   | B - Q T Z' | / ( |B| n ulp )
00068 *>
00069 *> (3)   | I - QQ' | / ( n ulp )
00070 *>
00071 *> (4)   | I - ZZ' | / ( n ulp )
00072 *>
00073 *> (5)   if A is in Schur form (i.e. quasi-triangular form)
00074 *>
00075 *> (6)   maximum over j of D(j)  where:
00076 *>
00077 *>       if alpha(j) is real:
00078 *>                     |alpha(j) - S(j,j)|        |beta(j) - T(j,j)|
00079 *>           D(j) = ------------------------ + -----------------------
00080 *>                  max(|alpha(j)|,|S(j,j)|)   max(|beta(j)|,|T(j,j)|)
00081 *>
00082 *>       if alpha(j) is complex:
00083 *>                                 | det( s S - w T ) |
00084 *>           D(j) = ---------------------------------------------------
00085 *>                  ulp max( s norm(S), |w| norm(T) )*norm( s S - w T )
00086 *>
00087 *>           and S and T are here the 2 x 2 diagonal blocks of S and T
00088 *>           corresponding to the j-th and j+1-th eigenvalues.
00089 *>
00090 *> (7)   if sorting worked and SDIM is the number of eigenvalues
00091 *>       which were selected.
00092 *>
00093 *> (8)   the estimated value DIF does not differ from the true values of
00094 *>       Difu and Difl more than a factor 10*THRESH. If the estimate DIF
00095 *>       equals zero the corresponding true values of Difu and Difl
00096 *>       should be less than EPS*norm(A, B). If the true value of Difu
00097 *>       and Difl equal zero, the estimate DIF should be less than
00098 *>       EPS*norm(A, B).
00099 *>
00100 *> (9)   If INFO = N+3 is returned by SGGESX, the reordering "failed"
00101 *>       and we check that DIF = PL = PR = 0 and that the true value of
00102 *>       Difu and Difl is < EPS*norm(A, B). We count the events when
00103 *>       INFO=N+3.
00104 *>
00105 *> For read-in test matrices, the above tests are run except that the
00106 *> exact value for DIF (and PL) is input data.  Additionally, there is
00107 *> one more test run for read-in test matrices:
00108 *>
00109 *> (10)  the estimated value PL does not differ from the true value of
00110 *>       PLTRU more than a factor THRESH. If the estimate PL equals
00111 *>       zero the corresponding true value of PLTRU should be less than
00112 *>       EPS*norm(A, B). If the true value of PLTRU equal zero, the
00113 *>       estimate PL should be less than EPS*norm(A, B).
00114 *>
00115 *> Note that for the built-in tests, a total of 10*NSIZE*(NSIZE-1)
00116 *> matrix pairs are generated and tested. NSIZE should be kept small.
00117 *>
00118 *> SVD (routine SGESVD) is used for computing the true value of DIF_u
00119 *> and DIF_l when testing the built-in test problems.
00120 *>
00121 *> Built-in Test Matrices
00122 *> ======================
00123 *>
00124 *> All built-in test matrices are the 2 by 2 block of triangular
00125 *> matrices
00126 *>
00127 *>          A = [ A11 A12 ]    and      B = [ B11 B12 ]
00128 *>              [     A22 ]                 [     B22 ]
00129 *>
00130 *> where for different type of A11 and A22 are given as the following.
00131 *> A12 and B12 are chosen so that the generalized Sylvester equation
00132 *>
00133 *>          A11*R - L*A22 = -A12
00134 *>          B11*R - L*B22 = -B12
00135 *>
00136 *> have prescribed solution R and L.
00137 *>
00138 *> Type 1:  A11 = J_m(1,-1) and A_22 = J_k(1-a,1).
00139 *>          B11 = I_m, B22 = I_k
00140 *>          where J_k(a,b) is the k-by-k Jordan block with ``a'' on
00141 *>          diagonal and ``b'' on superdiagonal.
00142 *>
00143 *> Type 2:  A11 = (a_ij) = ( 2(.5-sin(i)) ) and
00144 *>          B11 = (b_ij) = ( 2(.5-sin(ij)) ) for i=1,...,m, j=i,...,m
00145 *>          A22 = (a_ij) = ( 2(.5-sin(i+j)) ) and
00146 *>          B22 = (b_ij) = ( 2(.5-sin(ij)) ) for i=m+1,...,k, j=i,...,k
00147 *>
00148 *> Type 3:  A11, A22 and B11, B22 are chosen as for Type 2, but each
00149 *>          second diagonal block in A_11 and each third diagonal block
00150 *>          in A_22 are made as 2 by 2 blocks.
00151 *>
00152 *> Type 4:  A11 = ( 20(.5 - sin(ij)) ) and B22 = ( 2(.5 - sin(i+j)) )
00153 *>             for i=1,...,m,  j=1,...,m and
00154 *>          A22 = ( 20(.5 - sin(i+j)) ) and B22 = ( 2(.5 - sin(ij)) )
00155 *>             for i=m+1,...,k,  j=m+1,...,k
00156 *>
00157 *> Type 5:  (A,B) and have potentially close or common eigenvalues and
00158 *>          very large departure from block diagonality A_11 is chosen
00159 *>          as the m x m leading submatrix of A_1:
00160 *>                  |  1  b                            |
00161 *>                  | -b  1                            |
00162 *>                  |        1+d  b                    |
00163 *>                  |         -b 1+d                   |
00164 *>           A_1 =  |                  d  1            |
00165 *>                  |                 -1  d            |
00166 *>                  |                        -d  1     |
00167 *>                  |                        -1 -d     |
00168 *>                  |                               1  |
00169 *>          and A_22 is chosen as the k x k leading submatrix of A_2:
00170 *>                  | -1  b                            |
00171 *>                  | -b -1                            |
00172 *>                  |       1-d  b                     |
00173 *>                  |       -b  1-d                    |
00174 *>           A_2 =  |                 d 1+b            |
00175 *>                  |               -1-b d             |
00176 *>                  |                       -d  1+b    |
00177 *>                  |                      -1+b  -d    |
00178 *>                  |                              1-d |
00179 *>          and matrix B are chosen as identity matrices (see SLATM5).
00180 *>
00181 *> \endverbatim
00182 *
00183 *  Arguments:
00184 *  ==========
00185 *
00186 *> \param[in] NSIZE
00187 *> \verbatim
00188 *>          NSIZE is INTEGER
00189 *>          The maximum size of the matrices to use. NSIZE >= 0.
00190 *>          If NSIZE = 0, no built-in tests matrices are used, but
00191 *>          read-in test matrices are used to test SGGESX.
00192 *> \endverbatim
00193 *>
00194 *> \param[in] NCMAX
00195 *> \verbatim
00196 *>          NCMAX is INTEGER
00197 *>          Maximum allowable NMAX for generating Kroneker matrix
00198 *>          in call to SLAKF2
00199 *> \endverbatim
00200 *>
00201 *> \param[in] THRESH
00202 *> \verbatim
00203 *>          THRESH is REAL
00204 *>          A test will count as "failed" if the "error", computed as
00205 *>          described above, exceeds THRESH.  Note that the error
00206 *>          is scaled to be O(1), so THRESH should be a reasonably
00207 *>          small multiple of 1, e.g., 10 or 100.  In particular,
00208 *>          it should not depend on the precision (single vs. double)
00209 *>          or the size of the matrix.  THRESH >= 0.
00210 *> \endverbatim
00211 *>
00212 *> \param[in] NIN
00213 *> \verbatim
00214 *>          NIN is INTEGER
00215 *>          The FORTRAN unit number for reading in the data file of
00216 *>          problems to solve.
00217 *> \endverbatim
00218 *>
00219 *> \param[in] NOUT
00220 *> \verbatim
00221 *>          NOUT is INTEGER
00222 *>          The FORTRAN unit number for printing out error messages
00223 *>          (e.g., if a routine returns IINFO not equal to 0.)
00224 *> \endverbatim
00225 *>
00226 *> \param[out] A
00227 *> \verbatim
00228 *>          A is REAL array, dimension (LDA, NSIZE)
00229 *>          Used to store the matrix whose eigenvalues are to be
00230 *>          computed.  On exit, A contains the last matrix actually used.
00231 *> \endverbatim
00232 *>
00233 *> \param[in] LDA
00234 *> \verbatim
00235 *>          LDA is INTEGER
00236 *>          The leading dimension of A, B, AI, BI, Z and Q,
00237 *>          LDA >= max( 1, NSIZE ). For the read-in test,
00238 *>          LDA >= max( 1, N ), N is the size of the test matrices.
00239 *> \endverbatim
00240 *>
00241 *> \param[out] B
00242 *> \verbatim
00243 *>          B is REAL array, dimension (LDA, NSIZE)
00244 *>          Used to store the matrix whose eigenvalues are to be
00245 *>          computed.  On exit, B contains the last matrix actually used.
00246 *> \endverbatim
00247 *>
00248 *> \param[out] AI
00249 *> \verbatim
00250 *>          AI is REAL array, dimension (LDA, NSIZE)
00251 *>          Copy of A, modified by SGGESX.
00252 *> \endverbatim
00253 *>
00254 *> \param[out] BI
00255 *> \verbatim
00256 *>          BI is REAL array, dimension (LDA, NSIZE)
00257 *>          Copy of B, modified by SGGESX.
00258 *> \endverbatim
00259 *>
00260 *> \param[out] Z
00261 *> \verbatim
00262 *>          Z is REAL array, dimension (LDA, NSIZE)
00263 *>          Z holds the left Schur vectors computed by SGGESX.
00264 *> \endverbatim
00265 *>
00266 *> \param[out] Q
00267 *> \verbatim
00268 *>          Q is REAL array, dimension (LDA, NSIZE)
00269 *>          Q holds the right Schur vectors computed by SGGESX.
00270 *> \endverbatim
00271 *>
00272 *> \param[out] ALPHAR
00273 *> \verbatim
00274 *>          ALPHAR is REAL array, dimension (NSIZE)
00275 *> \endverbatim
00276 *>
00277 *> \param[out] ALPHAI
00278 *> \verbatim
00279 *>          ALPHAI is REAL array, dimension (NSIZE)
00280 *> \endverbatim
00281 *>
00282 *> \param[out] BETA
00283 *> \verbatim
00284 *>          BETA is REAL array, dimension (NSIZE)
00285 *> \verbatim
00286 *>          On exit, (ALPHAR + ALPHAI*i)/BETA are the eigenvalues.
00287 *> \endverbatim
00288 *>
00289 *> \param[out] C
00290 *> \verbatim
00291 *>          C is REAL array, dimension (LDC, LDC)
00292 *>          Store the matrix generated by subroutine SLAKF2, this is the
00293 *>          matrix formed by Kronecker products used for estimating
00294 *>          DIF.
00295 *> \endverbatim
00296 *>
00297 *> \param[in] LDC
00298 *> \verbatim
00299 *>          LDC is INTEGER
00300 *>          The leading dimension of C. LDC >= max(1, LDA*LDA/2 ).
00301 *> \endverbatim
00302 *>
00303 *> \param[out] S
00304 *> \verbatim
00305 *>          S is REAL array, dimension (LDC)
00306 *>          Singular values of C
00307 *> \endverbatim
00308 *>
00309 *> \param[out] WORK
00310 *> \verbatim
00311 *>          WORK is REAL array, dimension (LWORK)
00312 *> \endverbatim
00313 *>
00314 *> \param[in] LWORK
00315 *> \verbatim
00316 *>          LWORK is INTEGER
00317 *>          The dimension of the array WORK.
00318 *>          LWORK >= MAX( 5*NSIZE*NSIZE/2 - 2, 10*(NSIZE+1) )
00319 *> \endverbatim
00320 *>
00321 *> \param[out] IWORK
00322 *> \verbatim
00323 *>          IWORK is INTEGER array, dimension (LIWORK)
00324 *> \endverbatim
00325 *>
00326 *> \param[in] LIWORK
00327 *> \verbatim
00328 *>          LIWORK is INTEGER
00329 *>          The dimension of the array IWORK. LIWORK >= NSIZE + 6.
00330 *> \endverbatim
00331 *>
00332 *> \param[out] BWORK
00333 *> \verbatim
00334 *>          BWORK is LOGICAL array, dimension (LDA)
00335 *> \endverbatim
00336 *>
00337 *> \param[out] INFO
00338 *> \verbatim
00339 *>          INFO is INTEGER
00340 *>          = 0:  successful exit
00341 *>          < 0:  if INFO = -i, the i-th argument had an illegal value.
00342 *>          > 0:  A routine returned an error code.
00343 *> \endverbatim
00344 *
00345 *  Authors:
00346 *  ========
00347 *
00348 *> \author Univ. of Tennessee 
00349 *> \author Univ. of California Berkeley 
00350 *> \author Univ. of Colorado Denver 
00351 *> \author NAG Ltd. 
00352 *
00353 *> \date November 2011
00354 *
00355 *> \ingroup single_eig
00356 *
00357 *  =====================================================================
00358       SUBROUTINE SDRGSX( NSIZE, NCMAX, THRESH, NIN, NOUT, A, LDA, B,
00359      $                   AI, BI, Z, Q, ALPHAR, ALPHAI, BETA, C, LDC, S,
00360      $                   WORK, LWORK, IWORK, LIWORK, BWORK, INFO )
00361 *
00362 *  -- LAPACK test routine (version 3.4.0) --
00363 *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
00364 *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
00365 *     November 2011
00366 *
00367 *     .. Scalar Arguments ..
00368       INTEGER            INFO, LDA, LDC, LIWORK, LWORK, NCMAX, NIN,
00369      $                   NOUT, NSIZE
00370       REAL               THRESH
00371 *     ..
00372 *     .. Array Arguments ..
00373       LOGICAL            BWORK( * )
00374       INTEGER            IWORK( * )
00375       REAL               A( LDA, * ), AI( LDA, * ), ALPHAI( * ),
00376      $                   ALPHAR( * ), B( LDA, * ), BETA( * ),
00377      $                   BI( LDA, * ), C( LDC, * ), Q( LDA, * ), S( * ),
00378      $                   WORK( * ), Z( LDA, * )
00379 *     ..
00380 *
00381 *  =====================================================================
00382 *
00383 *     .. Parameters ..
00384       REAL               ZERO, ONE, TEN
00385       PARAMETER          ( ZERO = 0.0E+0, ONE = 1.0E+0, TEN = 1.0E+1 )
00386 *     ..
00387 *     .. Local Scalars ..
00388       LOGICAL            ILABAD
00389       CHARACTER          SENSE
00390       INTEGER            BDSPAC, I, I1, IFUNC, IINFO, J, LINFO, MAXWRK,
00391      $                   MINWRK, MM, MN2, NERRS, NPTKNT, NTEST, NTESTT,
00392      $                   PRTYPE, QBA, QBB
00393       REAL               ABNRM, BIGNUM, DIFTRU, PLTRU, SMLNUM, TEMP1,
00394      $                   TEMP2, THRSH2, ULP, ULPINV, WEIGHT
00395 *     ..
00396 *     .. Local Arrays ..
00397       REAL               DIFEST( 2 ), PL( 2 ), RESULT( 10 )
00398 *     ..
00399 *     .. External Functions ..
00400       LOGICAL            SLCTSX
00401       INTEGER            ILAENV
00402       REAL               SLAMCH, SLANGE
00403       EXTERNAL           SLCTSX, ILAENV, SLAMCH, SLANGE
00404 *     ..
00405 *     .. External Subroutines ..
00406       EXTERNAL           ALASVM, SGESVD, SGET51, SGET53, SGGESX, SLABAD,
00407      $                   SLACPY, SLAKF2, SLASET, SLATM5, XERBLA
00408 *     ..
00409 *     .. Intrinsic Functions ..
00410       INTRINSIC          ABS, MAX, SQRT
00411 *     ..
00412 *     .. Scalars in Common ..
00413       LOGICAL            FS
00414       INTEGER            K, M, MPLUSN, N
00415 *     ..
00416 *     .. Common blocks ..
00417       COMMON             / MN / M, N, MPLUSN, K, FS
00418 *     ..
00419 *     .. Executable Statements ..
00420 *
00421 *     Check for errors
00422 *
00423       IF( NSIZE.LT.0 ) THEN
00424          INFO = -1
00425       ELSE IF( THRESH.LT.ZERO ) THEN
00426          INFO = -2
00427       ELSE IF( NIN.LE.0 ) THEN
00428          INFO = -3
00429       ELSE IF( NOUT.LE.0 ) THEN
00430          INFO = -4
00431       ELSE IF( LDA.LT.1 .OR. LDA.LT.NSIZE ) THEN
00432          INFO = -6
00433       ELSE IF( LDC.LT.1 .OR. LDC.LT.NSIZE*NSIZE / 2 ) THEN
00434          INFO = -17
00435       ELSE IF( LIWORK.LT.NSIZE+6 ) THEN
00436          INFO = -21
00437       END IF
00438 *
00439 *     Compute workspace
00440 *      (Note: Comments in the code beginning "Workspace:" describe the
00441 *       minimal amount of workspace needed at that point in the code,
00442 *       as well as the preferred amount for good performance.
00443 *       NB refers to the optimal block size for the immediately
00444 *       following subroutine, as returned by ILAENV.)
00445 *
00446       MINWRK = 1
00447       IF( INFO.EQ.0 .AND. LWORK.GE.1 ) THEN
00448 c        MINWRK = MAX( 10*( NSIZE+1 ), 5*NSIZE*NSIZE / 2-2 )
00449          MINWRK = MAX( 10*( NSIZE+1 ), 5*NSIZE*NSIZE / 2 )
00450 *
00451 *        workspace for sggesx
00452 *
00453          MAXWRK = 9*( NSIZE+1 ) + NSIZE*
00454      $            ILAENV( 1, 'SGEQRF', ' ', NSIZE, 1, NSIZE, 0 )
00455          MAXWRK = MAX( MAXWRK, 9*( NSIZE+1 )+NSIZE*
00456      $            ILAENV( 1, 'SORGQR', ' ', NSIZE, 1, NSIZE, -1 ) )
00457 *
00458 *        workspace for sgesvd
00459 *
00460          BDSPAC = 5*NSIZE*NSIZE / 2
00461          MAXWRK = MAX( MAXWRK, 3*NSIZE*NSIZE / 2+NSIZE*NSIZE*
00462      $            ILAENV( 1, 'SGEBRD', ' ', NSIZE*NSIZE / 2,
00463      $            NSIZE*NSIZE / 2, -1, -1 ) )
00464          MAXWRK = MAX( MAXWRK, BDSPAC )
00465 *
00466          MAXWRK = MAX( MAXWRK, MINWRK )
00467 *
00468          WORK( 1 ) = MAXWRK
00469       END IF
00470 *
00471       IF( LWORK.LT.MINWRK )
00472      $   INFO = -19
00473 *
00474       IF( INFO.NE.0 ) THEN
00475          CALL XERBLA( 'SDRGSX', -INFO )
00476          RETURN
00477       END IF
00478 *
00479 *     Important constants
00480 *
00481       ULP = SLAMCH( 'P' )
00482       ULPINV = ONE / ULP
00483       SMLNUM = SLAMCH( 'S' ) / ULP
00484       BIGNUM = ONE / SMLNUM
00485       CALL SLABAD( SMLNUM, BIGNUM )
00486       THRSH2 = TEN*THRESH
00487       NTESTT = 0
00488       NERRS = 0
00489 *
00490 *     Go to the tests for read-in matrix pairs
00491 *
00492       IFUNC = 0
00493       IF( NSIZE.EQ.0 )
00494      $   GO TO 70
00495 *
00496 *     Test the built-in matrix pairs.
00497 *     Loop over different functions (IFUNC) of SGGESX, types (PRTYPE)
00498 *     of test matrices, different size (M+N)
00499 *
00500       PRTYPE = 0
00501       QBA = 3
00502       QBB = 4
00503       WEIGHT = SQRT( ULP )
00504 *
00505       DO 60 IFUNC = 0, 3
00506          DO 50 PRTYPE = 1, 5
00507             DO 40 M = 1, NSIZE - 1
00508                DO 30 N = 1, NSIZE - M
00509 *
00510                   WEIGHT = ONE / WEIGHT
00511                   MPLUSN = M + N
00512 *
00513 *                 Generate test matrices
00514 *
00515                   FS = .TRUE.
00516                   K = 0
00517 *
00518                   CALL SLASET( 'Full', MPLUSN, MPLUSN, ZERO, ZERO, AI,
00519      $                         LDA )
00520                   CALL SLASET( 'Full', MPLUSN, MPLUSN, ZERO, ZERO, BI,
00521      $                         LDA )
00522 *
00523                   CALL SLATM5( PRTYPE, M, N, AI, LDA, AI( M+1, M+1 ),
00524      $                         LDA, AI( 1, M+1 ), LDA, BI, LDA,
00525      $                         BI( M+1, M+1 ), LDA, BI( 1, M+1 ), LDA,
00526      $                         Q, LDA, Z, LDA, WEIGHT, QBA, QBB )
00527 *
00528 *                 Compute the Schur factorization and swapping the
00529 *                 m-by-m (1,1)-blocks with n-by-n (2,2)-blocks.
00530 *                 Swapping is accomplished via the function SLCTSX
00531 *                 which is supplied below.
00532 *
00533                   IF( IFUNC.EQ.0 ) THEN
00534                      SENSE = 'N'
00535                   ELSE IF( IFUNC.EQ.1 ) THEN
00536                      SENSE = 'E'
00537                   ELSE IF( IFUNC.EQ.2 ) THEN
00538                      SENSE = 'V'
00539                   ELSE IF( IFUNC.EQ.3 ) THEN
00540                      SENSE = 'B'
00541                   END IF
00542 *
00543                   CALL SLACPY( 'Full', MPLUSN, MPLUSN, AI, LDA, A, LDA )
00544                   CALL SLACPY( 'Full', MPLUSN, MPLUSN, BI, LDA, B, LDA )
00545 *
00546                   CALL SGGESX( 'V', 'V', 'S', SLCTSX, SENSE, MPLUSN, AI,
00547      $                         LDA, BI, LDA, MM, ALPHAR, ALPHAI, BETA,
00548      $                         Q, LDA, Z, LDA, PL, DIFEST, WORK, LWORK,
00549      $                         IWORK, LIWORK, BWORK, LINFO )
00550 *
00551                   IF( LINFO.NE.0 .AND. LINFO.NE.MPLUSN+2 ) THEN
00552                      RESULT( 1 ) = ULPINV
00553                      WRITE( NOUT, FMT = 9999 )'SGGESX', LINFO, MPLUSN,
00554      $                  PRTYPE
00555                      INFO = LINFO
00556                      GO TO 30
00557                   END IF
00558 *
00559 *                 Compute the norm(A, B)
00560 *
00561                   CALL SLACPY( 'Full', MPLUSN, MPLUSN, AI, LDA, WORK,
00562      $                         MPLUSN )
00563                   CALL SLACPY( 'Full', MPLUSN, MPLUSN, BI, LDA,
00564      $                         WORK( MPLUSN*MPLUSN+1 ), MPLUSN )
00565                   ABNRM = SLANGE( 'Fro', MPLUSN, 2*MPLUSN, WORK, MPLUSN,
00566      $                    WORK )
00567 *
00568 *                 Do tests (1) to (4)
00569 *
00570                   CALL SGET51( 1, MPLUSN, A, LDA, AI, LDA, Q, LDA, Z,
00571      $                         LDA, WORK, RESULT( 1 ) )
00572                   CALL SGET51( 1, MPLUSN, B, LDA, BI, LDA, Q, LDA, Z,
00573      $                         LDA, WORK, RESULT( 2 ) )
00574                   CALL SGET51( 3, MPLUSN, B, LDA, BI, LDA, Q, LDA, Q,
00575      $                         LDA, WORK, RESULT( 3 ) )
00576                   CALL SGET51( 3, MPLUSN, B, LDA, BI, LDA, Z, LDA, Z,
00577      $                         LDA, WORK, RESULT( 4 ) )
00578                   NTEST = 4
00579 *
00580 *                 Do tests (5) and (6): check Schur form of A and
00581 *                 compare eigenvalues with diagonals.
00582 *
00583                   TEMP1 = ZERO
00584                   RESULT( 5 ) = ZERO
00585                   RESULT( 6 ) = ZERO
00586 *
00587                   DO 10 J = 1, MPLUSN
00588                      ILABAD = .FALSE.
00589                      IF( ALPHAI( J ).EQ.ZERO ) THEN
00590                         TEMP2 = ( ABS( ALPHAR( J )-AI( J, J ) ) /
00591      $                          MAX( SMLNUM, ABS( ALPHAR( J ) ),
00592      $                          ABS( AI( J, J ) ) )+
00593      $                          ABS( BETA( J )-BI( J, J ) ) /
00594      $                          MAX( SMLNUM, ABS( BETA( J ) ),
00595      $                          ABS( BI( J, J ) ) ) ) / ULP
00596                         IF( J.LT.MPLUSN ) THEN
00597                            IF( AI( J+1, J ).NE.ZERO ) THEN
00598                               ILABAD = .TRUE.
00599                               RESULT( 5 ) = ULPINV
00600                            END IF
00601                         END IF
00602                         IF( J.GT.1 ) THEN
00603                            IF( AI( J, J-1 ).NE.ZERO ) THEN
00604                               ILABAD = .TRUE.
00605                               RESULT( 5 ) = ULPINV
00606                            END IF
00607                         END IF
00608                      ELSE
00609                         IF( ALPHAI( J ).GT.ZERO ) THEN
00610                            I1 = J
00611                         ELSE
00612                            I1 = J - 1
00613                         END IF
00614                         IF( I1.LE.0 .OR. I1.GE.MPLUSN ) THEN
00615                            ILABAD = .TRUE.
00616                         ELSE IF( I1.LT.MPLUSN-1 ) THEN
00617                            IF( AI( I1+2, I1+1 ).NE.ZERO ) THEN
00618                               ILABAD = .TRUE.
00619                               RESULT( 5 ) = ULPINV
00620                            END IF
00621                         ELSE IF( I1.GT.1 ) THEN
00622                            IF( AI( I1, I1-1 ).NE.ZERO ) THEN
00623                               ILABAD = .TRUE.
00624                               RESULT( 5 ) = ULPINV
00625                            END IF
00626                         END IF
00627                         IF( .NOT.ILABAD ) THEN
00628                            CALL SGET53( AI( I1, I1 ), LDA, BI( I1, I1 ),
00629      $                                  LDA, BETA( J ), ALPHAR( J ),
00630      $                                  ALPHAI( J ), TEMP2, IINFO )
00631                            IF( IINFO.GE.3 ) THEN
00632                               WRITE( NOUT, FMT = 9997 )IINFO, J,
00633      $                           MPLUSN, PRTYPE
00634                               INFO = ABS( IINFO )
00635                            END IF
00636                         ELSE
00637                            TEMP2 = ULPINV
00638                         END IF
00639                      END IF
00640                      TEMP1 = MAX( TEMP1, TEMP2 )
00641                      IF( ILABAD ) THEN
00642                         WRITE( NOUT, FMT = 9996 )J, MPLUSN, PRTYPE
00643                      END IF
00644    10             CONTINUE
00645                   RESULT( 6 ) = TEMP1
00646                   NTEST = NTEST + 2
00647 *
00648 *                 Test (7) (if sorting worked)
00649 *
00650                   RESULT( 7 ) = ZERO
00651                   IF( LINFO.EQ.MPLUSN+3 ) THEN
00652                      RESULT( 7 ) = ULPINV
00653                   ELSE IF( MM.NE.N ) THEN
00654                      RESULT( 7 ) = ULPINV
00655                   END IF
00656                   NTEST = NTEST + 1
00657 *
00658 *                 Test (8): compare the estimated value DIF and its
00659 *                 value. first, compute the exact DIF.
00660 *
00661                   RESULT( 8 ) = ZERO
00662                   MN2 = MM*( MPLUSN-MM )*2
00663                   IF( IFUNC.GE.2 .AND. MN2.LE.NCMAX*NCMAX ) THEN
00664 *
00665 *                    Note: for either following two causes, there are
00666 *                    almost same number of test cases fail the test.
00667 *
00668                      CALL SLAKF2( MM, MPLUSN-MM, AI, LDA,
00669      $                            AI( MM+1, MM+1 ), BI,
00670      $                            BI( MM+1, MM+1 ), C, LDC )
00671 *
00672                      CALL SGESVD( 'N', 'N', MN2, MN2, C, LDC, S, WORK,
00673      $                            1, WORK( 2 ), 1, WORK( 3 ), LWORK-2,
00674      $                            INFO )
00675                      DIFTRU = S( MN2 )
00676 *
00677                      IF( DIFEST( 2 ).EQ.ZERO ) THEN
00678                         IF( DIFTRU.GT.ABNRM*ULP )
00679      $                     RESULT( 8 ) = ULPINV
00680                      ELSE IF( DIFTRU.EQ.ZERO ) THEN
00681                         IF( DIFEST( 2 ).GT.ABNRM*ULP )
00682      $                     RESULT( 8 ) = ULPINV
00683                      ELSE IF( ( DIFTRU.GT.THRSH2*DIFEST( 2 ) ) .OR.
00684      $                        ( DIFTRU*THRSH2.LT.DIFEST( 2 ) ) ) THEN
00685                         RESULT( 8 ) = MAX( DIFTRU / DIFEST( 2 ),
00686      $                                DIFEST( 2 ) / DIFTRU )
00687                      END IF
00688                      NTEST = NTEST + 1
00689                   END IF
00690 *
00691 *                 Test (9)
00692 *
00693                   RESULT( 9 ) = ZERO
00694                   IF( LINFO.EQ.( MPLUSN+2 ) ) THEN
00695                      IF( DIFTRU.GT.ABNRM*ULP )
00696      $                  RESULT( 9 ) = ULPINV
00697                      IF( ( IFUNC.GT.1 ) .AND. ( DIFEST( 2 ).NE.ZERO ) )
00698      $                  RESULT( 9 ) = ULPINV
00699                      IF( ( IFUNC.EQ.1 ) .AND. ( PL( 1 ).NE.ZERO ) )
00700      $                  RESULT( 9 ) = ULPINV
00701                      NTEST = NTEST + 1
00702                   END IF
00703 *
00704                   NTESTT = NTESTT + NTEST
00705 *
00706 *                 Print out tests which fail.
00707 *
00708                   DO 20 J = 1, 9
00709                      IF( RESULT( J ).GE.THRESH ) THEN
00710 *
00711 *                       If this is the first test to fail,
00712 *                       print a header to the data file.
00713 *
00714                         IF( NERRS.EQ.0 ) THEN
00715                            WRITE( NOUT, FMT = 9995 )'SGX'
00716 *
00717 *                          Matrix types
00718 *
00719                            WRITE( NOUT, FMT = 9993 )
00720 *
00721 *                          Tests performed
00722 *
00723                            WRITE( NOUT, FMT = 9992 )'orthogonal', '''',
00724      $                        'transpose', ( '''', I = 1, 4 )
00725 *
00726                         END IF
00727                         NERRS = NERRS + 1
00728                         IF( RESULT( J ).LT.10000.0 ) THEN
00729                            WRITE( NOUT, FMT = 9991 )MPLUSN, PRTYPE,
00730      $                        WEIGHT, M, J, RESULT( J )
00731                         ELSE
00732                            WRITE( NOUT, FMT = 9990 )MPLUSN, PRTYPE,
00733      $                        WEIGHT, M, J, RESULT( J )
00734                         END IF
00735                      END IF
00736    20             CONTINUE
00737 *
00738    30          CONTINUE
00739    40       CONTINUE
00740    50    CONTINUE
00741    60 CONTINUE
00742 *
00743       GO TO 150
00744 *
00745    70 CONTINUE
00746 *
00747 *     Read in data from file to check accuracy of condition estimation
00748 *     Read input data until N=0
00749 *
00750       NPTKNT = 0
00751 *
00752    80 CONTINUE
00753       READ( NIN, FMT = *, END = 140 )MPLUSN
00754       IF( MPLUSN.EQ.0 )
00755      $   GO TO 140
00756       READ( NIN, FMT = *, END = 140 )N
00757       DO 90 I = 1, MPLUSN
00758          READ( NIN, FMT = * )( AI( I, J ), J = 1, MPLUSN )
00759    90 CONTINUE
00760       DO 100 I = 1, MPLUSN
00761          READ( NIN, FMT = * )( BI( I, J ), J = 1, MPLUSN )
00762   100 CONTINUE
00763       READ( NIN, FMT = * )PLTRU, DIFTRU
00764 *
00765       NPTKNT = NPTKNT + 1
00766       FS = .TRUE.
00767       K = 0
00768       M = MPLUSN - N
00769 *
00770       CALL SLACPY( 'Full', MPLUSN, MPLUSN, AI, LDA, A, LDA )
00771       CALL SLACPY( 'Full', MPLUSN, MPLUSN, BI, LDA, B, LDA )
00772 *
00773 *     Compute the Schur factorization while swaping the
00774 *     m-by-m (1,1)-blocks with n-by-n (2,2)-blocks.
00775 *
00776       CALL SGGESX( 'V', 'V', 'S', SLCTSX, 'B', MPLUSN, AI, LDA, BI, LDA,
00777      $             MM, ALPHAR, ALPHAI, BETA, Q, LDA, Z, LDA, PL, DIFEST,
00778      $             WORK, LWORK, IWORK, LIWORK, BWORK, LINFO )
00779 *
00780       IF( LINFO.NE.0 .AND. LINFO.NE.MPLUSN+2 ) THEN
00781          RESULT( 1 ) = ULPINV
00782          WRITE( NOUT, FMT = 9998 )'SGGESX', LINFO, MPLUSN, NPTKNT
00783          GO TO 130
00784       END IF
00785 *
00786 *     Compute the norm(A, B)
00787 *        (should this be norm of (A,B) or (AI,BI)?)
00788 *
00789       CALL SLACPY( 'Full', MPLUSN, MPLUSN, AI, LDA, WORK, MPLUSN )
00790       CALL SLACPY( 'Full', MPLUSN, MPLUSN, BI, LDA,
00791      $             WORK( MPLUSN*MPLUSN+1 ), MPLUSN )
00792       ABNRM = SLANGE( 'Fro', MPLUSN, 2*MPLUSN, WORK, MPLUSN, WORK )
00793 *
00794 *     Do tests (1) to (4)
00795 *
00796       CALL SGET51( 1, MPLUSN, A, LDA, AI, LDA, Q, LDA, Z, LDA, WORK,
00797      $             RESULT( 1 ) )
00798       CALL SGET51( 1, MPLUSN, B, LDA, BI, LDA, Q, LDA, Z, LDA, WORK,
00799      $             RESULT( 2 ) )
00800       CALL SGET51( 3, MPLUSN, B, LDA, BI, LDA, Q, LDA, Q, LDA, WORK,
00801      $             RESULT( 3 ) )
00802       CALL SGET51( 3, MPLUSN, B, LDA, BI, LDA, Z, LDA, Z, LDA, WORK,
00803      $             RESULT( 4 ) )
00804 *
00805 *     Do tests (5) and (6): check Schur form of A and compare
00806 *     eigenvalues with diagonals.
00807 *
00808       NTEST = 6
00809       TEMP1 = ZERO
00810       RESULT( 5 ) = ZERO
00811       RESULT( 6 ) = ZERO
00812 *
00813       DO 110 J = 1, MPLUSN
00814          ILABAD = .FALSE.
00815          IF( ALPHAI( J ).EQ.ZERO ) THEN
00816             TEMP2 = ( ABS( ALPHAR( J )-AI( J, J ) ) /
00817      $              MAX( SMLNUM, ABS( ALPHAR( J ) ), ABS( AI( J,
00818      $              J ) ) )+ABS( BETA( J )-BI( J, J ) ) /
00819      $              MAX( SMLNUM, ABS( BETA( J ) ), ABS( BI( J, J ) ) ) )
00820      $               / ULP
00821             IF( J.LT.MPLUSN ) THEN
00822                IF( AI( J+1, J ).NE.ZERO ) THEN
00823                   ILABAD = .TRUE.
00824                   RESULT( 5 ) = ULPINV
00825                END IF
00826             END IF
00827             IF( J.GT.1 ) THEN
00828                IF( AI( J, J-1 ).NE.ZERO ) THEN
00829                   ILABAD = .TRUE.
00830                   RESULT( 5 ) = ULPINV
00831                END IF
00832             END IF
00833          ELSE
00834             IF( ALPHAI( J ).GT.ZERO ) THEN
00835                I1 = J
00836             ELSE
00837                I1 = J - 1
00838             END IF
00839             IF( I1.LE.0 .OR. I1.GE.MPLUSN ) THEN
00840                ILABAD = .TRUE.
00841             ELSE IF( I1.LT.MPLUSN-1 ) THEN
00842                IF( AI( I1+2, I1+1 ).NE.ZERO ) THEN
00843                   ILABAD = .TRUE.
00844                   RESULT( 5 ) = ULPINV
00845                END IF
00846             ELSE IF( I1.GT.1 ) THEN
00847                IF( AI( I1, I1-1 ).NE.ZERO ) THEN
00848                   ILABAD = .TRUE.
00849                   RESULT( 5 ) = ULPINV
00850                END IF
00851             END IF
00852             IF( .NOT.ILABAD ) THEN
00853                CALL SGET53( AI( I1, I1 ), LDA, BI( I1, I1 ), LDA,
00854      $                      BETA( J ), ALPHAR( J ), ALPHAI( J ), TEMP2,
00855      $                      IINFO )
00856                IF( IINFO.GE.3 ) THEN
00857                   WRITE( NOUT, FMT = 9997 )IINFO, J, MPLUSN, NPTKNT
00858                   INFO = ABS( IINFO )
00859                END IF
00860             ELSE
00861                TEMP2 = ULPINV
00862             END IF
00863          END IF
00864          TEMP1 = MAX( TEMP1, TEMP2 )
00865          IF( ILABAD ) THEN
00866             WRITE( NOUT, FMT = 9996 )J, MPLUSN, NPTKNT
00867          END IF
00868   110 CONTINUE
00869       RESULT( 6 ) = TEMP1
00870 *
00871 *     Test (7) (if sorting worked)  <--------- need to be checked.
00872 *
00873       NTEST = 7
00874       RESULT( 7 ) = ZERO
00875       IF( LINFO.EQ.MPLUSN+3 )
00876      $   RESULT( 7 ) = ULPINV
00877 *
00878 *     Test (8): compare the estimated value of DIF and its true value.
00879 *
00880       NTEST = 8
00881       RESULT( 8 ) = ZERO
00882       IF( DIFEST( 2 ).EQ.ZERO ) THEN
00883          IF( DIFTRU.GT.ABNRM*ULP )
00884      $      RESULT( 8 ) = ULPINV
00885       ELSE IF( DIFTRU.EQ.ZERO ) THEN
00886          IF( DIFEST( 2 ).GT.ABNRM*ULP )
00887      $      RESULT( 8 ) = ULPINV
00888       ELSE IF( ( DIFTRU.GT.THRSH2*DIFEST( 2 ) ) .OR.
00889      $         ( DIFTRU*THRSH2.LT.DIFEST( 2 ) ) ) THEN
00890          RESULT( 8 ) = MAX( DIFTRU / DIFEST( 2 ), DIFEST( 2 ) / DIFTRU )
00891       END IF
00892 *
00893 *     Test (9)
00894 *
00895       NTEST = 9
00896       RESULT( 9 ) = ZERO
00897       IF( LINFO.EQ.( MPLUSN+2 ) ) THEN
00898          IF( DIFTRU.GT.ABNRM*ULP )
00899      $      RESULT( 9 ) = ULPINV
00900          IF( ( IFUNC.GT.1 ) .AND. ( DIFEST( 2 ).NE.ZERO ) )
00901      $      RESULT( 9 ) = ULPINV
00902          IF( ( IFUNC.EQ.1 ) .AND. ( PL( 1 ).NE.ZERO ) )
00903      $      RESULT( 9 ) = ULPINV
00904       END IF
00905 *
00906 *     Test (10): compare the estimated value of PL and it true value.
00907 *
00908       NTEST = 10
00909       RESULT( 10 ) = ZERO
00910       IF( PL( 1 ).EQ.ZERO ) THEN
00911          IF( PLTRU.GT.ABNRM*ULP )
00912      $      RESULT( 10 ) = ULPINV
00913       ELSE IF( PLTRU.EQ.ZERO ) THEN
00914          IF( PL( 1 ).GT.ABNRM*ULP )
00915      $      RESULT( 10 ) = ULPINV
00916       ELSE IF( ( PLTRU.GT.THRESH*PL( 1 ) ) .OR.
00917      $         ( PLTRU*THRESH.LT.PL( 1 ) ) ) THEN
00918          RESULT( 10 ) = ULPINV
00919       END IF
00920 *
00921       NTESTT = NTESTT + NTEST
00922 *
00923 *     Print out tests which fail.
00924 *
00925       DO 120 J = 1, NTEST
00926          IF( RESULT( J ).GE.THRESH ) THEN
00927 *
00928 *           If this is the first test to fail,
00929 *           print a header to the data file.
00930 *
00931             IF( NERRS.EQ.0 ) THEN
00932                WRITE( NOUT, FMT = 9995 )'SGX'
00933 *
00934 *              Matrix types
00935 *
00936                WRITE( NOUT, FMT = 9994 )
00937 *
00938 *              Tests performed
00939 *
00940                WRITE( NOUT, FMT = 9992 )'orthogonal', '''',
00941      $            'transpose', ( '''', I = 1, 4 )
00942 *
00943             END IF
00944             NERRS = NERRS + 1
00945             IF( RESULT( J ).LT.10000.0 ) THEN
00946                WRITE( NOUT, FMT = 9989 )NPTKNT, MPLUSN, J, RESULT( J )
00947             ELSE
00948                WRITE( NOUT, FMT = 9988 )NPTKNT, MPLUSN, J, RESULT( J )
00949             END IF
00950          END IF
00951 *
00952   120 CONTINUE
00953 *
00954   130 CONTINUE
00955       GO TO 80
00956   140 CONTINUE
00957 *
00958   150 CONTINUE
00959 *
00960 *     Summary
00961 *
00962       CALL ALASVM( 'SGX', NOUT, NERRS, NTESTT, 0 )
00963 *
00964       WORK( 1 ) = MAXWRK
00965 *
00966       RETURN
00967 *
00968  9999 FORMAT( ' SDRGSX: ', A, ' returned INFO=', I6, '.', / 9X, 'N=',
00969      $      I6, ', JTYPE=', I6, ')' )
00970 *
00971  9998 FORMAT( ' SDRGSX: ', A, ' returned INFO=', I6, '.', / 9X, 'N=',
00972      $      I6, ', Input Example #', I2, ')' )
00973 *
00974  9997 FORMAT( ' SDRGSX: SGET53 returned INFO=', I1, ' for eigenvalue ',
00975      $      I6, '.', / 9X, 'N=', I6, ', JTYPE=', I6, ')' )
00976 *
00977  9996 FORMAT( ' SDRGSX: S not in Schur form at eigenvalue ', I6, '.',
00978      $      / 9X, 'N=', I6, ', JTYPE=', I6, ')' )
00979 *
00980  9995 FORMAT( / 1X, A3, ' -- Real Expert Generalized Schur form',
00981      $      ' problem driver' )
00982 *
00983  9994 FORMAT( 'Input Example' )
00984 *
00985  9993 FORMAT( ' Matrix types: ', /
00986      $      '  1:  A is a block diagonal matrix of Jordan blocks ',
00987      $      'and B is the identity ', / '      matrix, ',
00988      $      / '  2:  A and B are upper triangular matrices, ',
00989      $      / '  3:  A and B are as type 2, but each second diagonal ',
00990      $      'block in A_11 and ', /
00991      $      '      each third diaongal block in A_22 are 2x2 blocks,',
00992      $      / '  4:  A and B are block diagonal matrices, ',
00993      $      / '  5:  (A,B) has potentially close or common ',
00994      $      'eigenvalues.', / )
00995 *
00996  9992 FORMAT( / ' Tests performed:  (S is Schur, T is triangular, ',
00997      $      'Q and Z are ', A, ',', / 19X,
00998      $      ' a is alpha, b is beta, and ', A, ' means ', A, '.)',
00999      $      / '  1 = | A - Q S Z', A,
01000      $      ' | / ( |A| n ulp )      2 = | B - Q T Z', A,
01001      $      ' | / ( |B| n ulp )', / '  3 = | I - QQ', A,
01002      $      ' | / ( n ulp )             4 = | I - ZZ', A,
01003      $      ' | / ( n ulp )', / '  5 = 1/ULP  if A is not in ',
01004      $      'Schur form S', / '  6 = difference between (alpha,beta)',
01005      $      ' and diagonals of (S,T)', /
01006      $      '  7 = 1/ULP  if SDIM is not the correct number of ',
01007      $      'selected eigenvalues', /
01008      $      '  8 = 1/ULP  if DIFEST/DIFTRU > 10*THRESH or ',
01009      $      'DIFTRU/DIFEST > 10*THRESH',
01010      $      / '  9 = 1/ULP  if DIFEST <> 0 or DIFTRU > ULP*norm(A,B) ',
01011      $      'when reordering fails', /
01012      $      ' 10 = 1/ULP  if PLEST/PLTRU > THRESH or ',
01013      $      'PLTRU/PLEST > THRESH', /
01014      $      '    ( Test 10 is only for input examples )', / )
01015  9991 FORMAT( ' Matrix order=', I2, ', type=', I2, ', a=', E10.3,
01016      $      ', order(A_11)=', I2, ', result ', I2, ' is ', 0P, F8.2 )
01017  9990 FORMAT( ' Matrix order=', I2, ', type=', I2, ', a=', E10.3,
01018      $      ', order(A_11)=', I2, ', result ', I2, ' is ', 0P, E10.3 )
01019  9989 FORMAT( ' Input example #', I2, ', matrix order=', I4, ',',
01020      $      ' result ', I2, ' is', 0P, F8.2 )
01021  9988 FORMAT( ' Input example #', I2, ', matrix order=', I4, ',',
01022      $      ' result ', I2, ' is', 1P, E10.3 )
01023 *
01024 *     End of SDRGSX
01025 *
01026       END
 All Files Functions