Télécharger dgeev.eso

Retour à la liste

Numérotation des lignes :

  1. C DGEEV SOURCE BP208322 20/09/18 21:15:48 10718
  2. *> \brief <b> DGEEV computes the eigenvalues and, optionally, the left and/or right eigenvectors for GE matrices</b>
  3. *
  4. * =========== DOCUMENTATION ===========
  5. *
  6. * Online html documentation available at
  7. * http://www.netlib.org/lapack/explore-html/
  8. *
  9. *> \htmlonly
  10. *> Download DGEEV + dependencies
  11. *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dgeev.f">
  12. *> [TGZ]</a>
  13. *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dgeev.f">
  14. *> [ZIP]</a>
  15. *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dgeev.f">
  16. *> [TXT]</a>
  17. *> \endhtmlonly
  18. *
  19. * Definition:
  20. * ===========
  21. *
  22. * SUBROUTINE DGEEV( JOBVL, JOBVR, N, A, LDA, WR, WI, VL, LDVL, VR,
  23. * LDVR, WORK, LWORK, INFO )
  24. *
  25. * .. Scalar Arguments ..
  26. * CHARACTER JOBVL, JOBVR
  27. * INTEGER INFO, LDA, LDVL, LDVR, LWORK, N
  28. * ..
  29. * .. Array Arguments ..
  30. * REAL*8 A( LDA, * ), VL( LDVL, * ), VR( LDVR, * ),
  31. * $ WI( * ), WORK( * ), WR( * )
  32. * ..
  33. *
  34. *
  35. *> \par Purpose:
  36. * =============
  37. *>
  38. *> \verbatim
  39. *>
  40. *> DGEEV computes for an N-by-N real nonsymmetric matrix A, the
  41. *> eigenvalues and, optionally, the left and/or right eigenvectors.
  42. *>
  43. *> The right eigenvector v(j) of A satisfies
  44. *> A * v(j) = lambda(j) * v(j)
  45. *> where lambda(j) is its eigenvalue.
  46. *> The left eigenvector u(j) of A satisfies
  47. *> u(j)**H * A = lambda(j) * u(j)**H
  48. *> where u(j)**H denotes the conjugate-transpose of u(j).
  49. *>
  50. *> The computed eigenvectors are normalized to have Euclidean norm
  51. *> equal to 1 and largest component real.
  52. *> \endverbatim
  53. *
  54. * Arguments:
  55. * ==========
  56. *
  57. *> \param[in] JOBVL
  58. *> \verbatim
  59. *> JOBVL is CHARACTER*1
  60. *> = 'N': left eigenvectors of A are not computed;
  61. *> = 'V': left eigenvectors of A are computed.
  62. *> \endverbatim
  63. *>
  64. *> \param[in] JOBVR
  65. *> \verbatim
  66. *> JOBVR is CHARACTER*1
  67. *> = 'N': right eigenvectors of A are not computed;
  68. *> = 'V': right eigenvectors of A are computed.
  69. *> \endverbatim
  70. *>
  71. *> \param[in] N
  72. *> \verbatim
  73. *> N is INTEGER
  74. *> The order of the matrix A. N >= 0.
  75. *> \endverbatim
  76. *>
  77. *> \param[in,out] A
  78. *> \verbatim
  79. *> A is REAL*8 array, dimension (LDA,N)
  80. *> On entry, the N-by-N matrix A.
  81. *> On exit, A has been overwritten.
  82. *> \endverbatim
  83. *>
  84. *> \param[in] LDA
  85. *> \verbatim
  86. *> LDA is INTEGER
  87. *> The leading dimension of the array A. LDA >= max(1,N).
  88. *> \endverbatim
  89. *>
  90. *> \param[out] WR
  91. *> \verbatim
  92. *> WR is REAL*8 array, dimension (N)
  93. *> \endverbatim
  94. *>
  95. *> \param[out] WI
  96. *> \verbatim
  97. *> WI is REAL*8 array, dimension (N)
  98. *> WR and WI contain the real and imaginary parts,
  99. *> respectively, of the computed eigenvalues. Complex
  100. *> conjugate pairs of eigenvalues appear consecutively
  101. *> with the eigenvalue having the positive imaginary part
  102. *> first.
  103. *> \endverbatim
  104. *>
  105. *> \param[out] VL
  106. *> \verbatim
  107. *> VL is REAL*8 array, dimension (LDVL,N)
  108. *> If JOBVL = 'V', the left eigenvectors u(j) are stored one
  109. *> after another in the columns of VL, in the same order
  110. *> as their eigenvalues.
  111. *> If JOBVL = 'N', VL is not referenced.
  112. *> If the j-th eigenvalue is real, then u(j) = VL(:,j),
  113. *> the j-th column of VL.
  114. *> If the j-th and (j+1)-st eigenvalues form a complex
  115. *> conjugate pair, then u(j) = VL(:,j) + i*VL(:,j+1) and
  116. *> u(j+1) = VL(:,j) - i*VL(:,j+1).
  117. *> \endverbatim
  118. *>
  119. *> \param[in] LDVL
  120. *> \verbatim
  121. *> LDVL is INTEGER
  122. *> The leading dimension of the array VL. LDVL >= 1; if
  123. *> JOBVL = 'V', LDVL >= N.
  124. *> \endverbatim
  125. *>
  126. *> \param[out] VR
  127. *> \verbatim
  128. *> VR is REAL*8 array, dimension (LDVR,N)
  129. *> If JOBVR = 'V', the right eigenvectors v(j) are stored one
  130. *> after another in the columns of VR, in the same order
  131. *> as their eigenvalues.
  132. *> If JOBVR = 'N', VR is not referenced.
  133. *> If the j-th eigenvalue is real, then v(j) = VR(:,j),
  134. *> the j-th column of VR.
  135. *> If the j-th and (j+1)-st eigenvalues form a complex
  136. *> conjugate pair, then v(j) = VR(:,j) + i*VR(:,j+1) and
  137. *> v(j+1) = VR(:,j) - i*VR(:,j+1).
  138. *> \endverbatim
  139. *>
  140. *> \param[in] LDVR
  141. *> \verbatim
  142. *> LDVR is INTEGER
  143. *> The leading dimension of the array VR. LDVR >= 1; if
  144. *> JOBVR = 'V', LDVR >= N.
  145. *> \endverbatim
  146. *>
  147. *> \param[out] WORK
  148. *> \verbatim
  149. *> WORK is REAL*8 array, dimension (MAX(1,LWORK))
  150. *> On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
  151. *> \endverbatim
  152. *>
  153. *> \param[in] LWORK
  154. *> \verbatim
  155. *> LWORK is INTEGER
  156. *> The dimension of the array WORK. LWORK >= max(1,3*N), and
  157. *> if JOBVL = 'V' or JOBVR = 'V', LWORK >= 4*N. For good
  158. *> performance, LWORK must generally be larger.
  159. *>
  160. *> If LWORK = -1, then a workspace query is assumed; the routine
  161. *> only calculates the optimal size of the WORK array, returns
  162. *> this value as the first entry of the WORK array, and no error
  163. *> message related to LWORK is issued by XERBLA.
  164. *> \endverbatim
  165. *>
  166. *> \param[out] INFO
  167. *> \verbatim
  168. *> INFO is INTEGER
  169. *> = 0: successful exit
  170. *> < 0: if INFO = -i, the i-th argument had an illegal value.
  171. *> > 0: if INFO = i, the QR algorithm failed to compute all the
  172. *> eigenvalues, and no eigenvectors have been computed;
  173. *> elements i+1:N of WR and WI contain eigenvalues which
  174. *> have converged.
  175. *> \endverbatim
  176. *
  177. * Authors:
  178. * ========
  179. *
  180. *> \author Univ. of Tennessee
  181. *> \author Univ. of California Berkeley
  182. *> \author Univ. of Colorado Denver
  183. *> \author NAG Ltd.
  184. *
  185. *> \date June 2016
  186. *
  187. * @precisions fortran d -> s
  188. *
  189. *> \ingroup doubleGEeigen
  190. *
  191. * =====================================================================
  192. SUBROUTINE DGEEV( JOBVL, JOBVR, N, A, LDA, WR, WI, VL, LDVL, VR,
  193. $ LDVR, WORK, LWORK, INFO )
  194. * implicit none
  195. IMPLICIT INTEGER(I-N)
  196. IMPLICIT REAL*8(A-H,O-Z)
  197. *
  198. * -- LAPACK driver routine (version 3.7.0) --
  199. * -- LAPACK is a software package provided by Univ. of Tennessee, --
  200. * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
  201. * June 2016
  202. *
  203. * .. Scalar Arguments ..
  204. CHARACTER JOBVL, JOBVR
  205. INTEGER INFO, LDA, LDVL, LDVR, LWORK, N
  206. * ..
  207. * .. Array Arguments ..
  208. REAL*8 A( LDA, * ), VL( LDVL, * ), VR( LDVR, * ),
  209. $ WI( * ), WORK( * ), WR( * )
  210. * ..
  211. *
  212. * =====================================================================
  213. *
  214. * .. Parameters ..
  215. REAL*8 ZERO, ONE
  216. PARAMETER ( ZERO = 0.0D0, ONE = 1.0D0 )
  217. * ..
  218. * .. Local Scalars ..
  219. LOGICAL LQUERY, SCALEA, WANTVL, WANTVR
  220. CHARACTER SIDE
  221. INTEGER HSWORK, I, IBAL, IERR, IHI, ILO, ITAU, IWRK, K,
  222. $ LWORK_TREVC, MAXWRK, MINWRK, NOUT
  223. REAL*8 ANRM, BIGNUM, CS, CSCALE, EPS, R, SCL, SMLNUM,
  224. $ SN
  225. * ..
  226. * .. Local Arrays ..
  227. LOGICAL SELECT( 1 )
  228. REAL*8 DUM( 1 )
  229. * ..
  230. * .. External Subroutines ..
  231. * EXTERNAL DGEBAK, DGEBAL, DGEHRD, DHSEQR, DLABAD, DLACPY,
  232. * $ DLARTG, DLASCL, DORGHR, DROT, DSCAL, DTREVC3,
  233. * $ XERBLA
  234. * ..
  235. * .. External Functions ..
  236. * LOGICAL LSAME
  237. * INTEGER IDAMAX, ILAENV
  238. * REAL*8 DLAMCH, DLANGE, DLAPY2, DNRM2
  239. * EXTERNAL LSAME, IDAMAX, ILAENV, DLAMCH, DLANGE, DLAPY2,
  240. * $ DNRM2
  241. * ..
  242. * .. Intrinsic Functions ..
  243. * INTRINSIC MAX, SQRT
  244. * ..
  245. * .. Executable Statements ..
  246. *
  247. * Test the input arguments
  248. *
  249. INFO = 0
  250. LQUERY = ( LWORK.EQ.-1 )
  251. WANTVL = ( JOBVL.EQ. 'V' )
  252. WANTVR = ( JOBVR.EQ. 'V' )
  253. IF( ( .NOT.WANTVL ) .AND. ( .NOT.( JOBVL.EQ. 'N' ) ) ) THEN
  254. INFO = -1
  255. ELSE IF( ( .NOT.WANTVR ) .AND. ( .NOT.( JOBVR.EQ. 'N' ) ) ) THEN
  256. INFO = -2
  257. ELSE IF( N.LT.0 ) THEN
  258. INFO = -3
  259. ELSE IF( LDA.LT.MAX( 1, N ) ) THEN
  260. INFO = -5
  261. ELSE IF( LDVL.LT.1 .OR. ( WANTVL .AND. LDVL.LT.N ) ) THEN
  262. INFO = -9
  263. ELSE IF( LDVR.LT.1 .OR. ( WANTVR .AND. LDVR.LT.N ) ) THEN
  264. INFO = -11
  265. END IF
  266. *
  267. * Compute workspace
  268. * (Note: Comments in the code beginning "Workspace:" describe the
  269. * minimal amount of workspace needed at that point in the code,
  270. * as well as the preferred amount for good performance.
  271. * NB refers to the optimal block size for the immediately
  272. * following subroutine, as returned by ILAENV.
  273. * HSWORK refers to the workspace preferred by DHSEQR, as
  274. * calculated below. HSWORK is computed assuming ILO=1 and IHI=N,
  275. * the worst case.)
  276. *
  277. IF( INFO.EQ.0 ) THEN
  278. IF( N.EQ.0 ) THEN
  279. MINWRK = 1
  280. MAXWRK = 1
  281. ELSE
  282. MAXWRK = 2*N + N*ILAENV( 1, 'DGEHRD', ' ', N, 1, N, 0 )
  283. IF( WANTVL ) THEN
  284. MINWRK = 4*N
  285. MAXWRK = MAX( MAXWRK, 2*N + ( N - 1 )*ILAENV( 1,
  286. $ 'DORGHR', ' ', N, 1, N, -1 ) )
  287. CALL DHSEQR( 'S', 'V', N, 1, N, A, LDA, WR, WI, VL, LDVL,
  288. $ WORK, -1, INFO )
  289. HSWORK = INT( WORK(1) )
  290. MAXWRK = MAX( MAXWRK, N + 1, N + HSWORK )
  291. CALL DTREVC3( 'L', 'B', SELECT, N, A, LDA,
  292. $ VL, LDVL, VR, LDVR, N, NOUT,
  293. $ WORK, -1, IERR )
  294. LWORK_TREVC = INT( WORK(1) )
  295. MAXWRK = MAX( MAXWRK, N + LWORK_TREVC )
  296. MAXWRK = MAX( MAXWRK, 4*N )
  297. ELSE IF( WANTVR ) THEN
  298. MINWRK = 4*N
  299. MAXWRK = MAX( MAXWRK, 2*N + ( N - 1 )*ILAENV( 1,
  300. $ 'DORGHR', ' ', N, 1, N, -1 ) )
  301. CALL DHSEQR( 'S', 'V', N, 1, N, A, LDA, WR, WI, VR, LDVR,
  302. $ WORK, -1, INFO )
  303. HSWORK = INT( WORK(1) )
  304. MAXWRK = MAX( MAXWRK, N + 1, N + HSWORK )
  305. CALL DTREVC3( 'R', 'B', SELECT, N, A, LDA,
  306. $ VL, LDVL, VR, LDVR, N, NOUT,
  307. $ WORK, -1, IERR )
  308. LWORK_TREVC = INT( WORK(1) )
  309. MAXWRK = MAX( MAXWRK, N + LWORK_TREVC )
  310. MAXWRK = MAX( MAXWRK, 4*N )
  311. ELSE
  312. MINWRK = 3*N
  313. CALL DHSEQR( 'E', 'N', N, 1, N, A, LDA, WR, WI, VR, LDVR,
  314. $ WORK, -1, INFO )
  315. HSWORK = INT( WORK(1) )
  316. MAXWRK = MAX( MAXWRK, N + 1, N + HSWORK )
  317. END IF
  318. MAXWRK = MAX( MAXWRK, MINWRK )
  319. END IF
  320. WORK( 1 ) = MAXWRK
  321. *
  322. IF( LWORK.LT.MINWRK .AND. .NOT.LQUERY ) THEN
  323. INFO = -13
  324. END IF
  325. END IF
  326. *
  327. IF( INFO.NE.0 ) THEN
  328. CALL XERBLA( 'DGEEV ', -INFO )
  329. RETURN
  330. ELSE IF( LQUERY ) THEN
  331. RETURN
  332. END IF
  333. *
  334. * Quick return if possible
  335. *
  336. IF( N.EQ.0 )
  337. $ RETURN
  338. *
  339. * Get machine constants
  340. *
  341. EPS = DLAMCH( 'P' )
  342. SMLNUM = DLAMCH( 'S' )
  343. BIGNUM = ONE / SMLNUM
  344. CALL DLABAD( SMLNUM, BIGNUM )
  345. SMLNUM = SQRT( SMLNUM ) / EPS
  346. BIGNUM = ONE / SMLNUM
  347. *
  348. * Scale A if max element outside range [SMLNUM,BIGNUM]
  349. *
  350. ANRM = DLANGE( 'M', N, N, A, LDA, DUM )
  351. SCALEA = .FALSE.
  352. IF( ANRM.GT.ZERO .AND. ANRM.LT.SMLNUM ) THEN
  353. SCALEA = .TRUE.
  354. CSCALE = SMLNUM
  355. ELSE IF( ANRM.GT.BIGNUM ) THEN
  356. SCALEA = .TRUE.
  357. CSCALE = BIGNUM
  358. END IF
  359. IF( SCALEA )
  360. $ CALL DLASCL( 'G', 0, 0, ANRM, CSCALE, N, N, A, LDA, IERR )
  361. *
  362. * Balance the matrix
  363. * (Workspace: need N)
  364. *
  365. IBAL = 1
  366. CALL DGEBAL( 'B', N, A, LDA, ILO, IHI, WORK( IBAL ), IERR )
  367. *
  368. * Reduce to upper Hessenberg form
  369. * (Workspace: need 3*N, prefer 2*N+N*NB)
  370. *
  371. ITAU = IBAL + N
  372. IWRK = ITAU + N
  373. CALL DGEHRD( N, ILO, IHI, A, LDA, WORK( ITAU ), WORK( IWRK ),
  374. $ LWORK-IWRK+1, IERR )
  375. *
  376. IF( WANTVL ) THEN
  377. *
  378. * Want left eigenvectors
  379. * Copy Householder vectors to VL
  380. *
  381. SIDE = 'L'
  382. CALL DLACPY( 'L', N, N, A, LDA, VL, LDVL )
  383. *
  384. * Generate orthogonal matrix in VL
  385. * (Workspace: need 3*N-1, prefer 2*N+(N-1)*NB)
  386. *
  387. CALL DORGHR( N, ILO, IHI, VL, LDVL, WORK( ITAU ), WORK( IWRK ),
  388. $ LWORK-IWRK+1, IERR )
  389. *
  390. * Perform QR iteration, accumulating Schur vectors in VL
  391. * (Workspace: need N+1, prefer N+HSWORK (see comments) )
  392. *
  393. IWRK = ITAU
  394. CALL DHSEQR( 'S', 'V', N, ILO, IHI, A, LDA, WR, WI, VL, LDVL,
  395. $ WORK( IWRK ), LWORK-IWRK+1, INFO )
  396. *
  397. IF( WANTVR ) THEN
  398. *
  399. * Want left and right eigenvectors
  400. * Copy Schur vectors to VR
  401. *
  402. SIDE = 'B'
  403. CALL DLACPY( 'F', N, N, VL, LDVL, VR, LDVR )
  404. END IF
  405. *
  406. ELSE IF( WANTVR ) THEN
  407. *
  408. * Want right eigenvectors
  409. * Copy Householder vectors to VR
  410. *
  411. SIDE = 'R'
  412. CALL DLACPY( 'L', N, N, A, LDA, VR, LDVR )
  413. *
  414. * Generate orthogonal matrix in VR
  415. * (Workspace: need 3*N-1, prefer 2*N+(N-1)*NB)
  416. *
  417. CALL DORGHR( N, ILO, IHI, VR, LDVR, WORK( ITAU ), WORK( IWRK ),
  418. $ LWORK-IWRK+1, IERR )
  419. *
  420. * Perform QR iteration, accumulating Schur vectors in VR
  421. * (Workspace: need N+1, prefer N+HSWORK (see comments) )
  422. *
  423. IWRK = ITAU
  424. CALL DHSEQR( 'S', 'V', N, ILO, IHI, A, LDA, WR, WI, VR, LDVR,
  425. $ WORK( IWRK ), LWORK-IWRK+1, INFO )
  426. *
  427. ELSE
  428. *
  429. * Compute eigenvalues only
  430. * (Workspace: need N+1, prefer N+HSWORK (see comments) )
  431. *
  432. IWRK = ITAU
  433. CALL DHSEQR( 'E', 'N', N, ILO, IHI, A, LDA, WR, WI, VR, LDVR,
  434. $ WORK( IWRK ), LWORK-IWRK+1, INFO )
  435. END IF
  436. *
  437. * If INFO .NE. 0 from DHSEQR, then quit
  438. *
  439. IF( INFO.NE.0 )
  440. $ GO TO 50
  441. *
  442. IF( WANTVL .OR. WANTVR ) THEN
  443. *
  444. * Compute left and/or right eigenvectors
  445. * (Workspace: need 4*N, prefer N + N + 2*N*NB)
  446. *
  447. CALL DTREVC3( SIDE, 'B', SELECT, N, A, LDA, VL, LDVL, VR, LDVR,
  448. $ N, NOUT, WORK( IWRK ), LWORK-IWRK+1, IERR )
  449. END IF
  450. *
  451. IF( WANTVL ) THEN
  452. *
  453. * Undo balancing of left eigenvectors
  454. * (Workspace: need N)
  455. *
  456. CALL DGEBAK( 'B', 'L', N, ILO, IHI, WORK( IBAL ), N, VL, LDVL,
  457. $ IERR )
  458. *
  459. * Normalize left eigenvectors and make largest component real
  460. *
  461. DO 20 I = 1, N
  462. IF( WI( I ).EQ.ZERO ) THEN
  463. SCL = ONE / DNRM2( N, VL( 1, I ), 1 )
  464. CALL DSCAL( N, SCL, VL( 1, I ), 1 )
  465. ELSE IF( WI( I ).GT.ZERO ) THEN
  466. SCL = ONE / DLAPY2( DNRM2( N, VL( 1, I ), 1 ),
  467. $ DNRM2( N, VL( 1, I+1 ), 1 ) )
  468. CALL DSCAL( N, SCL, VL( 1, I ), 1 )
  469. CALL DSCAL( N, SCL, VL( 1, I+1 ), 1 )
  470. DO 10 K = 1, N
  471. WORK( IWRK+K-1 ) = VL( K, I )**2 + VL( K, I+1 )**2
  472. 10 CONTINUE
  473. K = IDAMAX( N, WORK( IWRK ), 1 )
  474. CALL DLARTG( VL( K, I ), VL( K, I+1 ), CS, SN, R )
  475. CALL DROT( N, VL( 1, I ), 1, VL( 1, I+1 ), 1, CS, SN )
  476. VL( K, I+1 ) = ZERO
  477. END IF
  478. 20 CONTINUE
  479. END IF
  480. *
  481. IF( WANTVR ) THEN
  482. *
  483. * Undo balancing of right eigenvectors
  484. * (Workspace: need N)
  485. *
  486. CALL DGEBAK( 'B', 'R', N, ILO, IHI, WORK( IBAL ), N, VR, LDVR,
  487. $ IERR )
  488. *
  489. * Normalize right eigenvectors and make largest component real
  490. *
  491. DO 40 I = 1, N
  492. IF( WI( I ).EQ.ZERO ) THEN
  493. SCL = ONE / DNRM2( N, VR( 1, I ), 1 )
  494. CALL DSCAL( N, SCL, VR( 1, I ), 1 )
  495. ELSE IF( WI( I ).GT.ZERO ) THEN
  496. SCL = ONE / DLAPY2( DNRM2( N, VR( 1, I ), 1 ),
  497. $ DNRM2( N, VR( 1, I+1 ), 1 ) )
  498. CALL DSCAL( N, SCL, VR( 1, I ), 1 )
  499. CALL DSCAL( N, SCL, VR( 1, I+1 ), 1 )
  500. DO 30 K = 1, N
  501. WORK( IWRK+K-1 ) = VR( K, I )**2 + VR( K, I+1 )**2
  502. 30 CONTINUE
  503. K = IDAMAX( N, WORK( IWRK ), 1 )
  504. CALL DLARTG( VR( K, I ), VR( K, I+1 ), CS, SN, R )
  505. CALL DROT( N, VR( 1, I ), 1, VR( 1, I+1 ), 1, CS, SN )
  506. VR( K, I+1 ) = ZERO
  507. END IF
  508. 40 CONTINUE
  509. END IF
  510. *
  511. * Undo scaling if necessary
  512. *
  513. 50 CONTINUE
  514. IF( SCALEA ) THEN
  515. CALL DLASCL( 'G', 0, 0, CSCALE, ANRM, N-INFO, 1, WR( INFO+1 ),
  516. $ MAX( N-INFO, 1 ), IERR )
  517. CALL DLASCL( 'G', 0, 0, CSCALE, ANRM, N-INFO, 1, WI( INFO+1 ),
  518. $ MAX( N-INFO, 1 ), IERR )
  519. IF( INFO.GT.0 ) THEN
  520. CALL DLASCL( 'G', 0, 0, CSCALE, ANRM, ILO-1, 1, WR, N,
  521. $ IERR )
  522. CALL DLASCL( 'G', 0, 0, CSCALE, ANRM, ILO-1, 1, WI, N,
  523. $ IERR )
  524. END IF
  525. END IF
  526. *
  527. WORK( 1 ) = MAXWRK
  528. RETURN
  529. *
  530. * End of DGEEV
  531. *
  532. END
  533.  
  534.  
  535.  
  536.  

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