refstr=(index(controlcard,'REFSTR').gt.0)
write (iout,*) "REFSTR",refstr
pdbref=(index(controlcard,'PDBREF').gt.0)
+ dyn_ss=(index(controlcard,'DYN_SS').gt.0)
iscode=index(controlcard,'ONE_LETTER')
tree=(index(controlcard,'MAKE_TREE').gt.0)
min_var=(index(controlcard,'MINVAR').gt.0)
plot_tree=(index(controlcard,'PLOT_TREE').gt.0)
punch_dist=(index(controlcard,'PUNCH_DIST').gt.0)
- call readi(controlcard,'NCUT',ncut,1)
+ call readi(controlcard,'NCUT',ncut,0)
+ call readi(controlcard,'NCLUST',nclust,5)
call readi(controlcard,'NSTART',nstart,0)
call readi(controlcard,'NEND',nend,0)
call reada(controlcard,'ECUT',ecut,10.0d0)
lgrp=(index(controlcard,'LGRP').gt.0)
caonly=(index(controlcard,'CA_ONLY').gt.0)
print_dist=(index(controlcard,'PRINT_DIST').gt.0)
- call multreada(controlcard,'CUTOFF',rcutoff,ncut,-1.0d0)
+ if (ncut.gt.0)
+ & call multreada(controlcard,'CUTOFF',rcutoff,ncut,-1.0d0)
call readi(controlcard,'IOPT',iopt,2)
lside = index(controlcard,"SIDE").gt.0
efree = index(controlcard,"EFREE").gt.0
write (iout,*) 'beta_h',(beta_h(i),i=1,nT)
lprint_cart=index(controlcard,"PRINT_CART") .gt.0
lprint_int=index(controlcard,"PRINT_INT") .gt.0
+ with_dihed_constr = index(controlcard,"WITH_DIHED_CONSTR").gt.0
+ call readi(controlcard,'CONSTR_DIST',constr_dist,0)
+ write (iout,*) "with_dihed_constr ",with_dihed_constr,
+ & " CONSTR_DIST",constr_dist
+ call flush(iout)
if (min_var) iopt=1
return
end
include 'COMMON.CONTROL'
include 'COMMON.CONTACTS'
include 'COMMON.TIME1'
+ include 'COMMON.TORCNSTR'
#ifdef MPL
include 'COMMON.INFO'
#endif
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,'WSCCOR',wsccor,1.0D0)
+ 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)
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)
+ 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
+ write (iout,'(/a)') "Disulfide bridge parameters:"
+ write (iout,'(a,f10.2)') 'S-S bridge energy: ',ebr
+ write (iout,'(a,f10.2)') 'S-S depth: ',ss_depth
+ write (iout,'(2(a,f10.2))') 'd0cm:',d0cm,' akcm:',akcm
+ write (iout,'(2(a,f10.2))') 'akth:',akth,' akct:',akct
+ write (iout,'(3(a,f10.2))') 'v1ss:',v1ss,' v2ss:',v2ss,
+ & ' v3ss:',v3ss
+ write (iout,'(2(a,f10.2))') 'ht:',ht,' eps:', eps(1,1)
if (wcorr4.gt.0.0d0) wcorr=wcorr4
weights(1)=wsc
weights(2)=wscp
weights(18)=scal14
write (iout,10) wsc,wscp,welec,wvdwpp,wbond,wang,wscloc,wtor,
& wtor_d,wstrain,wel_loc,wcorr,wcorr5,wcorr6,wturn3,
- & wturn4,wturn6
+ & wturn4,wturn6,wsccor
10 format (/'Energy-term weights (unscaled):'//
& 'WSCC= ',f10.6,' (SC-SC)'/
& 'WSCP= ',f10.6,' (SC-p)'/
& 'WCORR6= ',f10.6,' (multi-body 6th order)'/
& 'WTURN3= ',f10.6,' (turns, 3rd order)'/
& 'WTURN4= ',f10.6,' (turns, 4th order)'/
- & 'WTURN6= ',f10.6,' (turns, 6th order)')
+ & 'WTURN6= ',f10.6,' (turns, 6th order)'/
+ & 'WSCCOR= ',f10.6,' (SC-backbone torsinal correalations)')
if (wcorr4.gt.0.0d0) then
write (iout,'(/2a/)') 'Local-electrostatic type correlation ',
& 'between contact pairs of peptide groups'
do i=1,nres
itype(i)=rescode(i,sequence(i),iscode)
enddo
+ if (itype(2).eq.10.and.itype(1).eq.ntyp1) then
+ write (iout,*)
+ & "Glycine is the first full residue, initial dummy deleted"
+ do i=1,nres
+ itype(i)=itype(i+1)
+ enddo
+ nres=nres-1
+ endif
+ if (itype(nres-1).eq.10.and.itype(nres).eq.ntyp1) then
+ write (iout,*)
+ & "Glycine is the last full residue, terminal dummy deleted"
+ nres=nres-1
+ endif
print *,nres
print '(20i4)',(itype(i),i=1,nres)
print *,'Call Read_Bridge.'
call read_bridge
+ if (with_dihed_constr) then
+
+ read (inp,*) ndih_constr
+ if (ndih_constr.gt.0) then
+ read (inp,*) ftors
+ write (iout,*) 'FTORS',ftors
+ read (inp,*) (idih_constr(i),phi0(i),drange(i),i=1,ndih_constr)
+ 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
+ do i=1,ndih_constr
+ phi0(i)=deg2rad*phi0(i)
+ drange(i)=deg2rad*drange(i)
+ enddo
+ endif
+
+ endif
nnt=1
nct=nres
print *,'NNT=',NNT,' NCT=',NCT
enddo
enddo
endif
- call contact(.true.,ncont_ref,icont_ref)
+ call contact(print_contact_map,ncont_ref,icont_ref)
+ endif
+c Read distance restraints
+ if (constr_dist.gt.0) then
+ call read_dist_constr
+ call hpb_partition
endif
return
end
enddo
write (iout,'(a,i3,a)') 'Pair',i,' contains unknown cystine.'
20 continue
- dhpb(i)=dbr
- forcon(i)=fbr
+c dhpb(i)=dbr
+c forcon(i)=fbr
enddo
do i=1,nss
ihpb(i)=ihpb(i)+nres
enddo
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.
+c write(iout,*) i,iss(i),dyn_ss_mask(iss(i)),"ATU"
+ enddo
+ endif
+ print *, "Leaving brigde read"
return
end
c----------------------------------------------------------------------------
return
end
c----------------------------------------------------------------------------
+ subroutine multreadi(rekord,lancuch,tablica,dim,default)
+ implicit none
+ integer dim,i
+ integer tablica(dim),default
+ character*(*) rekord,lancuch
+ character*80 aux
+ integer ilen,iread
+ 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
+c----------------------------------------------------------------------------
subroutine card_concat(card)
include 'DIMENSIONS'
include 'COMMON.IOUNITS'
#endif
return
end
+c-------------------------------------------------------------------------------
+ subroutine read_dist_constr
+ implicit real*8 (a-h,o-z)
+ include 'DIMENSIONS'
+ include 'COMMON.CONTROL'
+ include 'COMMON.CHAIN'
+ include 'COMMON.IOUNITS'
+ include 'COMMON.SBRIDGE'
+ integer ifrag_(2,100),ipair_(2,100)
+ double precision wfrag_(100),wpair_(100)
+ character*500 controlcard
+c write (iout,*) "Calling read_dist_constr"
+c write (iout,*) "nres",nres," nstart_sup",nstart_sup," nsup",nsup
+ call card_concat(controlcard)
+c call flush(iout)
+c write (iout,'(a)') controlcard
+ 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
+ call flush(iout)
+ if (.not.refstr .and. nfrag_.gt.0) then
+ write (iout,*)
+ & "ERROR: no reference structure to compute distance restraints"
+ write (iout,*)
+ & "Restraints must be specified explicitly (NDIST=number)"
+ stop
+ endif
+ if (nfrag_.lt.2 .and. npair_.gt.0) then
+ write (iout,*) "ERROR: Less than 2 fragments specified",
+ & " but distance restraints between pairs requested"
+ stop
+ endif
+ 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
+c 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
+ write (iout,'(a,3i5,f8.2,1pe12.2)') "+dist.constr ",
+ & nhpb,ihpb(nhpb),jhpb(nhpb),dhpb(nhpb),forcon(nhpb)
+ 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)
+ write (iout,'(a,3i5,f8.2,f10.1)') "+dist.constr ",
+ & nhpb,ihpb(nhpb),jhpb(nhpb),dhpb(nhpb),forcon(nhpb)
+ enddo
+ enddo
+ endif
+ enddo
+ do i=1,ndist_
+ read (inp,*) ihpb(nhpb+1),jhpb(nhpb+1),dhpb(i),dhpb1(i),
+ & ibecarb(i),forcon(nhpb+1)
+ if (forcon(nhpb+1).gt.0.0d0) then
+ nhpb=nhpb+1
+ 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
+ enddo
+ do i=1,nhpb
+ write (iout,'(a,3i5,2f8.2,i2,f10.1)') "+dist.constr ",
+ & i,ihpb(i),jhpb(i),dhpb(i),dhpb1(i),ibecarb(i),forcon(i)
+ enddo
+ call flush(iout)
+ return
+ end
projects / unres.git / blobdiff