Annotation of rpl/lapack/lapack/dlaln2.f, revision 1.9
1.9 ! bertrand 1: *> \brief \b DLALN2
! 2: *
! 3: * =========== DOCUMENTATION ===========
! 4: *
! 5: * Online html documentation available at
! 6: * http://www.netlib.org/lapack/explore-html/
! 7: *
! 8: *> \htmlonly
! 9: *> Download DLALN2 + dependencies
! 10: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dlaln2.f">
! 11: *> [TGZ]</a>
! 12: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dlaln2.f">
! 13: *> [ZIP]</a>
! 14: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dlaln2.f">
! 15: *> [TXT]</a>
! 16: *> \endhtmlonly
! 17: *
! 18: * Definition:
! 19: * ===========
! 20: *
! 21: * SUBROUTINE DLALN2( LTRANS, NA, NW, SMIN, CA, A, LDA, D1, D2, B,
! 22: * LDB, WR, WI, X, LDX, SCALE, XNORM, INFO )
! 23: *
! 24: * .. Scalar Arguments ..
! 25: * LOGICAL LTRANS
! 26: * INTEGER INFO, LDA, LDB, LDX, NA, NW
! 27: * DOUBLE PRECISION CA, D1, D2, SCALE, SMIN, WI, WR, XNORM
! 28: * ..
! 29: * .. Array Arguments ..
! 30: * DOUBLE PRECISION A( LDA, * ), B( LDB, * ), X( LDX, * )
! 31: * ..
! 32: *
! 33: *
! 34: *> \par Purpose:
! 35: * =============
! 36: *>
! 37: *> \verbatim
! 38: *>
! 39: *> DLALN2 solves a system of the form (ca A - w D ) X = s B
! 40: *> or (ca A**T - w D) X = s B with possible scaling ("s") and
! 41: *> perturbation of A. (A**T means A-transpose.)
! 42: *>
! 43: *> A is an NA x NA real matrix, ca is a real scalar, D is an NA x NA
! 44: *> real diagonal matrix, w is a real or complex value, and X and B are
! 45: *> NA x 1 matrices -- real if w is real, complex if w is complex. NA
! 46: *> may be 1 or 2.
! 47: *>
! 48: *> If w is complex, X and B are represented as NA x 2 matrices,
! 49: *> the first column of each being the real part and the second
! 50: *> being the imaginary part.
! 51: *>
! 52: *> "s" is a scaling factor (.LE. 1), computed by DLALN2, which is
! 53: *> so chosen that X can be computed without overflow. X is further
! 54: *> scaled if necessary to assure that norm(ca A - w D)*norm(X) is less
! 55: *> than overflow.
! 56: *>
! 57: *> If both singular values of (ca A - w D) are less than SMIN,
! 58: *> SMIN*identity will be used instead of (ca A - w D). If only one
! 59: *> singular value is less than SMIN, one element of (ca A - w D) will be
! 60: *> perturbed enough to make the smallest singular value roughly SMIN.
! 61: *> If both singular values are at least SMIN, (ca A - w D) will not be
! 62: *> perturbed. In any case, the perturbation will be at most some small
! 63: *> multiple of max( SMIN, ulp*norm(ca A - w D) ). The singular values
! 64: *> are computed by infinity-norm approximations, and thus will only be
! 65: *> correct to a factor of 2 or so.
! 66: *>
! 67: *> Note: all input quantities are assumed to be smaller than overflow
! 68: *> by a reasonable factor. (See BIGNUM.)
! 69: *> \endverbatim
! 70: *
! 71: * Arguments:
! 72: * ==========
! 73: *
! 74: *> \param[in] LTRANS
! 75: *> \verbatim
! 76: *> LTRANS is LOGICAL
! 77: *> =.TRUE.: A-transpose will be used.
! 78: *> =.FALSE.: A will be used (not transposed.)
! 79: *> \endverbatim
! 80: *>
! 81: *> \param[in] NA
! 82: *> \verbatim
! 83: *> NA is INTEGER
! 84: *> The size of the matrix A. It may (only) be 1 or 2.
! 85: *> \endverbatim
! 86: *>
! 87: *> \param[in] NW
! 88: *> \verbatim
! 89: *> NW is INTEGER
! 90: *> 1 if "w" is real, 2 if "w" is complex. It may only be 1
! 91: *> or 2.
! 92: *> \endverbatim
! 93: *>
! 94: *> \param[in] SMIN
! 95: *> \verbatim
! 96: *> SMIN is DOUBLE PRECISION
! 97: *> The desired lower bound on the singular values of A. This
! 98: *> should be a safe distance away from underflow or overflow,
! 99: *> say, between (underflow/machine precision) and (machine
! 100: *> precision * overflow ). (See BIGNUM and ULP.)
! 101: *> \endverbatim
! 102: *>
! 103: *> \param[in] CA
! 104: *> \verbatim
! 105: *> CA is DOUBLE PRECISION
! 106: *> The coefficient c, which A is multiplied by.
! 107: *> \endverbatim
! 108: *>
! 109: *> \param[in] A
! 110: *> \verbatim
! 111: *> A is DOUBLE PRECISION array, dimension (LDA,NA)
! 112: *> The NA x NA matrix A.
! 113: *> \endverbatim
! 114: *>
! 115: *> \param[in] LDA
! 116: *> \verbatim
! 117: *> LDA is INTEGER
! 118: *> The leading dimension of A. It must be at least NA.
! 119: *> \endverbatim
! 120: *>
! 121: *> \param[in] D1
! 122: *> \verbatim
! 123: *> D1 is DOUBLE PRECISION
! 124: *> The 1,1 element in the diagonal matrix D.
! 125: *> \endverbatim
! 126: *>
! 127: *> \param[in] D2
! 128: *> \verbatim
! 129: *> D2 is DOUBLE PRECISION
! 130: *> The 2,2 element in the diagonal matrix D. Not used if NW=1.
! 131: *> \endverbatim
! 132: *>
! 133: *> \param[in] B
! 134: *> \verbatim
! 135: *> B is DOUBLE PRECISION array, dimension (LDB,NW)
! 136: *> The NA x NW matrix B (right-hand side). If NW=2 ("w" is
! 137: *> complex), column 1 contains the real part of B and column 2
! 138: *> contains the imaginary part.
! 139: *> \endverbatim
! 140: *>
! 141: *> \param[in] LDB
! 142: *> \verbatim
! 143: *> LDB is INTEGER
! 144: *> The leading dimension of B. It must be at least NA.
! 145: *> \endverbatim
! 146: *>
! 147: *> \param[in] WR
! 148: *> \verbatim
! 149: *> WR is DOUBLE PRECISION
! 150: *> The real part of the scalar "w".
! 151: *> \endverbatim
! 152: *>
! 153: *> \param[in] WI
! 154: *> \verbatim
! 155: *> WI is DOUBLE PRECISION
! 156: *> The imaginary part of the scalar "w". Not used if NW=1.
! 157: *> \endverbatim
! 158: *>
! 159: *> \param[out] X
! 160: *> \verbatim
! 161: *> X is DOUBLE PRECISION array, dimension (LDX,NW)
! 162: *> The NA x NW matrix X (unknowns), as computed by DLALN2.
! 163: *> If NW=2 ("w" is complex), on exit, column 1 will contain
! 164: *> the real part of X and column 2 will contain the imaginary
! 165: *> part.
! 166: *> \endverbatim
! 167: *>
! 168: *> \param[in] LDX
! 169: *> \verbatim
! 170: *> LDX is INTEGER
! 171: *> The leading dimension of X. It must be at least NA.
! 172: *> \endverbatim
! 173: *>
! 174: *> \param[out] SCALE
! 175: *> \verbatim
! 176: *> SCALE is DOUBLE PRECISION
! 177: *> The scale factor that B must be multiplied by to insure
! 178: *> that overflow does not occur when computing X. Thus,
! 179: *> (ca A - w D) X will be SCALE*B, not B (ignoring
! 180: *> perturbations of A.) It will be at most 1.
! 181: *> \endverbatim
! 182: *>
! 183: *> \param[out] XNORM
! 184: *> \verbatim
! 185: *> XNORM is DOUBLE PRECISION
! 186: *> The infinity-norm of X, when X is regarded as an NA x NW
! 187: *> real matrix.
! 188: *> \endverbatim
! 189: *>
! 190: *> \param[out] INFO
! 191: *> \verbatim
! 192: *> INFO is INTEGER
! 193: *> An error flag. It will be set to zero if no error occurs,
! 194: *> a negative number if an argument is in error, or a positive
! 195: *> number if ca A - w D had to be perturbed.
! 196: *> The possible values are:
! 197: *> = 0: No error occurred, and (ca A - w D) did not have to be
! 198: *> perturbed.
! 199: *> = 1: (ca A - w D) had to be perturbed to make its smallest
! 200: *> (or only) singular value greater than SMIN.
! 201: *> NOTE: In the interests of speed, this routine does not
! 202: *> check the inputs for errors.
! 203: *> \endverbatim
! 204: *
! 205: * Authors:
! 206: * ========
! 207: *
! 208: *> \author Univ. of Tennessee
! 209: *> \author Univ. of California Berkeley
! 210: *> \author Univ. of Colorado Denver
! 211: *> \author NAG Ltd.
! 212: *
! 213: *> \date November 2011
! 214: *
! 215: *> \ingroup doubleOTHERauxiliary
! 216: *
! 217: * =====================================================================
1.1 bertrand 218: SUBROUTINE DLALN2( LTRANS, NA, NW, SMIN, CA, A, LDA, D1, D2, B,
219: $ LDB, WR, WI, X, LDX, SCALE, XNORM, INFO )
220: *
1.9 ! bertrand 221: * -- LAPACK auxiliary routine (version 3.4.0) --
1.1 bertrand 222: * -- LAPACK is a software package provided by Univ. of Tennessee, --
223: * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
1.9 ! bertrand 224: * November 2011
1.1 bertrand 225: *
226: * .. Scalar Arguments ..
227: LOGICAL LTRANS
228: INTEGER INFO, LDA, LDB, LDX, NA, NW
229: DOUBLE PRECISION CA, D1, D2, SCALE, SMIN, WI, WR, XNORM
230: * ..
231: * .. Array Arguments ..
232: DOUBLE PRECISION A( LDA, * ), B( LDB, * ), X( LDX, * )
233: * ..
234: *
235: * =====================================================================
236: *
237: * .. Parameters ..
238: DOUBLE PRECISION ZERO, ONE
239: PARAMETER ( ZERO = 0.0D0, ONE = 1.0D0 )
240: DOUBLE PRECISION TWO
241: PARAMETER ( TWO = 2.0D0 )
242: * ..
243: * .. Local Scalars ..
244: INTEGER ICMAX, J
245: DOUBLE PRECISION BBND, BI1, BI2, BIGNUM, BNORM, BR1, BR2, CI21,
246: $ CI22, CMAX, CNORM, CR21, CR22, CSI, CSR, LI21,
247: $ LR21, SMINI, SMLNUM, TEMP, U22ABS, UI11, UI11R,
248: $ UI12, UI12S, UI22, UR11, UR11R, UR12, UR12S,
249: $ UR22, XI1, XI2, XR1, XR2
250: * ..
251: * .. Local Arrays ..
252: LOGICAL RSWAP( 4 ), ZSWAP( 4 )
253: INTEGER IPIVOT( 4, 4 )
254: DOUBLE PRECISION CI( 2, 2 ), CIV( 4 ), CR( 2, 2 ), CRV( 4 )
255: * ..
256: * .. External Functions ..
257: DOUBLE PRECISION DLAMCH
258: EXTERNAL DLAMCH
259: * ..
260: * .. External Subroutines ..
261: EXTERNAL DLADIV
262: * ..
263: * .. Intrinsic Functions ..
264: INTRINSIC ABS, MAX
265: * ..
266: * .. Equivalences ..
267: EQUIVALENCE ( CI( 1, 1 ), CIV( 1 ) ),
268: $ ( CR( 1, 1 ), CRV( 1 ) )
269: * ..
270: * .. Data statements ..
271: DATA ZSWAP / .FALSE., .FALSE., .TRUE., .TRUE. /
272: DATA RSWAP / .FALSE., .TRUE., .FALSE., .TRUE. /
273: DATA IPIVOT / 1, 2, 3, 4, 2, 1, 4, 3, 3, 4, 1, 2, 4,
274: $ 3, 2, 1 /
275: * ..
276: * .. Executable Statements ..
277: *
278: * Compute BIGNUM
279: *
280: SMLNUM = TWO*DLAMCH( 'Safe minimum' )
281: BIGNUM = ONE / SMLNUM
282: SMINI = MAX( SMIN, SMLNUM )
283: *
284: * Don't check for input errors
285: *
286: INFO = 0
287: *
288: * Standard Initializations
289: *
290: SCALE = ONE
291: *
292: IF( NA.EQ.1 ) THEN
293: *
294: * 1 x 1 (i.e., scalar) system C X = B
295: *
296: IF( NW.EQ.1 ) THEN
297: *
298: * Real 1x1 system.
299: *
300: * C = ca A - w D
301: *
302: CSR = CA*A( 1, 1 ) - WR*D1
303: CNORM = ABS( CSR )
304: *
305: * If | C | < SMINI, use C = SMINI
306: *
307: IF( CNORM.LT.SMINI ) THEN
308: CSR = SMINI
309: CNORM = SMINI
310: INFO = 1
311: END IF
312: *
313: * Check scaling for X = B / C
314: *
315: BNORM = ABS( B( 1, 1 ) )
316: IF( CNORM.LT.ONE .AND. BNORM.GT.ONE ) THEN
317: IF( BNORM.GT.BIGNUM*CNORM )
318: $ SCALE = ONE / BNORM
319: END IF
320: *
321: * Compute X
322: *
323: X( 1, 1 ) = ( B( 1, 1 )*SCALE ) / CSR
324: XNORM = ABS( X( 1, 1 ) )
325: ELSE
326: *
327: * Complex 1x1 system (w is complex)
328: *
329: * C = ca A - w D
330: *
331: CSR = CA*A( 1, 1 ) - WR*D1
332: CSI = -WI*D1
333: CNORM = ABS( CSR ) + ABS( CSI )
334: *
335: * If | C | < SMINI, use C = SMINI
336: *
337: IF( CNORM.LT.SMINI ) THEN
338: CSR = SMINI
339: CSI = ZERO
340: CNORM = SMINI
341: INFO = 1
342: END IF
343: *
344: * Check scaling for X = B / C
345: *
346: BNORM = ABS( B( 1, 1 ) ) + ABS( B( 1, 2 ) )
347: IF( CNORM.LT.ONE .AND. BNORM.GT.ONE ) THEN
348: IF( BNORM.GT.BIGNUM*CNORM )
349: $ SCALE = ONE / BNORM
350: END IF
351: *
352: * Compute X
353: *
354: CALL DLADIV( SCALE*B( 1, 1 ), SCALE*B( 1, 2 ), CSR, CSI,
355: $ X( 1, 1 ), X( 1, 2 ) )
356: XNORM = ABS( X( 1, 1 ) ) + ABS( X( 1, 2 ) )
357: END IF
358: *
359: ELSE
360: *
361: * 2x2 System
362: *
1.8 bertrand 363: * Compute the real part of C = ca A - w D (or ca A**T - w D )
1.1 bertrand 364: *
365: CR( 1, 1 ) = CA*A( 1, 1 ) - WR*D1
366: CR( 2, 2 ) = CA*A( 2, 2 ) - WR*D2
367: IF( LTRANS ) THEN
368: CR( 1, 2 ) = CA*A( 2, 1 )
369: CR( 2, 1 ) = CA*A( 1, 2 )
370: ELSE
371: CR( 2, 1 ) = CA*A( 2, 1 )
372: CR( 1, 2 ) = CA*A( 1, 2 )
373: END IF
374: *
375: IF( NW.EQ.1 ) THEN
376: *
377: * Real 2x2 system (w is real)
378: *
379: * Find the largest element in C
380: *
381: CMAX = ZERO
382: ICMAX = 0
383: *
384: DO 10 J = 1, 4
385: IF( ABS( CRV( J ) ).GT.CMAX ) THEN
386: CMAX = ABS( CRV( J ) )
387: ICMAX = J
388: END IF
389: 10 CONTINUE
390: *
391: * If norm(C) < SMINI, use SMINI*identity.
392: *
393: IF( CMAX.LT.SMINI ) THEN
394: BNORM = MAX( ABS( B( 1, 1 ) ), ABS( B( 2, 1 ) ) )
395: IF( SMINI.LT.ONE .AND. BNORM.GT.ONE ) THEN
396: IF( BNORM.GT.BIGNUM*SMINI )
397: $ SCALE = ONE / BNORM
398: END IF
399: TEMP = SCALE / SMINI
400: X( 1, 1 ) = TEMP*B( 1, 1 )
401: X( 2, 1 ) = TEMP*B( 2, 1 )
402: XNORM = TEMP*BNORM
403: INFO = 1
404: RETURN
405: END IF
406: *
407: * Gaussian elimination with complete pivoting.
408: *
409: UR11 = CRV( ICMAX )
410: CR21 = CRV( IPIVOT( 2, ICMAX ) )
411: UR12 = CRV( IPIVOT( 3, ICMAX ) )
412: CR22 = CRV( IPIVOT( 4, ICMAX ) )
413: UR11R = ONE / UR11
414: LR21 = UR11R*CR21
415: UR22 = CR22 - UR12*LR21
416: *
417: * If smaller pivot < SMINI, use SMINI
418: *
419: IF( ABS( UR22 ).LT.SMINI ) THEN
420: UR22 = SMINI
421: INFO = 1
422: END IF
423: IF( RSWAP( ICMAX ) ) THEN
424: BR1 = B( 2, 1 )
425: BR2 = B( 1, 1 )
426: ELSE
427: BR1 = B( 1, 1 )
428: BR2 = B( 2, 1 )
429: END IF
430: BR2 = BR2 - LR21*BR1
431: BBND = MAX( ABS( BR1*( UR22*UR11R ) ), ABS( BR2 ) )
432: IF( BBND.GT.ONE .AND. ABS( UR22 ).LT.ONE ) THEN
433: IF( BBND.GE.BIGNUM*ABS( UR22 ) )
434: $ SCALE = ONE / BBND
435: END IF
436: *
437: XR2 = ( BR2*SCALE ) / UR22
438: XR1 = ( SCALE*BR1 )*UR11R - XR2*( UR11R*UR12 )
439: IF( ZSWAP( ICMAX ) ) THEN
440: X( 1, 1 ) = XR2
441: X( 2, 1 ) = XR1
442: ELSE
443: X( 1, 1 ) = XR1
444: X( 2, 1 ) = XR2
445: END IF
446: XNORM = MAX( ABS( XR1 ), ABS( XR2 ) )
447: *
448: * Further scaling if norm(A) norm(X) > overflow
449: *
450: IF( XNORM.GT.ONE .AND. CMAX.GT.ONE ) THEN
451: IF( XNORM.GT.BIGNUM / CMAX ) THEN
452: TEMP = CMAX / BIGNUM
453: X( 1, 1 ) = TEMP*X( 1, 1 )
454: X( 2, 1 ) = TEMP*X( 2, 1 )
455: XNORM = TEMP*XNORM
456: SCALE = TEMP*SCALE
457: END IF
458: END IF
459: ELSE
460: *
461: * Complex 2x2 system (w is complex)
462: *
463: * Find the largest element in C
464: *
465: CI( 1, 1 ) = -WI*D1
466: CI( 2, 1 ) = ZERO
467: CI( 1, 2 ) = ZERO
468: CI( 2, 2 ) = -WI*D2
469: CMAX = ZERO
470: ICMAX = 0
471: *
472: DO 20 J = 1, 4
473: IF( ABS( CRV( J ) )+ABS( CIV( J ) ).GT.CMAX ) THEN
474: CMAX = ABS( CRV( J ) ) + ABS( CIV( J ) )
475: ICMAX = J
476: END IF
477: 20 CONTINUE
478: *
479: * If norm(C) < SMINI, use SMINI*identity.
480: *
481: IF( CMAX.LT.SMINI ) THEN
482: BNORM = MAX( ABS( B( 1, 1 ) )+ABS( B( 1, 2 ) ),
483: $ ABS( B( 2, 1 ) )+ABS( B( 2, 2 ) ) )
484: IF( SMINI.LT.ONE .AND. BNORM.GT.ONE ) THEN
485: IF( BNORM.GT.BIGNUM*SMINI )
486: $ SCALE = ONE / BNORM
487: END IF
488: TEMP = SCALE / SMINI
489: X( 1, 1 ) = TEMP*B( 1, 1 )
490: X( 2, 1 ) = TEMP*B( 2, 1 )
491: X( 1, 2 ) = TEMP*B( 1, 2 )
492: X( 2, 2 ) = TEMP*B( 2, 2 )
493: XNORM = TEMP*BNORM
494: INFO = 1
495: RETURN
496: END IF
497: *
498: * Gaussian elimination with complete pivoting.
499: *
500: UR11 = CRV( ICMAX )
501: UI11 = CIV( ICMAX )
502: CR21 = CRV( IPIVOT( 2, ICMAX ) )
503: CI21 = CIV( IPIVOT( 2, ICMAX ) )
504: UR12 = CRV( IPIVOT( 3, ICMAX ) )
505: UI12 = CIV( IPIVOT( 3, ICMAX ) )
506: CR22 = CRV( IPIVOT( 4, ICMAX ) )
507: CI22 = CIV( IPIVOT( 4, ICMAX ) )
508: IF( ICMAX.EQ.1 .OR. ICMAX.EQ.4 ) THEN
509: *
510: * Code when off-diagonals of pivoted C are real
511: *
512: IF( ABS( UR11 ).GT.ABS( UI11 ) ) THEN
513: TEMP = UI11 / UR11
514: UR11R = ONE / ( UR11*( ONE+TEMP**2 ) )
515: UI11R = -TEMP*UR11R
516: ELSE
517: TEMP = UR11 / UI11
518: UI11R = -ONE / ( UI11*( ONE+TEMP**2 ) )
519: UR11R = -TEMP*UI11R
520: END IF
521: LR21 = CR21*UR11R
522: LI21 = CR21*UI11R
523: UR12S = UR12*UR11R
524: UI12S = UR12*UI11R
525: UR22 = CR22 - UR12*LR21
526: UI22 = CI22 - UR12*LI21
527: ELSE
528: *
529: * Code when diagonals of pivoted C are real
530: *
531: UR11R = ONE / UR11
532: UI11R = ZERO
533: LR21 = CR21*UR11R
534: LI21 = CI21*UR11R
535: UR12S = UR12*UR11R
536: UI12S = UI12*UR11R
537: UR22 = CR22 - UR12*LR21 + UI12*LI21
538: UI22 = -UR12*LI21 - UI12*LR21
539: END IF
540: U22ABS = ABS( UR22 ) + ABS( UI22 )
541: *
542: * If smaller pivot < SMINI, use SMINI
543: *
544: IF( U22ABS.LT.SMINI ) THEN
545: UR22 = SMINI
546: UI22 = ZERO
547: INFO = 1
548: END IF
549: IF( RSWAP( ICMAX ) ) THEN
550: BR2 = B( 1, 1 )
551: BR1 = B( 2, 1 )
552: BI2 = B( 1, 2 )
553: BI1 = B( 2, 2 )
554: ELSE
555: BR1 = B( 1, 1 )
556: BR2 = B( 2, 1 )
557: BI1 = B( 1, 2 )
558: BI2 = B( 2, 2 )
559: END IF
560: BR2 = BR2 - LR21*BR1 + LI21*BI1
561: BI2 = BI2 - LI21*BR1 - LR21*BI1
562: BBND = MAX( ( ABS( BR1 )+ABS( BI1 ) )*
563: $ ( U22ABS*( ABS( UR11R )+ABS( UI11R ) ) ),
564: $ ABS( BR2 )+ABS( BI2 ) )
565: IF( BBND.GT.ONE .AND. U22ABS.LT.ONE ) THEN
566: IF( BBND.GE.BIGNUM*U22ABS ) THEN
567: SCALE = ONE / BBND
568: BR1 = SCALE*BR1
569: BI1 = SCALE*BI1
570: BR2 = SCALE*BR2
571: BI2 = SCALE*BI2
572: END IF
573: END IF
574: *
575: CALL DLADIV( BR2, BI2, UR22, UI22, XR2, XI2 )
576: XR1 = UR11R*BR1 - UI11R*BI1 - UR12S*XR2 + UI12S*XI2
577: XI1 = UI11R*BR1 + UR11R*BI1 - UI12S*XR2 - UR12S*XI2
578: IF( ZSWAP( ICMAX ) ) THEN
579: X( 1, 1 ) = XR2
580: X( 2, 1 ) = XR1
581: X( 1, 2 ) = XI2
582: X( 2, 2 ) = XI1
583: ELSE
584: X( 1, 1 ) = XR1
585: X( 2, 1 ) = XR2
586: X( 1, 2 ) = XI1
587: X( 2, 2 ) = XI2
588: END IF
589: XNORM = MAX( ABS( XR1 )+ABS( XI1 ), ABS( XR2 )+ABS( XI2 ) )
590: *
591: * Further scaling if norm(A) norm(X) > overflow
592: *
593: IF( XNORM.GT.ONE .AND. CMAX.GT.ONE ) THEN
594: IF( XNORM.GT.BIGNUM / CMAX ) THEN
595: TEMP = CMAX / BIGNUM
596: X( 1, 1 ) = TEMP*X( 1, 1 )
597: X( 2, 1 ) = TEMP*X( 2, 1 )
598: X( 1, 2 ) = TEMP*X( 1, 2 )
599: X( 2, 2 ) = TEMP*X( 2, 2 )
600: XNORM = TEMP*XNORM
601: SCALE = TEMP*SCALE
602: END IF
603: END IF
604: END IF
605: END IF
606: *
607: RETURN
608: *
609: * End of DLALN2
610: *
611: END
CVSweb interface <joel.bertrand@systella.fr>