1* 2* $Id$ 3* 4 5 SUBROUTINE vib_RDINP( 6 & HESS,HESST,COORD,ATMASS,ZAN,OLDLBLAT,IZMAT,NELS, 7 & PROJEC,ZEROPE,HESOUT,INTERN, 8 & rtdb,first_pass) 9C 10C This routine reads the user input and the hessian matrix from 11C the appropriate place. The default is to use the "default" 12C masses and to read the analytic hessian from tape 10. 13C 14 IMPLICIT NONE ! REAL*8 (A-H,O-Z) 15#include "errquit.fh" 16#include "mafdecls.fh" 17#include "stdio.fh" 18#include "rtdb.fh" 19#include "geom.fh" 20#include "bas.fh" 21 INTEGER NELS 22 LOGICAL SETMAS,PROJEC,ZEROPE,HESOUT,INTERN 23 logical dum_log 24 integer natom, nat3, nhess, nhesst 25 logical first_pass 26 logical o_new_masses 27 DOUBLE PRECISION HESS(NAT3,NAT3) ! Hessian matrix 28 DOUBLE PRECISION HESST(NHESST) ! Lower triangle of hessian matrix 29 DOUBLE PRECISION COORD(3,NATOM) ! Atom x,y,z coordinates 30 DOUBLE PRECISION ATMASS(NATOM) ! Atoms mass 31 DOUBLE PRECISION ZAN(NATOM) ! Nuclear charge (i.e. atomic number) 32 INTEGER OLDLBLAT(2,NATOM) ! Atom character labels 33 INTEGER IZMAT(NELS) 34 integer num_new_masses 35 integer geom_index 36 double precision new_mass 37 character*32 mass_rtdb_id 38 character*16 mass_tag 39 character*16 tag 40 character*22 element 41 character*2 symbol 42 double precision q 43 double precision xyz(3) 44 COMMON /cvib_HESS/ NATOM,NAT3,NHESS,NHESST ! Hessian information 45 integer numans 46 double precision ams, wave 47 COMMON /cvib_SETCON/ AMS(36),WAVE,NUMANS ! setup parameters 48c 49 integer h_lblb, i_lblb, geom 50 integer il_cnt, iii, ijunk, jjj 51c 52 integer rtdb 53 HESOUT = .FALSE. ! do not write Hessian to tape 10 54 o_new_masses = .false. 55* 56 if (first_pass) WRITE(luout,*)' Vib: Default input used ' 57 PROJEC = .FALSE. 58 ZEROPE = .FALSE. 59 HESOUT = .FALSE. 60 INTERN = .FALSE. 61 if (rtdb_get(rtdb,'vib:project',MT_LOG,1,dum_log)) then 62 projec = dum_log 63 write(luout,*)' vib:project option set to',projec 64 endif 65 if (rtdb_get(rtdb,'vib:zero point energy',MT_LOG,1,dum_log)) then 66 zerope = dum_log 67 write(luout,*)' vib:zero point energy option set to',zerope 68 endif 69 if (zerope) projec = .true. 70* 71* over ride all input to do normal and projected analysis 72* 73 if (first_pass) then 74 projec = .false. 75 zerope = .false. 76 else 77 projec = .true. 78 zerope = .true. 79 endif 80*.. hessian read in vib_vib hess and hesst set there check it here 81 call vib_chkhess(hess,nat3,first_pass) 82* set labels, zan and coords information and mass information 83c 84 if (.not.ma_push_get(MT_BYTE,(2*natom),' labels for vib', 85 & h_lblb,i_lblb)) 86 & call errquit('vib_rdinp ma_get for labels failed',911, MA_ERR) 87 if (.not.geom_create(geom,'geometry')) call errquit 88 & ('vib_rdinp: error creating geometry',911, GEOM_ERR) 89 if (.not.geom_rtdb_load(rtdb,geom,'geometry')) call errquit 90 & ('vib_rdinp: error loading geometry',911, RTDB_ERR) 91c 92c see if ANY mass is overwritten by freq/vib input block 93c 94 if (rtdb_get 95 & (rtdb,'vib:remass:count',mt_int,1,num_new_masses)) then 96 o_new_masses = .true. 97 else 98 num_new_masses = 0 99 o_new_masses = .false. 100 endif 101 il_cnt = 0 102 do iii = 1,natom 103 if(.not.geom_cent_get(geom,iii,tag,xyz,q)) 104 & call errquit('vib_rdinp: geom_cent_get failed',911, 105 & GEOM_ERR) 106 zan(iii) = q 107 coord(1,iii) = xyz(1) 108 coord(2,iii) = xyz(2) 109 coord(3,iii) = xyz(3) 110 if(.not.geom_tag_to_element(tag,symbol,element,ijunk)) then 111 if (symbol.ne.'bq') 112 & call errquit('vib_rdinp: tag2elem fail',911, GEOM_ERR) 113 endif 114 byte_mb(il_cnt+i_lblb) = symbol(1:1) 115 il_cnt = il_cnt+1 116 byte_mb(il_cnt+i_lblb) = symbol(2:2) 117 il_cnt = il_cnt+1 118 enddo 119 120 121* 122* this loop structure separated out to do one pass over data in rtdb. 123* set the incore geometry with any mass information from the freq/vib 124* input. 125* 126 if (o_new_masses) then 127 do jjj = 1,num_new_masses 128 write(mass_rtdb_id,11111)jjj ! lexical index 129 if (rtdb_get(rtdb,mass_rtdb_id,mt_int,1,geom_index)) then 130 write(mass_rtdb_id,11113)jjj 131 if (.not.rtdb_get(rtdb,mass_rtdb_id,mt_dbl,1,new_mass)) 132 & call errquit 133 & ('vib_rdinp: rtdb get failed for mass',911, RTDB_ERR) 134 if (.not.geom_mass_set(geom,geom_index,new_mass)) 135 & call errquit 136 & ('vib_rdinp: geom_mass_set failed for mass',911, 137 & GEOM_ERR) 138 else ! check for tag input 139 write(mass_rtdb_id,11112)jjj 140 if (rtdb_cget(rtdb,mass_rtdb_id,1,mass_tag)) then 141 write(mass_rtdb_id,11113)jjj 142 if (.not.rtdb_get(rtdb,mass_rtdb_id,mt_dbl,1,new_mass)) 143 & call errquit 144 & ('vib_rdinp: rtdb get failed for mass',911, RTDB_ERR) 145* tag input found loop over all geometry atoms for a possible match 146 do iii = 1,natom 147 if(.not.geom_cent_get(geom,iii,tag,xyz,q)) 148 & call errquit('vib_rdinp: geom_cent_get failed',911, 149 & GEOM_ERR) 150 if (bas_do_tags_match(mass_tag,tag)) then 151 if (.not.geom_mass_set(geom,iii,new_mass)) 152 & call errquit 153 & ('vib_rdinp: geom_mass_set failed for mass',911, 154 & GEOM_ERR) 155 endif 156 enddo 157 else 158 if(.not.rtdb_print(rtdb,.true.)) call errquit('911',911, 159 & RTDB_ERR) 160 write(luout,*)' expected mass input for index ',jjj 161 write(luout,*)' It was neither index or tag based ' 162 call errquit 163 & ('vib_rdinp: fatal error for mass informaition',911, 164 & RTDB_ERR) 165 endif 166 endif 167 enddo 168 endif 169* 170* best place to set mass is in geometry object 171* modified in-core geometry object with new masses from 172* vibrational/freq input now reset's masses. 173* 174 setmas = geom_masses_get(geom,natom,atmass) 175 if (.not.setmas) 176 & call errquit 177 & ('vib_rdinp: geom_masses_get failed',911, GEOM_ERR) 178* 179 if (.not.geom_destroy(geom)) call errquit 180 & ('vib_rdinp: geom_destroy failed',911, GEOM_ERR) 181C 182C Write out atom information read here 183C 184 if (first_pass) then 185 WRITE(6,10000) 186C 187 il_cnt = 0 188 DO 00300 III = 1,NATOM 189 WRITE(6,10001) 190 & byte_mb(il_cnt+i_lblb), 191 & byte_mb(il_cnt+i_lblb+1), 192 & III, 193 & COORD(1,III),COORD(2,III),COORD(3,III), 194 & ATMASS(III) 195 il_cnt = il_cnt + 2 19600300 CONTINUE 197 WRITE(6,10002) 198 endif 199 if (.not.ma_pop_stack(h_lblb)) call errquit 200 & ('vib_rdinp ma_pop failed',911, MA_ERR) 201 RETURN ! leave routine 20210000 FORMAT(//,1X,28('-'),' Atom information ',28('-'),/, 203 +5X,'atom',4X,'#',8X,'X',14X,'Y',14X,'Z',12X,'mass',/, 204 +1X,74('-')) 20510001 FORMAT(4x,a1,a1,3x,I5,4(1PD15.7)) 20610002 FORMAT(1X,74('-')///) 20711111 format('vib:remass:',i4,':lexi:') 20811112 format('vib:remass:',i4,':tags:') 20911113 format('vib:remass:',i4,':mass:') 210C 211 END 212 213 214