implicit none
!-----------------------------------------------------------------------------
! Max. number of AA residues
- integer,parameter :: maxres=4000!1200
+ 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
!-----------------------------------------------------------------------------
!
!
write(iout,*) ' iss:',(iss(i),i=1,ns)
! Check whether the specified bridging residues are cystines.
do i=1,ns
+ write(iout,*) i,iss(i)
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 ',&
#endif
endif
enddo
+ if (dyn_ss) then
+ if(.not.allocated(ihpb)) allocate(ihpb(ns))
+ if(.not.allocated(jhpb)) allocate(jhpb(ns))
+ endif
! Read preformed bridges.
if (ns.gt.0) then
- read (inp,*) nss
+ read (inp,*) nss
if (nss.gt.0) then
if(.not.allocated(ihpb)) allocate(ihpb(nss))
if(.not.allocated(jhpb)) allocate(jhpb(nss))
enddo
endif
endif
-!write(iout,*) "end read_bridge"
+! write(iout,*) "end read_bridge"
return
end subroutine read_bridge
!-----------------------------------------------------------------------------
! include 'COMMON.CONTROL'
! include 'COMMON.LOCAL'
! include 'COMMON.INTERACT'
+! Read coordinates from input
!
!el local variables
integer :: l,k,j,i,kanal
-! Read coordinates from input
-!
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)
end subroutine read_threadbase
!-----------------------------------------------------------------------------
#ifdef WHAM_RUN
- subroutine read_angles(kanal,iscor,energ,iprot,*)
+!el subroutine read_angles(kanal,iscor,energ,iprot,*)
+ subroutine read_angles(kanal,*)
use geometry_data
use energy_data
subroutine read_angles(kanal,*)
use geometry_data
+ ! use energy
+ ! use control
! implicit real*8 (a-h,o-z)
! include 'DIMENSIONS'
! include 'COMMON.GEO'
card=card(:ilen(card)+1)//karta
return
end subroutine card_concat
+!----------------------------------------------------------------------------
+ subroutine read_afminp
+ use geometry_data
+ use energy_data
+ use control_data, only:out1file
+ use MPI_data
+ character*320 afmcard
+ integer i
+ print *, "wchodze"
+ call card_concat(afmcard,.true.)
+ call readi(afmcard,"BEG",afmbeg,0)
+ call readi(afmcard,"END",afmend,0)
+ call reada(afmcard,"FORCE",forceAFMconst,0.0d0)
+ call reada(afmcard,"VEL",velAFMconst,0.0d0)
+ print *,'FORCE=' ,forceAFMconst
+!------ NOW PROPERTIES FOR AFM
+ distafminit=0.0d0
+ do i=1,3
+ distafminit=(c(i,afmend)-c(i,afmbeg))**2+distafminit
+ enddo
+ distafminit=dsqrt(distafminit)
+ print *,'initdist',distafminit
+ return
+ end subroutine read_afminp
!-----------------------------------------------------------------------------
subroutine read_dist_constr
use MPI_data
+ ! use control
use geometry, only: dist
use geometry_data
use control_data
! write (iout,*) "Calling read_dist_constr"
! write (iout,*) "nres",nres," nstart_sup",nstart_sup," nsup",nsup
! call flush(iout)
+ if(.not.allocated(ihpb)) allocate(ihpb(maxdim))
+ if(.not.allocated(jhpb)) allocate(jhpb(maxdim))
+ if(.not.allocated(dhpb)) allocate(dhpb(maxdim))
+ if(.not.allocated(dhpb1)) allocate(dhpb1(maxdim))
+ if(.not.allocated(forcon)) allocate(forcon(maxdim))
+ if(.not.allocated(fordepth)) allocate(fordepth(maxdim))
+ if ((genconstr.gt.0).and.(constr_dist.eq.11)) then
+ call gen_dist_constr2
+ go to 1712
+ endif
call card_concat(controlcard,.true.)
call readi(controlcard,"NFRAG",nfrag_,0)
call readi(controlcard,"NPAIR",npair_,0)
! 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))
-
+! 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
endif
enddo
do i=1,ndist_
- read (inp,*) ihpb(nhpb+1),jhpb(nhpb+1),forcon(nhpb+1)
+! 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))
+ if (constr_dist.eq.11) then
+ read (inp,*) ihpb(nhpb+1),jhpb(nhpb+1),dhpb(i),dhpb1(i), &
+ ibecarb(i),forcon(nhpb+1),fordepth(nhpb+1)
+ fordepth(nhpb+1)=fordepth(nhpb+1)/forcon(nhpb+1)
+ else
+!C print *,"in else"
+ read (inp,*) ihpb(nhpb+1),jhpb(nhpb+1),dhpb(i),dhpb1(i), &
+ ibecarb(i),forcon(nhpb+1)
+ endif
if (forcon(nhpb+1).gt.0.0d0) then
nhpb=nhpb+1
- dhpb(nhpb)=dist(ihpb(nhpb),jhpb(nhpb))
+ if (ibecarb(i).gt.0) then
+ ihpb(i)=ihpb(i)+nres
+ jhpb(i)=jhpb(i)+nres
+ endif
+ if (dhpb(nhpb).eq.0.0d0) &
+ dhpb(nhpb)=dist(ihpb(nhpb),jhpb(nhpb))
+ endif
+
#ifdef MPI
if (.not.out1file .or. me.eq.king) &
write (iout,'(a,3i5,f8.2,f10.1)') "+dist.constr ",&
write (iout,'(a,3i5,f8.2,f10.1)') "+dist.constr ",&
nhpb,ihpb(nhpb),jhpb(nhpb),dhpb(nhpb),forcon(nhpb)
#endif
- endif
enddo
+ 1712 continue
call flush(iout)
return
end subroutine read_dist_constr
!-----------------------------------------------------------------------------
+ subroutine gen_dist_constr2
+ use MPI_data
+ ! use control
+ use geometry, only: dist
+ use geometry_data
+ use control_data
+ use energy_data
+ integer :: i,j
+ real(kind=8) :: distance
+ if (constr_dist.eq.11) then
+ do i=nstart_sup,nstart_sup+nsup-1
+ do j=i+2,nstart_sup+nsup-1
+ distance=dist(i,j)
+ if (distance.le.15.0) then
+ nhpb=nhpb+1
+ ihpb(nhpb)=i+nstart_seq-nstart_sup
+ jhpb(nhpb)=j+nstart_seq-nstart_sup
+ forcon(nhpb)=sqrt(0.04*distance)
+ fordepth(nhpb)=sqrt(40.0/distance)
+ dhpb(nhpb)=distance-0.1d0
+ dhpb1(nhpb)=distance+0.1d0
+
+#ifdef MPI
+ if (.not.out1file .or. me.eq.king) &
+ write (iout,'(a,3i5,f8.2,f10.2)') "+dist.constr ", &
+ nhpb,ihpb(nhpb),jhpb(nhpb),dhpb(nhpb),forcon(nhpb)
+#else
+ write (iout,'(a,3i5,f8.2,f10.2)') "+dist.constr ", &
+ nhpb,ihpb(nhpb),jhpb(nhpb),dhpb(nhpb),forcon(nhpb)
+#endif
+ endif
+ enddo
+ enddo
+ else
+ do i=nstart_sup,nstart_sup+nsup-1
+ do j=i+2,nstart_sup+nsup-1
+ nhpb=nhpb+1
+ ihpb(nhpb)=i+nstart_seq-nstart_sup
+ jhpb(nhpb)=j+nstart_seq-nstart_sup
+ forcon(nhpb)=weidis
+ dhpb(nhpb)=dist(i,j)
+ enddo
+ enddo
+ endif
+ return
+ end subroutine gen_dist_constr2
+
+!-----------------------------------------------------------------------------
#ifdef WINIFL
subroutine flush(iu)
return
use geometry_data, only: c,nres
use energy_data
+ ! use control
use compare_data
use MD_data
! implicit real*8 (a-h,o-z)
character*(*) :: tytul
character(len=1),dimension(10) :: chainid= (/'A','B','C','D','E','F','G','H','I','J'/)
integer,dimension(nres) :: ica !(maxres)
-
+ integer iti1
!el local variables
integer :: j,iti,itj,itk,itl,i,iatom,ichain,ires,iunit
real(kind=8) :: etot
ires=0
do i=nnt,nct
iti=itype(i)
+ iti1=itype(i+1)
+ if ((iti.eq.ntyp1).and.(iti1.eq.ntyp1)) cycle
if (iti.eq.ntyp1) then
ichain=ichain+1
- ires=0
+! ires=0
write (iunit,'(a)') 'TER'
else
ires=ires+1
! format.
use geometry_data, only: c
use energy_data
+ ! use control
! implicit real*8 (a-h,o-z)
! include 'DIMENSIONS'
! include 'COMMON.CHAIN'
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'
! print '(a,i5)',intname,igeom
!el local variables
integer :: i,it
- real(kind=8) :: ener
-#ifdef WHAM_RUN
- integer :: iii
-#endif
+ real(kind=8) :: ener,free
+! character(len=80) :: plik
+! integer :: iii
#if defined(AIX) || defined(PGI)
open (igeom,file=intname,position='append')
open (igeom,file=intname,access='append')
#endif
#ifdef WHAM_RUN
- iii=igeom
+! iii=igeom
igeom=iout
#endif
+ print *,nss
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)
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