include 'COMMON.FFIELD'
include 'COMMON.FREE'
include 'COMMON.INTERACT'
+ include "COMMON.SPLITELE"
+ include 'COMMON.SHIELD'
character*320 controlcard,ucase
#ifdef MPL
include 'COMMON.INFO'
#endif
- integer i
-
+ integer i,i1,i2,it1,it2
+ double precision pi
read (INP,'(a80)') titel
call card_concat(controlcard)
call readi(controlcard,'RESCALE',rescale_mode,2)
call reada(controlcard,'DISTCHAINMAX',distchainmax,50.0d0)
write (iout,*) "DISTCHAINMAX",distchainmax
+C Reading the dimensions of box in x,y,z coordinates
+ call reada(controlcard,'BOXX',boxxsize,100.0d0)
+ call reada(controlcard,'BOXY',boxysize,100.0d0)
+ call reada(controlcard,'BOXZ',boxzsize,100.0d0)
+c Cutoff range for interactions
+ call reada(controlcard,"R_CUT",r_cut,15.0d0)
+ call reada(controlcard,"LAMBDA",rlamb,0.3d0)
+ call reada(controlcard,"LIPTHICK",lipthick,0.0d0)
+ call reada(controlcard,"LIPAQBUF",lipbufthick,0.0d0)
+ if (lipthick.gt.0.0d0) then
+ bordliptop=(boxzsize+lipthick)/2.0
+ bordlipbot=bordliptop-lipthick
+C endif
+ if ((bordliptop.gt.boxzsize).or.(bordlipbot.lt.0.0))
+ & write(iout,*) "WARNING WRONG SIZE OF LIPIDIC PHASE"
+ buflipbot=bordlipbot+lipbufthick
+ bufliptop=bordliptop-lipbufthick
+ if ((lipbufthick*2.0d0).gt.lipthick)
+ &write(iout,*) "WARNING WRONG SIZE OF LIP AQ BUF"
+ endif
+ write(iout,*) "bordliptop=",bordliptop
+ write(iout,*) "bordlipbot=",bordlipbot
+ write(iout,*) "bufliptop=",bufliptop
+ write(iout,*) "buflipbot=",buflipbot
+C Shielding mode
+ call readi(controlcard,'SHIELD',shield_mode,0)
+ write (iout,*) "SHIELD MODE",shield_mode
+ if (shield_mode.gt.0) then
+ pi=3.141592d0
+C VSolvSphere the volume of solving sphere
+C print *,pi,"pi"
+C rpp(1,1) is the energy r0 for peptide group contact and will be used for it
+C there will be no distinction between proline peptide group and normal peptide
+C group in case of shielding parameters
+ VSolvSphere=4.0/3.0*pi*rpp(1,1)**3
+ VSolvSphere_div=VSolvSphere-4.0/3.0*pi*(rpp(1,1)/2.0)**3
+ write (iout,*) VSolvSphere,VSolvSphere_div
+C long axis of side chain
+ do i=1,ntyp
+ long_r_sidechain(i)=vbldsc0(1,i)
+ short_r_sidechain(i)=sigma0(i)
+ enddo
+ buff_shield=1.0d0
+ endif
call readi(controlcard,'PDBOUT',outpdb,0)
call readi(controlcard,'MOL2OUT',outmol2,0)
refstr=(index(controlcard,'REFSTR').gt.0)
pdbref=(index(controlcard,'PDBREF').gt.0)
iscode=index(controlcard,'ONE_LETTER')
tree=(index(controlcard,'MAKE_TREE').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
+ with_theta_constr = index(controlcard,"WITH_THETA_CONSTR").gt.0
+ write (iout,*) "with_theta_constr ",with_theta_constr
+ call flush(iout)
min_var=(index(controlcard,'MINVAR').gt.0)
plot_tree=(index(controlcard,'PLOT_TREE').gt.0)
punch_dist=(index(controlcard,'PUNCH_DIST').gt.0)
include 'COMMON.CONTROL'
include 'COMMON.CONTACTS'
include 'COMMON.TIME1'
+ include 'COMMON.TORCNSTR'
+ include 'COMMON.SHIELD'
#ifdef MPL
include 'COMMON.INFO'
#endif
double precision x(maxvar)
integer itype_pdb(maxres)
logical seq_comp
- integer i,j,kkk
+ integer i,j,kkk,i1,i2,it1,it2
C
C Body
C
C Read weights of the subsequent energy terms.
call card_concat(weightcard)
- call reada(weightcard,'WSC',wsc,1.0d0)
+ write(iout,*) weightcard
+C call reada(weightcard,'WSC',wsc,1.0d0)
+ write(iout,*) wsc
call reada(weightcard,'WLONG',wsc,wsc)
call reada(weightcard,'WSCP',wscp,1.0d0)
call reada(weightcard,'WELEC',welec,1.0D0)
call reada(weightcard,"V2SS",v2ss,7.61d0)
call reada(weightcard,"V3SS",v3ss,13.7d0)
call reada(weightcard,"EBR",ebr,-5.50D0)
+ call reada(weightcard,'WSHIELD',wshield,1.0d0)
+ write(iout,*) 'WSHIELD',wshield
+ call reada(weightcard,'WLT',wliptran,0.0D0)
call reada(weightcard,"ATRISS",atriss,0.301D0)
call reada(weightcard,"BTRISS",btriss,0.021D0)
call reada(weightcard,"CTRISS",ctriss,1.001D0)
print *,'Call Read_Bridge.'
call read_bridge
+C this fragment reads diheadral constrains
+ if (with_dihed_constr) then
+
+ read (inp,*) ndih_constr
+ if (ndih_constr.gt.0) then
+C read (inp,*) ftors
+C write (iout,*) 'FTORS',ftors
+C ftors is the force constant for torsional quartic constrains
+ read (inp,*) (idih_constr(i),phi0(i),drange(i),ftors(i),
+ & i=1,ndih_constr)
+ write (iout,*)
+ & 'There are',ndih_constr,' constraints on phi angles.'
+ do i=1,ndih_constr
+ write (iout,'(i5,3f8.3)') idih_constr(i),phi0(i),drange(i),
+ & ftors(i)
+ enddo
+ do i=1,ndih_constr
+ phi0(i)=deg2rad*phi0(i)
+ drange(i)=deg2rad*drange(i)
+ enddo
+ endif ! endif ndif_constr.gt.0
+ endif ! with_dihed_constr
+ if (with_theta_constr) then
+C with_theta_constr is keyword allowing for occurance of theta constrains
+ read (inp,*) ntheta_constr
+C ntheta_constr is the number of theta constrains
+ if (ntheta_constr.gt.0) then
+C read (inp,*) ftors
+ read (inp,*) (itheta_constr(i),theta_constr0(i),
+ & theta_drange(i),for_thet_constr(i),
+ & i=1,ntheta_constr)
+C the above code reads from 1 to ntheta_constr
+C itheta_constr(i) residue i for which is theta_constr
+C theta_constr0 the global minimum value
+C theta_drange is range for which there is no energy penalty
+C for_thet_constr is the force constant for quartic energy penalty
+C E=k*x**4
+C if(me.eq.king.or..not.out1file)then
+ write (iout,*)
+ & 'There are',ntheta_constr,' constraints on phi angles.'
+ do i=1,ntheta_constr
+ write (iout,'(i5,3f8.3)') itheta_constr(i),theta_constr0(i),
+ & theta_drange(i),
+ & for_thet_constr(i)
+ enddo
+C endif
+ do i=1,ntheta_constr
+ theta_constr0(i)=deg2rad*theta_constr0(i)
+ theta_drange(i)=deg2rad*theta_drange(i)
+ enddo
+C if(me.eq.king.or..not.out1file)
+C & write (iout,*) 'FTORS',ftors
+C do i=1,ntheta_constr
+C ii = itheta_constr(i)
+C thetabound(1,ii) = phi0(i)-drange(i)
+C thetabound(2,ii) = phi0(i)+drange(i)
+C enddo
+ endif ! ntheta_constr.gt.0
+ endif! with_theta_constr
+
nnt=1
nct=nres
print *,'NNT=',NNT,' NCT=',NCT
endif
call contact(.true.,ncont_ref,icont_ref)
endif
- if (ns.gt.0 .and. (me.eq.king .or. .not.out1file) ) then
- write (iout,'(/a,i3,a)')
- & 'The chain contains',ns,' disulfide-bridging cysteines.'
+ if (ns.gt.0) then
+C write (iout,'(/a,i3,a)')
+C & '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"
dyn_ss_mask(iss(i))=.true.
enddo
endif
+c Read distance restraints
+ if (constr_dist.gt.0) then
+ call read_dist_constr
+ call hpb_partition
+ endif
return
end
c-----------------------------------------------------------------------------
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
read (rekord(iread:),*) wartosc
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'
open (isidep1,file=sidepname,status="old")
call getenv('SCCORPAR',sccorname)
open (isccor,file=sccorname,status="old")
+ call getenv('LIPTRANPAR',liptranname)
+ open (iliptranpar,file=liptranname,status='old')
#ifndef OLDSCP
C
C 8/9/01 In the newest version SCp interaction constants are read from a file
#endif
return
end
+ subroutine read_dist_constr
+ implicit real*8 (a-h,o-z)
+ include 'DIMENSIONS'
+#ifdef MPI
+ include 'mpif.h'
+#endif
+ 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
+ logical lprn /.true./
+ write (iout,*) "Calling read_dist_constr"
+C write (iout,*) "nres",nres," nstart_sup",nstart_sup," nsup",nsup
+C call flush(iout)
+ write(iout,*) "TU sie wywalam?"
+ call card_concat(controlcard)
+ write (iout,*) controlcard
+ call flush(iout)
+ 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)
+ 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
+ if (lprn)
+ & 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_
+ 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)
+C write (iout,'(a,3i5,f8.2,f10.1)') "+dist.constr ",
+C & nhpb,ihpb(nhpb),jhpb(nhpb),dhpb(nhpb),forcon(nhpb)
+ else
+ read (inp,*) ihpb(nhpb+1),jhpb(nhpb+1),forcon(nhpb+1)
+ endif
+ 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))
+C dhpb(nhpb)=dist(ihpb(nhpb),jhpb(nhpb))
+ write (iout,'(a,3i5,f8.2,f10.1)') "+dist.constr ",
+ & nhpb,ihpb(nhpb),jhpb(nhpb),dhpb(nhpb),forcon(nhpb)
+ endif
+C endif
+ enddo
+ call hpb_partition
+ call flush(iout)
+ return
+ end
projects / unres.git / blobdiff