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

« back to all changes in this revision

Viewing changes to src/property/get_vecF1.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
* 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:
 
1
      subroutine get_vecF1(g_vecF1,    ! out:
 
2
     &                     g_vecF1_im, ! out:
 
3
     &                     g_rhs,      ! in :
 
4
     &                     g_rhs_im,   ! in :
 
5
     &                     vectors,    ! in : MO vectors
 
6
     &                     nbf,        ! in : nr. basis functions
 
7
     &                     nmo,        ! in : nr. MOs
 
8
     &                     ncomp,      ! in :
 
9
     &                     npol,       ! in : nr. polarizations
 
10
     &                     lifetime,   ! in : = (.true.,.false.) with/out damping
 
11
     &                     nocc,       ! in : nr. occupied MOs
 
12
     &                     nvirt,      ! in : nr. virtual  MOs
 
13
     &                     debug)      ! in : = .true. for debugging
 
14
c
 
15
c Author : Fredy W. Aquino
 
16
c Date   : 03-15-12
 
17
c Note.- Modified from original aoresponse source code
 
18
c        for extension to spin-unrestricted case
 
19
c        original aoresponse source code was written by 
 
20
c        J. Autschbach and appears on nwchem-devtrunk (date:03-02-12)
 
21
c         
 
22
c --> Experimental (not published yet)
 
23
 
 
24
       implicit none
 
25
#include "errquit.fh"
 
26
#include "global.fh"
 
27
#include "mafdecls.fh"
 
28
#include "msgids.fh"
 
29
#include "stdio.fh"
 
30
      integer ispin,npol,ncomp,nocc(npol),nvirt(npol)
 
31
c Note.- g_vecF1(npol,ncomp) does not work for two dimensions
 
32
      integer g_vecF1(2,2),g_vecF1_im(2,2),
 
33
     &        g_rhs(ncomp),g_rhs_im(ncomp),
 
34
     &        g_u
 
35
      integer vectors(npol)
 
36
      integer ipm,ifld,ndir,nbf,nmo,disp
 
37
      logical lifetime,debug
 
38
      integer alo(3), ahi(3), 
 
39
     &        blo(3), bhi(3),
 
40
     &        clo(3), chi(3),
 
41
     &        plo(3), phi(3),
 
42
     &        qlo(3), qhi(3)
 
43
      character*256 cstemp
 
44
      ndir=3 ! = nr directions (x,y,z)
 
45
c ---- Create (g_vecF1,g_vecF1_im) -------- START
 
46
      alo(1) = nbf
 
47
      alo(2) = -1
 
48
      alo(3) = -1
 
49
      ahi(1) = nbf
 
50
      ahi(2) = nbf
 
51
      ahi(3) = ndir
 
52
      do ispin=1,npol
 
53
       do ipm=1,ncomp
 
54
        write (cstemp,'(a,i1)') 'aor vecF1 ',ipm
 
55
        if (.not.nga_create(MT_DBL,3,ahi,trim(cstemp),
 
56
     &    alo,g_vecF1(ispin,ipm) ))
 
57
     &     call errquit('fiao_f1: nga_create failed vecF1',0,GA_ERR)
 
58
        call ga_zero(g_vecF1(ispin,ipm))
 
59
        if (lifetime) then
 
60
          write (cstemp,'(a,i1)') 'aor vecF1_Im ',ipm
 
61
          if (.not.nga_create(MT_DBL,3,ahi,trim(cstemp),
 
62
     &       alo,g_vecF1_im(ispin,ipm) ))
 
63
     &       call errquit('fiao_f1: nga_create failed E1_im',0,GA_ERR)
 
64
          call ga_zero(g_vecF1_im(ispin,ipm))
 
65
        end if
 
66
       enddo ! end-loop-ipm
 
67
      enddo ! end-loop-ispin
 
68
c ---- Create (g_vecF1,g_vecF1_im) -------- END
 
69
 
 
70
      do ispin=1,npol
 
71
       alo(1) = nbf
 
72
       alo(2) = -1
 
73
       alo(3) = -1
 
74
       ahi(1) = nbf  
 
75
       ahi(2) = nocc(ispin)
 
76
       ahi(3) = 1
 
77
       if (.not.nga_create(MT_DBL,3,ahi,'U matrix',alo,g_u)) call 
 
78
     &    errquit('giao_b1: nga_create failed g_u',0,GA_ERR)
 
79
c     C1 = C0 * U10
 
80
       alo(1) = 1
 
81
       ahi(1) = nmo
 
82
       alo(2) = 1
 
83
       ahi(2) = nocc(ispin)  
 
84
       alo(3) = 1
 
85
       ahi(3) = 1       
 
86
       blo(1) = 1
 
87
       bhi(1) = nbf
 
88
       blo(2) = 1
 
89
       bhi(2) = nmo 
 
90
       blo(3) = 1
 
91
       bhi(3) = 1
 
92
       clo(1) = 1
 
93
       chi(1) = nbf
 
94
       clo(2) = 1
 
95
       chi(2) = nocc(ispin)
 
96
       disp=nocc(1)*nvirt(1)*(ispin-1)
 
97
       plo(1) = disp+1
 
98
       phi(1) = disp+nocc(ispin)*nvirt(ispin)
 
99
       qlo(1) = nocc(ispin)+1
 
100
       qhi(1) = nmo
 
101
       qlo(2) = 1
 
102
       qhi(2) = nocc(ispin)
 
103
       qlo(3) = 1
 
104
       qhi(3) = 1
 
105
       do ifld = 1,ndir
 
106
        clo(3) = ifld
 
107
        chi(3) = ifld
 
108
c      Make C1           
 
109
        do ipm = 1,ncomp
 
110
c ======== Including g_u_ov (g_rhs --> g_u) ==== START
 
111
        plo(2) = ifld
 
112
        phi(2) = ifld
 
113
        call ga_zero(g_u)
 
114
        call nga_copy_patch('n',g_rhs(ipm),plo,phi,
 
115
     &                          g_u       ,qlo,qhi)
 
116
c ======== Including g_u_ov (g_rhs --> g_u) ==== END
 
117
         call nga_matmul_patch('n','n',1.0d0,0.0d0,
 
118
     &                         vectors(ispin)    ,blo,bhi,
 
119
     &                         g_u               ,alo,ahi,
 
120
     &                         g_vecF1(ispin,ipm),clo,chi)          
 
121
         if (lifetime) then
 
122
c ======== Including g_u_ov_im (g_rhs_im --> g_u_im) == START
 
123
         call ga_zero(g_u)
 
124
         call nga_copy_patch('n',g_rhs_im(ipm),plo,phi,
 
125
     &                           g_u          ,qlo,qhi)
 
126
c ======== Including g_u_ov_im (g_rhs_im --> g_u_im) == END
 
127
         call nga_matmul_patch('n','n',1.0d0,0.0d0,
 
128
     &                         vectors(ispin)       ,blo,bhi,
 
129
     &                         g_u                  ,alo,ahi,
 
130
     &                         g_vecF1_im(ispin,ipm),clo,chi)  
 
131
         end if ! end-if-lifetime (damping)     
 
132
        enddo ! end-loop-ipm
 
133
       enddo ! end-loop-ifld
 
134
      if (.not.ga_destroy(g_u)) call 
 
135
     &    errquit('fiao_b1: ga_destroy failed g_d1',0,GA_ERR)
 
136
      enddo ! end-loop-ispin
 
137
      return
 
138
      end
 
139
 
 
140
      subroutine update_rhs_dipole(
 
141
     &                     g_rhs,    ! in/out: 
 
142
     &                     vectors,  ! in : MO vectors
 
143
     &                     rtdb,     ! in : rtdb  handle
 
144
     &                     basis,    ! in : basis handle
 
145
     &                     lvelocity,! in : logical var
 
146
     &                     nat,      ! in : nr. atoms
 
147
     &                     npol,     ! in : nr of polarizations
 
148
     &                     nocc,     ! in : nr. occ  shells
 
149
     &                     nvirt,    ! in : nr. virt shells
 
150
     &                     nbf,      ! in : nr. basis functions
 
151
     &                     nmo,      ! in : nr. MOs   
 
152
     &                     ncomp,    ! in : nr components of ...
 
153
     &                     debug)    ! in : logical for debugging
 
154
c
 
155
c Author : Fredy W. Aquino
 
156
c Date   : 03-15-12
 
157
c Note.- Modified from original aoresponse source code
 
158
c        for extension to spin-unrestricted case
 
159
c        original aoresponse source code was written by 
 
160
c        J. Autschbach and appears on nwchem-devtrunk (date:03-02-12)
 
161
c         
 
162
c --> Experimental (not published yet)
 
163
 
 
164
       implicit none
 
165
#include "errquit.fh"
 
166
#include "global.fh"
 
167
#include "mafdecls.fh"
 
168
#include "msgids.fh"
 
169
#include "stdio.fh"
 
170
#include "geom.fh"
 
171
#include "rtdb.fh"
 
172
#include "bas.fh"
 
173
#include "prop.fh"
 
174
      integer npol,ncomp
 
175
      integer g_dipole,g_dipole1,g_rhs(ncomp)
 
176
      integer rtdb,basis
 
177
      integer vectors(npol)
 
178
      integer ifld,ndir,nbf,nmo,ipm,nat,
 
179
     &        nocc(npol),nvirt(npol)
 
180
      double precision origin(3)
 
181
      data origin/0d0,0d0,0d0/
 
182
      integer shift,disp,ispin
 
183
      integer alo(3), ahi(3),
 
184
     &        blo(3), bhi(3) 
 
185
      logical lvelocity,oskel,debug
 
186
      external giao_aotomo
 
187
      oskel = .false.
 
188
      ndir=3 ! nr directions (x,y,z)
 
189
c     Get dipole integrals in GA and transform to MO set (virt,occ)
 
190
      alo(1) = nbf
 
191
      alo(2) = -1
 
192
      alo(3) = -1
 
193
      ahi(1) = nbf
 
194
      ahi(2) = nbf
 
195
      ahi(3) = ndir
 
196
      if (.not.nga_create(MT_DBL,3,ahi,'dip matrix',alo,g_dipole1)) call 
 
197
     &    errquit('fiao_f1: nga_create failed g_dipole',0,GA_ERR)
 
198
      ahi(3) = ndir*npol
 
199
      if (.not.nga_create(MT_DBL,3,ahi,'dip matrix',alo,g_dipole)) call 
 
200
     &    errquit('fiao_f1: nga_create failed g_dipole',0,GA_ERR)
 
201
      if (debug) write (luout,*) 'fiao_f1: dipole matrix allocated'
 
202
c     Get H10 in GA, using g_dipole array
 
203
c     note: origin has been set to (0,0,0) for multipole integs.
 
204
      call ga_zero(g_dipole1)
 
205
      if (lvelocity) then 
 
206
        call int_giao_1ega(basis,basis,
 
207
     &                     g_dipole1,'velocity', ! out
 
208
     &                     origin,nat, oskel)
 
209
        call ga_scale (g_dipole1, -1d0)
 
210
      else
 
211
        call int_mpole_1ega(basis,basis,
 
212
     &                      g_dipole1,'dipole', ! out
 
213
     &                      origin,oskel)
 
214
      end if
 
215
      if (debug) write (luout,*) 'fiao_f1: AO integrals done'
 
216
c     ga_rhs(a,i) = ga_rhs(a,i) + H10(a,i)
 
217
c     Transform H10 to MO and add to g_rhs.
 
218
c     (the rhs is NOT divided by (e_a - e_i -/+ omega), this
 
219
c     will be considered in the CPKS solver, in the precon-
 
220
c     ditioner and the 1e part of the "product" routine)
 
221
c --------- g_dipole1 --> g_dipole --------- START
 
222
      call ga_zero(g_dipole)       
 
223
       blo(1) = 1
 
224
       bhi(1) = nbf
 
225
       blo(2) = 1
 
226
       bhi(2) = nbf
 
227
       blo(3) = 1
 
228
       bhi(3) = ndir
 
229
      do ispin=1,npol  
 
230
       disp=ndir*(ispin-1) 
 
231
       alo(1) = 1
 
232
       ahi(1) = nbf
 
233
       alo(2) = 1
 
234
       ahi(2) = nbf
 
235
       alo(3) = disp+1
 
236
       ahi(3) = disp+ndir  
 
237
       call nga_copy_patch('n',g_dipole1,blo,bhi,
 
238
     &                         g_dipole ,alo,ahi) 
 
239
      enddo ! end-loop-ispin
 
240
c --------- g_dipole1 --> g_dipole --------- END
 
241
      call giao_aotomo(g_dipole,vectors,nocc,nvirt,npol,ndir,nbf)
 
242
      do ispin=1,npol
 
243
       shift=ndir*(ispin-1)
 
244
       alo(1) = nocc(ispin)+1
 
245
       ahi(1) = nmo
 
246
       alo(2) = 1
 
247
       ahi(2) = nocc(ispin)
 
248
       alo(3) = shift+1
 
249
       ahi(3) = shift+ndir
 
250
       disp=nocc(1)*nvirt(1)*(ispin-1)
 
251
       blo(1) = disp+1
 
252
       bhi(1) = disp+nocc(ispin)*nvirt(ispin)
 
253
       blo(2) = 1
 
254
       bhi(2) = ndir
 
255
       blo(3) = 1
 
256
       bhi(3) = 1
 
257
       if (debug) then
 
258
        if (ga_nodeid().eq.0)
 
259
     &   write(*,*) '---- g_dipole-nw-------- START'
 
260
        call ga_print(g_dipole)
 
261
        if (ga_nodeid().eq.0)
 
262
     &   write(*,*) '---- g_dipole-nw-------- END'
 
263
       endif ! end-if-debug
 
264
       do ipm = 1,ncomp
 
265
        call nga_add_patch(1.0d0,g_rhs(ipm),blo,bhi,
 
266
     &                     1.0d0,g_dipole  ,alo,ahi,
 
267
     &                           g_rhs(ipm),blo,bhi)
 
268
       enddo ! end-loop-ipm
 
269
      enddo ! end-loop-ispin
 
270
      if (debug) write (luout,*) 'fiao_f1: dipole added to rhs'
 
271
c     Cleanup g_dipole as we do not need it right now
 
272
      if (.not.ga_destroy(g_dipole)) call 
 
273
     &    errquit('fiao_f1: ga_destroy failed g_dipole',0,GA_ERR)
 
274
      if (.not.ga_destroy(g_dipole1)) call 
 
275
     &    errquit('fiao_f1: ga_destroy failed g_dipole',0,GA_ERR)
 
276
      return
 
277
      end
 
278
 
 
279
      subroutine get_nocc(rtdb,   ! in : rtdb handle
 
280
     &                    nocc,   ! out: nr occupations
 
281
     &                    npol,   ! out: nr of polarization
 
282
     &                    nclosed,! in : nr closed shells
 
283
     &                    nopen,  ! in : nr open shells
 
284
     &                    nvirt,  ! in : nr virtual MOs
 
285
     &                    scftyp, ! in : string = UHF or RHF
 
286
     &                    ntot)   ! out: sum_{i,npol} nocc(i)*nvirt(i)
 
287
c
 
288
c Author : Fredy W. Aquino
 
289
c Date   : 03-15-12
 
290
c Note.- Modified from original aoresponse source code
 
291
c        for extension to spin-unrestricted case
 
292
c        original aoresponse source code was written by 
 
293
c        J. Autschbach and appears on nwchem-devtrunk (date:03-02-12)
 
294
c                
 
295
c --> Experimental (not published yet)
 
296
 
 
297
       implicit none
 
298
#include "errquit.fh"
 
299
#include "mafdecls.fh"
 
300
#include "msgids.fh"
 
301
#include "stdio.fh"
 
302
#include "rtdb.fh"
 
303
#include "prop.fh"
 
304
      character*3 scftyp
 
305
      integer rtdb,ispin,
 
306
     &        npol,nocc(2),nclosed(2),
 
307
     &        nopen(2),nvirt(2),ntot
 
308
        if      (scftyp .eq. 'UHF') then
 
309
         npol=2
 
310
         nocc(1)=nopen(1)   
 
311
         nocc(2)=nopen(2)    
 
312
c ------ Store nopen in rtdb so that CPHF routine is happy ---- START
 
313
c In scf_init(): to avoid error message: 
 
314
c          ===>  scf: no. of closed-shell electrons is not even!
 
315
          if (.not. rtdb_put(rtdb, 'scf:nopen', 
 
316
     &         MT_INT, 1, nocc(1)-nocc(2)))
 
317
     *         call errquit('get_nocc:rtdbput nopen failed',
 
318
     &         nocc(1)-nocc(2),
 
319
     &       RTDB_ERR)
 
320
c ------ Store nopen in rtdb so that CPHF routine is happy ---- END
 
321
        else if (scftyp .eq. 'RHF') then
 
322
         npol=1
 
323
         nocc(1)=nclosed(1)
 
324
         nocc(2)=0
 
325
        endif
 
326
        ntot=0
 
327
        do ispin=1,npol
 
328
         ntot=ntot+nocc(ispin)*nvirt(ispin)
 
329
        enddo
 
330
      return
 
331
      end