← Back to branch summary

~ubuntu-branches/ubuntu/trusty/nwchem/trusty-proposed

« back to all changes in this revision

Viewing changes to src/nwpw/band/kbpp/integrate_kbpp_e_band_local_new.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/nwpw/band/kbpp/integrate_kbpp_e_band_local_new.F
 
1
*
 
2
* $Id: integrate_kbpp_e_band_local_new.F 22503 2012-05-20 06:58:57Z d3y133 $
 
3
*
 
4
      subroutine integrate_kbpp_e_band_local_new(version,
 
5
     >                            nrho,drho,lmax,locp,nmax,
 
6
     >                            n_extra,n_expansion,zv,
 
7
     >                            vp,wp,rho,f,cs,sn,
 
8
     >                            nfft1,nfft2,nfft3,nprj,
 
9
     >                            G,vl,
 
10
     >                            n_prj,l_prj,m_prj,vnlnrm,
 
11
     >                            semicore,rho_sc_r,rho_sc_k,
 
12
     >                            nray,G_ray,vl_ray,vnl_ray,
 
13
     >                            rho_sc_k_ray,tmp_ray,
 
14
     >                            filter,
 
15
     >                            ierr)
 
16
      implicit none
 
17
      integer          version
 
18
      integer          nrho
 
19
      double precision drho
 
20
      integer          lmax
 
21
      integer          locp
 
22
      integer          nmax
 
23
      integer          n_extra,n_expansion(0:lmax)
 
24
      double precision zv
 
25
      double precision vp(nrho,0:lmax)
 
26
      double precision wp(nrho,0:(lmax+n_extra))
 
27
      double precision rho(nrho)
 
28
      double precision f(nrho)
 
29
      double precision cs(nrho)
 
30
      double precision sn(nrho)
 
31
 
 
32
      integer nfft1,nfft2,nfft3,nprj
 
33
      double precision G(nfft1,nfft2,nfft3,3)
 
34
      double precision vl(nfft1,nfft2,nfft3)
 
35
      integer          n_prj(nprj),l_prj(nprj),m_prj(nprj)
 
36
      double precision vnlnrm(nmax,nmax,0:lmax)
 
37
 
 
38
 
 
39
      logical semicore
 
40
      double precision rho_sc_r(nrho,2)
 
41
      double precision rho_sc_k(nfft1,nfft2,nfft3,4)
 
42
 
 
43
      integer nray
 
44
      double precision G_ray(nray)
 
45
      double precision vl_ray(nray,2)
 
46
      double precision vnl_ray(nray,0:(lmax+n_extra),2)
 
47
      double precision rho_sc_k_ray(nray,2,2)
 
48
      double precision tmp_ray(nray)
 
49
      logical filter
 
50
 
 
51
      integer ierr
 
52
 
 
53
#include "mafdecls.fh"
 
54
 
 
55
*     *** local variables ****
 
56
      integer np,taskid,MASTER
 
57
      parameter (MASTER=0)
 
58
 
 
59
      integer lcount,task_count,nfft3d
 
60
      integer k1,k2,k3,i,l,nx,n,nb,n1,n2
 
61
      double precision ecut,wcut
 
62
      double precision p,d
 
63
      double precision gx,gy,gz,a,q,dG,yp1
 
64
      integer indx(5,0:3)
 
65
 
 
66
*     **** external functions ****
 
67
      double precision dsum,simp,control_ecut,control_wcut,nwpw_splint
 
68
      external         dsum,simp,control_ecut,control_wcut,nwpw_splint
 
69
 
 
70
*     **** set up indx(n,l) --> to wp ****
 
71
      nb = lmax+1
 
72
      do l=0,lmax
 
73
         indx(1,l) = l
 
74
         do n=2,n_expansion(l)
 
75
            indx(n,l) = nb
 
76
            nb = nb+1
 
77
         end do
 
78
      end do
 
79
 
 
80
 
 
81
      call Parallel_np(np)
 
82
      call Parallel_taskid(taskid)
 
83
 
 
84
      nfft3d = (nfft1)*nfft2*nfft3
 
85
 
 
86
      if ((nrho/2)*2.eq.nrho) then
 
87
        write(*,*)"local pseudopotential not computed: nrho is not odd"
 
88
        ierr=2
 
89
        return
 
90
      end if
 
91
 
 
92
 
 
93
*::::::::::::::::::  Define non-local pseudopotential  ::::::::::::::::
 
94
      do l=0,lmax
 
95
        if (l.ne.locp) then
 
96
          do i=1,nrho
 
97
            vp(i,l)=vp(i,l)-vp(i,locp)
 
98
          end do
 
99
        end if
 
100
      end do
 
101
 
 
102
*:::::::::::::::::::::  Normarization constants  ::::::::::::::::::::::
 
103
      lcount = 0
 
104
      do l=0,lmax
 
105
        if (l.eq.locp) then
 
106
           do n2 = 1, n_expansion(l)
 
107
           do n1 = n2,n_expansion(l)
 
108
              vnlnrm(n1,n2,l) = 0.0d0
 
109
              vnlnrm(n2,n1,l) = 0.0d0
 
110
           end do
 
111
           end do
 
112
        else
 
113
           do n2 = 1, n_expansion(l)
 
114
              do i=1,nrho
 
115
                 f(i)=vp(i,l)*wp(i,indx(n2,l))**2
 
116
              end do
 
117
              a=simp(nrho,f,drho)
 
118
              vnlnrm(n2,n2,l) = a
 
119
              do n1 = n2+1,n_expansion(l)
 
120
                 do i=1,nrho
 
121
                    f(i)=vp(i,l)*wp(i,indx(n1,l))*wp(i,indx(n2,l))
 
122
                 end do
 
123
                 a=simp(nrho,f,drho)
 
124
                 vnlnrm(n1,n2,l) = a
 
125
                 vnlnrm(n2,n1,l) = a
 
126
              end do
 
127
           end do
 
128
           if (n_expansion(l).eq.1) then
 
129
              vnlnrm(1,1,l) = 1/a
 
130
           else if (n_expansion(l).eq.2) then
 
131
              d = vnlnrm(1,1,l)*vnlnrm(2,2,l)
 
132
     >          - vnlnrm(1,2,l)*vnlnrm(2,1,l)
 
133
              q = vnlnrm(1,1,l)
 
134
              vnlnrm(1,1,l) = vnlnrm(2,2,l)/d
 
135
              vnlnrm(2,2,l) = q/d
 
136
              vnlnrm(1,2,l) = -vnlnrm(1,2,l)/d
 
137
              vnlnrm(2,1,l) = -vnlnrm(2,1,l)/d
 
138
           else
 
139
              call nwpw_matrix_invert(n_expansion(l),vnlnrm(1,1,l),nmax)
 
140
           end if
 
141
        end if
 
142
      end do
 
143
 
 
144
*     ************* compute ray fourier transforms *********************
 
145
      call integrate_kbpp_e_band_ray(version,
 
146
     >                            nrho,drho,lmax,locp,nmax,
 
147
     >                            n_extra,n_expansion,zv,
 
148
     >                            vp,wp,rho,f,cs,sn,
 
149
     >                            nray,
 
150
     >                            G_ray,vl_ray,vnl_ray,
 
151
     >                            semicore,rho_sc_r,rho_sc_k_ray,
 
152
     >                            ierr)
 
153
 
 
154
*     **** filter the rays ****
 
155
      if (filter) then
 
156
         ecut = control_ecut()
 
157
         wcut = control_wcut()
 
158
         call kbpp_band_filter_ray(nray,G_ray,ecut,vl_ray)
 
159
         do l=0,lmax
 
160
            if (l.ne.locp) then
 
161
              do n=1,n_expansion(l)
 
162
                call kbpp_band_filter_ray(nray,G_ray,wcut,
 
163
     >                                    vnl_ray(1,indx(n,l),1))
 
164
              end do
 
165
            end if
 
166
         end do
 
167
         if (semicore) then
 
168
         call kbpp_band_filter_ray(nray,G_ray,ecut,rho_sc_k_ray(1,1,1))
 
169
         call kbpp_band_filter_ray(nray,G_ray,ecut,rho_sc_k_ray(1,2,1))
 
170
         end if
 
171
      end if
 
172
 
 
173
*     **** setup cubic bsplines ****
 
174
      dG = G_ray(3)-G_ray(2)
 
175
      !yp1 = (vl_ray(3,1)-vl_ray(2,1))/dG
 
176
      !**** five point formula ***
 
177
      yp1 = ( -50.0d0*vl_ray(2,1)
 
178
     >       + 96.0d0*vl_ray(3,1)
 
179
     >       - 72.0d0*vl_ray(4,1)
 
180
     >       + 32.0d0*vl_ray(5,1)
 
181
     >       -  6.0d0*vl_ray(6,1))/(24.0d0*dG)
 
182
      call nwpw_spline(G_ray(2),vl_ray(2,1),nray-1,yp1,0.0d0,
 
183
     >                          vl_ray(2,2),tmp_ray)
 
184
      do l=0,lmax
 
185
         if (l.ne.locp) then
 
186
            do n=1,n_expansion(l)
 
187
               call nwpw_spline(G_ray,vnl_ray(1,indx(n,l),1),nray,
 
188
     >                             0.0d0,0.0d0,
 
189
     >                             vnl_ray(1,indx(n,l),2),tmp_ray)
 
190
            end do
 
191
         end if
 
192
      end do
 
193
      if (semicore) then
 
194
         call nwpw_spline(G_ray,rho_sc_k_ray(1,1,1),nray,0.0d0,0.0d0,
 
195
     >                          rho_sc_k_ray(1,1,2),tmp_ray)
 
196
         call nwpw_spline(G_ray,rho_sc_k_ray(1,2,1),nray,0.0d0,0.0d0,
 
197
     >                          rho_sc_k_ray(1,2,2),tmp_ray)
 
198
      end if
 
199
 
 
200
 
 
201
*======================  Fourier transformation  ======================
 
202
      call dcopy(nfft3d,0.0d0,0,vl,1)
 
203
      call dcopy(4*nfft3d,0.0d0,0,rho_sc_k,1)
 
204
      task_count = -1
 
205
      DO 700 k3=1,nfft3
 
206
      DO 700 k2=1,nfft2
 
207
      DO 700 k1=1,nfft1
 
208
        task_count = task_count + 1
 
209
        if (mod(task_count,np).ne.taskid) go to 700
 
210
 
 
211
        Q=DSQRT(G(k1,k2,k3,1)**2
 
212
     >         +G(k1,k2,k3,2)**2
 
213
     >         +G(k1,k2,k3,3)**2)
 
214
        nx = (Q/dG) + 1.0d0
 
215
 
 
216
        if ((k1.eq.1).and.(k2.eq.1).and.(k3.eq.1)) go to 700
 
217
 
 
218
        
 
219
        GX=G(k1,k2,k3,1)/Q
 
220
        GY=G(k1,k2,k3,2)/Q
 
221
        GZ=G(k1,k2,k3,3)/Q
 
222
        do i=1,nrho
 
223
          cs(i)=DCOS(Q*rho(i))
 
224
          sn(i)=DSIN(Q*rho(i))
 
225
        end do
 
226
 
 
227
*::::::::::::::::::::::::::::::  local  :::::::::::::::::::::::::::::::
 
228
  600   CONTINUE
 
229
 
 
230
        vl(k1,k2,k3)= nwpw_splint(G_ray(2),vl_ray(2,1),
 
231
     >                                     vl_ray(2,2),nray-1,nx-1,Q)
 
232
 
 
233
 
 
234
*::::::::::::::::::::: semicore density :::::::::::::::::::::::::::::::
 
235
        if (semicore) then
 
236
           p = nwpw_splint(G_ray,rho_sc_k_ray(1,1,1),
 
237
     >                           rho_sc_k_ray(1,1,2),nray,nx,Q)
 
238
           rho_sc_k(k1,k2,k3,1) = p
 
239
 
 
240
           p = nwpw_splint(G_ray,rho_sc_k_ray(1,2,1),
 
241
     >                           rho_sc_k_ray(1,2,2),nray,nx,Q)
 
242
           rho_sc_k(k1,k2,k3,2)=p*GX
 
243
           rho_sc_k(k1,k2,k3,3)=p*GY
 
244
           rho_sc_k(k1,k2,k3,4)=p*GZ
 
245
 
 
246
        end if
 
247
    
 
248
  700 CONTINUE
 
249
      call Parallel_Vector_SumAll(4*nfft3d,rho_sc_k)
 
250
      call Parallel_Vector_SumAll(nfft3d,vl)
 
251
 
 
252
*:::::::::::::::::::::::::::::::  G=0  ::::::::::::::::::::::::::::::::      
 
253
 
 
254
*     **** local potential ****
 
255
      vl(1,1,1)=vl_ray(1,1)
 
256
 
 
257
*     **** semicore density ****
 
258
      if (semicore) then
 
259
         rho_sc_k(1,1,1,1) = rho_sc_k_ray(1,1,1)
 
260
         rho_sc_k(1,1,1,2) = 0.0d0
 
261
         rho_sc_k(1,1,1,3) = 0.0d0
 
262
         rho_sc_k(1,1,1,4) = 0.0d0
 
263
      end if
 
264
 
 
265
*     ********************************
 
266
*     **** define n_prj and l_prj ****
 
267
*     ********************************
 
268
      lcount = nprj+1
 
269
      GO TO (950,940,930,920), lmax+1
 
270
 
 
271
        !::::::  f-wave  :::::::
 
272
  920   CONTINUE
 
273
        if (locp.ne.3) then
 
274
          do n=1,n_expansion(3)
 
275
             lcount = lcount-1
 
276
             n_prj(lcount) = n
 
277
             l_prj(lcount) = 3
 
278
             m_prj(lcount) = -3
 
279
 
 
280
             lcount = lcount-1
 
281
             n_prj(lcount) = n
 
282
             l_prj(lcount) = 3
 
283
             m_prj(lcount) = -2
 
284
 
 
285
             lcount = lcount-1
 
286
             n_prj(lcount) = n
 
287
             l_prj(lcount) = 3
 
288
             m_prj(lcount) = -1
 
289
 
 
290
             lcount = lcount-1
 
291
             n_prj(lcount) = n
 
292
             l_prj(lcount) = 3
 
293
             m_prj(lcount) = 0
 
294
 
 
295
             lcount = lcount-1
 
296
             n_prj(lcount) = n
 
297
             l_prj(lcount) = 3
 
298
             m_prj(lcount) = 1
 
299
 
 
300
             lcount = lcount-1
 
301
             n_prj(lcount) = n
 
302
             l_prj(lcount) = 3
 
303
             m_prj(lcount) = 2
 
304
 
 
305
             lcount = lcount-1
 
306
             n_prj(lcount) = n
 
307
             l_prj(lcount) = 3
 
308
             m_prj(lcount) = 3
 
309
          end do
 
310
        end if
 
311
 
 
312
        !::::  d-wave  ::::
 
313
  930   CONTINUE
 
314
        if (locp.ne.2) then
 
315
          do n=1,n_expansion(2)
 
316
             lcount = lcount-1
 
317
             n_prj(lcount) = n
 
318
             l_prj(lcount) = 2
 
319
             m_prj(lcount) = -2
 
320
 
 
321
             lcount = lcount-1
 
322
             n_prj(lcount) = n
 
323
             l_prj(lcount) = 2
 
324
             m_prj(lcount) = -1
 
325
 
 
326
             lcount = lcount-1
 
327
             n_prj(lcount) = n
 
328
             l_prj(lcount) = 2
 
329
             m_prj(lcount) = 0
 
330
 
 
331
             lcount = lcount-1
 
332
             n_prj(lcount) = n
 
333
             l_prj(lcount) = 2
 
334
             m_prj(lcount) = 1
 
335
 
 
336
             lcount = lcount-1
 
337
             n_prj(lcount) = n
 
338
             l_prj(lcount) = 2
 
339
             m_prj(lcount) = 2
 
340
          end do
 
341
        end if
 
342
 
 
343
        !::::  p-wave  ::::
 
344
  940   CONTINUE
 
345
        if (locp.ne.1) then
 
346
          do n=1,n_expansion(1)
 
347
             lcount = lcount-1
 
348
             n_prj(lcount) = n
 
349
             l_prj(lcount) = 1
 
350
             m_prj(lcount) = -1
 
351
 
 
352
             lcount = lcount-1
 
353
             n_prj(lcount) = n
 
354
             l_prj(lcount) = 1
 
355
             m_prj(lcount) = 0
 
356
 
 
357
             lcount = lcount-1
 
358
             n_prj(lcount) = n
 
359
             l_prj(lcount) = 1
 
360
             m_prj(lcount) = 1
 
361
          end do
 
362
        end if
 
363
 
 
364
        !::::  s-wave  ::::
 
365
  950   CONTINUE
 
366
        if (locp.ne.0) then
 
367
          do n=1,n_expansion(0)
 
368
             lcount = lcount-1
 
369
             n_prj(lcount) = n
 
370
             l_prj(lcount) = 0
 
371
             m_prj(lcount) = 0
 
372
          end do
 
373
        end if
 
374
 
 
375
      ierr=0
 
376
      return
 
377
      end