← Back to branch summary

~ubuntu-branches/ubuntu/utopic/nwchem/utopic

« back to all changes in this revision

Viewing changes to src/lapack/double/dgeev.f

  • Committer: Package Import Robot
  • Author(s): Michael Banck, Daniel Leidert, Andreas Tille, Michael Banck
  • Date: 2013-07-04 12:14:55 UTC
  • mfrom: (1.1.2)
  • Revision ID: package-import@ubuntu.com-20130704121455-5tvsx2qabor3nrui
Tags: 6.3-1
http://bugs.debian.org/696361
* New upstream release.
* Fixes anisotropic properties (Closes: #696361).
* New features include:
  + Multi-reference coupled cluster (MRCC) approaches
  + Hybrid DFT calculations with short-range HF 
  + New density-functionals including Minnesota (M08, M11) and HSE hybrid
    functionals
  + X-ray absorption spectroscopy (XAS) with TDDFT
  + Analytical gradients for the COSMO solvation model
  + Transition densities from TDDFT 
  + DFT+U and Electron-Transfer (ET) methods for plane wave calculations
  + Exploitation of space group symmetry in plane wave geometry optimizations
  + Local density of states (LDOS) collective variable added to Metadynamics
  + Various new XC functionals added for plane wave calculations, including
    hybrid and range-corrected ones
  + Electric field gradients with relativistic corrections 
  + Nudged Elastic Band optimization method
  + Updated basis sets and ECPs 

[ Daniel Leidert ]
* debian/watch: Fixed.

[ Andreas Tille ]
* debian/upstream: References

[ Michael Banck ]
* debian/upstream (Name): New field.
* debian/patches/02_makefile_flags.patch: Refreshed.
* debian/patches/06_statfs_kfreebsd.patch: Likewise.
* debian/patches/07_ga_target_force_linux.patch: Likewise.
* debian/patches/05_avoid_inline_assembler.patch: Removed, no longer needed.
* debian/patches/09_backported_6.1.1_fixes.patch: Likewise.
* debian/control (Build-Depends): Added gfortran-4.7 and gcc-4.7.
* debian/patches/10_force_gcc-4.7.patch: New patch, explicitly sets
  gfortran-4.7 and gcc-4.7, fixes test suite hang with gcc-4.8 (Closes:
  #701328, #713262).
* debian/testsuite: Added tests for COSMO analytical gradients and MRCC.
* debian/rules (MRCC_METHODS): New variable, required to enable MRCC methods.

Show diffs side-by-side

added added

removed removed

Lines of Context:
src/lapack/double/dgeev.f
 
1
*> \brief <b> DGEEV computes the eigenvalues and, optionally, the left and/or right eigenvectors for GE matrices</b>
 
2
*
 
3
*  =========== DOCUMENTATION ===========
 
4
*
 
5
* Online html documentation available at 
 
6
*            http://www.netlib.org/lapack/explore-html/ 
 
7
*
 
8
*> \htmlonly
 
9
*> Download DGEEV + dependencies 
 
10
*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dgeev.f"> 
 
11
*> [TGZ]</a> 
 
12
*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dgeev.f"> 
 
13
*> [ZIP]</a> 
 
14
*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dgeev.f"> 
 
15
*> [TXT]</a>
 
16
*> \endhtmlonly 
 
17
*
 
18
*  Definition:
 
19
*  ===========
 
20
*
 
21
*       SUBROUTINE DGEEV( JOBVL, JOBVR, N, A, LDA, WR, WI, VL, LDVL, VR,
 
22
*                         LDVR, WORK, LWORK, INFO )
 
23
* 
 
24
*       .. Scalar Arguments ..
 
25
*       CHARACTER          JOBVL, JOBVR
 
26
*       INTEGER            INFO, LDA, LDVL, LDVR, LWORK, N
 
27
*       ..
 
28
*       .. Array Arguments ..
 
29
*       DOUBLE PRECISION   A( LDA, * ), VL( LDVL, * ), VR( LDVR, * ),
 
30
*      $                   WI( * ), WORK( * ), WR( * )
 
31
*       ..
 
32
*  
 
33
*
 
34
*> \par Purpose:
 
35
*  =============
 
36
*>
 
37
*> \verbatim
 
38
*>
 
39
*> DGEEV computes for an N-by-N real nonsymmetric matrix A, the
 
40
*> eigenvalues and, optionally, the left and/or right eigenvectors.
 
41
*>
 
42
*> The right eigenvector v(j) of A satisfies
 
43
*>                  A * v(j) = lambda(j) * v(j)
 
44
*> where lambda(j) is its eigenvalue.
 
45
*> The left eigenvector u(j) of A satisfies
 
46
*>               u(j)**H * A = lambda(j) * u(j)**H
 
47
*> where u(j)**H denotes the conjugate-transpose of u(j).
 
48
*>
 
49
*> The computed eigenvectors are normalized to have Euclidean norm
 
50
*> equal to 1 and largest component real.
 
51
*> \endverbatim
 
52
*
 
53
*  Arguments:
 
54
*  ==========
 
55
*
 
56
*> \param[in] JOBVL
 
57
*> \verbatim
 
58
*>          JOBVL is CHARACTER*1
 
59
*>          = 'N': left eigenvectors of A are not computed;
 
60
*>          = 'V': left eigenvectors of A are computed.
 
61
*> \endverbatim
 
62
*>
 
63
*> \param[in] JOBVR
 
64
*> \verbatim
 
65
*>          JOBVR is CHARACTER*1
 
66
*>          = 'N': right eigenvectors of A are not computed;
 
67
*>          = 'V': right eigenvectors of A are computed.
 
68
*> \endverbatim
 
69
*>
 
70
*> \param[in] N
 
71
*> \verbatim
 
72
*>          N is INTEGER
 
73
*>          The order of the matrix A. N >= 0.
 
74
*> \endverbatim
 
75
*>
 
76
*> \param[in,out] A
 
77
*> \verbatim
 
78
*>          A is DOUBLE PRECISION array, dimension (LDA,N)
 
79
*>          On entry, the N-by-N matrix A.
 
80
*>          On exit, A has been overwritten.
 
81
*> \endverbatim
 
82
*>
 
83
*> \param[in] LDA
 
84
*> \verbatim
 
85
*>          LDA is INTEGER
 
86
*>          The leading dimension of the array A.  LDA >= max(1,N).
 
87
*> \endverbatim
 
88
*>
 
89
*> \param[out] WR
 
90
*> \verbatim
 
91
*>          WR is DOUBLE PRECISION array, dimension (N)
 
92
*> \endverbatim
 
93
*>
 
94
*> \param[out] WI
 
95
*> \verbatim
 
96
*>          WI is DOUBLE PRECISION array, dimension (N)
 
97
*>          WR and WI contain the real and imaginary parts,
 
98
*>          respectively, of the computed eigenvalues.  Complex
 
99
*>          conjugate pairs of eigenvalues appear consecutively
 
100
*>          with the eigenvalue having the positive imaginary part
 
101
*>          first.
 
102
*> \endverbatim
 
103
*>
 
104
*> \param[out] VL
 
105
*> \verbatim
 
106
*>          VL is DOUBLE PRECISION array, dimension (LDVL,N)
 
107
*>          If JOBVL = 'V', the left eigenvectors u(j) are stored one
 
108
*>          after another in the columns of VL, in the same order
 
109
*>          as their eigenvalues.
 
110
*>          If JOBVL = 'N', VL is not referenced.
 
111
*>          If the j-th eigenvalue is real, then u(j) = VL(:,j),
 
112
*>          the j-th column of VL.
 
113
*>          If the j-th and (j+1)-st eigenvalues form a complex
 
114
*>          conjugate pair, then u(j) = VL(:,j) + i*VL(:,j+1) and
 
115
*>          u(j+1) = VL(:,j) - i*VL(:,j+1).
 
116
*> \endverbatim
 
117
*>
 
118
*> \param[in] LDVL
 
119
*> \verbatim
 
120
*>          LDVL is INTEGER
 
121
*>          The leading dimension of the array VL.  LDVL >= 1; if
 
122
*>          JOBVL = 'V', LDVL >= N.
 
123
*> \endverbatim
 
124
*>
 
125
*> \param[out] VR
 
126
*> \verbatim
 
127
*>          VR is DOUBLE PRECISION array, dimension (LDVR,N)
 
128
*>          If JOBVR = 'V', the right eigenvectors v(j) are stored one
 
129
*>          after another in the columns of VR, in the same order
 
130
*>          as their eigenvalues.
 
131
*>          If JOBVR = 'N', VR is not referenced.
 
132
*>          If the j-th eigenvalue is real, then v(j) = VR(:,j),
 
133
*>          the j-th column of VR.
 
134
*>          If the j-th and (j+1)-st eigenvalues form a complex
 
135
*>          conjugate pair, then v(j) = VR(:,j) + i*VR(:,j+1) and
 
136
*>          v(j+1) = VR(:,j) - i*VR(:,j+1).
 
137
*> \endverbatim
 
138
*>
 
139
*> \param[in] LDVR
 
140
*> \verbatim
 
141
*>          LDVR is INTEGER
 
142
*>          The leading dimension of the array VR.  LDVR >= 1; if
 
143
*>          JOBVR = 'V', LDVR >= N.
 
144
*> \endverbatim
 
145
*>
 
146
*> \param[out] WORK
 
147
*> \verbatim
 
148
*>          WORK is DOUBLE PRECISION array, dimension (MAX(1,LWORK))
 
149
*>          On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
 
150
*> \endverbatim
 
151
*>
 
152
*> \param[in] LWORK
 
153
*> \verbatim
 
154
*>          LWORK is INTEGER
 
155
*>          The dimension of the array WORK.  LWORK >= max(1,3*N), and
 
156
*>          if JOBVL = 'V' or JOBVR = 'V', LWORK >= 4*N.  For good
 
157
*>          performance, LWORK must generally be larger.
 
158
*>
 
159
*>          If LWORK = -1, then a workspace query is assumed; the routine
 
160
*>          only calculates the optimal size of the WORK array, returns
 
161
*>          this value as the first entry of the WORK array, and no error
 
162
*>          message related to LWORK is issued by XERBLA.
 
163
*> \endverbatim
 
164
*>
 
165
*> \param[out] INFO
 
166
*> \verbatim
 
167
*>          INFO is INTEGER
 
168
*>          = 0:  successful exit
 
169
*>          < 0:  if INFO = -i, the i-th argument had an illegal value.
 
170
*>          > 0:  if INFO = i, the QR algorithm failed to compute all the
 
171
*>                eigenvalues, and no eigenvectors have been computed;
 
172
*>                elements i+1:N of WR and WI contain eigenvalues which
 
173
*>                have converged.
 
174
*> \endverbatim
 
175
*
 
176
*  Authors:
 
177
*  ========
 
178
*
 
179
*> \author Univ. of Tennessee 
 
180
*> \author Univ. of California Berkeley 
 
181
*> \author Univ. of Colorado Denver 
 
182
*> \author NAG Ltd. 
 
183
*
 
184
*> \date September 2012
 
185
*
 
186
*> \ingroup doubleGEeigen
 
187
*
 
188
*  =====================================================================
1
189
      SUBROUTINE DGEEV( JOBVL, JOBVR, N, A, LDA, WR, WI, VL, LDVL, VR,
2
190
     $                  LDVR, WORK, LWORK, INFO )
3
 
c $Id: dgeev.f 19697 2010-10-29 16:57:34Z d3y133 $
4
191
*
5
 
*  -- LAPACK driver routine (version 2.0) --
6
 
*     Univ. of Tennessee, Univ. of California Berkeley, NAG Ltd.,
7
 
*     Courant Institute, Argonne National Lab, and Rice University
8
 
*     September 30, 1994
 
192
*  -- LAPACK driver routine (version 3.4.2) --
 
193
*  -- LAPACK is a software package provided by Univ. of Tennessee,    --
 
194
*  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
 
195
*     September 2012
9
196
*
10
197
*     .. Scalar Arguments ..
11
198
      CHARACTER          JOBVL, JOBVR
16
203
     $                   WI( * ), WORK( * ), WR( * )
17
204
*     ..
18
205
*
19
 
*  Purpose
20
 
*  =======
21
 
*
22
 
*  DGEEV computes for an N-by-N real nonsymmetric matrix A, the
23
 
*  eigenvalues and, optionally, the left and/or right eigenvectors.
24
 
*
25
 
*  The right eigenvector v(j) of A satisfies
26
 
*                   A * v(j) = lambda(j) * v(j)
27
 
*  where lambda(j) is its eigenvalue.
28
 
*  The left eigenvector u(j) of A satisfies
29
 
*                u(j)**H * A = lambda(j) * u(j)**H
30
 
*  where u(j)**H denotes the conjugate transpose of u(j).
31
 
*
32
 
*  The computed eigenvectors are normalized to have Euclidean norm
33
 
*  equal to 1 and largest component real.
34
 
*
35
 
*  Arguments
36
 
*  =========
37
 
*
38
 
*  JOBVL   (input) CHARACTER*1
39
 
*          = 'N': left eigenvectors of A are not computed;
40
 
*          = 'V': left eigenvectors of A are computed.
41
 
*
42
 
*  JOBVR   (input) CHARACTER*1
43
 
*          = 'N': right eigenvectors of A are not computed;
44
 
*          = 'V': right eigenvectors of A are computed.
45
 
*
46
 
*  N       (input) INTEGER
47
 
*          The order of the matrix A. N >= 0.
48
 
*
49
 
*  A       (input/output) DOUBLE PRECISION array, dimension (LDA,N)
50
 
*          On entry, the N-by-N matrix A.
51
 
*          On exit, A has been overwritten.
52
 
*
53
 
*  LDA     (input) INTEGER
54
 
*          The leading dimension of the array A.  LDA >= max(1,N).
55
 
*
56
 
*  WR      (output) DOUBLE PRECISION array, dimension (N)
57
 
*  WI      (output) DOUBLE PRECISION array, dimension (N)
58
 
*          WR and WI contain the real and imaginary parts,
59
 
*          respectively, of the computed eigenvalues.  Complex
60
 
*          conjugate pairs of eigenvalues appear consecutively
61
 
*          with the eigenvalue having the positive imaginary part
62
 
*          first.
63
 
*
64
 
*  VL      (output) DOUBLE PRECISION array, dimension (LDVL,N)
65
 
*          If JOBVL = 'V', the left eigenvectors u(j) are stored one
66
 
*          after another in the columns of VL, in the same order
67
 
*          as their eigenvalues.
68
 
*          If JOBVL = 'N', VL is not referenced.
69
 
*          If the j-th eigenvalue is real, then u(j) = VL(:,j),
70
 
*          the j-th column of VL.
71
 
*          If the j-th and (j+1)-st eigenvalues form a complex
72
 
*          conjugate pair, then u(j) = VL(:,j) + i*VL(:,j+1) and
73
 
*          u(j+1) = VL(:,j) - i*VL(:,j+1).
74
 
*
75
 
*  LDVL    (input) INTEGER
76
 
*          The leading dimension of the array VL.  LDVL >= 1; if
77
 
*          JOBVL = 'V', LDVL >= N.
78
 
*
79
 
*  VR      (output) DOUBLE PRECISION array, dimension (LDVR,N)
80
 
*          If JOBVR = 'V', the right eigenvectors v(j) are stored one
81
 
*          after another in the columns of VR, in the same order
82
 
*          as their eigenvalues.
83
 
*          If JOBVR = 'N', VR is not referenced.
84
 
*          If the j-th eigenvalue is real, then v(j) = VR(:,j),
85
 
*          the j-th column of VR.
86
 
*          If the j-th and (j+1)-st eigenvalues form a complex
87
 
*          conjugate pair, then v(j) = VR(:,j) + i*VR(:,j+1) and
88
 
*          v(j+1) = VR(:,j) - i*VR(:,j+1).
89
 
*
90
 
*  LDVR    (input) INTEGER
91
 
*          The leading dimension of the array VR.  LDVR >= 1; if
92
 
*          JOBVR = 'V', LDVR >= N.
93
 
*
94
 
*  WORK    (workspace/output) DOUBLE PRECISION array, dimension (LWORK)
95
 
*          On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
96
 
*
97
 
*  LWORK   (input) INTEGER
98
 
*          The dimension of the array WORK.  LWORK >= max(1,3*N), and
99
 
*          if JOBVL = 'V' or JOBVR = 'V', LWORK >= 4*N.  For good
100
 
*          performance, LWORK must generally be larger.
101
 
*
102
 
*  INFO    (output) INTEGER
103
 
*          = 0:  successful exit
104
 
*          < 0:  if INFO = -i, the i-th argument had an illegal value.
105
 
*          > 0:  if INFO = i, the QR algorithm failed to compute all the
106
 
*                eigenvalues, and no eigenvectors have been computed;
107
 
*                elements i+1:N of WR and WI contain eigenvalues which
108
 
*                have converged.
109
 
*
110
206
*  =====================================================================
111
207
*
112
208
*     .. Parameters ..
114
210
      PARAMETER          ( ZERO = 0.0D0, ONE = 1.0D0 )
115
211
*     ..
116
212
*     .. Local Scalars ..
117
 
      LOGICAL            SCALEA, WANTVL, WANTVR
 
213
      LOGICAL            LQUERY, SCALEA, WANTVL, WANTVR
118
214
      CHARACTER          SIDE
119
215
      INTEGER            HSWORK, I, IBAL, IERR, IHI, ILO, ITAU, IWRK, K,
120
 
     $                   MAXB, MAXWRK, MINWRK, NOUT
 
216
     $                   MAXWRK, MINWRK, NOUT
121
217
      DOUBLE PRECISION   ANRM, BIGNUM, CS, CSCALE, EPS, R, SCL, SMLNUM,
122
218
     $                   SN
123
219
*     ..
138
234
     $                   DNRM2
139
235
*     ..
140
236
*     .. Intrinsic Functions ..
141
 
      INTRINSIC          MAX, MIN, SQRT
 
237
      INTRINSIC          MAX, SQRT
142
238
*     ..
143
239
*     .. Executable Statements ..
144
240
*
145
241
*     Test the input arguments
146
242
*
147
243
      INFO = 0
 
244
      LQUERY = ( LWORK.EQ.-1 )
148
245
      WANTVL = LSAME( JOBVL, 'V' )
149
246
      WANTVR = LSAME( JOBVR, 'V' )
150
247
      IF( ( .NOT.WANTVL ) .AND. ( .NOT.LSAME( JOBVL, 'N' ) ) ) THEN
171
268
*       calculated below. HSWORK is computed assuming ILO=1 and IHI=N,
172
269
*       the worst case.)
173
270
*
174
 
      MINWRK = 1
175
 
      IF( INFO.EQ.0 .AND. LWORK.GE.1 ) THEN
176
 
         MAXWRK = 2*N + N*ILAENV( 1, 'DGEHRD', ' ', N, 1, N, 0 )
177
 
         IF( ( .NOT.WANTVL ) .AND. ( .NOT.WANTVR ) ) THEN
178
 
            MINWRK = MAX( 1, 3*N )
179
 
            MAXB = MAX( ILAENV( 8, 'DHSEQR', 'EN', N, 1, N, -1 ), 2 )
180
 
            K = MIN( MAXB, N, MAX( 2, ILAENV( 4, 'DHSEQR', 'EN', N, 1,
181
 
     $          N, -1 ) ) )
182
 
            HSWORK = MAX( K*( K+2 ), 2*N )
183
 
            MAXWRK = MAX( MAXWRK, N+1, N+HSWORK )
 
271
      IF( INFO.EQ.0 ) THEN
 
272
         IF( N.EQ.0 ) THEN
 
273
            MINWRK = 1
 
274
            MAXWRK = 1
184
275
         ELSE
185
 
            MINWRK = MAX( 1, 4*N )
186
 
            MAXWRK = MAX( MAXWRK, 2*N+( N-1 )*
187
 
     $               ILAENV( 1, 'DORGHR', ' ', N, 1, N, -1 ) )
188
 
            MAXB = MAX( ILAENV( 8, 'DHSEQR', 'SV', N, 1, N, -1 ), 2 )
189
 
            K = MIN( MAXB, N, MAX( 2, ILAENV( 4, 'DHSEQR', 'SV', N, 1,
190
 
     $          N, -1 ) ) )
191
 
            HSWORK = MAX( K*( K+2 ), 2*N )
192
 
            MAXWRK = MAX( MAXWRK, N+1, N+HSWORK )
193
 
            MAXWRK = MAX( MAXWRK, 4*N )
 
276
            MAXWRK = 2*N + N*ILAENV( 1, 'DGEHRD', ' ', N, 1, N, 0 )
 
277
            IF( WANTVL ) THEN
 
278
               MINWRK = 4*N
 
279
               MAXWRK = MAX( MAXWRK, 2*N + ( N - 1 )*ILAENV( 1,
 
280
     $                       'DORGHR', ' ', N, 1, N, -1 ) )
 
281
               CALL DHSEQR( 'S', 'V', N, 1, N, A, LDA, WR, WI, VL, LDVL,
 
282
     $                WORK, -1, INFO )
 
283
               HSWORK = WORK( 1 )
 
284
               MAXWRK = MAX( MAXWRK, N + 1, N + HSWORK )
 
285
               MAXWRK = MAX( MAXWRK, 4*N )
 
286
            ELSE IF( WANTVR ) THEN
 
287
               MINWRK = 4*N
 
288
               MAXWRK = MAX( MAXWRK, 2*N + ( N - 1 )*ILAENV( 1,
 
289
     $                       'DORGHR', ' ', N, 1, N, -1 ) )
 
290
               CALL DHSEQR( 'S', 'V', N, 1, N, A, LDA, WR, WI, VR, LDVR,
 
291
     $                WORK, -1, INFO )
 
292
               HSWORK = WORK( 1 )
 
293
               MAXWRK = MAX( MAXWRK, N + 1, N + HSWORK )
 
294
               MAXWRK = MAX( MAXWRK, 4*N )
 
295
            ELSE 
 
296
               MINWRK = 3*N
 
297
               CALL DHSEQR( 'E', 'N', N, 1, N, A, LDA, WR, WI, VR, LDVR,
 
298
     $                WORK, -1, INFO )
 
299
               HSWORK = WORK( 1 )
 
300
               MAXWRK = MAX( MAXWRK, N + 1, N + HSWORK )
 
301
            END IF
 
302
            MAXWRK = MAX( MAXWRK, MINWRK )
194
303
         END IF
195
304
         WORK( 1 ) = MAXWRK
196
 
      END IF
197
 
      IF( LWORK.LT.MINWRK ) THEN
198
 
         INFO = -13
199
 
      END IF
 
305
*
 
306
         IF( LWORK.LT.MINWRK .AND. .NOT.LQUERY ) THEN
 
307
            INFO = -13
 
308
         END IF
 
309
      END IF
 
310
*
200
311
      IF( INFO.NE.0 ) THEN
201
312
         CALL XERBLA( 'DGEEV ', -INFO )
202
313
         RETURN
 
314
      ELSE IF( LQUERY ) THEN
 
315
         RETURN
203
316
      END IF
204
317
*
205
318
*     Quick return if possible