From f6b5dd768a7aed80adf9cdfddabc0142a99cfdea Mon Sep 17 00:00:00 2001 From: Tjerk Straatsma Date: Mon, 20 Apr 1998 20:01:41 +0000 Subject: [PATCH] See ChangeLog --- src/prepar/ChangeLog | 2 + src/prepar/pre_atype.F | 133 +++++++++++++++++++++++++++++----------- src/prepar/pre_center.F | 123 +++++++++++++++++++++++-------------- src/prepar/pre_charge.F | 56 ++++++++++++++++- src/prepar/pre_dihimp.F | 30 +++++++-- src/prepar/pre_ffield.F | 38 ++++++++++-- src/prepar/pre_imp.F | 6 +- src/prepar/pre_mkfrg.F | 21 +++++-- src/prepar/pre_mklist.F | 4 +- src/prepar/pre_mkseq.F | 12 ++-- src/prepar/pre_mksgm.F | 3 +- src/prepar/pre_rdsgm.F | 10 +-- src/prepar/pre_wrttop.F | 12 ++-- src/prepar/prepar.F | 18 +++--- 14 files changed, 335 insertions(+), 133 deletions(-) diff --git a/src/prepar/ChangeLog b/src/prepar/ChangeLog index c0e0a3159b..0c68e9b76f 100644 --- a/src/prepar/ChangeLog +++ b/src/prepar/ChangeLog @@ -8,3 +8,5 @@ Apr 16 1998 Handle solvent Apr 17 1998 Generalize ring identification to 3,4,5 and 6 membered rings Modified code to retrieve improper dihedral parameters Allow force field atom type file in dir_s, dir_x, dir_u and dir_t +Apr 20 1998 Consolidate par and typ files + Corrections in improper dihedral center determination \ No newline at end of file diff --git a/src/prepar/pre_atype.F b/src/prepar/pre_atype.F index 782305412c..0de7be172e 100644 --- a/src/prepar/pre_atype.F +++ b/src/prepar/pre_atype.F @@ -1,11 +1,11 @@ - logical function pre_atype(lfnout,lprint,lfntyp,filtyp, + logical function pre_atype(lfnout,lprint,lfnpar,filpar, + dir_s,dir_x,dir_u,dir_t,latm,catm,matm,lbnd,mbnd,nbnd, + jlo,ilo,ihi,jhi,ltyp,mtyp,lring,aring,mring, + nring3,nring4,nring5,nring6) c implicit none c -c in : lfntyp = atom types file number +c in : lfnpar = atom types file number c filtyp = atom types file name c latm(2,matm) = atomic number c 3 = number of bonds @@ -31,9 +31,9 @@ c integer maty parameter (maty=1000) c - integer lfntyp,lfnout,lprint - character*255 filtyp - character*255 filpar,dir_s,dir_x,dir_u,dir_t + integer lfnpar,lfnout,lprint + character*255 filpar + character*255 dir_s,dir_x,dir_u,dir_t character*255 filnam integer matm,mtyp integer latm(5,matm) @@ -50,7 +50,7 @@ c character*4 atype(maty) integer itype(20,maty) character*100 card - integer iatnum,nhydr,nonh + integer iatnum,nhydr,nonh,irng logical skipa(5),skipt(5),skipaa(5),skiptt(5) c c setup typ array with latm data @@ -190,78 +190,97 @@ c read the atom type data c ----------------------- c ntype=0 - lt=index(filtyp,' ')-1 + lt=index(filpar,' ')-1 ld=index(dir_s,' ')-1 l=ld+lt - filnam=dir_s(1:ld)//filtyp(1:lt) - open(unit=lfntyp,file=filnam(1:l),status='old',form='formatted', + filnam=dir_s(1:ld)//filpar(1:lt) + open(unit=lfnpar,file=filnam(1:l),status='old',form='formatted', + err=901) - write(*,'(a,a)') 'TYPES FROM ',filnam(1:l) + 18 continue + read(lfnpar,1000,end=17,err=9999) card + if(card(1:10).ne.'Atom types') goto 18 + if(lprint.ge.2) + + write(*,'(a,a)') 'SATISFYING ATOM TYPES USING ',filnam(1:l) 16 continue - read(lfntyp,1000,end=17,err=9999) card + read(lfnpar,1000,end=17,err=9999) card 1000 format(a) if(card(1:1).eq.'#') goto 16 ntype=ntype+1 if(ntype.gt.maty) call errquit('increase maty',9999) - read(card,1001) atype(ntype),(itype(i,ntype),i=1,20) - 1001 format(a4,i4,i3,2i4,i3,3(i4,i3,3i4)) + read(card,1001) atype(ntype),(itype(i,ntype),i=1,10) + 1001 format(a4,i7,i3,2i5,i3,i7,i3,3i7) + read(lfnpar,1002) (itype(i,ntype),i=11,15) + read(lfnpar,1002) (itype(i,ntype),i=16,20) + 1002 format(27x,i7,i3,3i7) goto 16 17 continue - close(unit=lfntyp) + close(unit=lfnpar) c 901 continue - lt=index(filtyp,' ')-1 + lt=index(filpar,' ')-1 ld=index(dir_x,' ')-1 l=ld+lt - filnam=dir_x(1:ld)//filtyp(1:lt) - open(unit=lfntyp,file=filnam(1:l),status='old',form='formatted', + filnam=dir_x(1:ld)//filpar(1:lt) + open(unit=lfnpar,file=filnam(1:l),status='old',form='formatted', + err=902) - write(*,'(a,a)') 'TYPES FROM ',filnam(1:l) + 118 continue + read(lfnpar,1000,end=117,err=9999) card + if(card(1:10).ne.'Atom types') goto 118 + if(lprint.ge.2) + + write(*,'(a,a)') 'SATISFYING ATOM TYPES USING ',filnam(1:l) 116 continue - read(lfntyp,1000,end=117,err=9999) card + read(lfnpar,1000,end=117,err=9999) card if(card(1:1).eq.'#') goto 116 ntype=ntype+1 if(ntype.gt.maty) call errquit('increase maty',9999) read(card,1001) atype(ntype),(itype(i,ntype),i=1,20) goto 116 117 continue - close(unit=lfntyp) + close(unit=lfnpar) c 902 continue - lt=index(filtyp,' ')-1 + lt=index(filpar,' ')-1 ld=index(dir_u,' ')-1 l=ld+lt - filnam=dir_u(1:ld)//filtyp(1:lt) - open(unit=lfntyp,file=filnam(1:l),status='old',form='formatted', + filnam=dir_u(1:ld)//filpar(1:lt) + open(unit=lfnpar,file=filnam(1:l),status='old',form='formatted', + err=903) - write(*,'(a,a)') 'TYPES FROM ',filnam(1:l) + 218 continue + read(lfnpar,1000,end=217,err=9999) card + if(card(1:10).ne.'Atom types') goto 218 + if(lprint.ge.2) + + write(*,'(a,a)') 'SATISFYING ATOM TYPES USING ',filnam(1:l) 216 continue - read(lfntyp,1000,end=217,err=9999) card + read(lfnpar,1000,end=217,err=9999) card if(card(1:1).eq.'#') goto 216 ntype=ntype+1 if(ntype.gt.maty) call errquit('increase maty',9999) read(card,1001) atype(ntype),(itype(i,ntype),i=1,20) goto 216 217 continue - close(unit=lfntyp) + close(unit=lfnpar) c 903 continue - lt=index(filtyp,' ')-1 + lt=index(filpar,' ')-1 ld=index(dir_t,' ')-1 l=ld+lt - filnam=dir_t(1:ld)//filtyp(1:lt) - open(unit=lfntyp,file=filnam(1:l),status='old',form='formatted', + filnam=dir_t(1:ld)//filpar(1:lt) + open(unit=lfnpar,file=filnam(1:l),status='old',form='formatted', + err=904) - write(*,'(a,a)') 'TYPES FROM ',filnam(1:l) + 318 continue + read(lfnpar,1000,end=317,err=9999) card + if(card(1:10).ne.'Atom types') goto 318 + if(lprint.ge.2) + + write(*,'(a,a)') 'SATISFYING ATOM TYPES USING ',filnam(1:l) 316 continue - read(lfntyp,1000,end=317,err=9999) card + read(lfnpar,1000,end=317,err=9999) card if(card(1:1).eq.'#') goto 316 ntype=ntype+1 if(ntype.gt.maty) call errquit('increase maty',9999) read(card,1001) atype(ntype),(itype(i,ntype),i=1,20) goto 316 317 continue - close(unit=lfntyp) + close(unit=lfnpar) 904 continue c c determine the atom types @@ -270,12 +289,21 @@ c if(lprint.ge.5) +write(lfnout,'(a4,5i5)') (catm(2,i),(ltyp(j,i),j=1,5),i=ilo,ihi) c - do 18 i=ilo,ihi + do 38 i=ilo,ihi if(lprint.ge.5) + write(lfnout,'(a,a4,5i5)') 'Atom ',catm(2,i),(ltyp(j,i),j=1,5) if(latm(2,i).gt.0) then do 19 j=1,ntype iatnum=itype(1,j) + irng=0 + if(iatnum.ge.60000) then + irng=6 + iatnum=iatnum-60000 + endif + if(iatnum.ge.50000) then + irng=5 + iatnum=iatnum-50000 + endif nonh=0 do 20 k=1,4 if(iatnum.ge.1000) then @@ -299,7 +327,7 @@ c 4 match aromatic ring c 5 match explicit number of neighbors c if(iatnum.eq.latm(2,i).and. - + (nhydr.eq.0.or.(nhydr.eq.5.and.ltyp(4,j).eq.ltyp(5,j)).or. + + (nhydr.eq.0.or.(nhydr.eq.5.and.ltyp(4,i).eq.0).or. + nhydr.eq.ltyp(4,j)).and. + (itype(2,j).eq.0.or.itype(2,j).eq.ltyp(1,i)).and. + (itype(3,j).eq.0.or.itype(3,j).eq.ltyp(2,i).or. @@ -338,6 +366,15 @@ c + write(lfnout,'(5x,a,i5)') 'neighbor atom ', + latm(2,ltyp(5+k,i)) iatnum=itype(1+l*5,j) + irng=0 + if(iatnum.ge.60000) then + irng=6 + iatnum=iatnum-60000 + endif + if(iatnum.ge.50000) then + irng=5 + iatnum=iatnum-50000 + endif nonh=0 do 26 m=1,4 if(iatnum.ge.1000) then @@ -382,6 +419,13 @@ c ------------------------- c if(itype(2+l*5,j).gt.0.and. + itype(2+l*5,j).ne.ltyp(5,ltyp(5+k,i))) goto 25 +c +c check if neighbor in ring +c ------------------------- +c + if(irng.gt.0) then + if(irng.ne.ltyp(3,ltyp(5+k,i))) goto 25 + endif c if(lprint.ge.5) + write(lfnout,'(5x,a,2i5)') 'neighbor accepted ', @@ -416,6 +460,15 @@ c + latm(2,ltyp(5+kk,ltyp(5+k,i))) c iatnum=itype(2+l*5+ll,j) + irng=0 + if(iatnum.ge.60000) then + irng=6 + iatnum=iatnum-60000 + endif + if(iatnum.ge.50000) then + irng=5 + iatnum=iatnum-50000 + endif nonh=0 do 32 m=1,4 if(iatnum.ge.1000) then @@ -464,6 +517,14 @@ c if(nonh.ne.ltyp(5,ltyp(5+k,ltyp(5+k,i)))- + ltyp(4,ltyp(5+k,ltyp(5+k,i)))) goto 31 endif +c +c check if neighbor of neighbor in ring +c ------------------------------------- +c + if(irng.gt.0) then + if(irng.ne.ltyp(3,ltyp(5+kk,ltyp(5+k,i)))) goto 31 + endif +c if(lprint.ge.5) + write(lfnout,'(10x,a,4i5)') 'neighbor accepted ', + itype(2+l*5+ll,j),latm(2,ltyp(5+kk,ltyp(5+k,i))),kk,ll @@ -511,9 +572,9 @@ c endif 19 continue endif - 18 continue + 38 continue c - if(lprint.ge.5) then + if(lprint.ge.50) then do 100 i=jlo,jhi write(lfnout,'(16i5)') i,(ltyp(j,i),j=1,15) 100 continue diff --git a/src/prepar/pre_center.F b/src/prepar/pre_center.F index cfc007b914..47d0118cb9 100644 --- a/src/prepar/pre_center.F +++ b/src/prepar/pre_center.F @@ -40,55 +40,14 @@ c integer nring3,nring4,nring5,nring6 logical valid c - integer i,j,k,l,li,lj + integer i,j,k,l,li,lj,nh,ic real*8 angle c pre_center=.false. c -c peptide bonds -c ------------- -c - do 1 i=jlo,jhi - if(latm(2,i).eq.7.and.latm(3,i).eq.3) then - do 2 li=1,nbnd - j=0 - if(lbnd(1,li).eq.i) j=lbnd(2,li) - if(lbnd(2,li).eq.i) j=lbnd(1,li) - if(j.ge.jlo.and.j.le.jhi) then - if(latm(2,j).eq.6.and.latm(3,j).eq.3) then - do 3 lj=1,nbnd - k=0 - if(lbnd(1,lj).eq.j.and.lbnd(2,lj).ne.i) k=lbnd(2,lj) - if(lbnd(2,lj).eq.j.and.lbnd(1,lj).ne.i) k=lbnd(1,lj) - if(k.ge.jlo.and.k.le.jhi) then - if(latm(2,k).eq.8.and.latm(3,k).eq.1) then - latm(4,i)=1 - latm(4,j)=1 - goto 1 - endif - endif - 3 continue - endif - endif - 2 continue - endif - 1 continue -c -c planar C=O -c ---------- -c - do 26 i=jlo,jhi - if(latm(2,i).eq.6.and.latm(3,i).eq.3) then - do 27 li=1,nbnd - j=0 - if(lbnd(1,li).eq.i) j=lbnd(2,li) - if(lbnd(2,li).eq.i) j=lbnd(1,li) - if(j.ge.jlo.and.j.le.jhi) then - if(latm(2,j).eq.8.and.latm(3,j).eq.1) latm(4,i)=1 - endif - 27 continue - endif - 26 continue + do 53 i=jlo,jhi + latm(4,i)=0 + 53 continue c c 3 membered rings c ---------------- @@ -322,6 +281,80 @@ c 25 continue endif 24 continue +c +c peptide bonds +c ------------- +c + do 1 i=jlo,jhi + if(latm(2,i).eq.7.and.latm(3,i).eq.3) then + do 2 li=1,nbnd + j=0 + if(lbnd(1,li).eq.i) j=lbnd(2,li) + if(lbnd(2,li).eq.i) j=lbnd(1,li) + if(j.ge.jlo.and.j.le.jhi) then + if(latm(2,j).eq.6.and.latm(3,j).eq.3) then + do 3 lj=1,nbnd + k=0 + if(lbnd(1,lj).eq.j.and.lbnd(2,lj).ne.i) k=lbnd(2,lj) + if(lbnd(2,lj).eq.j.and.lbnd(1,lj).ne.i) k=lbnd(1,lj) + if(k.ge.jlo.and.k.le.jhi) then + if(latm(2,k).eq.8.and.latm(3,k).eq.1) then + latm(4,i)=1 + latm(4,j)=1 + goto 1 + endif + endif + 3 continue + endif + endif + 2 continue + endif + 1 continue +c +c planar C=O +c ---------- +c + do 26 i=jlo,jhi + if(latm(2,i).eq.6.and.latm(3,i).eq.3.and.latm(4,i).eq.0) then + do 27 li=1,nbnd + j=0 + if(lbnd(1,li).eq.i) j=lbnd(2,li) + if(lbnd(2,li).eq.i) j=lbnd(1,li) + if(j.ge.jlo.and.j.le.jhi) then + if(latm(2,j).eq.8.and.latm(3,j).eq.1) latm(4,i)=1 + endif + 27 continue + endif + 26 continue +c +c planar NH2 in C(sp2)-NH2 +c ------------------------ +c + do 45 i=jlo,jhi + if(latm(2,i).eq.7.and.latm(3,i).eq.3.and.latm(4,i).eq.0) then + nh=0 + ic=0 + do 46 li=1,nbnd + j=0 + if(lbnd(1,li).eq.i) j=lbnd(2,li) + if(lbnd(2,li).eq.i) j=lbnd(1,li) + if(j.ge.jlo.and.j.le.jhi) then + if(latm(2,j).eq.1) nh=nh+1 + if(latm(2,j).eq.6.and.latm(3,j).eq.3) ic=j + endif + 46 continue + if((nh.eq.1.or.nh.eq.2).and.ic.gt.0) latm(4,i)=1 + endif + 45 continue +c +c sp2 carbon +c ---------- +c + do 47 i=jlo,jhi + if(latm(2,i).eq.6.and.latm(3,i).eq.3.and.latm(4,i).eq.0) then + latm(4,i)=1 + endif + 47 continue c pre_center=.true. return diff --git a/src/prepar/pre_charge.F b/src/prepar/pre_charge.F index 9159bb0f86..8dc8711c40 100644 --- a/src/prepar/pre_charge.F +++ b/src/prepar/pre_charge.F @@ -1,4 +1,4 @@ - logical function pre_charge(latm,qatm,matm,lbnd,mbnd,nbnd) + logical function pre_charge(latm,qatm,matm,natm,lbnd,mbnd,nbnd) c implicit none c @@ -10,13 +10,13 @@ c nbnd = length bond list c c out : qatm(matm) = guestimated atomic partial charges c - integer matm + integer matm,natm integer latm(5,matm) real*8 qatm(matm) integer mbnd,nbnd integer lbnd(2,mbnd) c - integer i,ia,ja,it,jt,itmp + integer i,j,ia,ja,it,jt,itmp,iab(20),nb,nox,nhx c do 1 i=1,matm qatm(i)=0.0d0 @@ -91,6 +91,56 @@ c c endif 2 continue +c + do 3 i=1,natm + nb=0 + do 4 j=1,nbnd + ia=lbnd(1,j) + ja=lbnd(2,j) + if(ia.eq.i) then + nb=nb+1 + iab(nb)=ja + else + if(ja.eq.i) then + nb=nb+1 + iab(nb)=ia + endif + endif + 4 continue +c +c -COO (-) +c + if(latm(2,i).eq.6.and.nb.eq.3) then + nox=0 + do 5 j=1,nb + if(latm(2,iab(j)).eq.8.and.latm(3,iab(j)).eq.1) nox=nox+1 + 5 continue + if(nox.eq.2) then + do 6 j=1,nb + if(latm(2,iab(j)).eq.8.and.latm(3,iab(j)).eq.1) + + qatm(iab(j))=qatm(iab(j))-0.23 + 6 continue + qatm(i)=qatm(i)-0.54 + endif + endif +c +c -NH3 (+) +c + if(latm(2,i).eq.7.and.nb.eq.4) then + nhx=0 + do 7 j=1,nb + if(latm(2,iab(j)).eq.1.and.latm(3,iab(j)).eq.1) nhx=nhx+1 + 7 continue + if(nhx.eq.3) then + do 8 j=1,nb + if(latm(2,iab(j)).eq.1.and.latm(3,iab(j)).eq.1) + + qatm(iab(j))=qatm(iab(j))+0.07 + 8 continue + qatm(i)=qatm(i)+0.79 + endif + endif +c + 3 continue c pre_charge=.true. return diff --git a/src/prepar/pre_dihimp.F b/src/prepar/pre_dihimp.F index 6b1f6c0ce1..4917d06d3c 100644 --- a/src/prepar/pre_dihimp.F +++ b/src/prepar/pre_dihimp.F @@ -1,20 +1,30 @@ - logical function pre_dihimp(ffield,lang,mang,nang,ldih,mdih,ndih, - + limp,mimp,nimp) + logical function pre_dihimp(ffield,latm,matm,natm,lang,mang,nang, + + ldih,mdih,ndih,limp,mimp,nimp) c implicit none c - integer mang,nang,mdih,ndih,mimp,nimp - integer lang(3,mang),ldih(4,mdih),limp(4,mimp) + integer matm,natm,mang,nang,mdih,ndih,mimp,nimp + integer latm(5,matm),lang(3,mang),ldih(4,mdih),limp(4,mimp) character*80 ffield c - integer i,j,k,no + integer i,j,k,no,nd c if(ffield(1:6).eq.'gromos'.or.ffield(1:6).eq.'charmm') then c no=nimp +c +c dihedrals in aromatic ring are treated as improper by gromos +c both proper and improper by charmm c do 1 i=1,ndih + if(latm(2,ldih(1,i)).lt.4) goto 1 + if(latm(2,ldih(2,i)).lt.4) goto 1 + if(latm(2,ldih(3,i)).lt.4) goto 1 + if(latm(2,ldih(4,i)).lt.4) goto 1 do 2 j=1,nang + if(latm(2,lang(1,j)).lt.4) goto 2 + if(latm(2,lang(2,j)).lt.4) goto 2 + if(latm(2,lang(3,j)).lt.4) goto 2 if((ldih(1,i).eq.lang(1,j).and.ldih(4,i).eq.lang(3,j)).or. + (ldih(1,i).eq.lang(3,j).and.ldih(4,i).eq.lang(2,j))) then do 3 k=no+1,nimp @@ -27,7 +37,15 @@ c 1 continue c do 4 i=1,ndih-1 + if(latm(2,ldih(1,i)).lt.4) goto 4 + if(latm(2,ldih(2,i)).lt.4) goto 4 + if(latm(2,ldih(3,i)).lt.4) goto 4 + if(latm(2,ldih(4,i)).lt.4) goto 4 do 5 j=i+1,ndih + if(latm(2,ldih(1,j)).lt.4) goto 5 + if(latm(2,ldih(2,j)).lt.4) goto 5 + if(latm(2,ldih(3,j)).lt.4) goto 5 + if(latm(2,ldih(4,j)).lt.4) goto 5 if((ldih(1,i).eq.ldih(1,j).and.ldih(4,i).eq.ldih(4,j)).or. + (ldih(1,i).eq.ldih(4,j).and.ldih(4,i).eq.ldih(4,j))) then do 6 k=no+1,nimp @@ -50,7 +68,7 @@ c do 10 j=1,4 limp(j,i)=ldih(j,k) 10 continue - ldih(1,k)=0 + if(ffield(1:6).eq.'gromos') ldih(1,k)=0 9 continue c no=ndih diff --git a/src/prepar/pre_ffield.F b/src/prepar/pre_ffield.F index dfa7ce6fff..2b3dd5a447 100644 --- a/src/prepar/pre_ffield.F +++ b/src/prepar/pre_ffield.F @@ -98,6 +98,10 @@ c read(lfnpar,1000,end=999) card if(card(1:5).eq.'Cross') goto 2 if(card(1:5).eq.'Bonds') goto 2 + if(card(1:6).eq.'Angles') goto 2 + if(card(1:16).eq.'Proper dihedrals') goto 2 + if(card(1:18).eq.'Improper dihedrals') goto 2 + if(card(1:10).eq.'Atom types') goto 2 read(card,1001,err=2) atyp1(1:5),rmass,p1,p2,p3,p4,nt,ng 1001 format(a5,f10.5,4e12.5,i5,1x,i10) if(atyp1.eq.' ') goto 1 @@ -158,7 +162,12 @@ c atyp2(6:6)=' ' 6 continue read(lfnpar,1000,end=999) card + if(card(1:5).eq.'Atoms') goto 2 if(card(1:5).eq.'Bonds') goto 2 + if(card(1:6).eq.'Angles') goto 2 + if(card(1:16).eq.'Proper dihedrals') goto 2 + if(card(1:18).eq.'Improper dihedrals') goto 2 + if(card(1:10).eq.'Atom types') goto 2 read(card,1003,err=2) atyp1(1:5),atyp2(1:5),p1,p2 1003 format(a5,1x,a5,2x,2e12.5) if(atyp1.eq.' ') goto 1 @@ -204,7 +213,12 @@ c btyp2(6:6)=' ' 9 continue read(lfnpar,1000,end=999) card + if(card(1:5).eq.'Atoms') goto 2 + if(card(1:5).eq.'Cross') goto 2 if(card(1:6).eq.'Angles') goto 2 + if(card(1:16).eq.'Proper dihedrals') goto 2 + if(card(1:18).eq.'Improper dihedrals') goto 2 + if(card(1:10).eq.'Atom types') goto 2 read(card,1005,err=2) atyp1(1:5),atyp2(1:5),p1,p2 1005 format(a5,1x,a5,f10.5,e12.5) if(atyp1.eq.' ') goto 1 @@ -238,7 +252,12 @@ c btyp3(6:6)=' ' 12 continue read(lfnpar,1000,end=999) card - if(card(1:6).eq.'Proper') goto 2 + if(card(1:5).eq.'Atoms') goto 2 + if(card(1:5).eq.'Cross') goto 2 + if(card(1:5).eq.'Bonds') goto 2 + if(card(1:16).eq.'Proper dihedrals') goto 2 + if(card(1:18).eq.'Improper dihedrals') goto 2 + if(card(1:10).eq.'Atom types') goto 2 read(card,1006,err=2) atyp1(1:5),atyp2(1:5),atyp3(1:5),p1,p2 1006 format(a5,1x,a5,1x,a5,f10.5,e12.5) if(atyp1.eq.' ') goto 1 @@ -265,7 +284,7 @@ c c torsions c -------- c - if(card(1:6).eq.'Proper') then + if(card(1:16).eq.'Proper dihedrals') then atyp1(6:6)=' ' atyp2(6:6)=' ' atyp3(6:6)=' ' @@ -276,7 +295,12 @@ c btyp4(6:6)=' ' 15 continue read(lfnpar,1000,end=999) card - if(card(1:8).eq.'Improper') goto 2 + if(card(1:5).eq.'Atoms') goto 2 + if(card(1:5).eq.'Cross') goto 2 + if(card(1:5).eq.'Bonds') goto 2 + if(card(1:6).eq.'Angles') goto 2 + if(card(1:18).eq.'Improper dihedrals') goto 2 + if(card(1:10).eq.'Atom types') goto 2 read(card,1007,err=2) atyp1(1:5),atyp2(1:5),atyp3(1:5),atyp4(1:5), + p1,p2,mult 1007 format(a5,1x,a5,1x,a5,1x,a5,f10.5,e12.5,i5) @@ -315,7 +339,7 @@ c c impropers c --------- c - if(card(1:8).eq.'Improper') then + if(card(1:18).eq.'Improper dihedrals') then atyp1(6:6)=' ' atyp2(6:6)=' ' atyp3(6:6)=' ' @@ -326,6 +350,12 @@ c btyp4(6:6)=' ' 18 continue read(lfnpar,1000,end=999) card + if(card(1:5).eq.'Atoms') goto 2 + if(card(1:5).eq.'Cross') goto 2 + if(card(1:5).eq.'Bonds') goto 2 + if(card(1:6).eq.'Angles') goto 2 + if(card(1:16).eq.'Proper dihedrals') goto 2 + if(card(1:10).eq.'Atom types') goto 2 if(ffield(1:5).eq.'amber') then read(card,1008,err=2) atyp2(1:5),atyp3(1:5),atyp1(1:5),atyp4(1:5), + p1,p2,mult diff --git a/src/prepar/pre_imp.F b/src/prepar/pre_imp.F index 4a97657dc5..eb9dfc52d3 100644 --- a/src/prepar/pre_imp.F +++ b/src/prepar/pre_imp.F @@ -78,8 +78,9 @@ c c c aromatic ring c - if(ffield(1:5).eq.'amber'.or.ffield(1:5).eq.'charm') then - if(latm(2,i).ge.4.and.na.eq.3) then + if(latm(2,i).ge.4.and.na.eq.3.and. + + ((ffield(1:5).eq.'amber'.and.latm(2,i).le.6).or. + + (ffield(1:5).eq.'charm'))) then nimp=nimp+1 limp(1,nimp)=i limp(2,nimp)=ia(1) @@ -87,7 +88,6 @@ c limp(4,nimp)=ia(3) c write(*,'(6i5)') nimp,4,(limp(j,nimp),j=1,4) endif - endif c endif 1 continue diff --git a/src/prepar/pre_mkfrg.F b/src/prepar/pre_mkfrg.F index c065efd99a..e4a236a32e 100644 --- a/src/prepar/pre_mkfrg.F +++ b/src/prepar/pre_mkfrg.F @@ -1,4 +1,4 @@ - logical function pre_mkfrg(iunit,lfntyp,filtyp, + logical function pre_mkfrg(iunit,lfnpar,filpar, + dir_s,dir_x,dir_u,dir_t, + lfnout,lprint, + lseq,cseq,mseq,nseq,lsgm,csgm,msgm,nsgm,latm,catm,xatm,qatm,matm, @@ -36,7 +36,7 @@ c logical pre_atype,pre_charge external pre_atype,pre_charge c - integer iunit,lfntyp,lfnout + integer iunit,lfnpar,lfnout integer mseq,msgm,matm integer nseq,nsgm,natm integer mbnd,mang,mdih,mimp @@ -45,8 +45,8 @@ c character*6 catm(3,matm) real*8 xatm(3,matm),qatm(matm) integer lbnd(2,mbnd),lang(3,mang),ldih(4,mdih),limp(4,mimp) - character*255 filtyp,filnam - character*255 filpar,dir_s,dir_x,dir_u,dir_t + character*255 filpar,filnam + character*255 dir_s,dir_x,dir_u,dir_t integer lprint c integer len_t,length,len_f @@ -60,6 +60,7 @@ c integer link1,link2,link3,link4,link5,link integer ibnd,nbnd,nang,ndih,nconn real*8 rlink,dist + logical found c if(lprint.ge.5) then write(lfnout,2000) @@ -163,7 +164,7 @@ c c determine atom types c -------------------- c - if(.not.pre_atype(lfnout,lprint,lfntyp,filtyp, + if(.not.pre_atype(lfnout,lprint,lfnpar,filpar, + dir_s,dir_x,dir_u,dir_t,latm,catm,matm,lbnd,mbnd,nbnd, + jlo,ilo,ihi,jhi,int_mb(i_ltyp),mtyp,lring,aring,mring, + nring3,nring4,nring5,nring6)) @@ -179,7 +180,7 @@ c c guestimate partial charges c -------------------------- c - if(.not.pre_charge(latm,qatm,matm,lbnd,mbnd,nbnd)) + if(.not.pre_charge(latm,qatm,matm,natm,lbnd,mbnd,nbnd)) + call errquit('pre_charge failed',9999) c c write the fragment file @@ -206,6 +207,7 @@ c 3002 format(/,'FRAGMENT ',a,/, + ' num name type link cntr grp pgrp charge polarizab', + /) + found=.true. do 6 j=ilo,ihi if(latm(2,j).gt.0) then c @@ -214,6 +216,7 @@ c write(iunit,2003) j+1-ilo,catm(2,j),catm(3,j), + latm(5,j),latm(4,j),1,1,qatm(j),0.0 2003 format(i5,a4,2x,a6,4i5,2f12.6) + if(catm(3,j)(1:1).eq.' ') found=.false. endif 6 continue c @@ -238,6 +241,12 @@ c c write(lfnout,3005) filnam(1:len_f) 3005 format(/,'CREATED FRAGMENT FILE ',a) +c + if(.not.found) then + write(lfnout,3006) csgm(i)(1:length) + 3006 format(/,'ATOM TYPES COULD NOT BE DETERMINED FOR FRAGMENT ',a) + goto 9999 + endif c lsgm(2,i)=-4 c diff --git a/src/prepar/pre_mklist.F b/src/prepar/pre_mklist.F index d2d7fafeeb..be8fe6115c 100644 --- a/src/prepar/pre_mklist.F +++ b/src/prepar/pre_mklist.F @@ -368,8 +368,8 @@ c 20 continue endif c -c switch propers to impropers -c --------------------------- +c switch selected propers to impropers +c ------------------------------------ c if(ffield(1:6).eq.'gromos'.or.ffield(1:6).eq.'charmm') then no=nimp diff --git a/src/prepar/pre_mkseq.F b/src/prepar/pre_mkseq.F index d636e19b1b..9f01ea2ec0 100644 --- a/src/prepar/pre_mkseq.F +++ b/src/prepar/pre_mkseq.F @@ -1,5 +1,5 @@ logical function pre_mkseq(lfnout,lprint,ffield, - + lfnpdb,filpdb,lfnseq,filseq,lfntyp,filtyp, + + lfnpdb,filpdb,lfnseq,filseq,lfnpar,filpar, + lfnfrg,lfnsgm,dir_s,dir_x,dir_u,dir_t,slvnam,slvmdl) c c in : integer lfnout = logical file number output file @@ -8,8 +8,8 @@ c integer lfnpdb = logical file number for pdb file c char*80 filpdb = file name of pdb file c integer lfnseq = logical file number for seq file c char*80 filseq = file name of seq file -c integer lfntyp = logical file number for typ file -c char*80 filtyp = file name of typ file +c integer lfnpar = logical file number for typ file +c char*80 filpar = file name of typ file c integer lfnfrg = logical file number for frg file c char*80 dir_s = directory standard ffield files c char*80 dir_x = directory extensions ffield files @@ -36,13 +36,13 @@ c integer l_latm,i_latm,l_catm,i_catm,l_xatm,i_xatm,l_qatm,i_qatm integer l_bnd,i_bnd,l_ang,i_ang,l_dih,i_dih,l_imp,i_imp c - character*255 filpdb,filtyp,filseq + character*255 filpdb,filpar,filseq character*255 dir_s,dir_x,dir_u,dir_t character*80 ffield character*3 slvnam character*10 slvmdl c - integer lfnpdb,lfnout,lfnfrg,lfntyp,lfnseq,lfnsgm + integer lfnpdb,lfnout,lfnfrg,lfnpar,lfnseq,lfnsgm c if(lprint.ge.5) write(lfnout,1000) 1000 format('SEQUENCE GENERATION') @@ -259,7 +259,7 @@ c c c generate fragment files for missing segments in dir_t c - if(.not.pre_mkfrg(lfnfrg,lfntyp,filtyp, + if(.not.pre_mkfrg(lfnfrg,lfnpar,filpar, + dir_s,dir_x,dir_u,dir_t,lfnout,lprint, + int_mb(i_lseq),byte_mb(i_cseq),mseq,nseq, + int_mb(i_lsgm),byte_mb(i_csgm),msgm,nsgm,int_mb(i_latm), diff --git a/src/prepar/pre_mksgm.F b/src/prepar/pre_mksgm.F index ef226144d4..e29dcc1f8f 100644 --- a/src/prepar/pre_mksgm.F +++ b/src/prepar/pre_mksgm.F @@ -133,7 +133,8 @@ c c switch propers to impropers c --------------------------- c - if(.not.pre_dihimp(ffield,int_mb(i_ang),mang,nang, + if(.not.pre_dihimp(ffield,int_mb(i_latm),matm,natm, + + int_mb(i_ang),mang,nang, + int_mb(i_dih),mdih,ndih,int_mb(i_imp),mimp,nimp)) + call errquit('pre_dihsel failed',9999) c diff --git a/src/prepar/pre_rdsgm.F b/src/prepar/pre_rdsgm.F index 068795b93a..aab292120f 100644 --- a/src/prepar/pre_rdsgm.F +++ b/src/prepar/pre_rdsgm.F @@ -62,7 +62,7 @@ c read(lfnsgm,1003) (qatm(j,natm),j=1,6) 1003 format(6f12.6) c - if(lprint.ge.5) then + if(lprint.ge.50) then write(lfnout,1002) (catm(j,natm),j=1,4),(latm(j,natm),j=1,4) write(lfnout,1003) (qatm(j,natm),j=1,6) endif @@ -93,7 +93,7 @@ c read(lfnsgm,1005) (rbnd(j,nbnd),j=1,6) 1005 format(3(f12.6,e12.5)) c - if(lprint.ge.5) then + if(lprint.ge.50) then write(lfnout,1004) (lbnd(j,nbnd),j=1,4) write(lfnout,1005) (rbnd(j,nbnd),j=1,6) endif @@ -114,7 +114,7 @@ c lang(1,nang)=lang(1,nang)+na lang(2,nang)=lang(2,nang)+na lang(3,nang)=lang(3,nang)+na - if(lprint.ge.5) then + if(lprint.ge.50) then write(lfnout,1006) (lang(j,nang),j=1,5) write(lfnout,1007) (rang(j,nang),j=1,6) endif @@ -129,7 +129,7 @@ c read(lfnsgm,1009) (md(j),rd(2*j-1),rd(2*j),j=1,3) 1009 format(3(i5,f12.6,e12.5)) c - if(lprint.ge.5) then + if(lprint.ge.50) then write(lfnout,1008) (ld(j),j=1,6) write(lfnout,1009) (md(j),rd(2*j-1),rd(2*j),j=1,3) endif @@ -174,7 +174,7 @@ c + rimp(2*j,nimp),j=1,3) 1011 format(3(i5,f12.6,e12.5)) c - if(lprint.ge.5) then + if(lprint.ge.50) then write(lfnout,1010) (limp(j,nimp),j=1,6) write(lfnout,1011) (limp(6+j,nimp),rimp(2*j-1,nimp), + rimp(2*j,nimp),j=1,3) diff --git a/src/prepar/pre_wrttop.F b/src/prepar/pre_wrttop.F index 2308eda94c..d60b621611 100644 --- a/src/prepar/pre_wrttop.F +++ b/src/prepar/pre_wrttop.F @@ -110,7 +110,7 @@ c + 0,i-nang 1114 format(5i7) write(lfntop,1115) (rang(j,i),j=1,6) - 1115 format(3(0pf12.6,1pe12.5)) + 1115 format(3(0pf10.6,1pe12.5)) 21 continue endif c @@ -122,7 +122,7 @@ c 1116 format(6i7) write(lfntop,1117) + (kdih(k,j,i),rdih(k,2*j-1,i),rdih(k,2*j,i),j=1,3) - 1117 format(3(i3,0pf12.6,1pe12.5)) + 1117 format(3(i3,0pf10.6,1pe12.5)) 23 continue 22 continue endif @@ -137,7 +137,7 @@ c write(lfntop,1118) (limp(j,i)-natm,j=1,4),0,i-nimp 1118 format(6i7) write(lfntop,1119) (0,rimp(2*j-1,i),rimp(2*j,i),j=1,3) - 1119 format(3(i3,0pf12.6,1pe12.5)) + 1119 format(3(i3,0pf10.6,1pe12.5)) endif 24 continue endif @@ -183,7 +183,7 @@ c write(lfntop,1014) lang(1,i),lang(2,i),lang(3,i),0,i 1014 format(5i7) write(lfntop,1015) (rang(j,i),j=1,6) - 1015 format(3(0pf12.6,1pe12.5)) + 1015 format(3(0pf10.6,1pe12.5)) 12 continue c do 13 i=1,ndih @@ -193,7 +193,7 @@ c 1016 format(6i7) write(lfntop,1017) + (kdih(k,j,i),rdih(k,2*j-1,i),rdih(k,2*j,i),j=1,3) - 1017 format(3(i3,0pf12.6,1pe12.5)) + 1017 format(3(i3,0pf10.6,1pe12.5)) 14 continue 13 continue c @@ -205,7 +205,7 @@ c write(lfntop,1018) (limp(j,i),j=1,4),0,i 1018 format(6i7) write(lfntop,1019) (0,rimp(2*j-1,i),rimp(2*j,i),j=1,3) - 1019 format(3(i3,0pf12.6,1pe12.5)) + 1019 format(3(i3,0pf10.6,1pe12.5)) endif 15 continue c diff --git a/src/prepar/prepar.F b/src/prepar/prepar.F index 0c139a82ac..850f9700ec 100644 --- a/src/prepar/prepar.F +++ b/src/prepar/prepar.F @@ -1,6 +1,6 @@ logical function prepar(irtdb0) c -c $Id: prepar.F,v 1.7 1998-04-17 16:32:02 d3j191 Exp $ +c $Id: prepar.F,v 1.8 1998-04-20 20:01:41 d3j191 Exp $ c c ******************************************************** c ******************************************************** @@ -37,13 +37,13 @@ c c integer irtdb0,irtdb c - character*255 filpdb,filtyp,filseq,filtop,filpar + character*255 filpdb,filseq,filtop,filpar character*255 dir_s,dir_x,dir_u,dir_t character*80 sysnam,ffield character*10 slvmdl character*3 slvnam c - integer lfnpdb,lfnout,lfnfrg,lfntyp,lfnseq,lfnsgm,lfntop,lfnpar + integer lfnpdb,lfnout,lfnfrg,lfnseq,lfnsgm,lfntop,lfnpar integer len,lend,lprint integer newtop,newseq c @@ -80,7 +80,6 @@ c -------------------- c len=index(ffield,' ')-1 lend=index(dir_s,' ')-1 - filtyp=ffield(1:len)//'.typ ' filpar=ffield(1:len)//'.par ' c len=index(sysnam,' ')-1 @@ -95,11 +94,10 @@ c lfnout=6 lfnpdb=12 lfnfrg=13 - lfntyp=14 - lfnseq=15 - lfnsgm=16 - lfntop=17 - lfnpar=18 + lfnseq=14 + lfnsgm=15 + lfntop=16 + lfnpar=17 c lprint=10 if(util_print('sequence',print_none)) lprint=0 @@ -148,7 +146,7 @@ c generate sequence & segments from coordinates c --------------------------------------------- c if(.not.pre_mkseq(lfnout,lprint,ffield, - + lfnpdb,filpdb,lfnseq,filseq,lfntyp,filtyp, + + lfnpdb,filpdb,lfnseq,filseq,lfnpar,filpar, + lfnfrg,lfnsgm,dir_s,dir_x,dir_u,dir_t,slvnam,slvmdl)) + call errquit('pre_mkseq failed',9999) write(lfnout,1002) filseq(1:index(filseq,' ')-1)