module io_base !----------------------------------------------------------------------- use names use io_units implicit none !----------------------------------------------------------------------------- ! Max. number of AA residues integer,parameter :: maxres=6000!1200 ! Appr. max. number of interaction sites integer,parameter :: maxres2=2*maxres ! parameter (maxres6=6*maxres) ! parameter (mmaxres2=(maxres2*(maxres2+1)/2)) !----------------------------------------------------------------------------- ! Max. number of S-S bridges ! integer,parameter :: maxss=20 !----------------------------------------------------------------------------- ! Max. number of derivatives of virtual-bond and side-chain vectors in theta ! or phi. ! integer,parameter :: maxdim=(maxres-1)*(maxres-2)/2 !----------------------------------------------------------------------------- ! ! !----------------------------------------------------------------------------- contains !----------------------------------------------------------------------------- ! readrtns_CSA.F !----------------------------------------------------------------------------- subroutine read_bridge ! Read information about disulfide bridges. use geometry_data, only: nres use energy_data use control_data, only:out1file use MPI_data ! implicit real*8 (a-h,o-z) ! include 'DIMENSIONS' #ifdef MPI include 'mpif.h' #endif ! include 'COMMON.IOUNITS' ! include 'COMMON.GEO' ! include 'COMMON.VAR' ! include 'COMMON.INTERACT' ! include 'COMMON.LOCAL' ! include 'COMMON.NAMES' ! include 'COMMON.CHAIN' ! include 'COMMON.FFIELD' ! include 'COMMON.SBRIDGE' ! include 'COMMON.HEADER' ! include 'COMMON.CONTROL' ! include 'COMMON.DBASE' ! include 'COMMON.THREAD' ! include 'COMMON.TIME1' ! include 'COMMON.SETUP' !el local variables integer :: i,j,ierror ! Read bridging residues. read (inp,*) ns if (ns.gt.0) then allocate(iss(ns)) read (inp,*) (iss(i),i=1,ns) endif ! print *,'ns=',ns if(me.eq.king.or..not.out1file) & write (iout,*) 'ns=',ns if (ns.gt.0) & write(iout,*) ' iss:',(iss(i),i=1,ns) ! Check whether the specified bridging residues are cystines. do i=1,ns if (itype(iss(i)).ne.1) then if (me.eq.king.or..not.out1file) write (iout,'(2a,i3,a)') & 'Do you REALLY think that the residue ',& restyp(itype(iss(i))),i,& ' can form a disulfide bridge?!!!' write (*,'(2a,i3,a)') & 'Do you REALLY think that the residue ',& restyp(itype(iss(i))),i,& ' can form a disulfide bridge?!!!' #ifdef MPI call MPI_Finalize(MPI_COMM_WORLD,ierror) stop #endif endif enddo ! Read preformed bridges. if (ns.gt.0) then read (inp,*) nss if (nss.gt.0) then if(.not.allocated(ihpb)) allocate(ihpb(nss)) if(.not.allocated(jhpb)) allocate(jhpb(nss)) read (inp,*) (ihpb(i),jhpb(i),i=1,nss) if(.not.allocated(dhpb)) allocate(dhpb(nss)) if(.not.allocated(forcon)) allocate(forcon(nss))!(maxdim) !el maxdim=(maxres-1)*(maxres-2)/2 if(.not.allocated(dhpb1)) allocate(dhpb1(nss)) if(.not.allocated(ibecarb)) allocate(ibecarb(nss)) ! el Initialize the bridge array do i=1,nss dhpb(i)=0.0D0 enddo endif !-------------------- if(fg_rank.eq.0) & write(iout,*)'nss=',nss !el,' ihpb,jhpb: ',(ihpb(i),jhpb(i),i=1,nss) if (nss.gt.0) then write(iout,*)'ihpb,jhpb: ',(ihpb(i),jhpb(i),i=1,nss) nhpb=nss ! Check if the residues involved in bridges are in the specified list of ! bridging residues. do i=1,nss do j=1,i-1 if (ihpb(i).eq.ihpb(j).or.ihpb(i).eq.jhpb(j) & .or.jhpb(i).eq.ihpb(j).or.jhpb(i).eq.jhpb(j)) then write (iout,'(a,i3,a)') 'Disulfide pair',i,& ' contains residues present in other pairs.' write (*,'(a,i3,a)') 'Disulfide pair',i,& ' contains residues present in other pairs.' #ifdef MPI call MPI_Finalize(MPI_COMM_WORLD,ierror) stop #endif endif enddo do j=1,ns if (ihpb(i).eq.iss(j)) goto 10 enddo write (iout,'(a,i3,a)') 'Pair',i,' contains unknown cystine.' 10 continue do j=1,ns if (jhpb(i).eq.iss(j)) goto 20 enddo write (iout,'(a,i3,a)') 'Pair',i,' contains unknown cystine.' 20 continue dhpb(i)=d0cm forcon(i)=akcm enddo do i=1,nss ihpb(i)=ihpb(i)+nres jhpb(i)=jhpb(i)+nres enddo endif endif ! write(iout,*) "end read_bridge" return end subroutine read_bridge !----------------------------------------------------------------------------- subroutine read_x(kanal,*) use geometry_data use energy_data use geometry, only:int_from_cart1 ! implicit real*8 (a-h,o-z) ! include 'DIMENSIONS' ! include 'COMMON.GEO' ! include 'COMMON.VAR' ! include 'COMMON.CHAIN' ! include 'COMMON.IOUNITS' ! include 'COMMON.CONTROL' ! include 'COMMON.LOCAL' ! include 'COMMON.INTERACT' ! Read coordinates from input ! !el local variables integer :: l,k,j,i,kanal read(kanal,'(8f10.5)',end=10,err=10) & ((c(l,k),l=1,3),k=1,nres),& ((c(l,k+nres),l=1,3),k=nnt,nct) do j=1,3 c(j,nres+1)=c(j,1) c(j,2*nres)=c(j,nres) enddo call int_from_cart1(.false.) do i=1,nres-1 do j=1,3 dc(j,i)=c(j,i+1)-c(j,i) dc_norm(j,i)=dc(j,i)*vbld_inv(i+1) enddo enddo do i=nnt,nct if (itype(i).ne.10 .and. itype(i).ne.ntyp1) then do j=1,3 dc(j,i+nres)=c(j,i+nres)-c(j,i) dc_norm(j,i+nres)=dc(j,i+nres)*vbld_inv(i+nres) enddo endif enddo return 10 return 1 end subroutine read_x !----------------------------------------------------------------------------- subroutine read_threadbase use geometry_data use energy_data use compare_data ! implicit real*8 (a-h,o-z) ! include 'DIMENSIONS' ! include 'COMMON.IOUNITS' ! include 'COMMON.GEO' ! include 'COMMON.VAR' ! include 'COMMON.INTERACT' ! include 'COMMON.LOCAL' ! include 'COMMON.NAMES' ! include 'COMMON.CHAIN' ! include 'COMMON.FFIELD' ! include 'COMMON.SBRIDGE' ! include 'COMMON.HEADER' ! include 'COMMON.CONTROL' ! include 'COMMON.DBASE' ! include 'COMMON.THREAD' ! include 'COMMON.TIME1' !el local variables integer :: k,j,i ! Read pattern database for threading. read (icbase,*) nseq allocate(cart_base(3,maxres_base,nseq)) !(3,maxres_base,maxseq) allocate(nres_base(3,nseq)) !(3,maxseq) allocate(str_nam(nseq)) !(maxseq) do i=1,nseq read (icbase,'(i5,2x,a8,2i4)') nres_base(1,i),str_nam(i),& nres_base(2,i),nres_base(3,i) read (icbase,'(9f8.3)') ((cart_base(k,j,i),k=1,3),j=1,& nres_base(1,i)) ! write (iout,'(i5,2x,a8,2i4)') nres_base(1,i),str_nam(i), ! & nres_base(2,i),nres_base(3,i) ! write (iout,'(9f8.3)') ((cart_base(k,j,i),k=1,3),j=1, ! & nres_base(1,i)) enddo close (icbase) if (weidis.eq.0.0D0) weidis=0.1D0 do i=nnt,nct do j=i+2,nct nhpb=nhpb+1 ihpb(nhpb)=i jhpb(nhpb)=j forcon(nhpb)=weidis enddo enddo read (inp,*) nexcl,(iexam(1,i),iexam(2,i),i=1,nexcl) write (iout,'(a,i5)') 'nexcl: ',nexcl write (iout,'(2i5)') (iexam(1,i),iexam(2,i),i=1,nexcl) return end subroutine read_threadbase !----------------------------------------------------------------------------- #ifdef WHAM_RUN !el subroutine read_angles(kanal,iscor,energ,iprot,*) subroutine read_angles(kanal,*) use geometry_data use energy_data ! implicit real*8 (a-h,o-z) ! include 'DIMENSIONS' ! include 'DIMENSIONS.ZSCOPT' ! include 'COMMON.INTERACT' ! include 'COMMON.SBRIDGE' ! include 'COMMON.GEO' ! include 'COMMON.VAR' ! include 'COMMON.CHAIN' ! include 'COMMON.IOUNITS' character(len=80) :: lineh integer :: iscor,iprot,ic real(kind=8) :: energ #else subroutine read_angles(kanal,*) use geometry_data ! use energy ! use control ! implicit real*8 (a-h,o-z) ! include 'DIMENSIONS' ! include 'COMMON.GEO' ! include 'COMMON.VAR' ! include 'COMMON.CHAIN' ! include 'COMMON.IOUNITS' ! include 'COMMON.CONTROL' #endif ! Read angles from input ! !el local variables integer :: i,kanal #ifdef WHAM_RUN read(kanal,'(a80)',end=10,err=10) lineh read(lineh(:5),*,err=8) ic read(lineh(6:),*,err=8) energ goto 9 8 ic=1 print *,'error, assuming e=1d10',lineh energ=1d10 nss=0 9 continue read(lineh(18:),*,end=10,err=10) nss IF (NSS.LT.9) THEN read (lineh(20:),*,end=10,err=10) & (IHPB(I),JHPB(I),I=1,NSS),iscor ELSE read (lineh(20:),*,end=10,err=10) (IHPB(I),JHPB(I),I=1,8) read (kanal,*,end=10,err=10) (IHPB(I),JHPB(I),& I=9,NSS),iscor ENDIF ! print *,"energy",energ," iscor",iscor #endif read (kanal,*,err=10,end=10) (theta(i),i=3,nres) read (kanal,*,err=10,end=10) (phi(i),i=4,nres) read (kanal,*,err=10,end=10) (alph(i),i=2,nres-1) read (kanal,*,err=10,end=10) (omeg(i),i=2,nres-1) do i=1,nres ! 9/7/01 avoid 180 deg valence angle if (theta(i).gt.179.99d0) theta(i)=179.99d0 ! theta(i)=deg2rad*theta(i) phi(i)=deg2rad*phi(i) alph(i)=deg2rad*alph(i) omeg(i)=deg2rad*omeg(i) enddo return 10 return 1 end subroutine read_angles !----------------------------------------------------------------------------- subroutine reada(rekord,lancuch,wartosc,default) ! implicit none character*(*) :: rekord,lancuch real(kind=8) :: wartosc,default integer :: iread !,ilen !el external ilen iread=index(rekord,lancuch) if (iread.eq.0) then wartosc=default return endif iread=iread+ilen(lancuch)+1 read (rekord(iread:),*,err=10,end=10) wartosc return 10 wartosc=default return end subroutine reada !----------------------------------------------------------------------------- subroutine readi(rekord,lancuch,wartosc,default) ! implicit none character*(*) :: rekord,lancuch integer :: wartosc,default integer :: iread !,ilen !el external ilen iread=index(rekord,lancuch) if (iread.eq.0) then wartosc=default return endif iread=iread+ilen(lancuch)+1 read (rekord(iread:),*,err=10,end=10) wartosc return 10 wartosc=default return end subroutine readi !----------------------------------------------------------------------------- subroutine multreadi(rekord,lancuch,tablica,dim,default) ! implicit none integer :: dim,i integer :: tablica(dim),default character*(*) :: rekord,lancuch character(len=80) :: aux integer :: iread !,ilen !el external ilen do i=1,dim tablica(i)=default enddo iread=index(rekord,lancuch(:ilen(lancuch))//"=") if (iread.eq.0) return iread=iread+ilen(lancuch)+1 read (rekord(iread:),*,end=10,err=10) (tablica(i),i=1,dim) 10 return end subroutine multreadi !----------------------------------------------------------------------------- subroutine multreada(rekord,lancuch,tablica,dim,default) ! implicit none integer :: dim,i real(kind=8) :: tablica(dim),default character*(*) :: rekord,lancuch character(len=80) :: aux integer :: iread !,ilen !el external ilen do i=1,dim tablica(i)=default enddo iread=index(rekord,lancuch(:ilen(lancuch))//"=") if (iread.eq.0) return iread=iread+ilen(lancuch)+1 read (rekord(iread:),*,end=10,err=10) (tablica(i),i=1,dim) 10 return end subroutine multreada !----------------------------------------------------------------------------- subroutine card_concat(card,to_upper) ! dla UNRESA to_upper jest zawsze .true. ! implicit real*8 (a-h,o-z) ! include 'DIMENSIONS' ! include 'COMMON.IOUNITS' character(*) :: card character(len=80) :: karta !,ucase logical :: to_upper !el external ilen read (inp,'(a)') karta if (to_upper) karta=ucase(karta) card=' ' do while (karta(80:80).eq.'&') card=card(:ilen(card)+1)//karta(:79) read (inp,'(a)') karta if (to_upper) karta=ucase(karta) enddo card=card(:ilen(card)+1)//karta return end subroutine card_concat !----------------------------------------------------------------------------- subroutine read_dist_constr use MPI_data ! use control use geometry, only: dist use geometry_data use control_data use energy_data ! implicit real*8 (a-h,o-z) ! include 'DIMENSIONS' #ifdef MPI include 'mpif.h' #endif ! include 'COMMON.SETUP' ! include 'COMMON.CONTROL' ! include 'COMMON.CHAIN' ! include 'COMMON.IOUNITS' ! include 'COMMON.SBRIDGE' integer,dimension(2,100) :: ifrag_,ipair_ real(kind=8),dimension(100) :: wfrag_,wpair_ character(len=640) :: controlcard !el local variables integer :: i,k,j,ddjk,ii,jj,itemp integer :: nfrag_,npair_,ndist_ real(kind=8) :: dist_cut ! write (iout,*) "Calling read_dist_constr" ! write (iout,*) "nres",nres," nstart_sup",nstart_sup," nsup",nsup ! call flush(iout) call card_concat(controlcard,.true.) call readi(controlcard,"NFRAG",nfrag_,0) call readi(controlcard,"NPAIR",npair_,0) call readi(controlcard,"NDIST",ndist_,0) call reada(controlcard,'DIST_CUT',dist_cut,5.0d0) call multreadi(controlcard,"IFRAG",ifrag_(1,1),2*nfrag_,0) call multreadi(controlcard,"IPAIR",ipair_(1,1),2*npair_,0) call multreada(controlcard,"WFRAG",wfrag_(1),nfrag_,0.0d0) call multreada(controlcard,"WPAIR",wpair_(1),npair_,0.0d0) ! write (iout,*) "NFRAG",nfrag_," NPAIR",npair_," NDIST",ndist_ ! write (iout,*) "IFRAG" ! do i=1,nfrag_ ! write (iout,*) i,ifrag_(1,i),ifrag_(2,i),wfrag_(i) ! enddo ! write (iout,*) "IPAIR" ! do i=1,npair_ ! write (iout,*) i,ipair_(1,i),ipair_(2,i),wpair_(i) ! enddo if(.not.allocated(ihpb)) allocate(ihpb(maxdim)) if(.not.allocated(jhpb)) allocate(jhpb(maxdim)) if(.not.allocated(dhpb)) allocate(dhpb(maxdim)) if(.not.allocated(forcon)) allocate(forcon(maxdim)) call flush(iout) do i=1,nfrag_ if (ifrag_(1,i).lt.nstart_sup) ifrag_(1,i)=nstart_sup if (ifrag_(2,i).gt.nstart_sup+nsup-1) & ifrag_(2,i)=nstart_sup+nsup-1 ! write (iout,*) i,ifrag_(1,i),ifrag_(2,i),wfrag_(i) call flush(iout) if (wfrag_(i).gt.0.0d0) then do j=ifrag_(1,i),ifrag_(2,i)-1 do k=j+1,ifrag_(2,i) ! write (iout,*) "j",j," k",k ddjk=dist(j,k) if (constr_dist.eq.1) then nhpb=nhpb+1 ihpb(nhpb)=j jhpb(nhpb)=k dhpb(nhpb)=ddjk forcon(nhpb)=wfrag_(i) else if (constr_dist.eq.2) then if (ddjk.le.dist_cut) then nhpb=nhpb+1 ihpb(nhpb)=j jhpb(nhpb)=k dhpb(nhpb)=ddjk forcon(nhpb)=wfrag_(i) endif else nhpb=nhpb+1 ihpb(nhpb)=j jhpb(nhpb)=k dhpb(nhpb)=ddjk forcon(nhpb)=wfrag_(i)*dexp(-0.5d0*(ddjk/dist_cut)**2) endif #ifdef MPI if (.not.out1file .or. me.eq.king) & write (iout,'(a,3i5,f8.2,1pe12.2)') "+dist.constr ",& nhpb,ihpb(nhpb),jhpb(nhpb),dhpb(nhpb),forcon(nhpb) #else write (iout,'(a,3i5,f8.2,1pe12.2)') "+dist.constr ",& nhpb,ihpb(nhpb),jhpb(nhpb),dhpb(nhpb),forcon(nhpb) #endif enddo enddo endif enddo do i=1,npair_ if (wpair_(i).gt.0.0d0) then ii = ipair_(1,i) jj = ipair_(2,i) if (ii.gt.jj) then itemp=ii ii=jj jj=itemp endif do j=ifrag_(1,ii),ifrag_(2,ii) do k=ifrag_(1,jj),ifrag_(2,jj) nhpb=nhpb+1 ihpb(nhpb)=j jhpb(nhpb)=k forcon(nhpb)=wpair_(i) dhpb(nhpb)=dist(j,k) #ifdef MPI if (.not.out1file .or. me.eq.king) & write (iout,'(a,3i5,f8.2,f10.1)') "+dist.constr ",& nhpb,ihpb(nhpb),jhpb(nhpb),dhpb(nhpb),forcon(nhpb) #else write (iout,'(a,3i5,f8.2,f10.1)') "+dist.constr ",& nhpb,ihpb(nhpb),jhpb(nhpb),dhpb(nhpb),forcon(nhpb) #endif enddo enddo endif enddo do i=1,ndist_ read (inp,*) ihpb(nhpb+1),jhpb(nhpb+1),forcon(nhpb+1) if (forcon(nhpb+1).gt.0.0d0) then nhpb=nhpb+1 dhpb(nhpb)=dist(ihpb(nhpb),jhpb(nhpb)) #ifdef MPI if (.not.out1file .or. me.eq.king) & write (iout,'(a,3i5,f8.2,f10.1)') "+dist.constr ",& nhpb,ihpb(nhpb),jhpb(nhpb),dhpb(nhpb),forcon(nhpb) #else write (iout,'(a,3i5,f8.2,f10.1)') "+dist.constr ",& nhpb,ihpb(nhpb),jhpb(nhpb),dhpb(nhpb),forcon(nhpb) #endif endif enddo call flush(iout) return end subroutine read_dist_constr !----------------------------------------------------------------------------- #ifdef WINIFL subroutine flush(iu) return end subroutine flush #endif #ifdef AIX subroutine flush(iu) call flush_(iu) return end subroutine flush #endif !----------------------------------------------------------------------------- subroutine copy_to_tmp(source) ! include "DIMENSIONS" ! include "COMMON.IOUNITS" character*(*) :: source character(len=256) :: tmpfile ! integer ilen !el external ilen logical :: ex tmpfile=curdir(:ilen(curdir))//"/"//source(:ilen(source)) inquire(file=tmpfile,exist=ex) if (ex) then write (*,*) "Copying ",tmpfile(:ilen(tmpfile)),& " to temporary directory..." write (*,*) "/bin/cp "//tmpfile(:ilen(tmpfile))//" "//tmpdir call system("/bin/cp "//tmpfile(:ilen(tmpfile))//" "//tmpdir) endif return end subroutine copy_to_tmp !----------------------------------------------------------------------------- subroutine move_from_tmp(source) ! include "DIMENSIONS" ! include "COMMON.IOUNITS" character*(*) :: source ! integer ilen !el external ilen write (*,*) "Moving ",source(:ilen(source)),& " from temporary directory to working directory" write (*,*) "/bin/mv "//source(:ilen(source))//" "//curdir call system("/bin/mv "//source(:ilen(source))//" "//curdir) return end subroutine move_from_tmp !----------------------------------------------------------------------------- ! misc.f !----------------------------------------------------------------------------- ! $Date: 1994/10/12 17:24:21 $ ! $Revision: 2.5 $ logical function find_arg(ipos,line,errflag) integer, parameter :: maxlen = 80 character(len=80) :: line character(len=1) :: empty=' ',equal='=' logical :: errflag integer :: ipos ! This function returns .TRUE., if an argument follows keyword keywd; if so ! IPOS will point to the first non-blank character of the argument. Returns ! .FALSE., if no argument follows the keyword; in this case IPOS points ! to the first non-blank character of the next keyword. do while (line(ipos:ipos) .eq. empty .and. ipos.le.maxlen) ipos=ipos+1 enddo errflag=.false. if (line(ipos:ipos).eq.equal) then find_arg=.true. ipos=ipos+1 do while (line(ipos:ipos) .eq. empty .and. ipos.le.maxlen) ipos=ipos+1 enddo if (ipos.gt.maxlen) errflag=.true. else find_arg=.false. endif return end function find_arg !----------------------------------------------------------------------------- logical function find_group(iunit,jout,key1) character*(*) :: key1 character(len=80) :: karta !,ucase integer :: iunit,jout integer :: ll !,ilen !EL external ilen !EL logical lcom rewind (iunit) karta=' ' ll=ilen(key1) do while (index(ucase(karta),key1(1:ll)).eq.0.or.lcom(1,karta)) read (iunit,'(a)',end=10) karta enddo write (jout,'(2a)') '> ',karta(1:78) find_group=.true. return 10 find_group=.false. return end function find_group !----------------------------------------------------------------------------- logical function iblnk(charc) character(len=1) :: charc integer :: n n = ichar(charc) iblnk = (n.eq.9) .or. (n.eq.10) .or. (charc.eq.' ') return end function iblnk !----------------------------------------------------------------------------- integer function ilen(string) character*(*) :: string !EL logical :: iblnk ilen = len(string) 1 if ( ilen .gt. 0 ) then if ( iblnk( string(ilen:ilen) ) ) then ilen = ilen - 1 goto 1 endif endif return end function ilen !----------------------------------------------------------------------------- integer function in_keywd_set(nkey,ikey,narg,keywd,keywdset) integer :: nkey,i,ikey,narg character(len=16) :: keywd,keywdset(1:nkey,0:nkey) ! character(len=16) :: ucase do i=1,narg if (ucase(keywd).eq.keywdset(i,ikey)) then ! Match found in_keywd_set=i return endif enddo ! No match to the allowed set of keywords if this point is reached. in_keywd_set=0 return end function in_keywd_set !----------------------------------------------------------------------------- character function lcase(string) integer :: i, k, idiff character*(*) :: string character(len=1) :: c character(len=40) :: chtmp ! i = len(lcase) k = len(string) if (i .lt. k) then k = i if (string(k+1:) .ne. ' ') then chtmp = string endif endif idiff = ichar('a') - ichar('A') lcase = string do 99 i = 1, k c = string(i:i) if (lge(c,'A') .and. lle(c,'Z')) then lcase(i:i) = char(ichar(c) + idiff) endif 99 continue return end function lcase !----------------------------------------------------------------------------- logical function lcom(ipos,karta) character(len=80) :: karta character :: koment(2) = (/'!','#'/) integer :: ipos,i lcom=.false. do i=1,2 if (karta(ipos:ipos).eq.koment(i)) lcom=.true. enddo return end function lcom !----------------------------------------------------------------------------- logical function lower_case(ch) character*(*) :: ch lower_case=(ch.ge.'a' .and. ch.le.'z') return end function lower_case !----------------------------------------------------------------------------- subroutine mykey(line,keywd,ipos,blankline,errflag) ! This subroutine seeks a non-empty substring keywd in the string LINE. ! The substring begins with the first character different from blank and ! "=" encountered right to the pointer IPOS (inclusively) and terminates ! at the character left to the first blank or "=". When the subroutine is ! exited, the pointer IPOS is moved to the position of the terminator in LINE. ! The logical variable BLANKLINE is set at .TRUE., if LINE(IPOS:) contains ! only separators or the maximum length of the data line (80) has been reached. ! The logical variable ERRFLAG is set at .TRUE. if the string ! consists only from a "=". integer, parameter :: maxlen=80 character(len=1) :: empty=' ',equal='=',comma=',' character*(*) :: keywd character(len=80) :: line logical :: blankline,errflag !EL,lcom integer :: ipos,istart,iend errflag=.false. do while (line(ipos:ipos).eq.empty .and. (ipos.le.maxlen)) ipos=ipos+1 enddo if (ipos.gt.maxlen .or. lcom(ipos,line) ) then ! At this point the rest of the input line turned out to contain only blanks ! or to be commented out. blankline=.true. return endif blankline=.false. istart=ipos ! Checks whether the current char is a separator. do while (line(ipos:ipos).ne.empty .and. line(ipos:ipos).ne.equal & .and. line(ipos:ipos).ne.comma .and. ipos.le.maxlen) ipos=ipos+1 enddo iend=ipos-1 ! Error flag set to .true., if the length of the keyword was found less than 1. if (iend.lt.istart) then errflag=.true. return endif keywd=line(istart:iend) return end subroutine mykey !----------------------------------------------------------------------------- subroutine numstr(inum,numm) character(len=10) :: huj='0123456789' character*(*) :: numm integer :: inumm,inum,inum1,inum2 inumm=inum inum1=inumm/10 inum2=inumm-10*inum1 inumm=inum1 numm(3:3)=huj(inum2+1:inum2+1) inum1=inumm/10 inum2=inumm-10*inum1 inumm=inum1 numm(2:2)=huj(inum2+1:inum2+1) inum1=inumm/10 inum2=inumm-10*inum1 inumm=inum1 numm(1:1)=huj(inum2+1:inum2+1) return end subroutine numstr !----------------------------------------------------------------------------- function ucase(string) integer :: i, k, idiff character(*) :: string character(len=len(string)) :: ucase character(len=1) :: c character(len=40) :: chtmp ! i = len(ucase) k = len(string) if (i .lt. k) then k = i if (string(k+1:) .ne. ' ') then chtmp = string endif endif idiff = ichar('a') - ichar('A') ucase = string do 99 i = 1, k c = string(i:i) if (lge(c,'a') .and. lle(c,'z')) then ucase(i:i) = char(ichar(c) - idiff) endif 99 continue return end function ucase !----------------------------------------------------------------------------- ! geomout.F !----------------------------------------------------------------------------- subroutine pdbout(etot,tytul,iunit) use geometry_data, only: c,nres use energy_data ! use control use compare_data use MD_data ! implicit real*8 (a-h,o-z) ! include 'DIMENSIONS' ! include 'COMMON.CHAIN' ! include 'COMMON.INTERACT' ! include 'COMMON.NAMES' ! include 'COMMON.IOUNITS' ! include 'COMMON.HEADER' ! include 'COMMON.SBRIDGE' ! include 'COMMON.DISTFIT' ! include 'COMMON.MD' !el character(len=50) :: tytul character*(*) :: tytul character(len=1),dimension(10) :: chainid= (/'A','B','C','D','E','F','G','H','I','J'/) integer,dimension(nres) :: ica !(maxres) !el local variables integer :: j,iti,itj,itk,itl,i,iatom,ichain,ires,iunit real(kind=8) :: etot integer :: nres2 nres2=2*nres if(.not.allocated(vtot)) allocate(vtot(nres2)) !(maxres2) write (iunit,'(3a,1pe15.5)') 'REMARK ',tytul,' ENERGY ',etot !model write (iunit,'(a5,i6)') 'MODEL',1 if (nhfrag.gt.0) then do j=1,nhfrag iti=itype(hfrag(1,j)) itj=itype(hfrag(2,j)) if (j.lt.10) then write (iunit,'(a5,i5,1x,a1,i1,2x,a3,i7,2x,a3,i7,i3,t76,i5)') & 'HELIX',j,'H',j,& restyp(iti),hfrag(1,j)-1,& restyp(itj),hfrag(2,j)-1,1,hfrag(2,j)-hfrag(1,j) else write (iunit,'(a5,i5,1x,a1,i2,1x,a3,i7,2x,a3,i7,i3)') & 'HELIX',j,'H',j,& restyp(iti),hfrag(1,j)-1,& restyp(itj),hfrag(2,j)-1,1,hfrag(2,j)-hfrag(1,j) endif enddo endif if (nbfrag.gt.0) then do j=1,nbfrag iti=itype(bfrag(1,j)) itj=itype(bfrag(2,j)-1) write (iunit,'(a5,i5,1x,a1,i1,i3,1x,a3,i6,2x,a3,i6,i3)') & 'SHEET',1,'B',j,2,& restyp(iti),bfrag(1,j)-1,& restyp(itj),bfrag(2,j)-2,0 if (bfrag(3,j).gt.bfrag(4,j)) then itk=itype(bfrag(3,j)) itl=itype(bfrag(4,j)+1) write (iunit,'(a5,i5,1x,a1,i1,i3,1x,a3,i6,2x,a3,i6,i3,2x,a1,2x,a3,i6,3x,a1,2x,a3,i6)') & 'SHEET',2,'B',j,2,& restyp(itl),bfrag(4,j),& restyp(itk),bfrag(3,j)-1,-1,& "N",restyp(itk),bfrag(3,j)-1,& "O",restyp(iti),bfrag(1,j)-1 else itk=itype(bfrag(3,j)) itl=itype(bfrag(4,j)-1) write (iunit,'(a5,i5,1x,a1,i1,i3,1x,a3,i6,2x,a3,i6,i3,2x,a1,2x,a3,i6,3x,a1,2x,a3,i6)') & 'SHEET',2,'B',j,2,& restyp(itk),bfrag(3,j)-1,& restyp(itl),bfrag(4,j)-2,1,& "N",restyp(itk),bfrag(3,j)-1,& "O",restyp(iti),bfrag(1,j)-1 endif enddo endif if (nss.gt.0) then do i=1,nss if (dyn_ss) then write(iunit,'(a6,i4,1x,a3,i7,4x,a3,i7)') & 'SSBOND',i,'CYS',idssb(i)-nnt+1,& 'CYS',jdssb(i)-nnt+1 else write(iunit,'(a6,i4,1x,a3,i7,4x,a3,i7)') & 'SSBOND',i,'CYS',ihpb(i)-nnt+1-nres,& 'CYS',jhpb(i)-nnt+1-nres endif enddo endif iatom=0 ichain=1 ires=0 do i=nnt,nct iti=itype(i) if (iti.eq.ntyp1) then ichain=ichain+1 ires=0 write (iunit,'(a)') 'TER' else ires=ires+1 iatom=iatom+1 ica(i)=iatom write (iunit,10) iatom,restyp(iti),chainid(ichain),& ires,(c(j,i),j=1,3),vtot(i) if (iti.ne.10) then iatom=iatom+1 write (iunit,20) iatom,restyp(iti),chainid(ichain),& ires,(c(j,nres+i),j=1,3),& vtot(i+nres) endif endif enddo write (iunit,'(a)') 'TER' do i=nnt,nct-1 if (itype(i).eq.ntyp1) cycle if (itype(i).eq.10 .and. itype(i+1).ne.ntyp1) then write (iunit,30) ica(i),ica(i+1) else if (itype(i).ne.10 .and. itype(i+1).ne.ntyp1) then write (iunit,30) ica(i),ica(i+1),ica(i)+1 else if (itype(i).ne.10 .and. itype(i+1).eq.ntyp1) then write (iunit,30) ica(i),ica(i)+1 endif enddo if (itype(nct).ne.10) then write (iunit,30) ica(nct),ica(nct)+1 endif do i=1,nss if (dyn_ss) then write (iunit,30) ica(idssb(i))+1,ica(jdssb(i))+1 else write (iunit,30) ica(ihpb(i)-nres)+1,ica(jhpb(i)-nres)+1 endif enddo write (iunit,'(a6)') 'ENDMDL' 10 FORMAT ('ATOM',I7,' CA ',A3,1X,A1,I4,4X,3F8.3,f15.3) 20 FORMAT ('ATOM',I7,' CB ',A3,1X,A1,I4,4X,3F8.3,f15.3) 30 FORMAT ('CONECT',8I5) return end subroutine pdbout !----------------------------------------------------------------------------- subroutine MOL2out(etot,tytul) ! Prints the Cartesian coordinates of the alpha-carbons in the Tripos mol2 ! format. use geometry_data, only: c use energy_data ! use control ! implicit real*8 (a-h,o-z) ! include 'DIMENSIONS' ! include 'COMMON.CHAIN' ! include 'COMMON.INTERACT' ! include 'COMMON.NAMES' ! include 'COMMON.IOUNITS' ! include 'COMMON.HEADER' ! include 'COMMON.SBRIDGE' character(len=32) :: tytul,fd character(len=3) :: zahl character(len=6) :: res_num,pom !,ucase !el local variables integer :: i,j real(kind=8) :: etot #ifdef AIX call fdate_(fd) #elif (defined CRAY) call date(fd) #else call fdate(fd) #endif write (imol2,'(a)') '#' write (imol2,'(a)') & '# Creating user name: unres' write (imol2,'(2a)') '# Creation time: ',& fd write (imol2,'(/a)') '\@MOLECULE' write (imol2,'(a)') tytul write (imol2,'(5i5)') nct-nnt+1,nct-nnt+nss+1,nct-nnt+nss+1,0,0 write (imol2,'(a)') 'SMALL' write (imol2,'(a)') 'USER_CHARGES' write (imol2,'(a)') '\@ATOM' do i=nnt,nct write (zahl,'(i3)') i pom=ucase(restyp(itype(i))) res_num = pom(:3)//zahl(2:) write (imol2,10) i,(c(j,i),j=1,3),i,res_num,0.0 enddo write (imol2,'(a)') '\@BOND' do i=nnt,nct-1 write (imol2,'(i5,2i6,i2)') i-nnt+1,i-nnt+1,i-nnt+2,1 enddo do i=1,nss write (imol2,'(i5,2i6,i2)') nct-nnt+i,ihpb(i),jhpb(i),1 enddo write (imol2,'(a)') '\@SUBSTRUCTURE' do i=nnt,nct write (zahl,'(i3)') i pom = ucase(restyp(itype(i))) res_num = pom(:3)//zahl(2:) write (imol2,30) i-nnt+1,res_num,i-nnt+1,0 enddo 10 FORMAT (I7,' CA ',3F10.4,' C.3',I8,1X,A,F11.4,' ****') 30 FORMAT (I7,1x,A,I14,' RESIDUE',I13,' **** ****') return end subroutine MOL2out !----------------------------------------------------------------------------- subroutine intout use geometry_data use energy_data, only: itype ! use control ! implicit real*8 (a-h,o-z) ! include 'DIMENSIONS' ! include 'COMMON.IOUNITS' ! include 'COMMON.CHAIN' ! include 'COMMON.VAR' ! include 'COMMON.LOCAL' ! include 'COMMON.INTERACT' ! include 'COMMON.NAMES' ! include 'COMMON.GEO' ! include 'COMMON.TORSION' !el local variables integer :: i,iti write (iout,'(/a)') 'Geometry of the virtual chain.' write (iout,'(7a)') ' Res ',' d',' Theta',& ' Phi',' Dsc',' Alpha',' Omega' do i=1,nres iti=itype(i) write (iout,'(a3,i4,6f10.3)') restyp(iti),i,vbld(i),& rad2deg*theta(i),rad2deg*phi(i),vbld(nres+i),rad2deg*alph(i),& rad2deg*omeg(i) enddo return end subroutine intout !----------------------------------------------------------------------------- #ifdef CLUSTER subroutine briefout(it,ener,free)!,plik) #else subroutine briefout(it,ener) #endif use geometry_data use energy_data ! use control ! implicit real*8 (a-h,o-z) ! include 'DIMENSIONS' ! include 'COMMON.IOUNITS' ! include 'COMMON.CHAIN' ! include 'COMMON.VAR' ! include 'COMMON.LOCAL' ! include 'COMMON.INTERACT' ! include 'COMMON.NAMES' ! include 'COMMON.GEO' ! include 'COMMON.SBRIDGE' ! print '(a,i5)',intname,igeom !el local variables integer :: i,it real(kind=8) :: ener,free ! character(len=80) :: plik ! integer :: iii #if defined(AIX) || defined(PGI) open (igeom,file=intname,position='append') #else open (igeom,file=intname,access='append') #endif #ifdef WHAM_RUN ! iii=igeom igeom=iout #endif IF (NSS.LE.9) THEN #ifdef CLUSTER WRITE (igeom,180) IT,ENER,free,NSS,(IHPB(I),JHPB(I),I=1,NSS) ELSE WRITE (igeom,180) IT,ENER,free,NSS,(IHPB(I),JHPB(I),I=1,8) #else WRITE (igeom,180) IT,ENER,NSS,(IHPB(I),JHPB(I),I=1,NSS) ELSE WRITE (igeom,180) IT,ENER,NSS,(IHPB(I),JHPB(I),I=1,9) #endif WRITE (igeom,190) (IHPB(I),JHPB(I),I=10,NSS) ENDIF ! IF (nvar.gt.nphi) WRITE (igeom,200) (RAD2DEG*THETA(I),I=3,NRES) WRITE (igeom,200) (RAD2DEG*THETA(I),I=3,NRES) WRITE (igeom,200) (RAD2DEG*PHI(I),I=4,NRES) ! if (nvar.gt.nphi+ntheta) then write (igeom,200) (rad2deg*alph(i),i=2,nres-1) write (igeom,200) (rad2deg*omeg(i),i=2,nres-1) ! endif close(igeom) 180 format (I5,F12.3,I2,9(1X,2I3)) 190 format (3X,11(1X,2I3)) 200 format (8F10.4) return end subroutine briefout !----------------------------------------------------------------------------- #ifdef WINIFL subroutine fdate(fd) character(len=32) :: fd write(fd,'(32x)') return end subroutine fdate #endif !----------------------------------------------------------------------------- #ifdef WHAM_RUN real(kind=8) function gyrate(jcon) #else real(kind=8) function gyrate() #endif use geometry_data, only: c ! use geometry use energy_data ! implicit real*8 (a-h,o-z) ! include 'DIMENSIONS' ! include 'COMMON.INTERACT' ! include 'COMMON.CHAIN' real(kind=8) :: rg real(kind=8),dimension(3) :: cen !el local variables integer :: i,j,jcon do j=1,3 cen(j)=0.0d0 enddo do i=nnt,nct do j=1,3 cen(j)=cen(j)+c(j,i) enddo enddo do j=1,3 cen(j)=cen(j)/dble(nct-nnt+1) enddo rg = 0.0d0 do i = nnt, nct do j=1,3 rg = rg + (c(j,i)-cen(j))**2 enddo end do #ifdef WHAM_RUN gyrate = dsqrt(rg/dble(nct-nnt+1)) #else gyrate = sqrt(rg/dble(nct-nnt+1)) #endif return end function gyrate #ifdef WHAM_RUN !----------------------------------------------------------------------------- ! readrtns.F WHAM subroutine reads(rekord,lancuch,wartosc,default) ! implicit none character*(*) :: rekord,lancuch,wartosc,default character(len=80) :: aux integer :: lenlan,lenrec,iread,ireade !el external ilen !el logical iblnk !el external iblnk lenlan=ilen(lancuch) lenrec=ilen(rekord) iread=index(rekord,lancuch(:lenlan)//"=") ! print *,"rekord",rekord," lancuch",lancuch ! print *,"iread",iread," lenlan",lenlan," lenrec",lenrec if (iread.eq.0) then wartosc=default return endif iread=iread+lenlan+1 ! print *,"iread",iread ! print *,"|",rekord(iread:iread),"|",iblnk(rekord(iread:iread)) do while (iread.le.lenrec .and. iblnk(rekord(iread:iread))) iread=iread+1 ! print *,"|",rekord(iread:iread),"|",iblnk(rekord(iread:iread)) enddo ! print *,"iread",iread if (iread.gt.lenrec) then wartosc=default return endif ireade=iread+1 ! print *,"ireade",ireade do while (ireade.lt.lenrec .and. & .not.iblnk(rekord(ireade:ireade))) ireade=ireade+1 enddo wartosc=rekord(iread:ireade) return end subroutine reads #endif !----------------------------------------------------------------------------- ! permut.F !----------------------------------------------------------------------------- subroutine permut(isym) use geometry_data, only: tabperm ! implicit real*8 (a-h,o-z) ! include 'DIMENSIONS' ! include 'COMMON.LOCAL' ! include 'COMMON.VAR' ! include 'COMMON.CHAIN' ! include 'COMMON.INTERACT' ! include 'COMMON.IOUNITS' ! include 'COMMON.GEO' ! include 'COMMON.NAMES' ! include 'COMMON.CONTROL' integer :: n,isym ! logical nextp !el external nextp integer,dimension(isym) :: a ! parameter(n=symetr) !el local variables integer :: kkk,i n=isym if (n.eq.1) then tabperm(1,1)=1 return endif kkk=0 do i=1,n a(i)=i enddo 10 print *,(a(i),i=1,n) kkk=kkk+1 do i=1,n tabperm(kkk,i)=a(i) ! write (iout,*) "tututu", kkk enddo if(nextp(n,a)) go to 10 return end subroutine permut !----------------------------------------------------------------------------- logical function nextp(n,a) integer :: n,i,j,k,t ! logical :: nextp integer,dimension(n) :: a i=n-1 10 if(a(i).lt.a(i+1)) go to 20 i=i-1 if(i.eq.0) go to 20 go to 10 20 j=i+1 k=n 30 t=a(j) a(j)=a(k) a(k)=t j=j+1 k=k-1 if(j.lt.k) go to 30 j=i if(j.ne.0) go to 40 nextp=.false. return 40 j=j+1 if(a(j).lt.a(i)) go to 40 t=a(i) a(i)=a(j) a(j)=t nextp=.true. return end function nextp !----------------------------------------------------------------------------- !----------------------------------------------------------------------------- end module io_base