![]() |
LAPACK
3.4.1
LAPACK: Linear Algebra PACKage
|
00001 *> \brief \b STRSYL 00002 * 00003 * =========== DOCUMENTATION =========== 00004 * 00005 * Online html documentation available at 00006 * http://www.netlib.org/lapack/explore-html/ 00007 * 00008 *> \htmlonly 00009 *> Download STRSYL + dependencies 00010 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/strsyl.f"> 00011 *> [TGZ]</a> 00012 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/strsyl.f"> 00013 *> [ZIP]</a> 00014 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/strsyl.f"> 00015 *> [TXT]</a> 00016 *> \endhtmlonly 00017 * 00018 * Definition: 00019 * =========== 00020 * 00021 * SUBROUTINE STRSYL( TRANA, TRANB, ISGN, M, N, A, LDA, B, LDB, C, 00022 * LDC, SCALE, INFO ) 00023 * 00024 * .. Scalar Arguments .. 00025 * CHARACTER TRANA, TRANB 00026 * INTEGER INFO, ISGN, LDA, LDB, LDC, M, N 00027 * REAL SCALE 00028 * .. 00029 * .. Array Arguments .. 00030 * REAL A( LDA, * ), B( LDB, * ), C( LDC, * ) 00031 * .. 00032 * 00033 * 00034 *> \par Purpose: 00035 * ============= 00036 *> 00037 *> \verbatim 00038 *> 00039 *> STRSYL solves the real Sylvester matrix equation: 00040 *> 00041 *> op(A)*X + X*op(B) = scale*C or 00042 *> op(A)*X - X*op(B) = scale*C, 00043 *> 00044 *> where op(A) = A or A**T, and A and B are both upper quasi- 00045 *> triangular. A is M-by-M and B is N-by-N; the right hand side C and 00046 *> the solution X are M-by-N; and scale is an output scale factor, set 00047 *> <= 1 to avoid overflow in X. 00048 *> 00049 *> A and B must be in Schur canonical form (as returned by SHSEQR), that 00050 *> is, block upper triangular with 1-by-1 and 2-by-2 diagonal blocks; 00051 *> each 2-by-2 diagonal block has its diagonal elements equal and its 00052 *> off-diagonal elements of opposite sign. 00053 *> \endverbatim 00054 * 00055 * Arguments: 00056 * ========== 00057 * 00058 *> \param[in] TRANA 00059 *> \verbatim 00060 *> TRANA is CHARACTER*1 00061 *> Specifies the option op(A): 00062 *> = 'N': op(A) = A (No transpose) 00063 *> = 'T': op(A) = A**T (Transpose) 00064 *> = 'C': op(A) = A**H (Conjugate transpose = Transpose) 00065 *> \endverbatim 00066 *> 00067 *> \param[in] TRANB 00068 *> \verbatim 00069 *> TRANB is CHARACTER*1 00070 *> Specifies the option op(B): 00071 *> = 'N': op(B) = B (No transpose) 00072 *> = 'T': op(B) = B**T (Transpose) 00073 *> = 'C': op(B) = B**H (Conjugate transpose = Transpose) 00074 *> \endverbatim 00075 *> 00076 *> \param[in] ISGN 00077 *> \verbatim 00078 *> ISGN is INTEGER 00079 *> Specifies the sign in the equation: 00080 *> = +1: solve op(A)*X + X*op(B) = scale*C 00081 *> = -1: solve op(A)*X - X*op(B) = scale*C 00082 *> \endverbatim 00083 *> 00084 *> \param[in] M 00085 *> \verbatim 00086 *> M is INTEGER 00087 *> The order of the matrix A, and the number of rows in the 00088 *> matrices X and C. M >= 0. 00089 *> \endverbatim 00090 *> 00091 *> \param[in] N 00092 *> \verbatim 00093 *> N is INTEGER 00094 *> The order of the matrix B, and the number of columns in the 00095 *> matrices X and C. N >= 0. 00096 *> \endverbatim 00097 *> 00098 *> \param[in] A 00099 *> \verbatim 00100 *> A is REAL array, dimension (LDA,M) 00101 *> The upper quasi-triangular matrix A, in Schur canonical form. 00102 *> \endverbatim 00103 *> 00104 *> \param[in] LDA 00105 *> \verbatim 00106 *> LDA is INTEGER 00107 *> The leading dimension of the array A. LDA >= max(1,M). 00108 *> \endverbatim 00109 *> 00110 *> \param[in] B 00111 *> \verbatim 00112 *> B is REAL array, dimension (LDB,N) 00113 *> The upper quasi-triangular matrix B, in Schur canonical form. 00114 *> \endverbatim 00115 *> 00116 *> \param[in] LDB 00117 *> \verbatim 00118 *> LDB is INTEGER 00119 *> The leading dimension of the array B. LDB >= max(1,N). 00120 *> \endverbatim 00121 *> 00122 *> \param[in,out] C 00123 *> \verbatim 00124 *> C is REAL array, dimension (LDC,N) 00125 *> On entry, the M-by-N right hand side matrix C. 00126 *> On exit, C is overwritten by the solution matrix X. 00127 *> \endverbatim 00128 *> 00129 *> \param[in] LDC 00130 *> \verbatim 00131 *> LDC is INTEGER 00132 *> The leading dimension of the array C. LDC >= max(1,M) 00133 *> \endverbatim 00134 *> 00135 *> \param[out] SCALE 00136 *> \verbatim 00137 *> SCALE is REAL 00138 *> The scale factor, scale, set <= 1 to avoid overflow in X. 00139 *> \endverbatim 00140 *> 00141 *> \param[out] INFO 00142 *> \verbatim 00143 *> INFO is INTEGER 00144 *> = 0: successful exit 00145 *> < 0: if INFO = -i, the i-th argument had an illegal value 00146 *> = 1: A and B have common or very close eigenvalues; perturbed 00147 *> values were used to solve the equation (but the matrices 00148 *> A and B are unchanged). 00149 *> \endverbatim 00150 * 00151 * Authors: 00152 * ======== 00153 * 00154 *> \author Univ. of Tennessee 00155 *> \author Univ. of California Berkeley 00156 *> \author Univ. of Colorado Denver 00157 *> \author NAG Ltd. 00158 * 00159 *> \date November 2011 00160 * 00161 *> \ingroup realSYcomputational 00162 * 00163 * ===================================================================== 00164 SUBROUTINE STRSYL( TRANA, TRANB, ISGN, M, N, A, LDA, B, LDB, C, 00165 $ LDC, SCALE, INFO ) 00166 * 00167 * -- LAPACK computational routine (version 3.4.0) -- 00168 * -- LAPACK is a software package provided by Univ. of Tennessee, -- 00169 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 00170 * November 2011 00171 * 00172 * .. Scalar Arguments .. 00173 CHARACTER TRANA, TRANB 00174 INTEGER INFO, ISGN, LDA, LDB, LDC, M, N 00175 REAL SCALE 00176 * .. 00177 * .. Array Arguments .. 00178 REAL A( LDA, * ), B( LDB, * ), C( LDC, * ) 00179 * .. 00180 * 00181 * ===================================================================== 00182 * 00183 * .. Parameters .. 00184 REAL ZERO, ONE 00185 PARAMETER ( ZERO = 0.0E+0, ONE = 1.0E+0 ) 00186 * .. 00187 * .. Local Scalars .. 00188 LOGICAL NOTRNA, NOTRNB 00189 INTEGER IERR, J, K, K1, K2, KNEXT, L, L1, L2, LNEXT 00190 REAL A11, BIGNUM, DA11, DB, EPS, SCALOC, SGN, SMIN, 00191 $ SMLNUM, SUML, SUMR, XNORM 00192 * .. 00193 * .. Local Arrays .. 00194 REAL DUM( 1 ), VEC( 2, 2 ), X( 2, 2 ) 00195 * .. 00196 * .. External Functions .. 00197 LOGICAL LSAME 00198 REAL SDOT, SLAMCH, SLANGE 00199 EXTERNAL LSAME, SDOT, SLAMCH, SLANGE 00200 * .. 00201 * .. External Subroutines .. 00202 EXTERNAL SLABAD, SLALN2, SLASY2, SSCAL, XERBLA 00203 * .. 00204 * .. Intrinsic Functions .. 00205 INTRINSIC ABS, MAX, MIN, REAL 00206 * .. 00207 * .. Executable Statements .. 00208 * 00209 * Decode and Test input parameters 00210 * 00211 NOTRNA = LSAME( TRANA, 'N' ) 00212 NOTRNB = LSAME( TRANB, 'N' ) 00213 * 00214 INFO = 0 00215 IF( .NOT.NOTRNA .AND. .NOT.LSAME( TRANA, 'T' ) .AND. .NOT. 00216 $ LSAME( TRANA, 'C' ) ) THEN 00217 INFO = -1 00218 ELSE IF( .NOT.NOTRNB .AND. .NOT.LSAME( TRANB, 'T' ) .AND. .NOT. 00219 $ LSAME( TRANB, 'C' ) ) THEN 00220 INFO = -2 00221 ELSE IF( ISGN.NE.1 .AND. ISGN.NE.-1 ) THEN 00222 INFO = -3 00223 ELSE IF( M.LT.0 ) THEN 00224 INFO = -4 00225 ELSE IF( N.LT.0 ) THEN 00226 INFO = -5 00227 ELSE IF( LDA.LT.MAX( 1, M ) ) THEN 00228 INFO = -7 00229 ELSE IF( LDB.LT.MAX( 1, N ) ) THEN 00230 INFO = -9 00231 ELSE IF( LDC.LT.MAX( 1, M ) ) THEN 00232 INFO = -11 00233 END IF 00234 IF( INFO.NE.0 ) THEN 00235 CALL XERBLA( 'STRSYL', -INFO ) 00236 RETURN 00237 END IF 00238 * 00239 * Quick return if possible 00240 * 00241 SCALE = ONE 00242 IF( M.EQ.0 .OR. N.EQ.0 ) 00243 $ RETURN 00244 * 00245 * Set constants to control overflow 00246 * 00247 EPS = SLAMCH( 'P' ) 00248 SMLNUM = SLAMCH( 'S' ) 00249 BIGNUM = ONE / SMLNUM 00250 CALL SLABAD( SMLNUM, BIGNUM ) 00251 SMLNUM = SMLNUM*REAL( M*N ) / EPS 00252 BIGNUM = ONE / SMLNUM 00253 * 00254 SMIN = MAX( SMLNUM, EPS*SLANGE( 'M', M, M, A, LDA, DUM ), 00255 $ EPS*SLANGE( 'M', N, N, B, LDB, DUM ) ) 00256 * 00257 SGN = ISGN 00258 * 00259 IF( NOTRNA .AND. NOTRNB ) THEN 00260 * 00261 * Solve A*X + ISGN*X*B = scale*C. 00262 * 00263 * The (K,L)th block of X is determined starting from 00264 * bottom-left corner column by column by 00265 * 00266 * A(K,K)*X(K,L) + ISGN*X(K,L)*B(L,L) = C(K,L) - R(K,L) 00267 * 00268 * Where 00269 * M L-1 00270 * R(K,L) = SUM [A(K,I)*X(I,L)] + ISGN*SUM [X(K,J)*B(J,L)]. 00271 * I=K+1 J=1 00272 * 00273 * Start column loop (index = L) 00274 * L1 (L2) : column index of the first (first) row of X(K,L). 00275 * 00276 LNEXT = 1 00277 DO 70 L = 1, N 00278 IF( L.LT.LNEXT ) 00279 $ GO TO 70 00280 IF( L.EQ.N ) THEN 00281 L1 = L 00282 L2 = L 00283 ELSE 00284 IF( B( L+1, L ).NE.ZERO ) THEN 00285 L1 = L 00286 L2 = L + 1 00287 LNEXT = L + 2 00288 ELSE 00289 L1 = L 00290 L2 = L 00291 LNEXT = L + 1 00292 END IF 00293 END IF 00294 * 00295 * Start row loop (index = K) 00296 * K1 (K2): row index of the first (last) row of X(K,L). 00297 * 00298 KNEXT = M 00299 DO 60 K = M, 1, -1 00300 IF( K.GT.KNEXT ) 00301 $ GO TO 60 00302 IF( K.EQ.1 ) THEN 00303 K1 = K 00304 K2 = K 00305 ELSE 00306 IF( A( K, K-1 ).NE.ZERO ) THEN 00307 K1 = K - 1 00308 K2 = K 00309 KNEXT = K - 2 00310 ELSE 00311 K1 = K 00312 K2 = K 00313 KNEXT = K - 1 00314 END IF 00315 END IF 00316 * 00317 IF( L1.EQ.L2 .AND. K1.EQ.K2 ) THEN 00318 SUML = SDOT( M-K1, A( K1, MIN( K1+1, M ) ), LDA, 00319 $ C( MIN( K1+1, M ), L1 ), 1 ) 00320 SUMR = SDOT( L1-1, C( K1, 1 ), LDC, B( 1, L1 ), 1 ) 00321 VEC( 1, 1 ) = C( K1, L1 ) - ( SUML+SGN*SUMR ) 00322 SCALOC = ONE 00323 * 00324 A11 = A( K1, K1 ) + SGN*B( L1, L1 ) 00325 DA11 = ABS( A11 ) 00326 IF( DA11.LE.SMIN ) THEN 00327 A11 = SMIN 00328 DA11 = SMIN 00329 INFO = 1 00330 END IF 00331 DB = ABS( VEC( 1, 1 ) ) 00332 IF( DA11.LT.ONE .AND. DB.GT.ONE ) THEN 00333 IF( DB.GT.BIGNUM*DA11 ) 00334 $ SCALOC = ONE / DB 00335 END IF 00336 X( 1, 1 ) = ( VEC( 1, 1 )*SCALOC ) / A11 00337 * 00338 IF( SCALOC.NE.ONE ) THEN 00339 DO 10 J = 1, N 00340 CALL SSCAL( M, SCALOC, C( 1, J ), 1 ) 00341 10 CONTINUE 00342 SCALE = SCALE*SCALOC 00343 END IF 00344 C( K1, L1 ) = X( 1, 1 ) 00345 * 00346 ELSE IF( L1.EQ.L2 .AND. K1.NE.K2 ) THEN 00347 * 00348 SUML = SDOT( M-K2, A( K1, MIN( K2+1, M ) ), LDA, 00349 $ C( MIN( K2+1, M ), L1 ), 1 ) 00350 SUMR = SDOT( L1-1, C( K1, 1 ), LDC, B( 1, L1 ), 1 ) 00351 VEC( 1, 1 ) = C( K1, L1 ) - ( SUML+SGN*SUMR ) 00352 * 00353 SUML = SDOT( M-K2, A( K2, MIN( K2+1, M ) ), LDA, 00354 $ C( MIN( K2+1, M ), L1 ), 1 ) 00355 SUMR = SDOT( L1-1, C( K2, 1 ), LDC, B( 1, L1 ), 1 ) 00356 VEC( 2, 1 ) = C( K2, L1 ) - ( SUML+SGN*SUMR ) 00357 * 00358 CALL SLALN2( .FALSE., 2, 1, SMIN, ONE, A( K1, K1 ), 00359 $ LDA, ONE, ONE, VEC, 2, -SGN*B( L1, L1 ), 00360 $ ZERO, X, 2, SCALOC, XNORM, IERR ) 00361 IF( IERR.NE.0 ) 00362 $ INFO = 1 00363 * 00364 IF( SCALOC.NE.ONE ) THEN 00365 DO 20 J = 1, N 00366 CALL SSCAL( M, SCALOC, C( 1, J ), 1 ) 00367 20 CONTINUE 00368 SCALE = SCALE*SCALOC 00369 END IF 00370 C( K1, L1 ) = X( 1, 1 ) 00371 C( K2, L1 ) = X( 2, 1 ) 00372 * 00373 ELSE IF( L1.NE.L2 .AND. K1.EQ.K2 ) THEN 00374 * 00375 SUML = SDOT( M-K1, A( K1, MIN( K1+1, M ) ), LDA, 00376 $ C( MIN( K1+1, M ), L1 ), 1 ) 00377 SUMR = SDOT( L1-1, C( K1, 1 ), LDC, B( 1, L1 ), 1 ) 00378 VEC( 1, 1 ) = SGN*( C( K1, L1 )-( SUML+SGN*SUMR ) ) 00379 * 00380 SUML = SDOT( M-K1, A( K1, MIN( K1+1, M ) ), LDA, 00381 $ C( MIN( K1+1, M ), L2 ), 1 ) 00382 SUMR = SDOT( L1-1, C( K1, 1 ), LDC, B( 1, L2 ), 1 ) 00383 VEC( 2, 1 ) = SGN*( C( K1, L2 )-( SUML+SGN*SUMR ) ) 00384 * 00385 CALL SLALN2( .TRUE., 2, 1, SMIN, ONE, B( L1, L1 ), 00386 $ LDB, ONE, ONE, VEC, 2, -SGN*A( K1, K1 ), 00387 $ ZERO, X, 2, SCALOC, XNORM, IERR ) 00388 IF( IERR.NE.0 ) 00389 $ INFO = 1 00390 * 00391 IF( SCALOC.NE.ONE ) THEN 00392 DO 40 J = 1, N 00393 CALL SSCAL( M, SCALOC, C( 1, J ), 1 ) 00394 40 CONTINUE 00395 SCALE = SCALE*SCALOC 00396 END IF 00397 C( K1, L1 ) = X( 1, 1 ) 00398 C( K1, L2 ) = X( 2, 1 ) 00399 * 00400 ELSE IF( L1.NE.L2 .AND. K1.NE.K2 ) THEN 00401 * 00402 SUML = SDOT( M-K2, A( K1, MIN( K2+1, M ) ), LDA, 00403 $ C( MIN( K2+1, M ), L1 ), 1 ) 00404 SUMR = SDOT( L1-1, C( K1, 1 ), LDC, B( 1, L1 ), 1 ) 00405 VEC( 1, 1 ) = C( K1, L1 ) - ( SUML+SGN*SUMR ) 00406 * 00407 SUML = SDOT( M-K2, A( K1, MIN( K2+1, M ) ), LDA, 00408 $ C( MIN( K2+1, M ), L2 ), 1 ) 00409 SUMR = SDOT( L1-1, C( K1, 1 ), LDC, B( 1, L2 ), 1 ) 00410 VEC( 1, 2 ) = C( K1, L2 ) - ( SUML+SGN*SUMR ) 00411 * 00412 SUML = SDOT( M-K2, A( K2, MIN( K2+1, M ) ), LDA, 00413 $ C( MIN( K2+1, M ), L1 ), 1 ) 00414 SUMR = SDOT( L1-1, C( K2, 1 ), LDC, B( 1, L1 ), 1 ) 00415 VEC( 2, 1 ) = C( K2, L1 ) - ( SUML+SGN*SUMR ) 00416 * 00417 SUML = SDOT( M-K2, A( K2, MIN( K2+1, M ) ), LDA, 00418 $ C( MIN( K2+1, M ), L2 ), 1 ) 00419 SUMR = SDOT( L1-1, C( K2, 1 ), LDC, B( 1, L2 ), 1 ) 00420 VEC( 2, 2 ) = C( K2, L2 ) - ( SUML+SGN*SUMR ) 00421 * 00422 CALL SLASY2( .FALSE., .FALSE., ISGN, 2, 2, 00423 $ A( K1, K1 ), LDA, B( L1, L1 ), LDB, VEC, 00424 $ 2, SCALOC, X, 2, XNORM, IERR ) 00425 IF( IERR.NE.0 ) 00426 $ INFO = 1 00427 * 00428 IF( SCALOC.NE.ONE ) THEN 00429 DO 50 J = 1, N 00430 CALL SSCAL( M, SCALOC, C( 1, J ), 1 ) 00431 50 CONTINUE 00432 SCALE = SCALE*SCALOC 00433 END IF 00434 C( K1, L1 ) = X( 1, 1 ) 00435 C( K1, L2 ) = X( 1, 2 ) 00436 C( K2, L1 ) = X( 2, 1 ) 00437 C( K2, L2 ) = X( 2, 2 ) 00438 END IF 00439 * 00440 60 CONTINUE 00441 * 00442 70 CONTINUE 00443 * 00444 ELSE IF( .NOT.NOTRNA .AND. NOTRNB ) THEN 00445 * 00446 * Solve A**T *X + ISGN*X*B = scale*C. 00447 * 00448 * The (K,L)th block of X is determined starting from 00449 * upper-left corner column by column by 00450 * 00451 * A(K,K)**T*X(K,L) + ISGN*X(K,L)*B(L,L) = C(K,L) - R(K,L) 00452 * 00453 * Where 00454 * K-1 L-1 00455 * R(K,L) = SUM [A(I,K)**T*X(I,L)] +ISGN*SUM [X(K,J)*B(J,L)] 00456 * I=1 J=1 00457 * 00458 * Start column loop (index = L) 00459 * L1 (L2): column index of the first (last) row of X(K,L) 00460 * 00461 LNEXT = 1 00462 DO 130 L = 1, N 00463 IF( L.LT.LNEXT ) 00464 $ GO TO 130 00465 IF( L.EQ.N ) THEN 00466 L1 = L 00467 L2 = L 00468 ELSE 00469 IF( B( L+1, L ).NE.ZERO ) THEN 00470 L1 = L 00471 L2 = L + 1 00472 LNEXT = L + 2 00473 ELSE 00474 L1 = L 00475 L2 = L 00476 LNEXT = L + 1 00477 END IF 00478 END IF 00479 * 00480 * Start row loop (index = K) 00481 * K1 (K2): row index of the first (last) row of X(K,L) 00482 * 00483 KNEXT = 1 00484 DO 120 K = 1, M 00485 IF( K.LT.KNEXT ) 00486 $ GO TO 120 00487 IF( K.EQ.M ) THEN 00488 K1 = K 00489 K2 = K 00490 ELSE 00491 IF( A( K+1, K ).NE.ZERO ) THEN 00492 K1 = K 00493 K2 = K + 1 00494 KNEXT = K + 2 00495 ELSE 00496 K1 = K 00497 K2 = K 00498 KNEXT = K + 1 00499 END IF 00500 END IF 00501 * 00502 IF( L1.EQ.L2 .AND. K1.EQ.K2 ) THEN 00503 SUML = SDOT( K1-1, A( 1, K1 ), 1, C( 1, L1 ), 1 ) 00504 SUMR = SDOT( L1-1, C( K1, 1 ), LDC, B( 1, L1 ), 1 ) 00505 VEC( 1, 1 ) = C( K1, L1 ) - ( SUML+SGN*SUMR ) 00506 SCALOC = ONE 00507 * 00508 A11 = A( K1, K1 ) + SGN*B( L1, L1 ) 00509 DA11 = ABS( A11 ) 00510 IF( DA11.LE.SMIN ) THEN 00511 A11 = SMIN 00512 DA11 = SMIN 00513 INFO = 1 00514 END IF 00515 DB = ABS( VEC( 1, 1 ) ) 00516 IF( DA11.LT.ONE .AND. DB.GT.ONE ) THEN 00517 IF( DB.GT.BIGNUM*DA11 ) 00518 $ SCALOC = ONE / DB 00519 END IF 00520 X( 1, 1 ) = ( VEC( 1, 1 )*SCALOC ) / A11 00521 * 00522 IF( SCALOC.NE.ONE ) THEN 00523 DO 80 J = 1, N 00524 CALL SSCAL( M, SCALOC, C( 1, J ), 1 ) 00525 80 CONTINUE 00526 SCALE = SCALE*SCALOC 00527 END IF 00528 C( K1, L1 ) = X( 1, 1 ) 00529 * 00530 ELSE IF( L1.EQ.L2 .AND. K1.NE.K2 ) THEN 00531 * 00532 SUML = SDOT( K1-1, A( 1, K1 ), 1, C( 1, L1 ), 1 ) 00533 SUMR = SDOT( L1-1, C( K1, 1 ), LDC, B( 1, L1 ), 1 ) 00534 VEC( 1, 1 ) = C( K1, L1 ) - ( SUML+SGN*SUMR ) 00535 * 00536 SUML = SDOT( K1-1, A( 1, K2 ), 1, C( 1, L1 ), 1 ) 00537 SUMR = SDOT( L1-1, C( K2, 1 ), LDC, B( 1, L1 ), 1 ) 00538 VEC( 2, 1 ) = C( K2, L1 ) - ( SUML+SGN*SUMR ) 00539 * 00540 CALL SLALN2( .TRUE., 2, 1, SMIN, ONE, A( K1, K1 ), 00541 $ LDA, ONE, ONE, VEC, 2, -SGN*B( L1, L1 ), 00542 $ ZERO, X, 2, SCALOC, XNORM, IERR ) 00543 IF( IERR.NE.0 ) 00544 $ INFO = 1 00545 * 00546 IF( SCALOC.NE.ONE ) THEN 00547 DO 90 J = 1, N 00548 CALL SSCAL( M, SCALOC, C( 1, J ), 1 ) 00549 90 CONTINUE 00550 SCALE = SCALE*SCALOC 00551 END IF 00552 C( K1, L1 ) = X( 1, 1 ) 00553 C( K2, L1 ) = X( 2, 1 ) 00554 * 00555 ELSE IF( L1.NE.L2 .AND. K1.EQ.K2 ) THEN 00556 * 00557 SUML = SDOT( K1-1, A( 1, K1 ), 1, C( 1, L1 ), 1 ) 00558 SUMR = SDOT( L1-1, C( K1, 1 ), LDC, B( 1, L1 ), 1 ) 00559 VEC( 1, 1 ) = SGN*( C( K1, L1 )-( SUML+SGN*SUMR ) ) 00560 * 00561 SUML = SDOT( K1-1, A( 1, K1 ), 1, C( 1, L2 ), 1 ) 00562 SUMR = SDOT( L1-1, C( K1, 1 ), LDC, B( 1, L2 ), 1 ) 00563 VEC( 2, 1 ) = SGN*( C( K1, L2 )-( SUML+SGN*SUMR ) ) 00564 * 00565 CALL SLALN2( .TRUE., 2, 1, SMIN, ONE, B( L1, L1 ), 00566 $ LDB, ONE, ONE, VEC, 2, -SGN*A( K1, K1 ), 00567 $ ZERO, X, 2, SCALOC, XNORM, IERR ) 00568 IF( IERR.NE.0 ) 00569 $ INFO = 1 00570 * 00571 IF( SCALOC.NE.ONE ) THEN 00572 DO 100 J = 1, N 00573 CALL SSCAL( M, SCALOC, C( 1, J ), 1 ) 00574 100 CONTINUE 00575 SCALE = SCALE*SCALOC 00576 END IF 00577 C( K1, L1 ) = X( 1, 1 ) 00578 C( K1, L2 ) = X( 2, 1 ) 00579 * 00580 ELSE IF( L1.NE.L2 .AND. K1.NE.K2 ) THEN 00581 * 00582 SUML = SDOT( K1-1, A( 1, K1 ), 1, C( 1, L1 ), 1 ) 00583 SUMR = SDOT( L1-1, C( K1, 1 ), LDC, B( 1, L1 ), 1 ) 00584 VEC( 1, 1 ) = C( K1, L1 ) - ( SUML+SGN*SUMR ) 00585 * 00586 SUML = SDOT( K1-1, A( 1, K1 ), 1, C( 1, L2 ), 1 ) 00587 SUMR = SDOT( L1-1, C( K1, 1 ), LDC, B( 1, L2 ), 1 ) 00588 VEC( 1, 2 ) = C( K1, L2 ) - ( SUML+SGN*SUMR ) 00589 * 00590 SUML = SDOT( K1-1, A( 1, K2 ), 1, C( 1, L1 ), 1 ) 00591 SUMR = SDOT( L1-1, C( K2, 1 ), LDC, B( 1, L1 ), 1 ) 00592 VEC( 2, 1 ) = C( K2, L1 ) - ( SUML+SGN*SUMR ) 00593 * 00594 SUML = SDOT( K1-1, A( 1, K2 ), 1, C( 1, L2 ), 1 ) 00595 SUMR = SDOT( L1-1, C( K2, 1 ), LDC, B( 1, L2 ), 1 ) 00596 VEC( 2, 2 ) = C( K2, L2 ) - ( SUML+SGN*SUMR ) 00597 * 00598 CALL SLASY2( .TRUE., .FALSE., ISGN, 2, 2, A( K1, K1 ), 00599 $ LDA, B( L1, L1 ), LDB, VEC, 2, SCALOC, X, 00600 $ 2, XNORM, IERR ) 00601 IF( IERR.NE.0 ) 00602 $ INFO = 1 00603 * 00604 IF( SCALOC.NE.ONE ) THEN 00605 DO 110 J = 1, N 00606 CALL SSCAL( M, SCALOC, C( 1, J ), 1 ) 00607 110 CONTINUE 00608 SCALE = SCALE*SCALOC 00609 END IF 00610 C( K1, L1 ) = X( 1, 1 ) 00611 C( K1, L2 ) = X( 1, 2 ) 00612 C( K2, L1 ) = X( 2, 1 ) 00613 C( K2, L2 ) = X( 2, 2 ) 00614 END IF 00615 * 00616 120 CONTINUE 00617 130 CONTINUE 00618 * 00619 ELSE IF( .NOT.NOTRNA .AND. .NOT.NOTRNB ) THEN 00620 * 00621 * Solve A**T*X + ISGN*X*B**T = scale*C. 00622 * 00623 * The (K,L)th block of X is determined starting from 00624 * top-right corner column by column by 00625 * 00626 * A(K,K)**T*X(K,L) + ISGN*X(K,L)*B(L,L)**T = C(K,L) - R(K,L) 00627 * 00628 * Where 00629 * K-1 N 00630 * R(K,L) = SUM [A(I,K)**T*X(I,L)] + ISGN*SUM [X(K,J)*B(L,J)**T]. 00631 * I=1 J=L+1 00632 * 00633 * Start column loop (index = L) 00634 * L1 (L2): column index of the first (last) row of X(K,L) 00635 * 00636 LNEXT = N 00637 DO 190 L = N, 1, -1 00638 IF( L.GT.LNEXT ) 00639 $ GO TO 190 00640 IF( L.EQ.1 ) THEN 00641 L1 = L 00642 L2 = L 00643 ELSE 00644 IF( B( L, L-1 ).NE.ZERO ) THEN 00645 L1 = L - 1 00646 L2 = L 00647 LNEXT = L - 2 00648 ELSE 00649 L1 = L 00650 L2 = L 00651 LNEXT = L - 1 00652 END IF 00653 END IF 00654 * 00655 * Start row loop (index = K) 00656 * K1 (K2): row index of the first (last) row of X(K,L) 00657 * 00658 KNEXT = 1 00659 DO 180 K = 1, M 00660 IF( K.LT.KNEXT ) 00661 $ GO TO 180 00662 IF( K.EQ.M ) THEN 00663 K1 = K 00664 K2 = K 00665 ELSE 00666 IF( A( K+1, K ).NE.ZERO ) THEN 00667 K1 = K 00668 K2 = K + 1 00669 KNEXT = K + 2 00670 ELSE 00671 K1 = K 00672 K2 = K 00673 KNEXT = K + 1 00674 END IF 00675 END IF 00676 * 00677 IF( L1.EQ.L2 .AND. K1.EQ.K2 ) THEN 00678 SUML = SDOT( K1-1, A( 1, K1 ), 1, C( 1, L1 ), 1 ) 00679 SUMR = SDOT( N-L1, C( K1, MIN( L1+1, N ) ), LDC, 00680 $ B( L1, MIN( L1+1, N ) ), LDB ) 00681 VEC( 1, 1 ) = C( K1, L1 ) - ( SUML+SGN*SUMR ) 00682 SCALOC = ONE 00683 * 00684 A11 = A( K1, K1 ) + SGN*B( L1, L1 ) 00685 DA11 = ABS( A11 ) 00686 IF( DA11.LE.SMIN ) THEN 00687 A11 = SMIN 00688 DA11 = SMIN 00689 INFO = 1 00690 END IF 00691 DB = ABS( VEC( 1, 1 ) ) 00692 IF( DA11.LT.ONE .AND. DB.GT.ONE ) THEN 00693 IF( DB.GT.BIGNUM*DA11 ) 00694 $ SCALOC = ONE / DB 00695 END IF 00696 X( 1, 1 ) = ( VEC( 1, 1 )*SCALOC ) / A11 00697 * 00698 IF( SCALOC.NE.ONE ) THEN 00699 DO 140 J = 1, N 00700 CALL SSCAL( M, SCALOC, C( 1, J ), 1 ) 00701 140 CONTINUE 00702 SCALE = SCALE*SCALOC 00703 END IF 00704 C( K1, L1 ) = X( 1, 1 ) 00705 * 00706 ELSE IF( L1.EQ.L2 .AND. K1.NE.K2 ) THEN 00707 * 00708 SUML = SDOT( K1-1, A( 1, K1 ), 1, C( 1, L1 ), 1 ) 00709 SUMR = SDOT( N-L2, C( K1, MIN( L2+1, N ) ), LDC, 00710 $ B( L1, MIN( L2+1, N ) ), LDB ) 00711 VEC( 1, 1 ) = C( K1, L1 ) - ( SUML+SGN*SUMR ) 00712 * 00713 SUML = SDOT( K1-1, A( 1, K2 ), 1, C( 1, L1 ), 1 ) 00714 SUMR = SDOT( N-L2, C( K2, MIN( L2+1, N ) ), LDC, 00715 $ B( L1, MIN( L2+1, N ) ), LDB ) 00716 VEC( 2, 1 ) = C( K2, L1 ) - ( SUML+SGN*SUMR ) 00717 * 00718 CALL SLALN2( .TRUE., 2, 1, SMIN, ONE, A( K1, K1 ), 00719 $ LDA, ONE, ONE, VEC, 2, -SGN*B( L1, L1 ), 00720 $ ZERO, X, 2, SCALOC, XNORM, IERR ) 00721 IF( IERR.NE.0 ) 00722 $ INFO = 1 00723 * 00724 IF( SCALOC.NE.ONE ) THEN 00725 DO 150 J = 1, N 00726 CALL SSCAL( M, SCALOC, C( 1, J ), 1 ) 00727 150 CONTINUE 00728 SCALE = SCALE*SCALOC 00729 END IF 00730 C( K1, L1 ) = X( 1, 1 ) 00731 C( K2, L1 ) = X( 2, 1 ) 00732 * 00733 ELSE IF( L1.NE.L2 .AND. K1.EQ.K2 ) THEN 00734 * 00735 SUML = SDOT( K1-1, A( 1, K1 ), 1, C( 1, L1 ), 1 ) 00736 SUMR = SDOT( N-L2, C( K1, MIN( L2+1, N ) ), LDC, 00737 $ B( L1, MIN( L2+1, N ) ), LDB ) 00738 VEC( 1, 1 ) = SGN*( C( K1, L1 )-( SUML+SGN*SUMR ) ) 00739 * 00740 SUML = SDOT( K1-1, A( 1, K1 ), 1, C( 1, L2 ), 1 ) 00741 SUMR = SDOT( N-L2, C( K1, MIN( L2+1, N ) ), LDC, 00742 $ B( L2, MIN( L2+1, N ) ), LDB ) 00743 VEC( 2, 1 ) = SGN*( C( K1, L2 )-( SUML+SGN*SUMR ) ) 00744 * 00745 CALL SLALN2( .FALSE., 2, 1, SMIN, ONE, B( L1, L1 ), 00746 $ LDB, ONE, ONE, VEC, 2, -SGN*A( K1, K1 ), 00747 $ ZERO, X, 2, SCALOC, XNORM, IERR ) 00748 IF( IERR.NE.0 ) 00749 $ INFO = 1 00750 * 00751 IF( SCALOC.NE.ONE ) THEN 00752 DO 160 J = 1, N 00753 CALL SSCAL( M, SCALOC, C( 1, J ), 1 ) 00754 160 CONTINUE 00755 SCALE = SCALE*SCALOC 00756 END IF 00757 C( K1, L1 ) = X( 1, 1 ) 00758 C( K1, L2 ) = X( 2, 1 ) 00759 * 00760 ELSE IF( L1.NE.L2 .AND. K1.NE.K2 ) THEN 00761 * 00762 SUML = SDOT( K1-1, A( 1, K1 ), 1, C( 1, L1 ), 1 ) 00763 SUMR = SDOT( N-L2, C( K1, MIN( L2+1, N ) ), LDC, 00764 $ B( L1, MIN( L2+1, N ) ), LDB ) 00765 VEC( 1, 1 ) = C( K1, L1 ) - ( SUML+SGN*SUMR ) 00766 * 00767 SUML = SDOT( K1-1, A( 1, K1 ), 1, C( 1, L2 ), 1 ) 00768 SUMR = SDOT( N-L2, C( K1, MIN( L2+1, N ) ), LDC, 00769 $ B( L2, MIN( L2+1, N ) ), LDB ) 00770 VEC( 1, 2 ) = C( K1, L2 ) - ( SUML+SGN*SUMR ) 00771 * 00772 SUML = SDOT( K1-1, A( 1, K2 ), 1, C( 1, L1 ), 1 ) 00773 SUMR = SDOT( N-L2, C( K2, MIN( L2+1, N ) ), LDC, 00774 $ B( L1, MIN( L2+1, N ) ), LDB ) 00775 VEC( 2, 1 ) = C( K2, L1 ) - ( SUML+SGN*SUMR ) 00776 * 00777 SUML = SDOT( K1-1, A( 1, K2 ), 1, C( 1, L2 ), 1 ) 00778 SUMR = SDOT( N-L2, C( K2, MIN( L2+1, N ) ), LDC, 00779 $ B( L2, MIN(L2+1, N ) ), LDB ) 00780 VEC( 2, 2 ) = C( K2, L2 ) - ( SUML+SGN*SUMR ) 00781 * 00782 CALL SLASY2( .TRUE., .TRUE., ISGN, 2, 2, A( K1, K1 ), 00783 $ LDA, B( L1, L1 ), LDB, VEC, 2, SCALOC, X, 00784 $ 2, XNORM, IERR ) 00785 IF( IERR.NE.0 ) 00786 $ INFO = 1 00787 * 00788 IF( SCALOC.NE.ONE ) THEN 00789 DO 170 J = 1, N 00790 CALL SSCAL( M, SCALOC, C( 1, J ), 1 ) 00791 170 CONTINUE 00792 SCALE = SCALE*SCALOC 00793 END IF 00794 C( K1, L1 ) = X( 1, 1 ) 00795 C( K1, L2 ) = X( 1, 2 ) 00796 C( K2, L1 ) = X( 2, 1 ) 00797 C( K2, L2 ) = X( 2, 2 ) 00798 END IF 00799 * 00800 180 CONTINUE 00801 190 CONTINUE 00802 * 00803 ELSE IF( NOTRNA .AND. .NOT.NOTRNB ) THEN 00804 * 00805 * Solve A*X + ISGN*X*B**T = scale*C. 00806 * 00807 * The (K,L)th block of X is determined starting from 00808 * bottom-right corner column by column by 00809 * 00810 * A(K,K)*X(K,L) + ISGN*X(K,L)*B(L,L)**T = C(K,L) - R(K,L) 00811 * 00812 * Where 00813 * M N 00814 * R(K,L) = SUM [A(K,I)*X(I,L)] + ISGN*SUM [X(K,J)*B(L,J)**T]. 00815 * I=K+1 J=L+1 00816 * 00817 * Start column loop (index = L) 00818 * L1 (L2): column index of the first (last) row of X(K,L) 00819 * 00820 LNEXT = N 00821 DO 250 L = N, 1, -1 00822 IF( L.GT.LNEXT ) 00823 $ GO TO 250 00824 IF( L.EQ.1 ) THEN 00825 L1 = L 00826 L2 = L 00827 ELSE 00828 IF( B( L, L-1 ).NE.ZERO ) THEN 00829 L1 = L - 1 00830 L2 = L 00831 LNEXT = L - 2 00832 ELSE 00833 L1 = L 00834 L2 = L 00835 LNEXT = L - 1 00836 END IF 00837 END IF 00838 * 00839 * Start row loop (index = K) 00840 * K1 (K2): row index of the first (last) row of X(K,L) 00841 * 00842 KNEXT = M 00843 DO 240 K = M, 1, -1 00844 IF( K.GT.KNEXT ) 00845 $ GO TO 240 00846 IF( K.EQ.1 ) THEN 00847 K1 = K 00848 K2 = K 00849 ELSE 00850 IF( A( K, K-1 ).NE.ZERO ) THEN 00851 K1 = K - 1 00852 K2 = K 00853 KNEXT = K - 2 00854 ELSE 00855 K1 = K 00856 K2 = K 00857 KNEXT = K - 1 00858 END IF 00859 END IF 00860 * 00861 IF( L1.EQ.L2 .AND. K1.EQ.K2 ) THEN 00862 SUML = SDOT( M-K1, A( K1, MIN(K1+1, M ) ), LDA, 00863 $ C( MIN( K1+1, M ), L1 ), 1 ) 00864 SUMR = SDOT( N-L1, C( K1, MIN( L1+1, N ) ), LDC, 00865 $ B( L1, MIN( L1+1, N ) ), LDB ) 00866 VEC( 1, 1 ) = C( K1, L1 ) - ( SUML+SGN*SUMR ) 00867 SCALOC = ONE 00868 * 00869 A11 = A( K1, K1 ) + SGN*B( L1, L1 ) 00870 DA11 = ABS( A11 ) 00871 IF( DA11.LE.SMIN ) THEN 00872 A11 = SMIN 00873 DA11 = SMIN 00874 INFO = 1 00875 END IF 00876 DB = ABS( VEC( 1, 1 ) ) 00877 IF( DA11.LT.ONE .AND. DB.GT.ONE ) THEN 00878 IF( DB.GT.BIGNUM*DA11 ) 00879 $ SCALOC = ONE / DB 00880 END IF 00881 X( 1, 1 ) = ( VEC( 1, 1 )*SCALOC ) / A11 00882 * 00883 IF( SCALOC.NE.ONE ) THEN 00884 DO 200 J = 1, N 00885 CALL SSCAL( M, SCALOC, C( 1, J ), 1 ) 00886 200 CONTINUE 00887 SCALE = SCALE*SCALOC 00888 END IF 00889 C( K1, L1 ) = X( 1, 1 ) 00890 * 00891 ELSE IF( L1.EQ.L2 .AND. K1.NE.K2 ) THEN 00892 * 00893 SUML = SDOT( M-K2, A( K1, MIN( K2+1, M ) ), LDA, 00894 $ C( MIN( K2+1, M ), L1 ), 1 ) 00895 SUMR = SDOT( N-L2, C( K1, MIN( L2+1, N ) ), LDC, 00896 $ B( L1, MIN( L2+1, N ) ), LDB ) 00897 VEC( 1, 1 ) = C( K1, L1 ) - ( SUML+SGN*SUMR ) 00898 * 00899 SUML = SDOT( M-K2, A( K2, MIN( K2+1, M ) ), LDA, 00900 $ C( MIN( K2+1, M ), L1 ), 1 ) 00901 SUMR = SDOT( N-L2, C( K2, MIN( L2+1, N ) ), LDC, 00902 $ B( L1, MIN( L2+1, N ) ), LDB ) 00903 VEC( 2, 1 ) = C( K2, L1 ) - ( SUML+SGN*SUMR ) 00904 * 00905 CALL SLALN2( .FALSE., 2, 1, SMIN, ONE, A( K1, K1 ), 00906 $ LDA, ONE, ONE, VEC, 2, -SGN*B( L1, L1 ), 00907 $ ZERO, X, 2, SCALOC, XNORM, IERR ) 00908 IF( IERR.NE.0 ) 00909 $ INFO = 1 00910 * 00911 IF( SCALOC.NE.ONE ) THEN 00912 DO 210 J = 1, N 00913 CALL SSCAL( M, SCALOC, C( 1, J ), 1 ) 00914 210 CONTINUE 00915 SCALE = SCALE*SCALOC 00916 END IF 00917 C( K1, L1 ) = X( 1, 1 ) 00918 C( K2, L1 ) = X( 2, 1 ) 00919 * 00920 ELSE IF( L1.NE.L2 .AND. K1.EQ.K2 ) THEN 00921 * 00922 SUML = SDOT( M-K1, A( K1, MIN( K1+1, M ) ), LDA, 00923 $ C( MIN( K1+1, M ), L1 ), 1 ) 00924 SUMR = SDOT( N-L2, C( K1, MIN( L2+1, N ) ), LDC, 00925 $ B( L1, MIN( L2+1, N ) ), LDB ) 00926 VEC( 1, 1 ) = SGN*( C( K1, L1 )-( SUML+SGN*SUMR ) ) 00927 * 00928 SUML = SDOT( M-K1, A( K1, MIN( K1+1, M ) ), LDA, 00929 $ C( MIN( K1+1, M ), L2 ), 1 ) 00930 SUMR = SDOT( N-L2, C( K1, MIN( L2+1, N ) ), LDC, 00931 $ B( L2, MIN( L2+1, N ) ), LDB ) 00932 VEC( 2, 1 ) = SGN*( C( K1, L2 )-( SUML+SGN*SUMR ) ) 00933 * 00934 CALL SLALN2( .FALSE., 2, 1, SMIN, ONE, B( L1, L1 ), 00935 $ LDB, ONE, ONE, VEC, 2, -SGN*A( K1, K1 ), 00936 $ ZERO, X, 2, SCALOC, XNORM, IERR ) 00937 IF( IERR.NE.0 ) 00938 $ INFO = 1 00939 * 00940 IF( SCALOC.NE.ONE ) THEN 00941 DO 220 J = 1, N 00942 CALL SSCAL( M, SCALOC, C( 1, J ), 1 ) 00943 220 CONTINUE 00944 SCALE = SCALE*SCALOC 00945 END IF 00946 C( K1, L1 ) = X( 1, 1 ) 00947 C( K1, L2 ) = X( 2, 1 ) 00948 * 00949 ELSE IF( L1.NE.L2 .AND. K1.NE.K2 ) THEN 00950 * 00951 SUML = SDOT( M-K2, A( K1, MIN( K2+1, M ) ), LDA, 00952 $ C( MIN( K2+1, M ), L1 ), 1 ) 00953 SUMR = SDOT( N-L2, C( K1, MIN( L2+1, N ) ), LDC, 00954 $ B( L1, MIN( L2+1, N ) ), LDB ) 00955 VEC( 1, 1 ) = C( K1, L1 ) - ( SUML+SGN*SUMR ) 00956 * 00957 SUML = SDOT( M-K2, A( K1, MIN( K2+1, M ) ), LDA, 00958 $ C( MIN( K2+1, M ), L2 ), 1 ) 00959 SUMR = SDOT( N-L2, C( K1, MIN( L2+1, N ) ), LDC, 00960 $ B( L2, MIN( L2+1, N ) ), LDB ) 00961 VEC( 1, 2 ) = C( K1, L2 ) - ( SUML+SGN*SUMR ) 00962 * 00963 SUML = SDOT( M-K2, A( K2, MIN( K2+1, M ) ), LDA, 00964 $ C( MIN( K2+1, M ), L1 ), 1 ) 00965 SUMR = SDOT( N-L2, C( K2, MIN( L2+1, N ) ), LDC, 00966 $ B( L1, MIN( L2+1, N ) ), LDB ) 00967 VEC( 2, 1 ) = C( K2, L1 ) - ( SUML+SGN*SUMR ) 00968 * 00969 SUML = SDOT( M-K2, A( K2, MIN( K2+1, M ) ), LDA, 00970 $ C( MIN( K2+1, M ), L2 ), 1 ) 00971 SUMR = SDOT( N-L2, C( K2, MIN( L2+1, N ) ), LDC, 00972 $ B( L2, MIN( L2+1, N ) ), LDB ) 00973 VEC( 2, 2 ) = C( K2, L2 ) - ( SUML+SGN*SUMR ) 00974 * 00975 CALL SLASY2( .FALSE., .TRUE., ISGN, 2, 2, A( K1, K1 ), 00976 $ LDA, B( L1, L1 ), LDB, VEC, 2, SCALOC, X, 00977 $ 2, XNORM, IERR ) 00978 IF( IERR.NE.0 ) 00979 $ INFO = 1 00980 * 00981 IF( SCALOC.NE.ONE ) THEN 00982 DO 230 J = 1, N 00983 CALL SSCAL( M, SCALOC, C( 1, J ), 1 ) 00984 230 CONTINUE 00985 SCALE = SCALE*SCALOC 00986 END IF 00987 C( K1, L1 ) = X( 1, 1 ) 00988 C( K1, L2 ) = X( 1, 2 ) 00989 C( K2, L1 ) = X( 2, 1 ) 00990 C( K2, L2 ) = X( 2, 2 ) 00991 END IF 00992 * 00993 240 CONTINUE 00994 250 CONTINUE 00995 * 00996 END IF 00997 * 00998 RETURN 00999 * 01000 * End of STRSYL 01001 * 01002 END