Télécharger dlaln2.eso

Retour à la liste

Numérotation des lignes :

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

© Cast3M 2003 - Tous droits réservés.
Mentions légales