+++ /dev/null
- module io
-!-----------------------------------------------------------------------
- use io_units
- use names
- use io_base
- use io_config
- implicit none
-!-----------------------------------------------------------------------------
-!
-!
-!-----------------------------------------------------------------------------
- contains
-!-----------------------------------------------------------------------------
-! bank.F io_csa
-!-----------------------------------------------------------------------------
- subroutine write_csa_pdb(var,ene,nft,ik,iw_pdb)
-
- use csa_data
- use geometry_data, only:nres,nvar
- use geometry, only:var_to_geom,chainbuild
- use compare, only:secondary2
-! implicit real*8 (a-h,o-z)
-! include 'DIMENSIONS'
-! include 'COMMON.CSA'
-! include 'COMMON.BANK'
-! include 'COMMON.VAR'
-! include 'COMMON.IOUNITS'
-! include 'COMMON.MINIM'
-! include 'COMMON.SETUP'
-! include 'COMMON.GEO'
-! include 'COMMON.CHAIN'
-! include 'COMMON.LOCAL'
-! include 'COMMON.INTERACT'
-! include 'COMMON.NAMES'
-! include 'COMMON.SBRIDGE'
- integer :: lenpre,lenpot !,ilen
-!el external ilen
- real(kind=8),dimension(nvar) :: var !(maxvar) (maxvar=6*maxres)
- character(len=50) :: titelloc
- character(len=3) :: zahl
- real(kind=8),dimension(mxch*(mxch+1)/2+1) :: ene
-!el local variables
- integer :: nft,ik,iw_pdb
-
- nmin_csa=nmin_csa+1
- if(ene(1).lt.eglob_csa) then
- eglob_csa=ene(1)
- nglob_csa=nglob_csa+1
- call numstr(nglob_csa,zahl)
-
- call var_to_geom(nvar,var)
- call chainbuild
- call secondary2(.false.)
-
- lenpre=ilen(prefix)
- open(icsa_pdb,file=prefix(:lenpre)//'@'//zahl//'.pdb')
-
- if (iw_pdb.eq.1) then
- write(titelloc,'(a2,i3,a3,i9,a3,i6)') &
- 'GM',nglob_csa,' e ',nft,' m ',nmin_csa
- else
- write(titelloc,'(a2,i3,a3,i9,a3,i6,a5,f5.2,a5,f5.1)') &
- 'GM',nglob_csa,' e ',nft,' m ',nmin_csa,' rms ',&
- rmsn(ik),' %NC ',pncn(ik)*100
- endif
- call pdbout(eglob_csa,titelloc,icsa_pdb)
- close(icsa_pdb)
- endif
-
- return
- end subroutine write_csa_pdb
-!-----------------------------------------------------------------------------
-! csa.f io_csa
-!-----------------------------------------------------------------------------
- subroutine from_pdb(n,idum)
-! This subroutine stores the UNRES int variables generated from
-! subroutine readpdb into the 1st conformation of in dihang_in.
-! Subsequent n-1 conformations of dihang_in have identical values
-! of theta and phi as the 1st conformation but random values for
-! alph and omeg.
-! The array cref (also generated from subroutine readpdb) is stored
-! to crefjlee to be used for rmsd calculation in CSA, if necessary.
-
- use csa_data
- use geometry_data
- use random, only: ran1
-! implicit real*8 (a-h,o-z)
-! include 'DIMENSIONS'
-! include 'COMMON.IOUNITS'
-! include 'COMMON.CHAIN'
-! include 'COMMON.VAR'
-! include 'COMMON.BANK'
-! include 'COMMON.GEO'
-!el local variables
- integer :: n,idum,m,i,j,k,kk,kkk
- real(kind=8) :: e
-
- m=1
- do j=2,nres-1
- dihang_in(1,j,1,m)=theta(j+1)
- dihang_in(2,j,1,m)=phi(j+2)
- dihang_in(3,j,1,m)=alph(j)
- dihang_in(4,j,1,m)=omeg(j)
- enddo
- dihang_in(2,nres-1,1,k)=0.0d0
-
- do m=2,n
- do k=2,nres-1
- dihang_in(1,k,1,m)=dihang_in(1,k,1,1)
- dihang_in(2,k,1,m)=dihang_in(2,k,1,1)
- if(dabs(dihang_in(3,k,1,1)).gt.1.d-6) then
- dihang_in(3,k,1,m)=90.d0*ran1(idum)+90.d0
- dihang_in(3,k,1,m)=dihang_in(3,k,1,m)*deg2rad
- endif
- if(dabs(dihang_in(4,k,1,1)).gt.1.d-6) then
- dihang_in(4,k,1,m)=360.d0*ran1(idum)-180.d0
- dihang_in(4,k,1,m)=dihang_in(4,k,1,m)*deg2rad
- endif
- enddo
- enddo
-
-! Store cref to crefjlee (they are in COMMON.CHAIN).
- do k=1,2*nres
- do kk=1,3
- kkk=1
- crefjlee(kk,k)=cref(kk,k,kkk)
- enddo
- enddo
-
- open(icsa_native_int,file=csa_native_int,status="old")
- do m=1,n
- write(icsa_native_int,*) m,e
- write(icsa_native_int,200) &
- (dihang_in(1,k,1,m)*rad2deg,k=2,nres-1)
- write(icsa_native_int,200) &
- (dihang_in(2,k,1,m)*rad2deg,k=2,nres-2)
- write(icsa_native_int,200) &
- (dihang_in(3,k,1,m)*rad2deg,k=2,nres-1)
- write(icsa_native_int,200) &
- (dihang_in(4,k,1,m)*rad2deg,k=2,nres-1)
- enddo
-
- do k=1,nres
- write(icsa_native_int,200) (crefjlee(i,k),i=1,3)
- enddo
- close(icsa_native_int)
-
- 200 format (8f10.4)
-
- return
- end subroutine from_pdb
-!-----------------------------------------------------------------------------
- subroutine from_int(n,mm,idum)
-
- use csa_data
- use geometry_data
- use energy_data
- use geometry, only:chainbuild,gen_side
- use energy, only:etotal
- use compare
-! implicit real*8 (a-h,o-z)
-! include 'DIMENSIONS'
-! include 'COMMON.IOUNITS'
-! include 'COMMON.CHAIN'
-! include 'COMMON.VAR'
-! include 'COMMON.INTERACT'
-! include 'COMMON.BANK'
-! include 'COMMON.GEO'
-! include 'COMMON.CONTACTS'
-! integer ilen
-!el external ilen
- logical :: fail
- real(kind=8),dimension(0:n_ene) :: energia
-!el local variables
- integer :: n,mm,idum,i,ii,j,m,k,kk,maxcount_fail,icount_fail,maxsi
- real(kind=8) :: co
-
- open(icsa_native_int,file=csa_native_int,status="old")
- read (icsa_native_int,*)
- call read_angles(icsa_native_int,*10)
- goto 11
- 10 write (iout,'(2a)') "CHUJ NASTAPIL - error in ",&
- csa_native_int(:ilen(csa_native_int))
- 11 continue
- call intout
- do j=2,nres-1
- dihang_in(1,j,1,1)=theta(j+1)
- dihang_in(2,j,1,1)=phi(j+2)
- dihang_in(3,j,1,1)=alph(j)
- dihang_in(4,j,1,1)=omeg(j)
- enddo
- dihang_in(2,nres-1,1,1)=0.0d0
-
-! read(icsa_native_int,*) ind,e
-! read(icsa_native_int,200) (dihang_in(1,k,1,1),k=2,nres-1)
-! read(icsa_native_int,200) (dihang_in(2,k,1,1),k=2,nres-2)
-! read(icsa_native_int,200) (dihang_in(3,k,1,1),k=2,nres-1)
-! read(icsa_native_int,200) (dihang_in(4,k,1,1),k=2,nres-1)
-! dihang_in(2,nres-1,1,1)=0.d0
-
- maxsi=100
- maxcount_fail=100
-
- do m=mm+2,n
-! do k=2,nres-1
-! dihang_in(1,k,1,m)=dihang_in(1,k,1,1)
-! dihang_in(2,k,1,m)=dihang_in(2,k,1,1)
-! if(abs(dihang_in(3,k,1,1)).gt.1.d-3) then
-! dihang_in(3,k,1,m)=90.d0*ran1(idum)+90.d0
-! endif
-! if(abs(dihang_in(4,k,1,1)).gt.1.d-3) then
-! dihang_in(4,k,1,m)=360.d0*ran1(idum)-180.d0
-! endif
-! enddo
-! call intout
- fail=.true.
-
- icount_fail=0
-
- DO WHILE (FAIL .AND. ICOUNT_FAIL .LE. MAXCOUNT_FAIL)
-
- do i=nnt,nct
- if (itype(i).ne.10) then
-!d print *,'i=',i,' itype=',itype(i),' theta=',theta(i+1)
- fail=.true.
- ii=0
- do while (fail .and. ii .le. maxsi)
- call gen_side(itype(i),theta(i+1),alph(i),omeg(i),fail)
- ii = ii+1
- enddo
- endif
- enddo
- call chainbuild
- call etotal(energia)
- fail = (energia(0).ge.1.0d20)
- icount_fail=icount_fail+1
-
- ENDDO
-
- if (icount_fail.gt.maxcount_fail) then
- write (iout,*) &
- 'Failed to generate non-overlaping near-native conf.',&
- m
- endif
-
- do j=2,nres-1
- dihang_in(1,j,1,m)=theta(j+1)
- dihang_in(2,j,1,m)=phi(j+2)
- dihang_in(3,j,1,m)=alph(j)
- dihang_in(4,j,1,m)=omeg(j)
- enddo
- dihang_in(2,nres-1,1,m)=0.0d0
- enddo
-
-! do m=1,n
-! write(icsa_native_int,*) m,e
-! write(icsa_native_int,200) (dihang_in(1,k,1,m),k=2,nres-1)
-! write(icsa_native_int,200) (dihang_in(2,k,1,m),k=2,nres-2)
-! write(icsa_native_int,200) (dihang_in(3,k,1,m),k=2,nres-1)
-! write(icsa_native_int,200) (dihang_in(4,k,1,m),k=2,nres-1)
-! enddo
-! close(icsa_native_int)
-
-! do m=mm+2,n
-! do i=1,4
-! do j=2,nres-1
-! dihang_in(i,j,1,m)=dihang_in(i,j,1,m)*deg2rad
-! enddo
-! enddo
-! enddo
-
- call dihang_to_c(dihang_in(1,1,1,1))
-
-! Store c to cref (they are in COMMON.CHAIN).
- do k=1,2*nres
- do kk=1,3
- crefjlee(kk,k)=c(kk,k)
- enddo
- enddo
-
- call contact(.true.,ncont_ref,icont_ref,co)
-
-! do k=1,nres
-! write(icsa_native_int,200) (crefjlee(i,k),i=1,3)
-! enddo
- close(icsa_native_int)
-
- 200 format (8f10.4)
-
- return
- end subroutine from_int
-!-----------------------------------------------------------------------------
- subroutine dihang_to_c(aarray)
-
- use geometry_data
- use csa_data
- use geometry, only:chainbuild
-! implicit real*8 (a-h,o-z)
-! include 'DIMENSIONS'
-! include 'COMMON.CSA'
-! include 'COMMON.BANK'
-! include 'COMMON.CHAIN'
-! include 'COMMON.GEO'
-! include 'COMMON.VAR'
- integer :: i
- real(kind=8),dimension(mxang,nres,mxch) :: aarray !(mxang,maxres,mxch)
-
-! do i=4,nres
-! phi(i)=dihang_in(1,i-2,1,1)
-! enddo
- do i=2,nres-1
- theta(i+1)=aarray(1,i,1)
- phi(i+2)=aarray(2,i,1)
- alph(i)=aarray(3,i,1)
- omeg(i)=aarray(4,i,1)
- enddo
-
- call chainbuild
-
- return
- end subroutine dihang_to_c
-!-----------------------------------------------------------------------------
-! geomout.F
-!-----------------------------------------------------------------------------
-#ifdef NOXDR
- subroutine cartout(time)
-#else
- subroutine cartoutx(time)
-#endif
- use geometry_data, only: c,nres
- use energy_data
- use MD_data, only: potE,t_bath
-! 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'
- real(kind=8) :: time
-!el local variables
- integer :: j,k,i
-
-#if defined(AIX) || defined(PGI)
- open(icart,file=cartname,position="append")
-#else
- open(icart,file=cartname,access="append")
-#endif
- write (icart,'(e15.8,2e15.5,f12.5,$)') time,potE,uconst,t_bath
- if (dyn_ss) then
- write (icart,'(i4,$)') &
- nss,(idssb(j)+nres,jdssb(j)+nres,j=1,nss)
- else
- write (icart,'(i4,$)') &
- nss,(ihpb(j),jhpb(j),j=1,nss)
- endif
- write (icart,'(i4,20f7.4)') nfrag+npair+3*nfrag_back,&
- (qfrag(i),i=1,nfrag),(qpair(i),i=1,npair),&
- (utheta(i),ugamma(i),uscdiff(i),i=1,nfrag_back)
- write (icart,'(8f10.5)') &
- ((c(k,j),k=1,3),j=1,nres),&
- ((c(k,j+nres),k=1,3),j=nnt,nct)
- close(icart)
- return
-
-#ifdef NOXDR
- end subroutine cartout
-#else
- end subroutine cartoutx
-#endif
-!-----------------------------------------------------------------------------
-#ifndef NOXDR
- subroutine cartout(time)
-! implicit real*8 (a-h,o-z)
-! include 'DIMENSIONS'
- use geometry_data, only: c,nres
- use energy_data
- use MD_data, only: potE,t_bath
-#ifdef MPI
- use MPI_data
- include 'mpif.h'
-! include 'COMMON.SETUP'
-#else
- integer,parameter :: me=0
-#endif
-! include 'COMMON.CHAIN'
-! include 'COMMON.INTERACT'
-! include 'COMMON.NAMES'
-! include 'COMMON.IOUNITS'
-! include 'COMMON.HEADER'
-! include 'COMMON.SBRIDGE'
-! include 'COMMON.DISTFIT'
-! include 'COMMON.MD'
- real(kind=8) :: time
- integer :: iret,itmp
- real(kind=4) :: prec
- real(kind=4),dimension(3,2*nres+2) :: xcoord !(3,maxres2+2) (maxres2=2*maxres
-!el local variables
- integer :: j,i,ixdrf
-
-#ifdef AIX
- call xdrfopen_(ixdrf,cartname, "a", iret)
- call xdrffloat_(ixdrf, real(time), iret)
- call xdrffloat_(ixdrf, real(potE), iret)
- call xdrffloat_(ixdrf, real(uconst), iret)
- call xdrffloat_(ixdrf, real(uconst_back), iret)
- call xdrffloat_(ixdrf, real(t_bath), iret)
- call xdrfint_(ixdrf, nss, iret)
- do j=1,nss
- if (dyn_ss) then
- call xdrfint_(ixdrf, idssb(j)+nres, iret)
- call xdrfint_(ixdrf, jdssb(j)+nres, iret)
- else
- call xdrfint_(ixdrf, ihpb(j), iret)
- call xdrfint_(ixdrf, jhpb(j), iret)
- endif
- enddo
- call xdrfint_(ixdrf, nfrag+npair+3*nfrag_back, iret)
- do i=1,nfrag
- call xdrffloat_(ixdrf, real(qfrag(i)), iret)
- enddo
- do i=1,npair
- call xdrffloat_(ixdrf, real(qpair(i)), iret)
- enddo
- do i=1,nfrag_back
- call xdrffloat_(ixdrf, real(utheta(i)), iret)
- call xdrffloat_(ixdrf, real(ugamma(i)), iret)
- call xdrffloat_(ixdrf, real(uscdiff(i)), iret)
- enddo
-#else
- call xdrfopen(ixdrf,cartname, "a", iret)
- call xdrffloat(ixdrf, real(time), iret)
- call xdrffloat(ixdrf, real(potE), iret)
- call xdrffloat(ixdrf, real(uconst), iret)
- call xdrffloat(ixdrf, real(uconst_back), iret)
- call xdrffloat(ixdrf, real(t_bath), iret)
- call xdrfint(ixdrf, nss, iret)
- do j=1,nss
- if (dyn_ss) then
- call xdrfint(ixdrf, idssb(j)+nres, iret)
- call xdrfint(ixdrf, jdssb(j)+nres, iret)
- else
- call xdrfint(ixdrf, ihpb(j), iret)
- call xdrfint(ixdrf, jhpb(j), iret)
- endif
- enddo
- call xdrfint(ixdrf, nfrag+npair+3*nfrag_back, iret)
- do i=1,nfrag
- call xdrffloat(ixdrf, real(qfrag(i)), iret)
- enddo
- do i=1,npair
- call xdrffloat(ixdrf, real(qpair(i)), iret)
- enddo
- do i=1,nfrag_back
- call xdrffloat(ixdrf, real(utheta(i)), iret)
- call xdrffloat(ixdrf, real(ugamma(i)), iret)
- call xdrffloat(ixdrf, real(uscdiff(i)), iret)
- enddo
-#endif
- prec=10000.0
- do i=1,nres
- do j=1,3
- xcoord(j,i)=c(j,i)
- enddo
- enddo
- do i=nnt,nct
- do j=1,3
- xcoord(j,nres+i-nnt+1)=c(j,i+nres)
- enddo
- enddo
-
- itmp=nres+nct-nnt+1
-#ifdef AIX
- call xdrf3dfcoord_(ixdrf, xcoord, itmp, prec, iret)
- call xdrfclose_(ixdrf, iret)
-#else
- call xdrf3dfcoord(ixdrf, xcoord, itmp, prec, iret)
- call xdrfclose(ixdrf, iret)
-#endif
- return
- end subroutine cartout
-#endif
-!-----------------------------------------------------------------------------
- subroutine statout(itime)
-
- use energy_data
- use control_data, only:refstr
- use MD_data
- use MPI_data
- use compare, only:rms_nac_nnc
-! implicit real*8 (a-h,o-z)
-! include 'DIMENSIONS'
-! include 'COMMON.CONTROL'
-! include 'COMMON.CHAIN'
-! include 'COMMON.INTERACT'
-! include 'COMMON.NAMES'
-! include 'COMMON.IOUNITS'
-! include 'COMMON.HEADER'
-! include 'COMMON.SBRIDGE'
-! include 'COMMON.DISTFIT'
-! include 'COMMON.MD'
-! include 'COMMON.REMD'
-! include 'COMMON.SETUP'
- integer :: itime
- real(kind=8),dimension(0:n_ene) :: energia
-! double precision gyrate
-!el external gyrate
-!el common /gucio/ cm
- character(len=256) :: line1,line2
- character(len=4) :: format1,format2
- character(len=30) :: format
-!el local variables
- integer :: i,ii1,ii2
- real(kind=8) :: rms,frac,frac_nn,co
-
-#ifdef AIX
- if(itime.eq.0) then
- open(istat,file=statname,position="append")
- endif
-#else
-#ifdef PGI
- open(istat,file=statname,position="append")
-#else
- open(istat,file=statname,access="append")
-#endif
-#endif
- if (refstr) then
- call rms_nac_nnc(rms,frac,frac_nn,co,.false.)
- write (line1,'(i10,f15.2,3f12.3,f7.2,4f6.3,3f12.3,i5,$)') &
- itime,totT,EK,potE,totE,&
- rms,frac,frac_nn,co,amax,kinetic_T,t_bath,gyrate(),me
- format1="a133"
- else
- write (line1,'(i10,f15.2,7f12.3,i5,$)') &
- itime,totT,EK,potE,totE,&
- amax,kinetic_T,t_bath,gyrate(),me
- format1="a114"
- endif
- if(usampl.and.totT.gt.eq_time) then
- write(line2,'(i5,2f9.4,300f7.4)') iset,uconst,uconst_back,&
- (qfrag(ii1),ii1=1,nfrag),(qpair(ii2),ii2=1,npair),&
- (utheta(i),ugamma(i),uscdiff(i),i=1,nfrag_back)
- write(format2,'(a1,i3.3)') "a",23+7*nfrag+7*npair &
- +21*nfrag_back
- else
- format2="a001"
- line2=' '
- endif
- if (print_compon) then
- if(itime.eq.0) then
- write(format,'(a1,a4,a1,a4,a10)') "(",format1,",",format2,&
- ",20a12)"
- write (istat,format) "#","",&
- (ename(print_order(i)),i=1,nprint_ene)
- endif
- write(format,'(a1,a4,a1,a4,a10)') "(",format1,",",format2,&
- ",20f12.3)"
- write (istat,format) line1,line2,&
- (potEcomp(print_order(i)),i=1,nprint_ene)
- else
- write(format,'(a1,a4,a1,a4,a1)') "(",format1,",",format2,")"
- write (istat,format) line1,line2
- endif
-#if defined(AIX)
- call flush(istat)
-#else
- close(istat)
-#endif
- return
- end subroutine statout
-!-----------------------------------------------------------------------------
-! readrtns_CSA.F
-!-----------------------------------------------------------------------------
- subroutine readrtns
-
- use control_data
- use energy_data
- use MPI_data
- use muca_md, only:read_muca
-! implicit real*8 (a-h,o-z)
-! include 'DIMENSIONS'
-#ifdef MPI
- include 'mpif.h'
-#endif
-! include 'COMMON.SETUP'
-! include 'COMMON.CONTROL'
-! include 'COMMON.SBRIDGE'
-! include 'COMMON.IOUNITS'
- logical :: file_exist
- integer :: i
-! Read force-field parameters except weights
- call parmread
-! Read job setup parameters
- call read_control
-! Read control parameters for energy minimzation if required
- if (minim) call read_minim
-! Read MCM control parameters if required
- if (modecalc.eq.3 .or. modecalc.eq.6) call mcmread
-! Read MD control parameters if reqjuired
- if (modecalc.eq.12) call read_MDpar
-! Read MREMD control parameters if required
- if (modecalc.eq.14) then
- call read_MDpar
- call read_REMDpar
- endif
-! Read MUCA control parameters if required
- if (lmuca) call read_muca
-! Read CSA control parameters if required (from fort.40 if exists
-! otherwise from general input file)
- if (modecalc.eq.8) then
- inquire (file="fort.40",exist=file_exist)
- if (.not.file_exist) call csaread
- endif
-!fmc if (modecalc.eq.10) call mcmfread
-! Read molecule information, molecule geometry, energy-term weights, and
-! restraints if requested
- call molread
-! Print restraint information
-#ifdef MPI
- if (.not. out1file .or. me.eq.king) then
-#endif
- if (nhpb.gt.nss) &
- write (iout,'(a,i5,a)') "The following",nhpb-nss,&
- " distance constraints have been imposed"
- do i=nss+1,nhpb
- write (iout,'(3i6,f10.5)') i-nss,ihpb(i),jhpb(i),forcon(i)
- enddo
-#ifdef MPI
- endif
-#endif
-! print *,"Processor",myrank," leaves READRTNS"
-! write(iout,*) "end readrtns"
- return
- end subroutine readrtns
-!-----------------------------------------------------------------------------
- subroutine molread
-!
-! Read molecular data.
-!
-! use control, only: ilen
- use control_data
- use geometry_data
- use energy_data
- use energy
- use compare_data
- use MD_data, only: t_bath
- use MPI_data
- use compare, only:seq_comp,contact
- use control
-! implicit real*8 (a-h,o-z)
-! include 'DIMENSIONS'
-#ifdef MPI
- include 'mpif.h'
- integer :: error_msg,ierror,ierr,ierrcode
-#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.CONTACTS'
-! include 'COMMON.TORCNSTR'
-! include 'COMMON.TIME1'
-! include 'COMMON.BOUNDS'
-! include 'COMMON.MD'
-! include 'COMMON.SETUP'
- character(len=4),dimension(:),allocatable :: sequence !(maxres)
-! integer :: rescode
-! double precision x(maxvar)
- character(len=256) :: pdbfile
- character(len=320) :: weightcard
- character(len=80) :: weightcard_t!,ucase
-! integer,dimension(:),allocatable :: itype_pdb !(maxres)
-! common /pizda/ itype_pdb
- logical :: fail !seq_comp,
- real(kind=8) :: energia(0:n_ene)
-! integer ilen
-!el external ilen
-!el local varables
- integer :: i,j,l,k,kkk,ii,i1,i2,it1,it2
-
- real(kind=8),dimension(3,maxres2+2) :: c_alloc
- real(kind=8),dimension(3,0:maxres2) :: dc_alloc
- integer,dimension(maxres) :: itype_alloc
-
- integer :: iti,nsi,maxsi,itrial,itmp
- real(kind=8) :: wlong,scalscp,co
- allocate(weights(n_ene))
-!-----------------------------
- allocate(c(3,2*maxres+2)) !(3,maxres2+2) maxres2=2*maxres
- allocate(dc(3,0:2*maxres)) !(3,0:maxres2)
- allocate(itype(maxres)) !(maxres)
-!
-! Zero out tables.
-!
- c(:,:)=0.0D0
- dc(:,:)=0.0D0
- itype(:)=0
-!-----------------------------
-!
-! Body
-!
-! Read weights of the subsequent energy terms.
- call card_concat(weightcard,.true.)
- call reada(weightcard,'WLONG',wlong,1.0D0)
- call reada(weightcard,'WSC',wsc,wlong)
- call reada(weightcard,'WSCP',wscp,wlong)
- call reada(weightcard,'WELEC',welec,1.0D0)
- call reada(weightcard,'WVDWPP',wvdwpp,welec)
- call reada(weightcard,'WEL_LOC',wel_loc,1.0D0)
- call reada(weightcard,'WCORR4',wcorr4,0.0D0)
- call reada(weightcard,'WCORR5',wcorr5,0.0D0)
- call reada(weightcard,'WCORR6',wcorr6,0.0D0)
- call reada(weightcard,'WTURN3',wturn3,1.0D0)
- call reada(weightcard,'WTURN4',wturn4,1.0D0)
- call reada(weightcard,'WTURN6',wturn6,1.0D0)
- call reada(weightcard,'WSCCOR',wsccor,1.0D0)
- call reada(weightcard,'WSTRAIN',wstrain,1.0D0)
- call reada(weightcard,'WBOND',wbond,1.0D0)
- call reada(weightcard,'WTOR',wtor,1.0D0)
- call reada(weightcard,'WTORD',wtor_d,1.0D0)
- call reada(weightcard,'WANG',wang,1.0D0)
- call reada(weightcard,'WSCLOC',wscloc,1.0D0)
- call reada(weightcard,'SCAL14',scal14,0.4D0)
- call reada(weightcard,'SCALSCP',scalscp,1.0d0)
- call reada(weightcard,'CUTOFF',cutoff_corr,7.0d0)
- call reada(weightcard,'DELT_CORR',delt_corr,0.5d0)
- call reada(weightcard,'TEMP0',temp0,300.0d0)
- if (index(weightcard,'SOFT').gt.0) ipot=6
-! 12/1/95 Added weight for the multi-body term WCORR
- call reada(weightcard,'WCORRH',wcorr,1.0D0)
- if (wcorr4.gt.0.0d0) wcorr=wcorr4
- weights(1)=wsc
- weights(2)=wscp
- weights(3)=welec
- weights(4)=wcorr
- weights(5)=wcorr5
- weights(6)=wcorr6
- weights(7)=wel_loc
- weights(8)=wturn3
- weights(9)=wturn4
- weights(10)=wturn6
- weights(11)=wang
- weights(12)=wscloc
- weights(13)=wtor
- weights(14)=wtor_d
- weights(15)=wstrain
- weights(16)=wvdwpp
- weights(17)=wbond
- weights(18)=scal14
- weights(21)=wsccor
- if(me.eq.king.or..not.out1file) &
- write (iout,10) wsc,wscp,welec,wvdwpp,wbond,wang,wscloc,wtor,&
- wtor_d,wstrain,wel_loc,wcorr,wcorr5,wcorr6,wsccor,wturn3,&
- wturn4,wturn6
- 10 format (/'Energy-term weights (unscaled):'// &
- 'WSCC= ',f10.6,' (SC-SC)'/ &
- 'WSCP= ',f10.6,' (SC-p)'/ &
- 'WELEC= ',f10.6,' (p-p electr)'/ &
- 'WVDWPP= ',f10.6,' (p-p VDW)'/ &
- 'WBOND= ',f10.6,' (stretching)'/ &
- 'WANG= ',f10.6,' (bending)'/ &
- 'WSCLOC= ',f10.6,' (SC local)'/ &
- 'WTOR= ',f10.6,' (torsional)'/ &
- 'WTORD= ',f10.6,' (double torsional)'/ &
- 'WSTRAIN=',f10.6,' (SS bridges & dist. cnstr.)'/ &
- 'WEL_LOC=',f10.6,' (multi-body 3-rd order)'/ &
- 'WCORR4= ',f10.6,' (multi-body 4th order)'/ &
- 'WCORR5= ',f10.6,' (multi-body 5th order)'/ &
- 'WCORR6= ',f10.6,' (multi-body 6th order)'/ &
- 'WSCCOR= ',f10.6,' (back-scloc correlation)'/ &
- 'WTURN3= ',f10.6,' (turns, 3rd order)'/ &
- 'WTURN4= ',f10.6,' (turns, 4th order)'/ &
- 'WTURN6= ',f10.6,' (turns, 6th order)')
- if(me.eq.king.or..not.out1file)then
- if (wcorr4.gt.0.0d0) then
- write (iout,'(/2a/)') 'Local-electrostatic type correlation ',&
- 'between contact pairs of peptide groups'
- write (iout,'(2(a,f5.3/))') &
- 'Cutoff on 4-6th order correlation terms: ',cutoff_corr,&
- 'Range of quenching the correlation terms:',2*delt_corr
- else if (wcorr.gt.0.0d0) then
- write (iout,'(/2a/)') 'Hydrogen-bonding correlation ',&
- 'between contact pairs of peptide groups'
- endif
- write (iout,'(a,f8.3)') &
- 'Scaling factor of 1,4 SC-p interactions:',scal14
- write (iout,'(a,f8.3)') &
- 'General scaling factor of SC-p interactions:',scalscp
- endif
- r0_corr=cutoff_corr-delt_corr
- do i=1,ntyp
- aad(i,1)=scalscp*aad(i,1)
- aad(i,2)=scalscp*aad(i,2)
- bad(i,1)=scalscp*bad(i,1)
- bad(i,2)=scalscp*bad(i,2)
- enddo
- call rescale_weights(t_bath)
- if(me.eq.king.or..not.out1file) &
- write (iout,22) wsc,wscp,welec,wvdwpp,wbond,wang,wscloc,wtor,&
- wtor_d,wstrain,wel_loc,wcorr,wcorr5,wcorr6,wsccor,wturn3,&
- wturn4,wturn6
- 22 format (/'Energy-term weights (scaled):'// &
- 'WSCC= ',f10.6,' (SC-SC)'/ &
- 'WSCP= ',f10.6,' (SC-p)'/ &
- 'WELEC= ',f10.6,' (p-p electr)'/ &
- 'WVDWPP= ',f10.6,' (p-p VDW)'/ &
- 'WBOND= ',f10.6,' (stretching)'/ &
- 'WANG= ',f10.6,' (bending)'/ &
- 'WSCLOC= ',f10.6,' (SC local)'/ &
- 'WTOR= ',f10.6,' (torsional)'/ &
- 'WTORD= ',f10.6,' (double torsional)'/ &
- 'WSTRAIN=',f10.6,' (SS bridges & dist. cnstr.)'/ &
- 'WEL_LOC=',f10.6,' (multi-body 3-rd order)'/ &
- 'WCORR4= ',f10.6,' (multi-body 4th order)'/ &
- 'WCORR5= ',f10.6,' (multi-body 5th order)'/ &
- 'WCORR6= ',f10.6,' (multi-body 6th order)'/ &
- 'WSCCOR= ',f10.6,' (back-scloc correlatkion)'/ &
- 'WTURN3= ',f10.6,' (turns, 3rd order)'/ &
- 'WTURN4= ',f10.6,' (turns, 4th order)'/ &
- 'WTURN6= ',f10.6,' (turns, 6th order)')
- if(me.eq.king.or..not.out1file) &
- write (iout,*) "Reference temperature for weights calculation:",&
- temp0
- call reada(weightcard,"D0CM",d0cm,3.78d0)
- call reada(weightcard,"AKCM",akcm,15.1d0)
- call reada(weightcard,"AKTH",akth,11.0d0)
- call reada(weightcard,"AKCT",akct,12.0d0)
- call reada(weightcard,"V1SS",v1ss,-1.08d0)
- call reada(weightcard,"V2SS",v2ss,7.61d0)
- call reada(weightcard,"V3SS",v3ss,13.7d0)
- call reada(weightcard,"EBR",ebr,-5.50D0)
- dyn_ss=(index(weightcard,'DYN_SS').gt.0)
-
- call reada(weightcard,"HT",Ht,0.0D0)
- if (dyn_ss) then
- ss_depth=ebr/wsc-0.25*eps(1,1)
- Ht=Ht/wsc-0.25*eps(1,1)
- akcm=akcm*wstrain/wsc
- akth=akth*wstrain/wsc
- akct=akct*wstrain/wsc
- v1ss=v1ss*wstrain/wsc
- v2ss=v2ss*wstrain/wsc
- v3ss=v3ss*wstrain/wsc
- else
- ss_depth=ebr/wstrain-0.25*eps(1,1)*wsc/wstrain
- endif
-
- if(me.eq.king.or..not.out1file) then
- write (iout,*) "Parameters of the SS-bond potential:"
- write (iout,*) "D0CM",d0cm," AKCM",akcm," AKTH",akth,&
- " AKCT",akct
- write (iout,*) "V1SS",v1ss," V2SS",v2ss," V3SS",v3ss
- write (iout,*) "EBR",ebr," SS_DEPTH",ss_depth
- write (iout,*)" HT",Ht
- print *,'indpdb=',indpdb,' pdbref=',pdbref
- endif
- if (indpdb.gt.0 .or. pdbref) then
- read(inp,'(a)') pdbfile
- if(me.eq.king.or..not.out1file) &
- write (iout,'(2a)') 'PDB data will be read from file ',&
- pdbfile(:ilen(pdbfile))
- open(ipdbin,file=pdbfile,status='old',err=33)
- goto 34
- 33 write (iout,'(a)') 'Error opening PDB file.'
- stop
- 34 continue
-! print *,'Begin reading pdb data'
- call readpdb
-! print *,'Finished reading pdb data'
- if(me.eq.king.or..not.out1file) &
- write (iout,'(a,i3,a,i3)')'nsup=',nsup,&
- ' nstart_sup=',nstart_sup !,"ergwergewrgae"
-!el if(.not.allocated(itype_pdb))
- allocate(itype_pdb(nres))
- do i=1,nres
- itype_pdb(i)=itype(i)
- enddo
- close (ipdbin)
- nnt=nstart_sup
- nct=nstart_sup+nsup-1
-!el if(.not.allocated(icont_ref))
- allocate(icont_ref(2,12*nres)) ! maxcont=12*maxres
- call contact(.false.,ncont_ref,icont_ref,co)
-
- if (sideadd) then
- if(me.eq.king.or..not.out1file) &
- write(iout,*)'Adding sidechains'
- maxsi=1000
- do i=2,nres-1
- iti=itype(i)
- if (iti.ne.10 .and. itype(i).ne.ntyp1) then
- nsi=0
- fail=.true.
- do while (fail.and.nsi.le.maxsi)
- call gen_side(iti,theta(i+1),alph(i),omeg(i),fail)
- nsi=nsi+1
- enddo
- if(fail) write(iout,*)'Adding sidechain failed for res ',&
- i,' after ',nsi,' trials'
- endif
- enddo
- endif
- endif
-
- if (indpdb.eq.0) then
-! Read sequence if not taken from the pdb file.
- read (inp,*) nres
-! print *,'nres=',nres
- allocate(sequence(nres))
- if (iscode.gt.0) then
- read (inp,'(80a1)') (sequence(i)(1:1),i=1,nres)
- else
- read (inp,'(20(1x,a3))') (sequence(i),i=1,nres)
- endif
-! Convert sequence to numeric code
- do i=1,nres
- itype(i)=rescode(i,sequence(i),iscode)
- enddo
-! Assign initial virtual bond lengths
-!elwrite(iout,*) "test_alloc"
- if(.not.allocated(vbld)) allocate(vbld(2*nres))
-!elwrite(iout,*) "test_alloc"
- if(.not.allocated(vbld_inv)) allocate(vbld_inv(2*nres))
-!elwrite(iout,*) "test_alloc"
- do i=2,nres
- vbld(i)=vbl
- vbld_inv(i)=vblinv
- enddo
- do i=2,nres-1
- vbld(i+nres)=dsc(iabs(itype(i)))
- vbld_inv(i+nres)=dsc_inv(iabs(itype(i)))
-! write (iout,*) "i",i," itype",itype(i),
-! & " dsc",dsc(itype(i))," vbld",vbld(i),vbld(i+nres)
- enddo
- endif
-! print *,nres
-! print '(20i4)',(itype(i),i=1,nres)
-!----------------------------
-!el reallocate tables
-! do i=1,maxres2
-! do j=1,3
-! c_alloc(j,i)=c(j,i)
-! dc_alloc(j,i)=dc(j,i)
-! enddo
-! enddo
-! do i=1,maxres
-!elwrite(iout,*) "itype",i,itype(i)
-! itype_alloc(i)=itype(i)
-! enddo
-
-! deallocate(c)
-! deallocate(dc)
-! deallocate(itype)
-! allocate(c(3,2*nres+4))
-! allocate(dc(3,0:2*nres+2))
-! allocate(itype(nres+2))
- allocate(itel(nres+2))
- itel(:)=0
-
-! do i=1,2*nres+2
-! do j=1,3
-! c(j,i)=c_alloc(j,i)
-! dc(j,i)=dc_alloc(j,i)
-! enddo
-! enddo
-! do i=1,nres+2
-! itype(i)=itype_alloc(i)
-! itel(i)=0
-! enddo
-!--------------------------
- do i=1,nres
-#ifdef PROCOR
- if (itype(i).eq.ntyp1 .or. itype(i+1).eq.ntyp1) then
-#else
- if (itype(i).eq.ntyp1) then
-#endif
- itel(i)=0
-#ifdef PROCOR
- else if (iabs(itype(i+1)).ne.20) then
-#else
- else if (iabs(itype(i)).ne.20) then
-#endif
- itel(i)=1
- else
- itel(i)=2
- endif
- enddo
- if(me.eq.king.or..not.out1file)then
- write (iout,*) "ITEL"
- do i=1,nres-1
- write (iout,*) i,itype(i),itel(i)
- enddo
- print *,'Call Read_Bridge.'
- endif
- call read_bridge
-!--------------------------------
-! znamy nres oraz nss można zaalokowac potrzebne tablice
- call alloc_geo_arrays
- call alloc_ener_arrays
-!--------------------------------
-! 8/13/98 Set limits to generating the dihedral angles
- do i=1,nres
- phibound(1,i)=-pi
- phibound(2,i)=pi
- enddo
- read (inp,*) ndih_constr
- if (ndih_constr.gt.0) then
- allocate(idih_constr(ndih_constr),idih_nconstr(ndih_constr)) !(maxdih_constr)
- allocate(phi0(ndih_constr),drange(ndih_constr)) !(maxdih_constr)
- read (inp,*) ftors
- read (inp,*) (idih_constr(i),phi0(i),drange(i),i=1,ndih_constr)
- if(me.eq.king.or..not.out1file)then
- write (iout,*) &
- 'There are',ndih_constr,' constraints on phi angles.'
- do i=1,ndih_constr
- write (iout,'(i5,2f8.3)') idih_constr(i),phi0(i),drange(i)
- enddo
- endif
- do i=1,ndih_constr
- phi0(i)=deg2rad*phi0(i)
- drange(i)=deg2rad*drange(i)
- enddo
- if(me.eq.king.or..not.out1file) &
- write (iout,*) 'FTORS',ftors
- do i=1,ndih_constr
- ii = idih_constr(i)
- phibound(1,ii) = phi0(i)-drange(i)
- phibound(2,ii) = phi0(i)+drange(i)
- enddo
- endif
- nnt=1
-#ifdef MPI
- if (me.eq.king) then
-#endif
- write (iout,'(a)') 'Boundaries in phi angle sampling:'
- do i=1,nres
- write (iout,'(a3,i5,2f10.1)') &
- restyp(itype(i)),i,phibound(1,i)*rad2deg,phibound(2,i)*rad2deg
- enddo
-#ifdef MP
- endif
-#endif
- nct=nres
-!d print *,'NNT=',NNT,' NCT=',NCT
- if (itype(1).eq.ntyp1) nnt=2
- if (itype(nres).eq.ntyp1) nct=nct-1
- if (pdbref) then
- if(me.eq.king.or..not.out1file) &
- write (iout,'(a,i3)') 'nsup=',nsup
- nstart_seq=nnt
- if (nsup.le.(nct-nnt+1)) then
- do i=0,nct-nnt+1-nsup
- if (seq_comp(itype(nnt+i),itype_pdb(nstart_sup),nsup)) then
- nstart_seq=nnt+i
- goto 111
- endif
- enddo
- write (iout,'(a)') &
- 'Error - sequences to be superposed do not match.'
- stop
- else
- do i=0,nsup-(nct-nnt+1)
- if (seq_comp(itype(nnt),itype_pdb(nstart_sup+i),nct-nnt+1)) &
- then
- nstart_sup=nstart_sup+i
- nsup=nct-nnt+1
- goto 111
- endif
- enddo
- write (iout,'(a)') &
- 'Error - sequences to be superposed do not match.'
- endif
- 111 continue
- if (nsup.eq.0) nsup=nct-nnt
- if (nstart_sup.eq.0) nstart_sup=nnt
- if (nstart_seq.eq.0) nstart_seq=nnt
- if(me.eq.king.or..not.out1file) &
- write (iout,*) 'nsup=',nsup,' nstart_sup=',nstart_sup,&
- ' nstart_seq=',nstart_seq !,"242343453254"
- endif
-!--- Zscore rms -------
- if (nz_start.eq.0) nz_start=nnt
- if (nz_end.eq.0 .and. nsup.gt.0) then
- nz_end=nnt+nsup-1
- else if (nz_end.eq.0) then
- nz_end=nct
- endif
- if(me.eq.king.or..not.out1file)then
- write (iout,*) 'NZ_START=',nz_start,' NZ_END=',nz_end
- write (iout,*) 'IZ_SC=',iz_sc
- endif
-!----------------------
- call init_int_table
- if (refstr) then
- if (.not.pdbref) then
- call read_angles(inp,*38)
- goto 39
- 38 write (iout,'(a)') 'Error reading reference structure.'
-#ifdef MPI
- call MPI_Finalize(MPI_COMM_WORLD,IERROR)
- stop 'Error reading reference structure'
-#endif
- 39 call chainbuild
- call setup_var
-!zscore call geom_to_var(nvar,coord_exp_zs(1,1))
- nstart_sup=nnt
- nstart_seq=nnt
- nsup=nct-nnt+1
- kkk=1
- do i=1,2*nres
- do j=1,3
- cref(j,i,kkk)=c(j,i)
- enddo
- enddo
- call contact(.true.,ncont_ref,icont_ref,co)
- endif
-! write (iout,*) "constr_dist",constr_dist,nstart_sup,nsup
- call flush(iout)
- if (constr_dist.gt.0) call read_dist_constr
- write (iout,*) "After read_dist_constr nhpb",nhpb
- call hpb_partition
- if(me.eq.king.or..not.out1file) &
- write (iout,*) 'Contact order:',co
- if (pdbref) then
- if(me.eq.king.or..not.out1file) &
- write (2,*) 'Shifting contacts:',nstart_seq,nstart_sup
- do i=1,ncont_ref
- do j=1,2
- icont_ref(j,i)=icont_ref(j,i)+nstart_seq-nstart_sup
- enddo
- if(me.eq.king.or..not.out1file) &
- write (2,*) i,' ',restyp(itype(icont_ref(1,i))),' ',&
- icont_ref(1,i),' ',&
- restyp(itype(icont_ref(2,i))),' ',icont_ref(2,i)
- enddo
- endif
- endif
- if (indpdb.eq.0 .and. modecalc.ne.2 .and. modecalc.ne.4 &
- .and. modecalc.ne.8 .and. modecalc.ne.9 .and. &
- modecalc.ne.10) then
-! If input structure hasn't been supplied from the PDB file read or generate
-! initial geometry.
- if (iranconf.eq.0 .and. .not. extconf) then
- if(me.eq.king.or..not.out1file .and.fg_rank.eq.0) &
- write (iout,'(a)') 'Initial geometry will be read in.'
- if (read_cart) then
- read(inp,'(8f10.5)',end=36,err=36) &
- ((c(l,k),l=1,3),k=1,nres),&
- ((c(l,k+nres),l=1,3),k=nnt,nct)
- write (iout,*) "Exit READ_CART"
- write (iout,'(8f10.5)') &
- ((c(l,k),l=1,3),k=1,nres),&
- ((c(l,k+nres),l=1,3),k=nnt,nct)
- call int_from_cart1(.true.)
- write (iout,*) "Finish INT_TO_CART"
- do i=1,nres-1
- do j=1,3
- dc(j,i)=c(j,i+1)-c(j,i)
- dc_norm(j,i)=dc_norm(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_norm(j,i+nres)*vbld_inv(i+nres)
- enddo
- endif
- enddo
- return
- else
- call read_angles(inp,*36)
- endif
- goto 37
- 36 write (iout,'(a)') 'Error reading angle file.'
-#ifdef MPI
- call mpi_finalize( MPI_COMM_WORLD,IERR )
-#endif
- stop 'Error reading angle file.'
- 37 continue
- else if (extconf) then
- if(me.eq.king.or..not.out1file .and. fg_rank.eq.0) &
- write (iout,'(a)') 'Extended chain initial geometry.'
- do i=3,nres
- theta(i)=90d0*deg2rad
- enddo
- do i=4,nres
- phi(i)=180d0*deg2rad
- enddo
- do i=2,nres-1
- alph(i)=110d0*deg2rad
- enddo
-!elwrite (iout,*)"alph(i)*deg2rad",(alph(i), i=1,nres)
- do i=2,nres-1
- omeg(i)=-120d0*deg2rad
- if (itype(i).le.0) omeg(i)=-omeg(i)
- enddo
- else
- if(me.eq.king.or..not.out1file) &
- write (iout,'(a)') 'Random-generated initial geometry.'
-
-
-#ifdef MPI
- if (me.eq.king .or. fg_rank.eq.0 .and. &
- ( modecalc.eq.12 .or. modecalc.eq.14) ) then
-#endif
- do itrial=1,100
- itmp=1
- call gen_rand_conf(itmp,*30)
- goto 40
- 30 write (iout,*) 'Failed to generate random conformation',&
- ', itrial=',itrial
- write (*,*) 'Processor:',me,&
- ' Failed to generate random conformation',&
- ' itrial=',itrial
- call intout
-
-#ifdef AIX
- call flush_(iout)
-#else
- call flush(iout)
-#endif
- enddo
- write (iout,'(a,i3,a)') 'Processor:',me,&
- ' error in generating random conformation.'
- write (*,'(a,i3,a)') 'Processor:',me,&
- ' error in generating random conformation.'
- call flush(iout)
-#ifdef MPI
- call MPI_Abort(mpi_comm_world,error_msg,ierrcode)
- 40 continue
- endif
-#else
- do itrial=1,100
- itmp=1
- call gen_rand_conf(itmp,*335)
- goto 40
- 335 write (iout,*) 'Failed to generate random conformation',&
- ', itrial=',itrial
- write (*,*) 'Failed to generate random conformation',&
- ', itrial=',itrial
- enddo
- write (iout,'(a,i3,a)') 'Processor:',me,&
- ' error in generating random conformation.'
- write (*,'(a,i3,a)') 'Processor:',me,&
- ' error in generating random conformation.'
- stop
- 40 continue
-#endif
- endif
- elseif (modecalc.eq.4) then
- read (inp,'(a)') intinname
- open (intin,file=intinname,status='old',err=333)
- if (me.eq.king .or. .not.out1file.and.fg_rank.eq.0) &
- write (iout,'(a)') 'intinname',intinname
- write (*,'(a)') 'Processor',myrank,' intinname',intinname
- goto 334
- 333 write (iout,'(2a)') 'Error opening angle file ',intinname
-#ifdef MPI
- call MPI_Finalize(MPI_COMM_WORLD,IERR)
-#endif
- stop 'Error opening angle file.'
- 334 continue
-
- endif
-! Generate distance constraints, if the PDB structure is to be regularized.
- if (nthread.gt.0) then
- call read_threadbase
- endif
- call setup_var
-!elwrite (iout,*)"alph(i)*deg2rad",(alph(i), i=1,nres)
- if (me.eq.king .or. .not. out1file) &
- call intout
-!elwrite (iout,*)"alph(i)*rad2deg",(alph(i)*rad2deg, i=1,nres)
- if (ns.gt.0 .and. (me.eq.king .or. .not.out1file) ) then
- write (iout,'(/a,i3,a)') &
- 'The chain contains',ns,' disulfide-bridging cysteines.'
- write (iout,'(20i4)') (iss(i),i=1,ns)
- if (dyn_ss) then
- write(iout,*)"Running with dynamic disulfide-bond formation"
- else
- write (iout,'(/a/)') 'Pre-formed links are:'
- do i=1,nss
- i1=ihpb(i)-nres
- i2=jhpb(i)-nres
- it1=itype(i1)
- it2=itype(i2)
- if (me.eq.king.or..not.out1file) &
- write (iout,'(2a,i3,3a,i3,a,3f10.3)') &
- restyp(it1),'(',i1,') -- ',restyp(it2),'(',i2,')',dhpb(i),&
- ebr,forcon(i)
- enddo
- write (iout,'(a)')
- endif
- endif
- if (ns.gt.0.and.dyn_ss) then
- do i=nss+1,nhpb
- ihpb(i-nss)=ihpb(i)
- jhpb(i-nss)=jhpb(i)
- forcon(i-nss)=forcon(i)
- dhpb(i-nss)=dhpb(i)
- enddo
- nhpb=nhpb-nss
- nss=0
- call hpb_partition
- do i=1,ns
- dyn_ss_mask(iss(i))=.true.
- enddo
- endif
- if (i2ndstr.gt.0) call secstrp2dihc
-! call geom_to_var(nvar,x)
-! call etotal(energia(0))
-! call enerprint(energia(0))
-! call briefout(0,etot)
-! stop
-!d write (iout,'(2(a,i3))') 'NNT',NNT,' NCT',NCT
-!d write (iout,'(a)') 'Variable list:'
-!d write (iout,'(i4,f10.5)') (i,rad2deg*x(i),i=1,nvar)
-#ifdef MPI
- if (me.eq.king .or. (fg_rank.eq.0 .and. .not.out1file)) &
- write (iout,'(//80(1h*)/20x,a,i4,a/80(1h*)//)') &
- 'Processor',myrank,': end reading molecular data.'
-#endif
- return
- end subroutine molread
-!-----------------------------------------------------------------------------
-!-----------------------------------------------------------------------------
- end module io
projects / unres4.git / blobdiff