summaryrefslogtreecommitdiff
path: root/var/spack/repos/builtin/packages/nwchem/raman_displ.patch
blob: 7ff9b65ea5ebeda9e7ba213288cae26a79e2f817 (plain) (blame)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
Index: src/property/raman_input.F
===================================================================
--- src/property/raman_input.F	(revision 28032)
+++ src/property/raman_input.F	(revision 28033)
@@ -47,6 +47,7 @@
 c
 c  set some defaults
 c
+      field=' '
       plot = 'normal' ! normal or resonance
       line = 'lorentzian' !  lorentzian (l) or gaussian (g) lineshape
       width = 20.0D+00 ! full-width at half maximum (FWHM) in 1/cm
@@ -54,7 +55,6 @@
       hyperraman = .false. ! flag to calculate hyperaman terms
       vroa = .false. ! flag to calculate vibrational raman spec
       rmmodes = 0
-      first = 7
       last = 10000
       low = 0.0D+00
       high = 100000.0D+00
@@ -132,9 +132,9 @@
       else if(inp_compare(.false.,'first',test)) then
          if(.not. inp_i(first))
      $     call errquit(pname//'missing value for first',911, INPUT_ERR)
-         if (.not. rtdb_put(rtdb,'raman:first',mt_int,1,first))
-     $     call errquit(pname//'rtdb put failed',0, RTDB_ERR)
-c        --- determine first normal mode to use ---
+c        --- not setting default here, it will be set later after
+c            frequency calculation has been done so we know if we have
+c            a linear molecule or not
       else if(inp_compare(.false.,'last',test)) then
          if(.not. inp_i(last)) ! FA-06-16-12 bug-fixed (BEF: first AFT: last)
      $     call errquit(pname//'missing value for last',911, INPUT_ERR)
Index: src/property/task_raman.F
===================================================================
--- src/property/task_raman.F	(revision 28032)
+++ src/property/task_raman.F	(revision 28033)
@@ -59,6 +59,7 @@
 
       integer j,pos,first0 ! FA-06-15-12
       logical preraman     ! FA-06-18-12
+      logical linear
 
       character*32 pname
 
@@ -107,6 +108,12 @@
      $   call errquit(pname//'rtdb_put freq_done',911, RTDB_ERR)
       endif
 c
+c --------Figure out if molecule is linear------------    
+
+c     if vib module doesn't list molecule as linear, assume it is not
+      if (.not. rtdb_get(rtdb,'vib:linear',mt_log,1,linear))
+     $      linear=.false.
+c
 c --------Create/load reference geometry to get the number of atoms------------    
 
       if (.not.geom_create(geom,'geometry')) call errquit
@@ -116,7 +123,11 @@
       if (.not. geom_ncent(geom,nat)) 
      &  call errquit(pname//'geom_ncent failed?',3, GEOM_ERR)
       nc = nat*3
-      rmmodes = nc-6
+      if (linear) then
+        rmmodes = nc-5
+      else
+        rmmodes = nc-6
+      end if
 
 c      if (ga_nodeid().eq.0) then
 c       write(*,1) nat,nc,rmmodes
@@ -146,8 +157,13 @@
      $      low  = 0.0D+00 ! lowest wavenumber  normal mode to use
       if (.not. rtdb_get(rtdb,'raman:high',mt_dbl,1,high))
      $      high = 100000.0D+00 ! Highest wavenumber normal mode to use
-      if (.not. rtdb_get(rtdb,'raman:first',mt_int,1,first))
-     $      first = 7 ! first normal mode to use
+      if (.not. rtdb_get(rtdb,'raman:first',mt_int,1,first)) then
+            if (linear) then
+              first = 6 ! first normal mode to use
+            else
+              first = 7 ! first normal mode to use
+            end if
+      end if
       if (.not. rtdb_get(rtdb,'raman:last',mt_int,1,last))
      $      last = 10000 ! last normal mode to use
       if (.not. rtdb_get(rtdb,'raman:hyperraman',mt_log,1,hyperraman))
@@ -156,7 +172,11 @@
      $      vroa = .false. ! # flag to calculate vibrational 
       if (.not. rtdb_get(rtdb,'raman:preraman',mt_log,1,preraman))
      $      preraman = .false. ! # flag to do task_freq() and leave
-      first0=7 ! constant
+      if (linear) then
+        first0=6 ! constant
+      else
+        first0=7 ! constant
+      end if
 c ======== FA-debug =============== START
 c      if (ga_nodeid().eq.0) then
 c       write(*,2) plot,line,width,step_size,steps
@@ -172,8 +192,13 @@
          rmmodes = nc
 c
 c --- in case we want overide the defaults for modes to include ---
-         if (.not. rtdb_get(rtdb,'raman:first',mt_int,1,first))
-     $      first = 7 ! srtep size for displacement along modes
+         if (.not. rtdb_get(rtdb,'raman:first',mt_int,1,first)) then
+            if (linear) then
+              first = 6 ! srtep size for displacement along modes
+            else
+              first = 7 ! srtep size for displacement along modes
+            end if
+         end if
       endif 
 c
 c ----------alocate space for freq and normal modes----------------------------
@@ -294,7 +319,7 @@
 c ------------enough setup really do the calculation------------------------
       if (.not.preraman) then
       call task_raman_doit(rtdb,geom,nc,nat,
-     &                     first0, ! = 7 constant
+     &                     first0, ! = 6 or 7
      &                     first,last,rmmodes,
      &                     steps,nfreq,plot,line,width,
      &                     step_size,
@@ -336,7 +361,7 @@
 c
 c     == perform raman calculation ==
       subroutine task_raman_doit(rtdb,geom,nc,nat,
-     &                           first0, ! = 7 constant
+     &                           first0, ! = 6 or 7
      &                           first,last,rmmodes,
      &                           steps,nfreq,
      &                           plot,line,width,
@@ -495,7 +520,7 @@
      &             lbl_raman,  ! in: raman label
      &             begin,      ! in: 
      &             last,       ! in: 
-     &             first0,     ! in: = 7 constant
+     &             first0,     ! in: = 6 or 7
      &             eigenvecs,  ! in: hessian data (modes)
      &             eigenvals,  ! in: hessian data (frequencies)
      &             mass,       ! in: mass(i) i=1,nat
@@ -519,7 +544,7 @@
      &             lbl_raman,  ! in: raman label
      &             mode_ini,   ! in: 
      &             mode_end,   ! in: 
-     &             first0,     ! in: = 7 constant
+     &             first0,     ! in: = 6 or 7
      &             eigenvecs,  ! in: hessian data (modes)
      &             eigenvals,  ! in: hessian data (frequencies)
      &             mass,       ! in: mass(i) i=1,nat
@@ -541,7 +566,7 @@
      &             lbl_raman,  ! in: raman label
      &             begin,      ! in: starting mode
      &             last,       ! in: ending   mode
-     &             first0,     ! in: = 7 constant
+     &             first0,     ! in: = 6 or 7
      &             eigenvecs,  ! in: hessian data (modes)
      &             eigenvals,  ! in: hessian data (frequencies)
      &             mass,       ! in: mass(i) i=1,nat
@@ -596,7 +621,7 @@
      &             rmmodes,    ! in: total nr. modes
      &             rminfo,     ! in: stores raman info
      &             nc,nat,     ! in: (nc,nat)=(nr coord,nr atoms)
-     &             first0,     ! in: = 7 constant
+     &             first0,     ! in: = 6 or 7
      &             eigenvecs,  ! in: hessian data (modes)
      &             eigenvals,  ! in: hessian data (frequencies)
      &             mass,       ! in: mass(i) i=1,nat
@@ -757,7 +782,8 @@
      &                 step_size,
      &                 rminfo,
      &                 eigenvecs,
-     &                 mass)
+     &                 mass,
+     &                 first0)
 c ======== FA: Writing to file rminfo ========= START
 c         if (ga_nodeid().eq.0)
 c     &     write(*,*) 'BEF raman_write() ...'
@@ -783,7 +809,7 @@
      &             lbl_raman,  ! in: raman label
      &             begin,      ! in: starting mode
      &             last,       ! in: ending   mode
-     &             first0,     ! in: = 7 constant
+     &             first0,     ! in: = 6 or 7
      &             eigenvecs,  ! in: hessian data (modes)
      &             eigenvals,  ! in: hessian data (frequencies)
      &             mass,       ! in: mass(i) i=1,nat
@@ -890,7 +916,7 @@
      &             rmmodes,    ! in: total nr. modes
      &             rminfo,     ! in: stores raman info
      &             nc,nat,     ! in: (nc,nat)=(nr coord,nr atoms)
-     &             first0,     ! in: = 7 constant
+     &             first0,     ! in: = 6 or 7
      &             eigenvecs,  ! in: hessian data (modes)
      &             eigenvals,  ! in: hessian data (frequencies)
      &             mass,       ! in: mass(i) i=1,nat
@@ -915,7 +941,7 @@
      &             lbl_raman,  ! in: raman label
      &             mode_ini,   ! in: 
      &             mode_end,   ! in: 
-     &             first0,     ! in: = 7 constant
+     &             first0,     ! in: = 6 or 7
      &             eigenvecs,  ! in: hessian data (modes)
      &             eigenvals,  ! in: hessian data (frequencies)
      &             mass,       ! in: mass(i) i=1,nat
@@ -1036,7 +1062,7 @@
      &             rmmodes,    ! in: total nr. modes
      &             rminfo,     ! in: stores raman info
      &             nc,nat,     ! in: (nc,nat)=(nr coord,nr atoms)
-     &             first0,     ! in: = 7 constant
+     &             first0,     ! in: = 6 or 7
      &             eigenvecs,  ! in: hessian data (modes)
      &             eigenvals,  ! in: hessian data (frequencies)
      &             mass,       ! in: mass(i) i=1,nat
@@ -1058,7 +1084,7 @@
      &             lbl_raman,  ! in: raman label
      &             begin,      ! in: 
      &             last,       ! in: 
-     &             first0,     ! in: = 7 constant
+     &             first0,     ! in: = 6 or 7
      &             eigenvecs,  ! in: hessian data (modes)
      &             eigenvals,  ! in: hessian data (frequencies)
      &             mass,       ! in: mass(i) i=1,nat
@@ -1139,7 +1165,7 @@
      &             rmmodes,    ! in: total nr. modes
      &             rminfo,     ! in: stores raman info
      &             nc,nat,     ! in: (nc,nat)=(nr coord,nr atoms)
-     &             first0,     ! in: = 7 constant
+     &             first0,     ! in: = 6 or 7
      &             eigenvecs,  ! in: hessian data (modes)
      &             eigenvals,  ! in: hessian data (frequencies)
      &             mass,       ! in: mass(i) i=1,nat
Index: src/property/raman.F
===================================================================
--- src/property/raman.F	(revision 28032)
+++ src/property/raman.F	(revision 28033)
@@ -29,8 +29,8 @@
       integer rtdb    ! [input] rtdb handle
       integer natom   ! [input] number of atoms
       integer nat3    ! [input] 3*number of atoms
-      integer first   ! first mode to consider in aoresponse (default =7 ramana =1 hyperraman)
-      integer tmpmode ! set to fill rminfo from 1 ( not 7 for raman calc)
+      integer first   ! first mode to consider in aoresponse (default =6 or 7 raman =1 hyperraman)
+      integer tmpmode ! set to fill rminfo from 1 ( not 6 or 7 for raman calc)
       integer rmmodes ! # of raman active modes
 
       double precision rminfo(rmmodes,4) ! data for raman spec
@@ -41,6 +41,10 @@
       double precision ncoords(3,natom)    ! [scratch] coords after step
       double precision steps(3,natom)     ! [scratch] step generated by vector and scaled
 c
+      double precision length_of_step, scale
+      double precision ddot
+      external ddot
+c
       parameter (bohr2ang=0.52917724924D+00) ! CONVERSION OF BOHR to ANGSTROMS
 c -------------determine sign of the step---------------------------------
       if (iii.eq.1) then
@@ -57,13 +61,16 @@
 c     &       i4,',',i4,',',i4,',',i4,',',f15.8,')')
 c ======= FA-check rminfo(x,1) ======== END   
 c --------------------------------------------------------------------
-        ivec = 1
-        do iatom = 1,natom
-          do ixyz = 1,3
-            steps(ixyz,iatom)=sign*step_size*eigenvecs(ivec,imode)
-            ivec = ivec + 1
-          enddo ! ixyz
-        enddo ! iatom
+      ivec = 1
+      do iatom = 1,natom
+        do ixyz = 1,3
+          steps(ixyz,iatom)=eigenvecs(ivec,imode)
+          ivec = ivec + 1
+        enddo ! ixyz
+      enddo ! iatom
+      length_of_step = sqrt(ddot(nat3,steps,1,steps,1))
+      scale = sign*step_size/length_of_step
+      call dscal(nat3,scale,steps,1)
 
       call daxpy(nat3,1.0d00,steps,1,ncoords,1) ! mult coords 
       if (.not. geom_cart_coords_set(geom,ncoords))
@@ -85,7 +92,8 @@
      &                    step_size,! in : step of finite differencing
      &                    rminfo,   ! in : Raman data
      &                    eigenvecs,! in : normal modes eigenvectors (nat3,nat3)
-     &                    mass)     ! in : mass
+     &                    mass,     ! in : mass
+     &                    first0)   ! in : first nonzero mode (6 or 7)
 c
 c Authors: Jonathan Mullin, Northwestern University (ver 1: Jan. 2011)
 c          Fredy W. Aquino, Northwestern University (ver 2: Oct. 2012)
@@ -108,6 +116,7 @@
       integer imode   ! mode #
       integer natom   ! [input] number of atoms
       integer nat3    ! [input] 3*number of atoms
+      integer first0  ! [input] first nonzero mode (6 or 7)
 c
       double precision rminfo(rmmodes,4) ! raman data
       double precision step_size,stepsize ! [input] step of finite differencing
@@ -134,7 +143,7 @@
       call dfill(3*natom,0.0D+00,tmode,1) ! 
 c zero
       stepsize = zero
-      m = imode - 6 
+      m = imode - first0 + 1
       j=1
       i=1
       ar2   = zero ! alpha real