call readi(controlcard,'MAXGEN',maxgen,10000)
call readi(controlcard,'MAXOVERLAP',maxoverlap,1000)
call readi(controlcard,"KDIAG",kdiag,0)
- call readi(controlcard,"RESCALE_MODE",rescale_mode,1)
+ call readi(controlcard,"RESCALE_MODE",rescale_mode,2)
if(me.eq.king .or. .not. out1file .and. fg_rank.eq.0)
& write (iout,*) "RESCALE_MODE",rescale_mode
split_ene=index(controlcard,'SPLIT_ENE').gt.0
do i=1,nrep
iremd_m_total=iremd_m_total+remd_m(i)
enddo
- write (iout,*) 'Total number of replicas ',iremd_m_total
+ write (iout,*) 'Total number of replicas ',iremd_m_total
+ endif
endif
- endif
if(me.eq.king.or..not.out1file)
- & write (iout,'(/30(1h=),a,29(1h=)/)') " End of REMD run setup "
+ & write (iout,'(/30(1h=),a,29(1h=)/)') " End of REMD run setup "
return
end
c--------------------------------------------------------------------------
rest = index(controlcard,"REST").gt.0
tbf = index(controlcard,"TBF").gt.0
usampl = index(controlcard,"USAMPL").gt.0
+
mdpdb = index(controlcard,"MDPDB").gt.0
call reada(controlcard,"T_BATH",t_bath,300.0d0)
call reada(controlcard,"TAU_BATH",tau_bath,1.0d-1)
C Body
C
C Read weights of the subsequent energy terms.
- call card_concat(weightcard)
- 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
+ call card_concat(weightcard)
+ 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
C 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
+ 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,
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)
+ do i=1,maxres
+ dyn_ss_mask(i)=.false.
+ enddo
+ do i=1,maxres-1
+ do j=i+1,maxres
+ dyn_ssbond_ij(i,j)=1.0d300
+ enddo
+ enddo
+ 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
-c print *,'indpdb=',indpdb,' pdbref=',pdbref
+ 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
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)')
+ 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)
include 'COMMON.SETUP'
C Read bridging residues.
read (inp,*) ns,(iss(i),i=1,ns)
-c print *,'ns=',ns
+ print *,'ns=',ns
if(me.eq.king.or..not.out1file)
& write (iout,*) 'ns=',ns,' iss:',(iss(i),i=1,ns)
C 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(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
call readi(minimcard,'MINFUN',minfun,maxmin)
call reada(minimcard,'TOLF',tolf,1.0D-2)
call reada(minimcard,'RTOLF',rtolf,1.0D-4)
+ print_min_stat=min0(index(minimcard,'PRINT_MIN_STAT'),1)
+ print_min_res=min0(index(minimcard,'PRINT_MIN_RES'),1)
+ print_min_ini=min0(index(minimcard,'PRINT_MIN_INI'),1)
write (iout,'(/80(1h*)/20x,a/80(1h*))')
& 'Options in energy minimization:'
write (iout,'(4(a,i5),a,1pe14.5,a,1pe14.5)')
c 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
+ nhpb=nhpb+1
+ ihpb(nhpb)=j
+ jhpb(nhpb)=k
dhpb(nhpb)=ddjk
- forcon(nhpb)=wfrag_(i)
+ forcon(nhpb)=wfrag_(i)
else if (constr_dist.eq.2) then
if (ddjk.le.dist_cut) then
nhpb=nhpb+1
#endif
if (OKRandom) then
c r1 = prng_next(me)
- r1=ran_number(0.0D0,1.0D0)
+ r1=ran_number(0.0D0,1.0D0)
if(me.eq.king)
& write (iout,*) 'ran_num',r1
if (r1.lt.0.0d0) OKRandom=.false.