1 subroutine getmdp(natoms,heat,izo,ido,istatz,cc,ianc, 2 & bl,alph,bet,ibl,ialph,ibet,imap,ianz,iz, 3 & c,cz,alpha,beta,ian) 4 5c this is really getmop 6 7 implicit double precision (a-h,o-z), integer ( i-n) 8 integer getlin 9 character*137 string 10 character*2 catom,catomt,tolowf,iel 11 logical first,ottest, oerror,zmstar,getnxt,gnreal 12 13 common /zmfrst/ ihaszm, nz, mxzat 14 15 dimension ian(*),c(3,*),cz(3,*),alpha(*),beta(*) 16 dimension ianc(*),cc(3,*) 17 dimension bl(*),alph(*),bet(*),ibl(*),ialph(*),ibet(*), 18 & imap(*),ianz(*),iz(4,*) 19 dimension izt(4),iztt(3,3),zzz(3),idm(3),tmp(3) 20 dimension iel(100) 21 character*137 coplin 22 character*137 line 23 common /curlin/ line 24 common /rdwr/ iun1,iun2,iun3,iun4,iun5 25 data iel/'bq', 26 & 'h ', 'he', 27 & 'li', 'be', 'b ', 'c ', 'n ', 'o ', 'f ', 'ne', 28 & 'na', 'mg', 'al', 'si', 'p ', 's ', 'cl', 'ar', 29 & 'k ', 'ca', 30 & 'sc', 'ti', 'v ', 'cr', 'mn', 31 & 'fe', 'co', 'ni', 'cu', 'zn', 32 & 'ga', 'ge', 'as', 'se', 'br', 'kr', 33 & 'rb','sr','y ','zr','nb','mo','tc','ru','rh','pd','ag','cd', 34 & 'in','sn','sb','te','i ','xe','cs','ba','la','ce','pr','nd', 35 & 'pm','sm','eu','gd','tb','dy','ho','er','tm','yb','lu','hf', 36 & 'ta','w ','re','os','ir','pt','au','hg','tl','pb','bi','po', 37 & 'at','rn','fr','ra','ac','th','pa','u ','np','pu','am','cm', 38 & 'bk','cf','x '/ 39 40 istatz = 0 41 ifnd = 0 42 irwd = 0 43 toang = 0.52917706d0 44 455 continue 46 47 first = .true. 48 if (izo.eq.1) then 49 nzz = nz 50 else 51 nz = 0 52 endif 53 maxnz = mxzat 54 ottest=.true. 55 heat = 0.0d0 56 5710 continue 58 zmstar = .false. 59 nztt = 0 60 do while (getlin(0).eq.1) 61 if (icdex(line,'[AMB').ne.0) return 62 i1 = index(line,'(') 63 i2 = index(line,')') 64 if (i1.ne.0.and.i2.ne.0.and.i2.gt.i1) then 65 line = line(1:i1-1)//line(i2+1:) 66 endif 67 if (index(line,'VARIABLE FUNCTION').ne.0) then 68 nword = 3 69 if (index(line,'SECOND').ne.0) nword = 4 70 if (getlin(0).eq.0) goto 10 71 do ii=1,nword 72 ktype = nxtwrd(string,nstr,itype,rtype) 73 end do 74 if (ktype.eq.3) then 75 heat = rtype 76 endif 77 goto 10 78 endif 79 if (index(line,'HEAT OF FORMATION').ne.0) then 80 do while(.true.) 81 ktype = nxtwrd(string,nstr,itype,rtype) 82 if (ktype.eq.3) heat = rtype 83 if (ktype.eq.0.or.ktype.eq.3) goto 10 84 end do 85 endif 86 if (index(line,'POINT POTENTIAL').ne.0.or. 87 & index(line,'POINT POTENTIAL').ne.0) then 88 in = 2 89 if (index(line,'FEMTOSECONDS').ne.0) in = 3 90 if (getlin(0).eq.0) goto 10 91c do ii=1,4 92c Armando Juan Navarro Vazquez says total energy not interresting 93c Use potential energy 94c 95 do ii=1,in 96 ktype = nxtwrd(string,nstr,itype,rtype) 97 end do 98 if (ktype.eq.3) heat = rtype 99 goto 10 100 endif 101 if (index(line,'FINAL GEOMETRY OBTAINED').ne.0) then 102 call redel(line,4) 103 coplin = line 104 iws = 0 105 ktype = -1 106 do while (ktype.ne.0) 107 ktype = nxtwrd(string,nstr,itype,rtype) 108 if (ktype.ne.0) iws = iws + 1 109 end do 110 line = coplin 111 if (iws.eq.7) then 112 natoms = 0 113 do while(.true.) 114 ktype = nxtwrd(string,nstr,itype,rtype) 115 if (ktype.eq.1.and.(nstr.eq.1.or.nstr.eq.2)) then 116 if (nstr.eq.1) then 117 catomt(1:1) = string(1:1) 118 catomt(2:2) = ' ' 119 else 120 catomt = string(1:2) 121 endif 122 else 123 if (natoms.ne.0) then 124 istatz = 2 125 return 126 endif 127 goto 100 128 endif 129 130 iatmp = 0 131 catom = tolowf(catomt) 132 do i=1,100 133 if ( catom .eq. iel(i) ) iatmp = i - 1 134 end do 135 if (catom.eq.'xx') iatmp = 99 136 if (catom.eq.'tv') iatmp = 99 137 if (iatmp.eq.0.or.iatmp.gt.100) goto 100 138 139 natoms = natoms + 1 140 ianc(natoms) = iatmp 141 142 do ii=1,3 143 ktype = nxtwrd(string,nstr,itype,rtype) 144 if (ktype.eq.3) then 145 cc(ii,natoms) = rtype / toang 146 else 147 goto 100 148 endif 149 ktype = nxtwrd(string,nstr,itype,rtype) 150 if (ktype.ne.2) goto 100 151 end do 152 if (getlin(0).ne.1) return 153 end do 154 endif 155 endif 156 if (index(line,'scf done: ').ne.0) then 157 do while(.true.) 158 ktype = nxtwrd(string,nstr,itype,rtype) 159 if (ktype.eq.3) heat = rtype 160 if (ktype.eq.0.or.ktype.eq.3) goto 10 161 end do 162 endif 163 do i=1,3 164 idm(i) = 0 165 zzz(i) = 0.0d0 166 izt(i) = 0 167 end do 168 ktype = nxtwrd(string,nstr,itype,rtype) 169 if (ktype.eq.2) then 170 if (itype.eq.1.and.first) zmstar = .true. 171 if (zmstar) then 172 ktype = nxtwrd(string,nstr,itype,rtype) 173 else 174 goto 100 175 endif 176 endif 177c element label 178 if (ktype.eq.1.and.(nstr.eq.1.or.nstr.eq.2)) then 179 if (nstr.eq.1) then 180 catomt(1:1) = string(1:1) 181 catomt(2:2) = ' ' 182 else 183 catomt = string(1:2) 184 endif 185 else 186 goto 100 187 endif 188 189c check which zmstar type 190 191 if (zmstar.and.nztt.eq.0) then 192 izmtyp = 0 193 if (gnreal(tmp,3,.false.)) izmtyp = 1 194 endif 195 196 num = nztt 197 if (num.gt.3) num = 3 198 if (.not.zmstar.or.(zmstar.and.izmtyp.eq.1)) num = 3 199 200 getnxt = .true. 201 do i=1,num 202 if (getnxt) ktype = nxtwrd(string,nstr,itype,rtype) 203 204c bl,aplh,bet 205 206 getnxt = .true. 207 if (ktype.eq.2.or.ktype.eq.3) then 208 if (ktype.eq.2) then 209 zzz(i) = 1.0d0*itype 210 else 211 zzz(i) = rtype 212 endif 213 else 214 if (zmstar) then 215 if (ktype.eq.0) then 216c read in connectivity 217 goto 200 218 else 219 goto 100 220 endif 221 else 222 goto 100 223 endif 224 endif 225 226c optimize flag 1-3 227 228 ktype = nxtwrd(string,nstr,itype,rtype) 229 if (zmstar) then 230 if (ktype.eq.1) then 231 if (nstr.eq.1) then 232 if (string(1:1).eq.'*'.or. 233 & string(1:1).eq.'+') then 234 idm(i) = 1 235 else 236 goto 100 237 endif 238 else 239 goto 100 240 endif 241 else 242 getnxt = .false. 243 idm(i) = 0 244 endif 245 else 246 if (ktype.eq.2) then 247 idm(i) = abs(itype) 248 else 249 goto 100 250 endif 251 endif 252 end do 253200 continue 254c 255c connectivity 256c 257 num = nztt 258 if (num.gt.3) num = 3 259 if (.not.zmstar) num = 3 260 261 do i=1,num 262 if (getnxt) ktype = nxtwrd(string,nstr,itype,rtype) 263 getnxt = .true. 264 if (ktype.eq.2) then 265 izt(i) = itype 266 endif 267 end do 268 269c check first three atoms 270 271 if (izo.ne.1) then 272 ioke = 1 273 if (nztt.le.2) then 274 do i=1,3 275 iztt(i,nztt+1) = izt(i) 276 end do 277 endif 278 if (nztt.le.2) then 279 if (iztt(1,1).ne.0) ioke = 0 280 if (iztt(2,1).ne.0) ioke = 0 281 if (iztt(3,1).ne.0) ioke = 0 282 endif 283 if (nztt.ge.1.and.nztt.le.2) then 284 if (iztt(1,2).le.0) ioke = 0 285 if (iztt(2,2).ne.0) ioke = 0 286 if (iztt(3,2).ne.0) ioke = 0 287 endif 288 if (nztt.eq.2) then 289 if (iztt(1,3).le.0) ioke = 0 290 if (iztt(2,3).le.0) ioke = 0 291 if (iztt(3,3).ne.0) ioke = 0 292 endif 293 if (ioke.ne.1) then 294 istatz = 0 295 first = .true. 296 goto 5 297 endif 298 endif 299 300c check atom4 and upwards 301 302 iztot = 0 303 do i=1,3 304 iztot = iztot + izt(i) 305 end do 306 if (nztt.gt.4.and.iztot.eq.0) then 307 istatz = 0 308 first = .true. 309 goto 5 310 endif 311 312300 continue 313 if (zmstar.and.heat.eq.0.0d0.and.ifnd.eq.0) then 314 ifnd = 1 315 goto 100 316 endif 317c 318c copy temporary variables into z-mat 319c 320 iatmp = 0 321 catom = tolowf(catomt) 322 do i=1,100 323 if ( catom .eq. iel(i) ) iatmp = i - 1 324 end do 325 if (catom.eq.'xx') iatmp = 99 326 if (catom.eq.'tv') iatmp = 99 327 if (iatmp.eq.0.or.iatmp.gt.100) goto 100 328 first = .false. 329 if (nztt.ge.1) istatz = 1 330 nz = nz + 1 331 nztt = nztt + 1 332 ianz(nz) = iatmp 333 do i=1,3 334 iz(i,nz) = izt(i) 335 end do 336 337 iz(4,nz) = 0 338 bl(nz) = zzz(1) 339 ibl(nz) = idm(1) 340 alph(nz) = zzz(2) 341 ialph(nz) = idm(2) 342 bet(nz) = zzz(3) 343 ibet(nz) = idm(3) 344 end do 345 first = .false. 346 347100 if (first) goto 10 348 if (izo.eq.1) then 349 if (istatz.eq.0) nz = nzz 350 return 351 endif 352 if (istatz.eq.1) then 353c call prtzm 354 ihaszm = 1 355 call stocc(maxnz,nz,iz,bl,alph,bet,alpha,beta,ierror) 356 if (ierror.lt.5) then 357 call stoc(maxnz,nz,0,0,0,ianz,iz,bl,alph,bet, 358 & ottest,natoms,ian,c,cz,imap,alpha,beta, 359 & oerror,.true.,.true.) 360 if (oerror) then 361 istatz = 0 362 return 363 endif 364 else 365 istatz = 0 366 return 367 endif 368 do i=1,natoms 369 do j=1,3 370 cc(j,i) = c(j,i) 371 end do 372 ianc(i) = ian(i) 373 end do 374 else 375 if (ifnd.eq.1.and.irwd.eq.0.and.ido.eq.1) then 376 irwd = 1 377 call rewfil 378 goto 5 379 endif 380 endif 381 382 return 383 end 384