← Back to branch summary

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

« back to all changes in this revision

Viewing changes to .pc/09_backported_6.1.1_fixes.patch/src/property/raman.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:
.pc/09_backported_6.1.1_fixes.patch/src/property/raman.F
1
 
c
2
 
c ---------------------------------------------------------------------
3
 
c computes the xyz coords along the imode normal mode
4
 
c ---------------------------------------------------------------------
5
 
c
6
 
      subroutine raman_modestep(rtdb,nat3,natom,geom,rmmodes,imode,
7
 
     &    iii,first,eigenvecs,eigenvals,ncoords,rminfo,step_size)
8
 
c
9
 
      implicit none
10
 
c
11
 
#include "errquit.fh"
12
 
#include "rtdb.fh"
13
 
#include "mafdecls.fh"
14
 
#include "stdio.fh"
15
 
#include "geom.fh"
16
 
#include "nwc_const.fh"
17
 
#include "inp.fh"
18
 
#include "global.fh"
19
 
#include "bas.fh"
20
 
c
21
 
      integer iii ! step number which determines sign 
22
 
      integer imode ! mode cordinates are dipsplaced along
23
 
      integer iatom, ixyz, ivec ! counting index
24
 
      integer geom    ! [input] geom handle
25
 
      integer rtdb    ! [input] rtdb handle
26
 
      integer natom   ! [input] number of atoms
27
 
      integer nat3    ! [input] 3*number of atoms
28
 
      integer first   ! first mode to consider in aoresponse (default =7 ramana =1 hyperraman)
29
 
      integer tmpmode ! set to fill rminfo from 1 ( not 7 for raman calc)
30
 
      integer rmmodes ! # of raman active modes
31
 
c
32
 
      double precision rminfo(rmmodes,4) ! data for raman spec
33
 
      double precision step_size ! multiplictive factor for step along normal mode
34
 
      double precision sign,bohr2ang ! sign of the step , unit conversion
35
 
      double precision eigenvecs(nat3,nat3) ! [input](xyz&atom,mode)
36
 
      double precision eigenvals(nat3)      ! [input] (mode)
37
 
      double precision ncoords(3,natom)    ! [scratch] coords after step
38
 
      double precision steps(3,natom)     ! [scratch] step generated by vector and scaled
39
 
c
40
 
      parameter (bohr2ang=0.52917724924D+00) ! CONVERSION OF BOHR to ANGSTROMS
41
 
c -------------determine sign of the step---------------------------------
42
 
      if (iii.eq.1) then
43
 
        sign = 1.0D+00
44
 
      else
45
 
        sign = -1.0D+00
46
 
      endif
47
 
c  --- save imode values in rminfo  ---
48
 
      tmpmode=imode-first+1
49
 
      rminfo(tmpmode,1) = eigenvals(imode)
50
 
c --------------------------------------------------------------------
51
 
        ivec = 1
52
 
        do iatom = 1,natom
53
 
          do ixyz = 1,3
54
 
      steps(ixyz,iatom)=sign*step_size*eigenvecs(ivec,imode)
55
 
            ivec = ivec + 1
56
 
          enddo ! ixyz
57
 
        enddo ! iatom
58
 
c
59
 
      call daxpy(nat3,1.0d00,steps,1,ncoords,1) ! mult coords 
60
 
      if (.not. geom_cart_coords_set(geom,ncoords))
61
 
     $     call errquit('raman_modestep: bad geom',0, GEOM_ERR)
62
 
c
63
 
      return
64
 
      end
65
 
c
66
 
c --------------------------------------------------------------------
67
 
c    Finite differencing for Raman Scattering Calculations
68
 
c --------------------------------------------------------------------
69
 
c
70
 
      subroutine fd_raman(rtdb,imode,rmmodes,natom,nat3,alfarex2,
71
 
     &           alfarimx2,step_size,rminfo,eigenvecs,mass)
72
 
c
73
 
      implicit none
74
 
c
75
 
#include "errquit.fh"
76
 
#include "global.fh"
77
 
#include "mafdecls.fh"
78
 
#include "stdio.fh"
79
 
#include "rtdb.fh"
80
 
#include "util.fh"
81
 
c
82
 
      logical debug
83
 
c
84
 
      integer rmmodes !  # modes used in raman calculation (3N-6)
85
 
      integer  i,j,m,ii,jj,mm,ivec  ! counting indexes
86
 
      integer rtdb    ! [input] rtdb handle
87
 
      integer imode   ! mode #
88
 
      integer natom   ! [input] number of atoms
89
 
      integer nat3    ! [input] 3*number of atoms
90
 
c
91
 
      double precision rminfo(rmmodes,4) ! raman data
92
 
      double precision step_size,stepsize ! [input] step of finite differencing
93
 
c      double precision laser(nfreq) ! [input] frequency used in AORESPONSE 
94
 
      double precision eigenvecs(nat3,nat3) ! [input](xyz&atom,mode)
95
 
      double precision tmode(3,natom) ! [input](atom#,xyz)
96
 
      double precision mass(natom) ! [input](atom#)
97
 
      double precision norm ! constant
98
 
      double precision fdipol(3,3), fdrpol(3,3) ! finite difference 
99
 
      double precision alfarex2(6,3), alfarimx2(6,3) ! AORESPONSE date for plus and minus step
100
 
      double precision ar2, ai2, gr2, gi2 ! alpha and g for imode
101
 
      double precision zero, one, two, three, six, seven,
102
 
     &                  fourty_five, bohr2ang
103
 
      parameter (zero=0.0D+00, one=1.0D+00, two=2.0D+00, three=3.0D+00,
104
 
     &            six=6.0D+00, seven=7.0D+00, fourty_five=45.0D+00)
105
 
      parameter (bohr2ang=0.52917724924D+00) ! CONVERSION OF BOHR to ANGSTROMS
106
 
c
107
 
      debug =  ( .false. .and. ga_nodeid().eq.0 )
108
 
c --------------------------------------------------------------------  
109
 
      call dfill(9,0.0D+00,fdrpol,1)  ! real fdipole polariability for 
110
 
      call dfill(9,0.0D+00,fdipol,1) ! imaginary fdipole polariability for
111
 
      call dfill(3*natom,0.0D+00,tmode,1) ! 
112
 
c zero
113
 
      stepsize = zero
114
 
      m = imode - 6 
115
 
      j=1
116
 
      i=1
117
 
      ar2   = zero ! alpha real
118
 
      gr2   = zero ! gradient real
119
 
      ai2   = zero ! alpha imaginaruy
120
 
      gi2   = zero ! gradient imaginary
121
 
      MM=3
122
 
c -----------mass weight the normal coordinates--------------------------------
123
 
       norm= 0.0D+00
124
 
       ivec = 1
125
 
      do i=1,natom
126
 
       do j=1,3
127
 
         tmode(j,i) = eigenvecs(ivec,imode)*sqrt(mass(i))
128
 
          ivec = ivec + 1
129
 
       end do
130
 
      end do
131
 
c ---calculate the norm of the qmass---
132
 
      do i=1,natom
133
 
         norm = norm + tmode(1,i)*tmode(1,i) + 
134
 
     *    tmode(2,i)*tmode(2,i) + tmode(3,i)*tmode(3,i)
135
 
      end do
136
 
      norm = sqrt(norm)
137
 
c --- normalize ---
138
 
      do i=1,natom
139
 
       do j=1,3
140
 
         tmode(j,i)= tmode(j,i)/norm
141
 
       end do
142
 
      end do
143
 
c
144
 
      do i=1,natom
145
 
       do j=1,3
146
 
         tmode(j,i) = tmode(j,i)/sqrt(mass(i))
147
 
       end do
148
 
      end do
149
 
c
150
 
      norm = zero
151
 
      do i=1,natom
152
 
         norm = norm + tmode(1,i)*tmode(1,i) + tmode(2,i)*tmode(2,i)
153
 
     *          + tmode(3,i)*tmode(3,i)
154
 
      end do
155
 
c --- calculate stepsize --- 
156
 
      stepsize = step_size/sqrt(norm)
157
 
c
158
 
      if (ga_nodeid().eq.0)  write(luout,*) 'Stepsize :',stepsize
159
 
c --------------------------------------------------------------------
160
 
      j=1
161
 
      i=1
162
 
      DO ii=1,3
163
 
       DO jj=1,3
164
 
c       difference polaraizability (real)
165
 
      fdrpol(ii,jj)     = alfarex2(ii,jj) - alfarex2(ii+3,jj)
166
 
c        difference polaraizability (imaginary)
167
 
      fdipol(ii,jj)     = alfarimx2(ii,jj) - alfarimx2(ii+3,jj)
168
 
c        convert units
169
 
       fdrpol(ii,jj) = ( fdrpol(ii,jj)*(bohr2ang**2) )/(two*stepsize)
170
 
       fdipol(ii,jj) = ( fdipol(ii,jj)*(bohr2ang**2) )/(two*stepsize)
171
 
       enddo
172
 
      enddo
173
 
c ----------calculate a and g for real and impaginary--------
174
 
      ar2 = ( ( (fdrpol(i,j)+fdrpol(i+1,j+1)+fdrpol(i+2,j+2))/three) 
175
 
     *       *( (fdrpol(i,j)+fdrpol(i+1,j+1)+fdrpol(i+2,j+2))/three) ) 
176
 
c
177
 
      gr2 = ( ( ((fdrpol(i,j)    - fdrpol(i+1,j+1))**2)
178
 
     *      +  ((fdrpol(i+1,j+1) - fdrpol(i+2,j+2))**2) 
179
 
     *      +  ((fdrpol(i+2,j+2) - fdrpol(i,j    ))**2) )
180
 
     *      + ( six*((fdrpol(i,j+1)   * fdrpol(i,j+1))
181
 
     *      +        (fdrpol(i+2,j)   * fdrpol(i+2,j))
182
 
     *      +        (fdrpol(i+1,j+2) * fdrpol(i+1,j+2))) )  )/two
183
 
c imaginary
184
 
      ai2 = ( ( (fdipol(i,j)+fdipol(i+1,j+1)+fdipol(i+2,j+2))/three)   
185
 
     *       *( (fdipol(i,j)+fdipol(i+1,j+1)+fdipol(i+2,j+2))/three) ) 
186
 
c
187
 
      gi2 = (  ( ((fdipol(i,j)     - fdipol(i+1,j+1))**2)
188
 
     *      +    ((fdipol(i+1,j+1) - fdipol(i+2,j+2))**2)
189
 
     *      +    ((fdipol(i+2,j+2) - fdipol(i,j    ))**2) )
190
 
     *      + ( six*((fdipol(i,j+1)   * fdipol(i,j+1))
191
 
     *      +        (fdipol(i+2,j)   * fdipol(i+2,j))
192
 
     *      +        (fdipol(i+1,j+2) * fdipol(i+1,j+2))) )  )/two
193
 
c --------------------------------------------------------------------
194
 
      if (debug) then
195
 
          write(luout,*) 'alfarimx2'
196
 
          call output(alfarimx2,1,6,1,3,6,3,1)
197
 
          write(luout,*) 'alfarex2'
198
 
          call output(alfarex2,1,6,1,3,6,3,1)
199
 
          write(luout,*) 'fdrpol'
200
 
          call output(fdrpol,1,3,1,3,3,3,1)
201
 
          write(luout,*) 'fdipol'
202
 
          call output(fdipol,1,3,1,3,3,3,1)
203
 
          write(luout,*) 'ar2=',ar2
204
 
          write(luout,*) 'gr2=',gr2
205
 
          write(luout,*) 'ai2=',ai2
206
 
          write(luout,*) 'gi2=',gi2
207
 
          write(luout,*) 'stepsize=',stepsize
208
 
      endif !  for DEBUG
209
 
      rminfo(m,2) = fourty_five*ar2 + seven*gr2 ! real component  
210
 
      rminfo(m,3) = fourty_five*ai2 + seven*gi2 ! imaginary component
211
 
      rminfo(m,4) = fourty_five*ar2 + seven*gr2 + fourty_five*ai2 +
212
 
     &               seven*gi2 ! total fdipole_polarizability_derivative
213
 
c------------ print result from finite difference ---------------------
214
 
      if (ga_nodeid().eq.0) write(luout,9001) rminfo(m,1),rminfo(m,2)
215
 
      if (ga_nodeid().eq.0) write(luout,9002) rminfo(m,1),rminfo(m,3)
216
 
      if (ga_nodeid().eq.0) write(luout,9003) rminfo(m,1),rminfo(m,4)
217
 
c --------------- save the rminfo for restart -------------------------
218
 
      if (.not. rtdb_put(rtdb, 'raman:rminfo ', mt_dbl,
219
 
     &    rmmodes*4, rminfo))
220
 
     &  call errquit('FD_raman:failed to write rminfo', 0, RTDB_ERR)
221
 
c --------------------------------------------------------------------
222
 
 9001 FORMAT(/1X,F7.2,2X,'REAL     ',2X,F24.6/)
223
 
 9002 FORMAT(/1X,F7.2,2X,'IMAGINARY',2X,F24.6/)
224
 
 9003 FORMAT(/1X,F7.2,2X,'TOTAL    ',2X,F24.6/)
225
 
      return
226
 
      end 
227
 
c
228
 
c --------------------------------------------------------------------
229
 
c   calculate absolute differential raman scattering cross section.
230
 
c   Producing a raman scattering output file 
231
 
c --------------------------------------------------------------------
232
 
c
233
 
      subroutine raman_scattering(rtdb,start,last,rmmodes,nfreq,plot,
234
 
     &                            line,width,laser,rminfo)
235
 
c
236
 
      implicit none
237
 
c
238
 
#include "global.fh"
239
 
#include "stdio.fh"
240
 
#include "rtdb.fh"
241
 
#include "inp.fh"
242
 
c
243
 
      character*16 plot, line
244
 
c
245
 
      integer rmmodes ! number or raman active modes
246
 
      integer rtdb    ! [input] rtdb handle
247
 
      integer nfreq   ! [input] number of excitation freqencies
248
 
      integer i,j   ! counting indexes 
249
 
      integer numpts ! number of points in plot
250
 
      integer last  ! [input] last mode to calculate for aoresponse
251
 
      integer begin ! first=begin if not restarting, else modified
252
 
      integer start ! [input] first mode to use in calculation of plot
253
 
c
254
 
      integer ifreq
255
 
      Double precision laser(nfreq) ! [input] laser used in AORESPONSE calculation
256
 
      Double precision width ! [input] line width
257
 
      Double precision rminfo(rmmodes,4) ! raman data
258
 
c      Double precision freq(rmmodes) ! [input] number of excitation frequencies used in AORESPONSE
259
 
      Double precision tempfc, frq4th 
260
 
      Double precision crs_real(rmmodes) ! real differential cross section
261
 
      Double precision  crs_imag(rmmodes) ! imaginary differential cross section
262
 
      Double precision crs_tot(rmmodes) ! total differential cross section
263
 
      Double precision frq, frqmin, frqmax, dfrq ! current, minimum, maximum and step of frequency
264
 
      Double precision afac, bfac, factor ! lineshape factors.
265
 
      Double precision conv_tot, conv_real, conv_imag !
266
 
      Double precision lambda ! exciation inverse wavenumber (1/cm)
267
 
      Double precision widthfac, pi, planck, boltz, speed, epsilon0 
268
 
      Double precision avogadro, nm2icm, temp, scale, amu,au2nm
269
 
      Double precision exparg, conver, zero, two, one, forty_five
270
 
c
271
 
      character*255 filename
272
 
      integer unitno
273
 
      parameter (unitno = 77)
274
 
c
275
 
      parameter (widthfac = 2.35482D+00)
276
 
      parameter (zero = 0.0D+00, one = 1.0D+00, 
277
 
     *           two = 2.0D+00, forty_five = 45.0D+00)
278
 
      parameter (pi = 3.14159265358979323846D+00)
279
 
      parameter (planck = 6.6260755D-34) ! plank's constant (J*sec)
280
 
      parameter (boltz = 1.3806580D-23) ! bolztmann's constant  (J/K)
281
 
      parameter (speed = 2.99792458D+08) ! speed of light (m/sec)
282
 
      parameter (epsilon0 = 8.8541878D-12) ! free permitivity (J^-1C^2m^-1)
283
 
      parameter (avogadro = 6.02214199D+23) ! avagadro's number (unitless)
284
 
      parameter (amu = (one/avogadro)*1.0D-03 ) ! atomic mass unit 
285
 
      parameter (au2nm =45.563353D+00)  ! converion factor from au to nm
286
 
      parameter (nm2icm = 1.0D+07) !converion factor from nm to 1/cm
287
 
      parameter (temp = 3.0D+02) ! tempature (k)
288
 
      parameter (scale = 1.0D+34)
289
 
      parameter (exparg  = (planck*speed*100.0D+00)/boltz )
290
 
      parameter (conver  = ( (two*planck*pi*pi)/speed*1.0D-40)/amu ) ! conver = 2*planck_h*pi^2/c * 10^40/atomic_mass_unit               
291
 
      parameter (numpts  = 1000)
292
 
c
293
 
      logical debug
294
 
      debug      = .false. ! produce debug output
295
 
c
296
 
c ----------outline and units of the calculation:---------------------
297
 
c d(sigma)/d(omega) = h/(8*(epsilon0)^2 *c) *(lambda-lambda_z)^4/lambda_z
298
 
c                     /(1-exp(-h*c*lambda_z/(kT))*(scattering factor)/45                    
299
 
c differential cross sections:
300
 
c - constant pre-factor for differential scattering cross sections:
301
 
c        h/(8*(epsilon0)^2*c) = 3.524100053d-21 J^3*s^2*m/C^4
302
 
c - lambda_t is given in [1/m^3]
303
 
c - squared polarizability derivatives have to be given in
304
 
c   [C^2*m^2/V^2/kg], but they are given as squared derivatives of
305
 
c   polarizability volumes in [Angstrom^4/amu].
306
 
c   conversion factor: (4*pi*epsilon0)^2*10^-40{Angstrom^4/m^4}
307
 
c                      /atomic_mass_unit{amu/kg}
308
 
c --------------------------------------------------------------------
309
 
c
310
 
      if (debug) then
311
 
          write(6,*) "rmmodes",rmmodes
312
 
          do ifreq = 1,nfreq
313
 
            write(6,*) "laser",laser(ifreq)
314
 
          end do
315
 
          write(6,*) 'rminfo'
316
 
          call output(rminfo,1,rmmodes,1,4,rmmodes,4,1)
317
 
      endif
318
 
c
319
 
c     --- write raman scattering data to file ---
320
 
c
321
 
      call util_file_name(plot,.false.,.false.,filename)
322
 
      open(unitno, status='unknown', form='formatted',file=filename)
323
 
 
324
 
      write(unitno,9002) plot
325
 
      do ifreq = 1,nfreq
326
 
         write(unitno,9003) au2nm/laser(ifreq)
327
 
      end do
328
 
      write(unitno,9004) line 
329
 
      write(unitno,9005) width
330
 
      write(unitno,9010) rmmodes
331
 
c
332
 
      if (plot =='normal') then
333
 
         write(unitno,9011)
334
 
         do i=1,rmmodes             
335
 
            write(unitno,9012) i, rminfo(i,1), rminfo(i,4)
336
 
         end do
337
 
            write(unitno,9006)
338
 
      else if (plot =='resonance') then
339
 
         write(unitno,9013)
340
 
         do i=1,rmmodes
341
 
            write(unitno,9014) i, rminfo(i,1), rminfo(i,4), 
342
 
     &                        rminfo(i,2), rminfo(i,3)
343
 
         end do
344
 
            write(unitno,9006)
345
 
      end if ! normal
346
 
c
347
 
c ----------convert laser from au. -> nm -> 1/cm ----------------------------
348
 
c
349
 
        lambda = nm2icm/(au2nm/laser(1))
350
 
c
351
 
c ------------Calculate differential cross section-----------------------
352
 
c
353
 
      do i=start,rmmodes
354
 
         crs_real(i) = zero
355
 
         crs_imag(i) = zero
356
 
         crs_tot(i)  = zero
357
 
         tempfc=1.0D+06/(one-exp(-exparg*rminfo(i,1)/temp))
358
 
         frq4th = (lambda-rminfo(i,1))**4
359
 
         crs_tot(i) = tempfc*frq4th*conver*rminfo(i,4)/
360
 
     *                     (forty_five*rminfo(i,1))
361
 
         if (plot == 'resonance') then
362
 
            crs_real(i) = tempfc*frq4th*conver*rminfo(i,2)/
363
 
     *                     (forty_five*rminfo(i,1))
364
 
            crs_imag(i) = tempfc*frq4th*conver*rminfo(i,3)/
365
 
     *                     (forty_five*rminfo(i,1))
366
 
         end if
367
 
      enddo ! rmmodes
368
 
c
369
 
         frqmin = rminfo(1,1) - 60.0D+00  ! determine plot minimum
370
 
         frqmax = rminfo(rmmodes,1) + 60.0D+00 ! determine plot maximum
371
 
         dfrq   = ((frqmax-frqmin)) / (numpts-1) !energy steps determined by numpts
372
 
c
373
 
c ------calculate lineshape factors for the plot-----------------------
374
 
c
375
 
      if (line == 'lorentzian') then
376
 
         afac = width/(two*pi)
377
 
         bfac = width/two
378
 
      end if
379
 
c      
380
 
      if (line == 'gaussian') then
381
 
         width = width/widthfac
382
 
         afac = one/(sqrt(two*pi)*width)
383
 
         bfac = one/(two*width**two)
384
 
      end if
385
 
c
386
 
c ----------- make the plot---------------------------------------
387
 
c
388
 
      do i=1,numpts ! loop over number of points
389
 
         conv_tot = zero
390
 
         conv_real = zero
391
 
         conv_imag = zero
392
 
         frq = frqmin + i*dfrq
393
 
         do j=start,rmmodes ! loop over number of raman active modes
394
 
            if (line == 'lorentzian') then  
395
 
             factor = scale*afac/((frq-rminfo(j,1))**2 + bfac**2)
396
 
            else if (line == 'gaussian') then  
397
 
             factor = scale*afac*exp(-bfac* (frq-rminfo(j,1))**2)
398
 
            endif
399
 
            conv_tot = conv_tot + factor*crs_tot(j)
400
 
            if (plot == 'resonance') then
401
 
               conv_real = conv_real + factor*crs_real(j)
402
 
               conv_imag = conv_imag + factor*crs_imag(j)
403
 
            end if
404
 
         end do  ! end loop over rmmodes 
405
 
c
406
 
         if (plot == 'normal' ) then
407
 
            write (unitno,9017) frq, conv_tot
408
 
         else if (plot =='resonance') then
409
 
            write(unitno,9018) frq, conv_tot, conv_real, conv_imag
410
 
         endif  ! plot
411
 
      end do ! end loop over numpts
412
 
c
413
 
      close(unitno) 
414
 
c
415
 
      write(luout,22) filename(1:inp_strlen(filename))
416
 
 22   format(/' Raman scattering data written to ',a/)
417
 
      call util_flush(luout)
418
 
c
419
 
c --------------------------------------------------------------------
420
 
 9002 FORMAT(/1X,60(1H-)/5X,A16,'Raman Scattering Plot'/1X,60(1H-))
421
 
 9003 FORMAT(1X,'Excitation wavelength',F8.2,1X,'nm')
422
 
 9004 FORMAT(1X,'Convolution lineshape is ',A12)
423
 
 9005 FORMAT(1X,'FWHM',F8.2,1X,' 1/cm'/1X,60(1H-))
424
 
 9006 FORMAT(1X,60(1H-))
425
 
 9010 FORMAT(1X,60(1H-)/1X,'# of frequencies : ',I5)
426
 
 9011 FORMAT(1X,60(1H-)/1X,'#  freq [ 1/cm]  S [Ang**4/amu]'/
427
 
     *        1X,60(1H-))
428
 
 9012 FORMAT(1X,I3,F9.2,F14.4)
429
 
 9013 FORMAT(1X,60(1H-)/'#  freq [ 1/cm]  S-tot [a.u.]',  
430
 
     * 'S-real [a.u.]  S-imag [a.u.]'/1X,60(1H-))
431
 
 9014 FORMAT(1X,I3,F9.2,3F14.4)
432
 
 9015 FORMAT(1X,I4,F10.2,2E12.4)
433
 
 9016 FORMAT(1X,F10.2,3E12.4)
434
 
 9017 FORMAT(1X,F8.2,G14.4)
435
 
 9018 FORMAT(1X,F8.2,3G14.4)
436
 
c
437
 
      return
438
 
      end  
439
 
c ----------------end raman_scattering---------------------------------------
440
 
c
441
 
      subroutine raman_save(rtdb,ii,junk3,junk4)
442
 
c
443
 
      implicit none
444
 
#include "errquit.fh"
445
 
#include "rtdb.fh"
446
 
#include "mafdecls.fh"
447
 
#include "stdio.fh"
448
 
#include "geom.fh"
449
 
#include "nwc_const.fh" 
450
 
#include "inp.fh"
451
 
#include "global.fh"
452
 
c
453
 
      integer i,j,n,ii ! counting indexes
454
 
      integer rtdb     ! [input] rtdb handle
455
 
      double precision alfare(3,3),alfaim(3,3) ! [aoresponse] real and imaginary polarizanbility
456
 
      double precision junk3(6,3),junk4(6,3) ! temp. arrays for FD of polarizanbility
457
 
      logical debug
458
 
      logical status
459
 
c
460
 
      debug      = (.false. .and. ga_nodeid().eq.0) ! produce debug output
461
 
      status = rtdb_parallel(.true.)
462
 
c  -------------get raman ploarizability from rtdb ----------------------------
463
 
        if (.not. rtdb_get(rtdb, 'raman:alfare ', mt_dbl, 9, alfare))
464
 
     &   call errquit('raman_save:failed to get alfare', 1, RTDB_ERR)
465
 
        if (.not. rtdb_get(rtdb, 'raman:alfaim ', mt_dbl, 9, alfaim)) 
466
 
     &   call errquit('raman_save:failed to get alfaim', 2, RTDB_ERR)
467
 
c
468
 
      n=(ii-1)*3
469
 
      DO i=1,3
470
 
         n=n+1
471
 
        DO j=1,3
472
 
          junk3(n,j)= alfare(i,j) ! add AORESPONSE real data to array for finite_differencing
473
 
          junk4(n,j)= alfaim(i,j) ! add AORESPONSE imaginary data to array for finite_differencing
474
 
        enddo ! j loop
475
 
       enddo ! i loop
476
 
c  -------------delete raman ploarizability from rtdb ----------------------------
477
 
      IF (debug) THEN
478
 
          write(luout,*) 'junk3'
479
 
          call output(junk3,1,6,1,3,6,3,1)
480
 
          write(luout,*) 'junk4'
481
 
          call output(junk4,1,6,1,3,6,3,1)
482
 
          write(luout,*) 'alfare'
483
 
          call output(alfare,1,3,1,3,3,3,1)
484
 
          write(luout,*) 'alfaim'
485
 
          call output(alfaim,1,3,1,3,3,3,1)
486
 
      ENDIF !  for DEBUG
487
 
c
488
 
      return
489
 
      end
490
 
c ----------------end raman_save------------------------------------------------
491
 
c $Id: raman.F 21176 2011-10-10 06:35:49Z d3y133 $