C Set up the time limit (caution! The time must be input in minutes!)
read_cart=index(controlcard,'READ_CART').gt.0
call readi(controlcard,'CONSTR_DIST',constr_dist,0)
+ call readi(controlcard,'CONSTR_HOMOL',constr_homology,0)
call reada(controlcard,'TIMLIM',timlim,960.0D0) ! default 16 hours
unres_pdb = index(controlcard,'UNRES_PDB') .gt. 0
call reada(controlcard,'SAFETY',safety,30.0D0) ! default 30 minutes
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,'WDFAD',wdfa_dist,0.0d0)
+ call reada(weightcard,'WDFAT',wdfa_tor,0.0d0)
+ call reada(weightcard,'WDFAN',wdfa_nei,0.0d0)
+ call reada(weightcard,'WDFAB',wdfa_beta,0.0d0)
call reada(weightcard,'SCAL14',scal14,0.4D0)
call reada(weightcard,'SCALSCP',scalscp,1.0d0)
call reada(weightcard,'CUTOFF',cutoff_corr,7.0d0)
weights(18)=scal14
weights(21)=wsccor
endif
+ weights(25)=wdfa_dist
+ weights(26)=wdfa_tor
+ weights(27)=wdfa_nei
+ weights(28)=wdfa_beta
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
+ & wturn4,wturn6,
+ & wdfa_dist,wdfa_tor,wdfa_nei,wdfa_beta
10 format (/'Energy-term weights (unscaled):'//
& 'WSCC= ',f10.6,' (SC-SC)'/
& 'WSCP= ',f10.6,' (SC-p)'/
& 'WSCCOR= ',f10.6,' (back-scloc correlation)'/
& 'WTURN3= ',f10.6,' (turns, 3rd order)'/
& 'WTURN4= ',f10.6,' (turns, 4th order)'/
- & 'WTURN6= ',f10.6,' (turns, 6th order)')
+ & 'WTURN6= ',f10.6,' (turns, 6th order)'/
+ & 'WDFA_D= ',f10.6,' (DFA, distance)' /
+ & 'WDFA_T= ',f10.6,' (DFA, torsional)' /
+ & 'WDFA_N= ',f10.6,' (DFA, number of neighbor)' /
+ & 'WDFA_B= ',f10.6,' (DFA, beta formation)')
if(me.eq.king.or..not.out1file)then
if (wcorr4.gt.0.0d0) then
write (iout,'(/2a/)') 'Local-electrostatic type correlation ',
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
+ & wturn4,wturn6,
+ & wdfa_dist,wdfa_tor,wdfa_nei,wdfa_beta
22 format (/'Energy-term weights (scaled):'//
& 'WSCC= ',f10.6,' (SC-SC)'/
& 'WSCP= ',f10.6,' (SC-p)'/
& 'WSCCOR= ',f10.6,' (back-scloc correlatkion)'/
& 'WTURN3= ',f10.6,' (turns, 3rd order)'/
& 'WTURN4= ',f10.6,' (turns, 4th order)'/
- & 'WTURN6= ',f10.6,' (turns, 6th order)')
+ & 'WTURN6= ',f10.6,' (turns, 6th order)'/
+ & 'WDFA_D= ',f10.6,' (DFA, distance)' /
+ & 'WDFA_T= ',f10.6,' (DFA, torsional)' /
+ & 'WDFA_N= ',f10.6,' (DFA, number of neighbor)' /
+ & 'WDFA_B= ',f10.6,' (DFA, beta formation)')
if(me.eq.king.or..not.out1file)
& write (iout,*) "Reference temperature for weights calculation:",
& temp0
call chainbuild
endif
C Following 2 lines for diagnostics; comment out if not needed
- write (iout,*) "After sideadd"
- call intout
+c write (iout,*) "After sideadd"
+c call intout
endif
if (indpdb.eq.0) then
C Read sequence if not taken from the pdb file.
cd print *,'NNT=',NNT,' NCT=',NCT
if (itype(1).eq.21) nnt=2
if (itype(nres).eq.21) nct=nct-1
+
+C Bartek:READ init_vars
+C Initialize variables!
+C Juyong:READ read_info
+C READ fragment information!!
+C both routines should be in dfa.F file!!
+
+ if (.not. (wdfa_dist.eq.0.0 .and. wdfa_tor.eq.0.0 .and.
+ & wdfa_nei.eq.0.0 .and. wdfa_beta.eq.0.0)) then
+ call init_dfa_vars
+ print*, 'init_dfa_vars finished!'
+ call read_dfa_info
+ print*, 'read_dfa_info finished!'
+ endif
+C
if (pdbref) then
if(me.eq.king.or..not.out1file)
& write (iout,'(a,i3)') 'nsup=',nsup
c write (iout,*) "constr_dist",constr_dist,nstart_sup,nsup
if (constr_dist.gt.0) then
call read_dist_constr
- call hpb_partition
endif
+
+
+ if (constr_homology.gt.0) then
+ call read_constr_homology
+ else
+ homol_nset=0
+ endif
+
+
+ if (nhpb.gt.0) call hpb_partition
c write (iout,*) "After read_dist_constr nhpb",nhpb
c call flush(iout)
if (indpdb.eq.0 .and. modecalc.ne.2 .and. modecalc.ne.4
write (iout,'(20i4)') (iss(i),i=1,ns)
if (dyn_ss) then
write(iout,*)"Running with dynamic disulfide-bond formation"
- 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
else
write (iout,'(/a/)') 'Pre-formed links are:'
do i=1,nss
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)')
+ 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
c call geom_to_var(nvar,x)
c call etotal(energia(0))
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(iss(i)),i,
+ & '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(iss(i)),i,
+ & '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)
C Read preformed bridges.
if (ns.gt.0) then
read (inp,*) nss,(ihpb(i),jhpb(i),i=1,nss)
- write (iout,*) 'nss=',nss,' ihpb,jhpb: ',(ihpb(i),jhpb(i),i=1,nss)
+ if(fg_rank.eq.0)
+ & write(iout,*)'nss=',nss,' ihpb,jhpb: ',(ihpb(i),jhpb(i),i=1,nss)
if (nss.gt.0) then
nhpb=nss
C Check if the residues involved in bridges are in the specified list of
include 'COMMON.CHAIN'
include 'COMMON.IOUNITS'
include 'COMMON.MD'
+ include 'COMMON.CONTROL'
open(irest2,file=rest2name,status='unknown')
read(irest2,*) totT,EK,potE,totE,t_bath
do i=1,2*nres
do i=1,2*nres
read(irest2,'(3e15.5)') (dc(j,i),j=1,3)
enddo
- if(usampl) then
+ if(usampl.or.homol_nset.gt.1) then
read (irest2,*) iset
endif
close(irest2)
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
+c write (iout,*) "j",j," k",k
ddjk=dist(j,k)
if (constr_dist.eq.1) then
nhpb=nhpb+1
return
end
c-------------------------------------------------------------------------------
+
+ subroutine read_constr_homology
+
+ include 'DIMENSIONS'
+#ifdef MPI
+ include 'mpif.h'
+#endif
+ include 'COMMON.SETUP'
+ include 'COMMON.CONTROL'
+ include 'COMMON.CHAIN'
+ include 'COMMON.IOUNITS'
+ include 'COMMON.MD'
+ include 'COMMON.GEO'
+ include 'COMMON.INTERACT'
+ double precision odl_temp,sigma_odl_temp
+ common /przechowalnia/ odl_temp(maxres,maxres,max_template),
+ & sigma_odl_temp(maxres,maxres,max_template)
+ character*2 kic2
+ character*24 model_ki_dist, model_ki_angle
+ character*500 controlcard
+ character*3200 controlcard1
+ integer ki, i, j, k, l
+ logical lprn /.true./
+
+ call card_concat(controlcard)
+ call readi(controlcard,"HOMOL_NSET",homol_nset,1)
+ if (homol_nset.gt.1)then
+ call card_concat(controlcard)
+ read(controlcard,*) (waga_dist(i),i=1,homol_nset)
+ call card_concat(controlcard)
+ read(controlcard,*) (waga_angle(i),i=1,homol_nset)
+ write(iout,*) "iset distance_weight angle_weight"
+ do i=1,homol_nset
+ write(iout,*) i,waga_dist(i),waga_angle(i)
+ enddo
+ else
+ iset=1
+ call reada(controlcard,"HOMOL_DIST",waga_dist(1),1.0d0)
+ call reada(controlcard,"HOMOL_ANGLE",waga_angle(1),1.0d0)
+ endif
+
+ write (iout,*) "nnt",nnt," nct",nct
+ call flush(iout)
+ lim_odl=0
+ lim_dih=0
+ do i=1,nres
+ do j=i+2,nres
+ do ki=1,constr_homology
+ sigma_odl_temp(i,j,ki)=0.0d0
+ odl_temp(i,j,ki)=0.0d0
+ enddo
+ enddo
+ enddo
+ do i=1,nres-3
+ do ki=1,constr_homology
+ dih(ki,i)=0.0d0
+ sigma_dih(ki,i)=0.0d0
+ enddo
+ enddo
+ do ki=1,constr_homology
+ write(kic2,'(i2)') ki
+ if (ki.le.9) kic2="0"//kic2(2:2)
+
+ model_ki_dist="model"//kic2//".dist"
+ model_ki_angle="model"//kic2//".angle"
+ open (ientin,file=model_ki_dist,status='old')
+ do irec=1,maxdim !petla do czytania wiezow na odleglosc
+ read (ientin,*,end=1401) i, j, odl_temp(i+nnt-1,j+nnt-1,ki),
+ & sigma_odl_temp(i+nnt-1,j+nnt-1,ki)
+ odl_temp(j+nnt-1,i+nnt-1,ki)=odl_temp(i+nnt-1,j+nnt-1,ki)
+ sigma_odl_temp(j+nnt-1,i+nnt-1,ki)=
+ & sigma_odl_temp(i+nnt-1,j+nnt-1,ki)
+ enddo
+ 1401 continue
+ close (ientin)
+ open (ientin,file=model_ki_angle,status='old')
+ do irec=1,maxres-3 !petla do czytania wiezow na katy torsyjne
+ read (ientin,*,end=1402) i, j, k,l,dih(ki,i+nnt-1),
+ & sigma_dih(ki,i+nnt-1)
+ if (i+nnt-1.gt.lim_dih) lim_dih=i+nnt-1
+ sigma_dih(ki,i+nnt-1)=1.0d0/sigma_dih(ki,i+nnt-1)**2
+ enddo
+ 1402 continue
+ close (ientin)
+ enddo
+ ii=0
+ write (iout,*) "nnt",nnt," nct",nct
+ do i=nnt,nct-2
+ do j=i+2,nct
+ ki=1
+c write (iout,*) "i",i," j",j," constr_homology",constr_homology
+ do while (ki.le.constr_homology .and.
+ & sigma_odl_temp(i,j,ki).le.0.0d0)
+c write (iout,*) "ki",ki," sigma_odl",sigma_odl_temp(i,j,ki)
+ ki=ki+1
+ enddo
+c write (iout,*) "ki",ki
+ if (ki.gt.constr_homology) cycle
+ ii=ii+1
+ ires_homo(ii)=i
+ jres_homo(ii)=j
+ do ki=1,constr_homology
+ odl(ki,ii)=odl_temp(i,j,ki)
+ sigma_odl(ki,ii)=1.0d0/sigma_odl_temp(i,j,ki)**2
+ enddo
+ enddo
+ enddo
+ lim_odl=ii
+ if (constr_homology.gt.0) call homology_partition
+c Print restraints
+ if (.not.lprn) return
+ write (iout,*) "Distance restraints from templates"
+ do ii=1,lim_odl
+ write(iout,'(3i5,10(2f8.2,4x))') ii,ires_homo(ii),jres_homo(ii),
+ & (odl(ki,ii),1.0d0/dsqrt(sigma_odl(ki,ii)),ki=1,constr_homology)
+ enddo
+ write (iout,*) "Dihedral angle restraints from templates"
+ do i=nnt,lim_dih
+ write (iout,'(i5,10(2f8.2,4x))') i,(rad2deg*dih(ki,i),
+ & rad2deg/dsqrt(sigma_dih(ki,i)),ki=1,constr_homology)
+ enddo
+c write(iout,*) "TEST CZYTANIA1",odl(1,2,1),odl(1,3,1),odl(1,4,1)
+c write(iout,*) "TEST CZYTANIA2",dih(1,1),dih(2,1),dih(3,1)
+
+
+ return
+ end
+c----------------------------------------------------------------------
+
#ifdef WINIFL
subroutine flush(iu)
return
return
end
#endif
+
c------------------------------------------------------------------------------
subroutine copy_to_tmp(source)
include "DIMENSIONS"