subroutine etotal(energia)
- implicit real*8 (a-h,o-z)
+ implicit none
include 'DIMENSIONS'
#ifndef ISNAN
external proc_proc
#ifdef MPI
include "mpif.h"
double precision weights_(n_ene)
+ double precision time00
+ integer ierror,ierr
#endif
include 'COMMON.SETUP'
include 'COMMON.IOUNITS'
include 'COMMON.SBRIDGE'
include 'COMMON.CHAIN'
include 'COMMON.VAR'
- include 'COMMON.MD'
+c include 'COMMON.MD'
+ include 'COMMON.QRESTR'
include 'COMMON.CONTROL'
include 'COMMON.TIME1'
include 'COMMON.SPLITELE'
include 'COMMON.TORCNSTR'
+ include 'COMMON.SAXS'
+ double precision evdw,evdw1,evdw2,evdw2_14,ees,eel_loc,
+ & eello_turn3,eello_turn4,edfadis,estr,ehpb,ebe,ethetacnstr,
+ & escloc,etors,edihcnstr,etors_d,esccor,ecorr,ecorr5,ecorr6,eturn6,
+ & eliptran,Eafmforce,Etube,
+ & esaxs_constr,ehomology_constr,edfator,edfanei,edfabet
+ integer n_corr,n_corr1
#ifdef MPI
c print*,"ETOTAL Processor",fg_rank," absolute rank",myrank,
c & " nfgtasks",nfgtasks
weights_(17)=wbond
weights_(18)=scal14
weights_(21)=wsccor
- weights_(22)=wtube
+ weights_(22)=wliptran
+ weights_(25)=wtube
weights_(26)=wsaxs
weights_(28)=wdfa_dist
weights_(29)=wdfa_tor
wbond=weights(17)
scal14=weights(18)
wsccor=weights(21)
- wtube=weights(22)
+ wliptran=weights(22)
+ wtube=weights(25)
wsaxs=weights(26)
wdfa_dist=weights_(28)
wdfa_tor=weights_(29)
else
esccor=0.0d0
endif
+#ifdef FOURBODY
C print *,"PRZED MULIt"
c print *,"Processor",myrank," computed Usccorr"
C
c & n_corr1
c call flush(iout)
endif
+#endif
c print *,"Processor",myrank," computed Ucorr"
c write (iout,*) "nsaxs",nsaxs," saxs_mode",saxs_mode
if (nsaxs.gt.0 .and. saxs_mode.eq.0) then
C print *,"przed lipidami"
if (wliptran.gt.0) then
call Eliptransfer(eliptran)
+ else
+ eliptran=0.0d0
endif
C print *,"za lipidami"
if (AFMlog.gt.0) then
end
c-------------------------------------------------------------------------------
subroutine sum_energy(energia,reduce)
- implicit real*8 (a-h,o-z)
+ implicit none
include 'DIMENSIONS'
#ifndef ISNAN
external proc_proc
#endif
#ifdef MPI
include "mpif.h"
+ integer ierr
+ double precision time00
#endif
include 'COMMON.SETUP'
include 'COMMON.IOUNITS'
include 'COMMON.CONTROL'
include 'COMMON.TIME1'
logical reduce
+ integer i
+ double precision evdw,evdw1,evdw2,evdw2_14,ees,eel_loc,
+ & eello_turn3,eello_turn4,edfadis,estr,ehpb,ebe,ethetacnstr,
+ & escloc,etors,edihcnstr,etors_d,esccor,ecorr,ecorr5,ecorr6,eturn6,
+ & eliptran,Eafmforce,Etube,
+ & esaxs_constr,ehomology_constr,edfator,edfanei,edfabet
+ double precision Uconst,etot
#ifdef MPI
if (nfgtasks.gt.1 .and. reduce) then
#ifdef DEBUG
end
c-------------------------------------------------------------------------------
subroutine sum_gradient
- implicit real*8 (a-h,o-z)
+ implicit none
include 'DIMENSIONS'
#ifndef ISNAN
external proc_proc
#endif
#ifdef MPI
include 'mpif.h'
+ integer ierror,ierr
+ double precision time00,time01
#endif
double precision gradbufc(3,-1:maxres),gradbufx(3,-1:maxres),
& glocbuf(4*maxres),gradbufc_sum(3,-1:maxres)
include 'COMMON.TIME1'
include 'COMMON.MAXGRAD'
include 'COMMON.SCCOR'
- include 'COMMON.MD'
+c include 'COMMON.MD'
+ include 'COMMON.QRESTR'
+ integer i,j,k
+ double precision scalar
+ double precision gvdwc_norm,gvdwc_scp_norm,gelc_norm,gvdwpp_norm,
+ &gradb_norm,ghpbc_norm,gradcorr_norm,gel_loc_norm,gcorr3_turn_norm,
+ &gcorr4_turn_norm,gradcorr5_norm,gradcorr6_norm,
+ &gcorr6_turn_norm,gsccorrc_norm,gscloc_norm,gvdwx_norm,
+ &gradx_scp_norm,ghpbx_norm,gradxorr_norm,gsccorrx_norm,
+ &gsclocx_norm
#ifdef TIMING
time01=MPI_Wtime()
#endif
gradcorr5_max=0.0d0
gradcorr6_max=0.0d0
gcorr6_turn_max=0.0d0
- gsccorc_max=0.0d0
+ gsccorrc_max=0.0d0
gscloc_max=0.0d0
gvdwx_max=0.0d0
gradx_scp_max=0.0d0
ghpbx_max=0.0d0
gradxorr_max=0.0d0
- gsccorx_max=0.0d0
+ gsccorrx_max=0.0d0
gsclocx_max=0.0d0
do i=1,nct
gvdwc_norm=dsqrt(scalar(gvdwc(1,i),gvdwc(1,i)))
if (gradcorr5_norm.gt.gradcorr5_max)
& gradcorr5_max=gradcorr5_norm
gradcorr6_norm=dsqrt(scalar(gradcorr6(1,i),gradcorr6(1,i)))
- if (gradcorr6_norm.gt.gradcorr6_max) gcorr6_max=gradcorr6_norm
+ if (gradcorr6_norm.gt.gradcorr6_max)gradcorr6_max=gradcorr6_norm
gcorr6_turn_norm=dsqrt(scalar(gcorr6_turn(1,i),
& gcorr6_turn(1,i)))
if (gcorr6_turn_norm.gt.gcorr6_turn_max)
& gcorr6_turn_max=gcorr6_turn_norm
- gsccorr_norm=dsqrt(scalar(gsccorc(1,i),gsccorc(1,i)))
- if (gsccorr_norm.gt.gsccorr_max) gsccorr_max=gsccorr_norm
+ gsccorrc_norm=dsqrt(scalar(gsccorc(1,i),gsccorc(1,i)))
+ if (gsccorrc_norm.gt.gsccorrc_max) gsccorrc_max=gsccorrc_norm
gscloc_norm=dsqrt(scalar(gscloc(1,i),gscloc(1,i)))
if (gscloc_norm.gt.gscloc_max) gscloc_max=gscloc_norm
gvdwx_norm=dsqrt(scalar(gvdwx(1,i),gvdwx(1,i)))
write (istat,'(1h#,21f10.2)') gvdwc_max,gvdwc_scp_max,
& gelc_max,gvdwpp_max,gradb_max,ghpbc_max,
& gradcorr_max,gel_loc_max,gcorr3_turn_max,gcorr4_turn_max,
- & gradcorr5_max,gradcorr6_max,gcorr6_turn_max,gsccorc_max,
+ & gradcorr5_max,gradcorr6_max,gcorr6_turn_max,gsccorrc_max,
& gscloc_max,gvdwx_max,gradx_scp_max,ghpbx_max,gradxorr_max,
- & gsccorx_max,gsclocx_max
+ & gsccorrx_max,gsclocx_max
close(istat)
if (gvdwc_max.gt.1.0d4) then
write (iout,*) "gvdwc gvdwx gradb gradbx"
end
c-------------------------------------------------------------------------------
subroutine rescale_weights(t_bath)
- implicit real*8 (a-h,o-z)
+ implicit none
+#ifdef MPI
+ include 'mpif.h'
+ integer ierror
+#endif
include 'DIMENSIONS'
include 'COMMON.IOUNITS'
include 'COMMON.FFIELD'
include 'COMMON.SBRIDGE'
include 'COMMON.CONTROL'
+ double precision t_bath
+ double precision facT,facT2,facT3,facT4,facT5
double precision kfac /2.4d0/
double precision x,x2,x3,x4,x5,licznik /1.12692801104297249644/
c facT=temp0/t_bath
end
C------------------------------------------------------------------------
subroutine enerprint(energia)
- implicit real*8 (a-h,o-z)
+ implicit none
include 'DIMENSIONS'
include 'COMMON.IOUNITS'
include 'COMMON.FFIELD'
include 'COMMON.SBRIDGE'
- include 'COMMON.MD'
+ include 'COMMON.QRESTR'
double precision energia(0:n_ene)
+ double precision evdw,evdw1,evdw2,evdw2_14,ees,eel_loc,
+ & eello_turn3,eello_turn4,edfadis,estr,ehpb,ebe,ethetacnstr,
+ & escloc,etors,edihcnstr,etors_d,esccor,ecorr,ecorr5,ecorr6,
+ & eello_turn6,
+ & eliptran,Eafmforce,Etube,
+ & esaxs,ehomology_constr,edfator,edfanei,edfabet,etot
etot=energia(0)
evdw=energia(1)
evdw2=energia(2)
write (iout,10) evdw,wsc,evdw2,wscp,ees,welec,evdw1,wvdwpp,
& estr,wbond,ebe,wang,
& escloc,wscloc,etors,wtor,etors_d,wtor_d,ehpb,wstrain,
+#ifdef FOURBODY
& ecorr,wcorr,
- & ecorr5,wcorr5,ecorr6,wcorr6,eel_loc,wel_loc,eello_turn3,wturn3,
- & eello_turn4,wturn4,eello_turn6,wturn6,esccor,wsccor,edihcnstr,
- & ethetacnstr,ebr*nss,Uconst,wumb,eliptran,wliptran,Eafmforc,
+ & ecorr5,wcorr5,ecorr6,wcorr6,
+#endif
+ & eel_loc,wel_loc,eello_turn3,wturn3,
+ & eello_turn4,wturn4,
+#ifdef FOURBODY
+ & eello_turn6,wturn6,
+#endif
+ & esccor,wsccor,edihcnstr,
+ & ethetacnstr,ebr*nss,Uconst,wumb,eliptran,wliptran,Eafmforce,
& etube,wtube,esaxs,wsaxs,ehomology_constr,
& edfadis,wdfa_dist,edfator,wdfa_tor,edfanei,wdfa_nei,
& edfabet,wdfa_beta,
& 'ETORSD=',1pE16.6,' WEIGHT=',1pE16.6,' (double torsional)'/
& 'EHBP= ',1pE16.6,' WEIGHT=',1pE16.6,
& ' (SS bridges & dist. cnstr.)'/
+#ifdef FOURBODY
& 'ECORR4=',1pE16.6,' WEIGHT=',1pE16.6,' (multi-body)'/
& 'ECORR5=',1pE16.6,' WEIGHT=',1pE16.6,' (multi-body)'/
& 'ECORR6=',1pE16.6,' WEIGHT=',1pE16.6,' (multi-body)'/
+#endif
& 'EELLO= ',1pE16.6,' WEIGHT=',1pE16.6,' (electrostatic-local)'/
& 'ETURN3=',1pE16.6,' WEIGHT=',1pE16.6,' (turns, 3rd order)'/
& 'ETURN4=',1pE16.6,' WEIGHT=',1pE16.6,' (turns, 4th order)'/
+#ifdef FOURBODY
& 'ETURN6=',1pE16.6,' WEIGHT=',1pE16.6,' (turns, 6th order)'/
+#endif
& 'ESCCOR=',1pE16.6,' WEIGHT=',1pE16.6,' (backbone-rotamer corr)'/
& 'EDIHC= ',1pE16.6,' (virtual-bond dihedral angle restraints)'/
& 'ETHETC=',1pE16.6,' (virtual-bond angle restraints)'/
write (iout,10) evdw,wsc,evdw2,wscp,ees,welec,
& estr,wbond,ebe,wang,
& escloc,wscloc,etors,wtor,etors_d,wtor_d,ehpb,wstrain,
+#ifdef FOURBODY
& ecorr,wcorr,
- & ecorr5,wcorr5,ecorr6,wcorr6,eel_loc,wel_loc,eello_turn3,wturn3,
- & eello_turn4,wturn4,eello_turn6,wturn6,esccor,wsccor,edihcnstr,
+ & ecorr5,wcorr5,ecorr6,wcorr6,
+#endif
+ & eel_loc,wel_loc,eello_turn3,wturn3,
+ & eello_turn4,wturn4,
+#ifdef FOURBODY
+ & eello_turn6,wturn6,
+#endif
+ & esccor,wsccor,edihcnstr,
& ethetacnstr,ebr*nss,Uconst,wumb,eliptran,wliptran,Eafmforc,
& etube,wtube,esaxs,wsaxs,ehomology_constr,
& edfadis,wdfa_dist,edfator,wdfa_tor,edfanei,wdfa_nei,
& 'ETORSD=',1pE16.6,' WEIGHT=',1pE16.6,' (double torsional)'/
& 'EHBP= ',1pE16.6,' WEIGHT=',1pE16.6,
& ' (SS bridges & dist. restr.)'/
+#ifdef FOURBODY
& 'ECORR4=',1pE16.6,' WEIGHT=',1pE16.6,' (multi-body)'/
& 'ECORR5=',1pE16.6,' WEIGHT=',1pE16.6,' (multi-body)'/
& 'ECORR6=',1pE16.6,' WEIGHT=',1pE16.6,' (multi-body)'/
+#endif
& 'EELLO= ',1pE16.6,' WEIGHT=',1pE16.6,' (electrostatic-local)'/
& 'ETURN3=',1pE16.6,' WEIGHT=',1pE16.6,' (turns, 3rd order)'/
& 'ETURN4=',1pE16.6,' WEIGHT=',1pE16.6,' (turns, 4th order)'/
+#ifdef FOURBODY
& 'ETURN6=',1pE16.6,' WEIGHT=',1pE16.6,' (turns, 6th order)'/
+#endif
& 'ESCCOR=',1pE16.6,' WEIGHT=',1pE16.6,' (backbone-rotamer corr)'/
& 'EDIHC= ',1pE16.6,' (virtual-bond dihedral angle restraints)'/
& 'ETHETC=',1pE16.6,' (virtual-bond angle restraints)'/
C This subroutine calculates the interaction energy of nonbonded side chains
C assuming the LJ potential of interaction.
C
- implicit real*8 (a-h,o-z)
+ implicit none
+ double precision accur
include 'DIMENSIONS'
parameter (accur=1.0d-10)
include 'COMMON.GEO'
include 'COMMON.SBRIDGE'
include 'COMMON.NAMES'
include 'COMMON.IOUNITS'
+ include 'COMMON.SPLITELE'
+#ifdef FOURBODY
include 'COMMON.CONTACTS'
- dimension gg(3)
+ include 'COMMON.CONTMAT'
+#endif
+ double precision gg(3)
+ double precision evdw,evdwij
+ integer i,j,k,itypi,itypj,itypi1,num_conti,iint
+ double precision xi,yi,zi,xj,yj,zj,rij,eps0ij,fac,e1,e2,rrij,
+ & sigij,r0ij,rcut,sqrij,sss1,sssgrad1
+ double precision fcont,fprimcont
+ double precision sscale,sscagrad
c write(iout,*)'Entering ELJ nnt=',nnt,' nct=',nct,' expon=',expon
evdw=0.0D0
do i=iatsc_s,iatsc_e
C Change 12/1/95 to calculate four-body interactions
rij=xj*xj+yj*yj+zj*zj
rrij=1.0D0/rij
+ sqrij=dsqrt(rij)
+ sss1=sscale(sqrij,r_cut_int)
+ if (sss1.eq.0.0d0) cycle
+ sssgrad1=sscagrad(sqrij,r_cut_int)
+
c write (iout,*)'i=',i,' j=',j,' itypi=',itypi,' itypj=',itypj
eps0ij=eps(itypi,itypj)
fac=rrij**expon2
cd & restyp(itypi),i,restyp(itypj),j,a(itypi,itypj),
cd & bb(itypi,itypj),1.0D0/dsqrt(rrij),evdwij,epsi,sigm,
cd & (c(k,i),k=1,3),(c(k,j),k=1,3)
- evdw=evdw+evdwij
+ evdw=evdw+sss1*evdwij
C
C Calculate the components of the gradient in DC and X
C
- fac=-rrij*(e1+evdwij)
+ fac=-rrij*(e1+evdwij)*sss1
+ & +evdwij*sssgrad1/sqrij/expon
gg(1)=xj*fac
gg(2)=yj*fac
gg(3)=zj*fac
cgrad enddo
cgrad enddo
C
+#ifdef FOURBODY
C 12/1/95, revised on 5/20/97
C
C Calculate the contact function. The ith column of the array JCONT will
cd & i,j,(gacont(kk,num_conti,i),kk=1,3)
endif
endif
+#endif
enddo ! j
enddo ! iint
C Change 12/1/95
+#ifdef FOURBODY
num_cont(i)=num_conti
+#endif
enddo ! i
do i=1,nct
do j=1,3
C This subroutine calculates the interaction energy of nonbonded side chains
C assuming the LJK potential of interaction.
C
- implicit real*8 (a-h,o-z)
+ implicit none
include 'DIMENSIONS'
include 'COMMON.GEO'
include 'COMMON.VAR'
include 'COMMON.INTERACT'
include 'COMMON.IOUNITS'
include 'COMMON.NAMES'
- dimension gg(3)
+ include 'COMMON.SPLITELE'
+ double precision gg(3)
+ double precision evdw,evdwij
+ integer i,j,k,itypi,itypj,itypi1,iint
+ double precision xi,yi,zi,xj,yj,zj,rij,eps0ij,fac,e1,e2,rrij,
+ & fac_augm,e_augm,r_inv_ij,r_shift_inv,sss1,sssgrad1
logical scheck
+ double precision sscale,sscagrad
c print *,'Entering ELJK nnt=',nnt,' nct=',nct,' expon=',expon
evdw=0.0D0
do i=iatsc_s,iatsc_e
e_augm=augm(itypi,itypj)*fac_augm
r_inv_ij=dsqrt(rrij)
rij=1.0D0/r_inv_ij
+ sss1=sscale(rij,r_cut_int)
+ if (sss1.eq.0.0d0) cycle
+ sssgrad1=sscagrad(rij,r_cut_int)
r_shift_inv=1.0D0/(rij+r0(itypi,itypj)-sigma(itypi,itypj))
fac=r_shift_inv**expon
C have you changed here?
cd & bb(itypi,itypj),augm(itypi,itypj),epsi,sigm,
cd & sigma(itypi,itypj),1.0D0/dsqrt(rrij),evdwij,
cd & (c(k,i),k=1,3),(c(k,j),k=1,3)
- evdw=evdw+evdwij
+ evdw=evdw+evdwij*sss1
C
C Calculate the components of the gradient in DC and X
C
fac=-2.0D0*rrij*e_augm-r_inv_ij*r_shift_inv*(e1+e1+e2)
+ & +evdwij*sssgrad1*r_inv_ij/expon
gg(1)=xj*fac
gg(2)=yj*fac
gg(3)=zj*fac
C This subroutine calculates the interaction energy of nonbonded side chains
C assuming the Berne-Pechukas potential of interaction.
C
- implicit real*8 (a-h,o-z)
+ implicit none
include 'DIMENSIONS'
include 'COMMON.GEO'
include 'COMMON.VAR'
include 'COMMON.INTERACT'
include 'COMMON.IOUNITS'
include 'COMMON.CALC'
+ include 'COMMON.SPLITELE'
+ integer icall
common /srutu/ icall
+ double precision evdw
+ integer itypi,itypj,itypi1,iint,ind
+ double precision eps0ij,epsi,sigm,fac,e1,e2,rrij,xi,yi,zi,
+ & sss1,sssgrad1
+ double precision sscale,sscagrad
c double precision rrsave(maxdim)
logical lprn
evdw=0.0D0
cd rrij=rrsave(ind)
cd endif
rij=dsqrt(rrij)
+ sss1=sscale(1.0d0/rij,r_cut_int)
+ if (sss1.eq.0.0d0) cycle
+ sssgrad1=sscagrad(1.0d0/rij,r_cut_int)
C Calculate the angle-dependent terms of energy & contributions to derivatives.
call sc_angular
C Calculate whole angle-dependent part of epsilon and contributions
eps2der=evdwij*eps3rt
eps3der=evdwij*eps2rt
evdwij=evdwij*eps2rt*eps3rt
- evdw=evdw+evdwij
+ evdw=evdw+sss1*evdwij
if (lprn) then
sigm=dabs(aa/bb)**(1.0D0/6.0D0)
epsi=bb**2/aa
fac=-expon*(e1+evdwij)
sigder=fac/sigsq
fac=rrij*fac
+ & +evdwij*sssgrad1/sss1*rij
C Calculate radial part of the gradient
gg(1)=xj*fac
gg(2)=yj*fac
C This subroutine calculates the interaction energy of nonbonded side chains
C assuming the Gay-Berne potential of interaction.
C
- implicit real*8 (a-h,o-z)
+ implicit none
include 'DIMENSIONS'
include 'COMMON.GEO'
include 'COMMON.VAR'
include 'COMMON.SPLITELE'
include 'COMMON.SBRIDGE'
logical lprn
- integer xshift,yshift,zshift
-
+ integer xshift,yshift,zshift,subchap
+ double precision evdw
+ integer itypi,itypj,itypi1,iint,ind
+ double precision eps0ij,epsi,sigm,fac,e1,e2,rrij,xi,yi,zi
+ double precision fracinbuf,sslipi,evdwij_przed_tri,sig0ij,
+ & sslipj,ssgradlipj,ssgradlipi,dist_init,xj_safe,yj_safe,zj_safe,
+ & xj_temp,yj_temp,zj_temp,dist_temp,sig,rij_shift,faclip
+ double precision dist,sscale,sscagrad,sscagradlip,sscalelip
evdw=0.0D0
ccccc energy_dec=.false.
C print *,'Entering EGB nnt=',nnt,' nct=',nct,' expon=',expon
c & dc_norm(1,nres+j),dc_norm(2,nres+j),dc_norm(3,nres+j)
rrij=1.0D0/(xj*xj+yj*yj+zj*zj)
rij=dsqrt(rrij)
- sss=sscale((1.0d0/rij)/sigma(itypi,itypj))
- sssgrad=sscagrad((1.0d0/rij)/sigma(itypi,itypj))
-
+ sss=sscale(1.0d0/rij,r_cut_int)
c write (iout,'(a7,4f8.3)')
c & "ssscale",sss,((1.0d0/rij)/sigma(itypi,itypj)),r_cut,rlamb
- if (sss.gt.0.0d0) then
+ if (sss.eq.0.0d0) cycle
+ sssgrad=sscagrad(1.0d0/rij,r_cut_int)
C Calculate angle-dependent terms of energy and contributions to their
C derivatives.
call sc_angular
& evdwij
endif
- if (energy_dec) write (iout,'(a6,2i5,0pf7.3)')
- & 'evdw',i,j,evdwij
+ if (energy_dec) write (iout,'(a,2i5,3f10.5)')
+ & 'r sss evdw',i,j,rij,sss,evdwij
C Calculate gradient components.
e1=e1*eps1*eps2rt**2*eps3rt**2
fac=rij*fac
c print '(2i4,6f8.4)',i,j,sss,sssgrad*
c & evdwij,fac,sigma(itypi,itypj),expon
- fac=fac+evdwij/sss*sssgrad/sigma(itypi,itypj)*rij
+ fac=fac+evdwij*sssgrad/sss*rij
c fac=0.0d0
C Calculate the radial part of the gradient
gg_lipi(3)=eps1*(eps2rt*eps2rt)
- &*(eps3rt*eps3rt)*sss/2.0d0*(faclip*faclip*
- & (aa_lip(itypi,itypj)-aa_aq(itypi,itypj))
- &+faclip*(bb_lip(itypi,itypj)-bb_aq(itypi,itypj)))
+ & *(eps3rt*eps3rt)*sss/2.0d0*(faclip*faclip*
+ & (aa_lip(itypi,itypj)-aa_aq(itypi,itypj))
+ & +faclip*(bb_lip(itypi,itypj)-bb_aq(itypi,itypj)))
gg_lipj(3)=ssgradlipj*gg_lipi(3)
gg_lipi(3)=gg_lipi(3)*ssgradlipi
C gg_lipi(3)=0.0d0
gg(2)=yj*fac
gg(3)=zj*fac
C Calculate angular part of the gradient.
+c call sc_grad_scale(sss)
call sc_grad
- endif
ENDIF ! dyn_ss
enddo ! j
enddo ! iint
C This subroutine calculates the interaction energy of nonbonded side chains
C assuming the Gay-Berne-Vorobjev potential of interaction.
C
- implicit real*8 (a-h,o-z)
+ implicit none
include 'DIMENSIONS'
include 'COMMON.GEO'
include 'COMMON.VAR'
include 'COMMON.INTERACT'
include 'COMMON.IOUNITS'
include 'COMMON.CALC'
- integer xshift,yshift,zshift
+ include 'COMMON.SPLITELE'
+ integer xshift,yshift,zshift,subchap
+ integer icall
common /srutu/ icall
logical lprn
+ double precision evdw
+ integer itypi,itypj,itypi1,iint,ind
+ double precision eps0ij,epsi,sigm,fac,e1,e2,rrij,r0ij,
+ & xi,yi,zi,fac_augm,e_augm
+ double precision fracinbuf,sslipi,evdwij_przed_tri,sig0ij,
+ & sslipj,ssgradlipj,ssgradlipi,dist_init,xj_safe,yj_safe,zj_safe,
+ & xj_temp,yj_temp,zj_temp,dist_temp,sig,rij_shift,faclip,sssgrad1
+ double precision dist,sscale,sscagrad,sscagradlip,sscalelip
evdw=0.0D0
c print *,'Entering EGB nnt=',nnt,' nct=',nct,' expon=',expon
evdw=0.0D0
dzj=dc_norm(3,nres+j)
rrij=1.0D0/(xj*xj+yj*yj+zj*zj)
rij=dsqrt(rrij)
+ sss=sscale(1.0d0/rij,r_cut_int)
+ if (sss.eq.0.0d0) cycle
+ sssgrad=sscagrad(1.0d0/rij,r_cut_int)
C Calculate angle-dependent terms of energy and contributions to their
C derivatives.
call sc_angular
fac=-expon*(e1+evdwij)*rij_shift
sigder=fac*sigder
fac=rij*fac-2*expon*rrij*e_augm
- fac=fac+evdwij/sss*sssgrad/sigma(itypi,itypj)*rij
+ fac=fac+(evdwij+e_augm)*sssgrad/sss*rij
C Calculate the radial part of the gradient
gg(1)=xj*fac
gg(2)=yj*fac
gg(3)=zj*fac
C Calculate angular part of the gradient.
+c call sc_grad_scale(sss)
call sc_grad
enddo ! j
enddo ! iint
include 'COMMON.SBRIDGE'
include 'COMMON.NAMES'
include 'COMMON.IOUNITS'
- include 'COMMON.CONTACTS'
+c include 'COMMON.CONTACTS'
dimension gg(3)
cd print *,'Entering Esoft_sphere nnt=',nnt,' nct=',nct
evdw=0.0D0
include 'COMMON.CHAIN'
include 'COMMON.DERIV'
include 'COMMON.INTERACT'
- include 'COMMON.CONTACTS'
+c include 'COMMON.CONTACTS'
include 'COMMON.TORSION'
include 'COMMON.VECTORS'
include 'COMMON.FFIELD'
zj=zj_safe-zmedi
endif
rij=xj*xj+yj*yj+zj*zj
- sss=sscale(sqrt(rij))
- sssgrad=sscagrad(sqrt(rij))
+ sss=sscale(sqrt(rij),r_cut_int)
+ sssgrad=sscagrad(sqrt(rij),r_cut_int)
if (rij.lt.r0ijsq) then
evdw1ij=0.25d0*(rij-r0ijsq)**2
fac=rij-r0ijsq
#endif
return
end
-C-----------------------------------------------------------------------------
- subroutine check_vecgrad
- implicit real*8 (a-h,o-z)
- include 'DIMENSIONS'
- include 'COMMON.IOUNITS'
- include 'COMMON.GEO'
- include 'COMMON.VAR'
- include 'COMMON.LOCAL'
- include 'COMMON.CHAIN'
- include 'COMMON.VECTORS'
- dimension uygradt(3,3,2,maxres),uzgradt(3,3,2,maxres)
- dimension uyt(3,maxres),uzt(3,maxres)
- dimension uygradn(3,3,2),uzgradn(3,3,2),erij(3)
- double precision delta /1.0d-7/
- call vec_and_deriv
-cd do i=1,nres
-crc write(iout,'(2i5,2(3f10.5,5x))') i,1,dc_norm(:,i)
-crc write(iout,'(2i5,2(3f10.5,5x))') i,2,uy(:,i)
-crc write(iout,'(2i5,2(3f10.5,5x)/)')i,3,uz(:,i)
-cd write(iout,'(2i5,2(3f10.5,5x))') i,1,
-cd & (dc_norm(if90,i),if90=1,3)
-cd write(iout,'(2i5,2(3f10.5,5x))') i,2,(uy(if90,i),if90=1,3)
-cd write(iout,'(2i5,2(3f10.5,5x)/)')i,3,(uz(if90,i),if90=1,3)
-cd write(iout,'(a)')
-cd enddo
- do i=1,nres
- do j=1,2
- do k=1,3
- do l=1,3
- uygradt(l,k,j,i)=uygrad(l,k,j,i)
- uzgradt(l,k,j,i)=uzgrad(l,k,j,i)
- enddo
- enddo
- enddo
- enddo
- call vec_and_deriv
- do i=1,nres
- do j=1,3
- uyt(j,i)=uy(j,i)
- uzt(j,i)=uz(j,i)
- enddo
- enddo
- do i=1,nres
-cd write (iout,*) 'i=',i
- do k=1,3
- erij(k)=dc_norm(k,i)
- enddo
- do j=1,3
- do k=1,3
- dc_norm(k,i)=erij(k)
- enddo
- dc_norm(j,i)=dc_norm(j,i)+delta
-c fac=dsqrt(scalar(dc_norm(1,i),dc_norm(1,i)))
-c do k=1,3
-c dc_norm(k,i)=dc_norm(k,i)/fac
-c enddo
-c write (iout,*) (dc_norm(k,i),k=1,3)
-c write (iout,*) (erij(k),k=1,3)
- call vec_and_deriv
- do k=1,3
- uygradn(k,j,1)=(uy(k,i)-uyt(k,i))/delta
- uygradn(k,j,2)=(uy(k,i-1)-uyt(k,i-1))/delta
- uzgradn(k,j,1)=(uz(k,i)-uzt(k,i))/delta
- uzgradn(k,j,2)=(uz(k,i-1)-uzt(k,i-1))/delta
- enddo
-c write (iout,'(i5,3f8.5,3x,3f8.5,5x,3f8.5,3x,3f8.5)')
-c & j,(uzgradt(k,j,1,i),k=1,3),(uzgradn(k,j,1),k=1,3),
-c & (uzgradt(k,j,2,i-1),k=1,3),(uzgradn(k,j,2),k=1,3)
- enddo
- do k=1,3
- dc_norm(k,i)=erij(k)
- enddo
-cd do k=1,3
-cd write (iout,'(i5,3f8.5,3x,3f8.5,5x,3f8.5,3x,3f8.5)')
-cd & k,(uygradt(k,l,1,i),l=1,3),(uygradn(k,l,1),l=1,3),
-cd & (uygradt(k,l,2,i-1),l=1,3),(uygradn(k,l,2),l=1,3)
-cd write (iout,'(i5,3f8.5,3x,3f8.5,5x,3f8.5,3x,3f8.5)')
-cd & k,(uzgradt(k,l,1,i),l=1,3),(uzgradn(k,l,1),l=1,3),
-cd & (uzgradt(k,l,2,i-1),l=1,3),(uzgradn(k,l,2),l=1,3)
-cd write (iout,'(a)')
-cd enddo
- enddo
- return
- end
C--------------------------------------------------------------------------
subroutine set_matrices
implicit real*8 (a-h,o-z)
include 'COMMON.CHAIN'
include 'COMMON.DERIV'
include 'COMMON.INTERACT'
- include 'COMMON.CONTACTS'
+ include 'COMMON.CORRMAT'
include 'COMMON.TORSION'
include 'COMMON.VECTORS'
include 'COMMON.FFIELD'
#else
do i=3,nres+1
#endif
- if (i.gt. nnt+2 .and. i.lt.nct+2) then
+ ii=ireschain(i-2)
+c write (iout,*) "i",i,i-2," ii",ii
+ if (ii.eq.0) cycle
+ innt=chain_border(1,ii)
+ inct=chain_border(2,ii)
+c write (iout,*) "i",i,i-2," ii",ii," innt",innt," inct",inct
+c if (i.gt. nnt+2 .and. i.lt.nct+2) then
+ if (i.gt. innt+2 .and. i.lt.inct+2) then
iti = itype2loc(itype(i-2))
else
iti=nloctyp
endif
c if (i.gt. iatel_s+1 .and. i.lt.iatel_e+4) then
- if (i.gt. nnt+1 .and. i.lt.nct+1) then
+ if (i.gt. innt+1 .and. i.lt.inct+1) then
iti1 = itype2loc(itype(i-1))
else
iti1=nloctyp
endif
-c write(iout,*),i
+c write(iout,*),"i",i,i-2," iti",itype(i-2),iti,
+c & " iti1",itype(i-1),iti1
#ifdef NEWCORR
cost1=dcos(theta(i-1))
sint1=dsin(theta(i-1))
write (iout,*) 'theta=', theta(i-1)
#endif
#else
- if (i.gt. nnt+2 .and. i.lt.nct+2) then
+ if (i.gt. innt+2 .and. i.lt.inct+2) then
+c if (i.gt. nnt+2 .and. i.lt.nct+2) then
iti = itype2loc(itype(i-2))
else
iti=nloctyp
write(iout,*) 'b2=',(b2(k,i-2),k=1,2)
#endif
enddo
+ mu=0.0d0
#ifdef PARMAT
do i=ivec_start+2,ivec_end+2
#else
do i=3,nres+1
#endif
- if (i .lt. nres+1) then
+c if (itype(i-1).eq.ntyp1 .and. itype(i).eq.ntyp1) cycle
+ if (i .lt. nres+1 .and. itype(i-1).lt.ntyp1) then
sin1=dsin(phi(i))
cos1=dcos(phi(i))
sintab(i-2)=sin1
Ug2(2,1,i-2)=0.0d0
Ug2(2,2,i-2)=0.0d0
endif
- if (i .gt. 3 .and. i .lt. nres+1) then
+ if (i .gt. 3) then
obrot_der(1,i-2)=-sin1
obrot_der(2,i-2)= cos1
Ugder(1,1,i-2)= sin1
Ug2der(2,2,i-2)=0.0d0
endif
c if (i.gt. iatel_s+2 .and. i.lt.iatel_e+5) then
- if (i.gt. nnt+2 .and. i.lt.nct+2) then
+c if (i.gt. nnt+2 .and. i.lt.nct+2) then
+ if (i.gt.nnt+2 .and.i.lt.nct+2) then
iti = itype2loc(itype(i-2))
else
iti=nloctyp
c write(iout,*) "Macierz EUG",
c & eug(1,1,i-2),eug(1,2,i-2),eug(2,1,i-2),
c & eug(2,2,i-2)
+#ifdef FOURBODY
if (wcorr4.gt.0.0d0 .or. wcorr5.gt.0.0d0 .or. wcorr6.gt.0.0d0)
& then
call matmat2(CC(1,1,i-2),Ug(1,1,i-2),CUg(1,1,i-2))
call matvec2(Ctilde(1,1,i-1),obrot(1,i-2),Ctobr(1,i-2))
call matvec2(Dtilde(1,1,i-2),obrot2(1,i-2),Dtobr2(1,i-2))
endif
+#endif
else
do k=1,2
Ub2(k,i-2)=0.0d0
cd write (iout,*) 'mu1',mu1(:,i-2)
cd write (iout,*) 'mu2',mu2(:,i-2)
cd write (iout,*) 'mu',i-2,mu(:,i-2)
+#ifdef FOURBODY
if (wcorr4.gt.0.0d0 .or. wcorr5.gt.0.0d0 .or.wcorr6.gt.0.0d0)
& then
call matmat2(CC(1,1,i-1),Ugder(1,1,i-2),CUgder(1,1,i-2))
call matmat2(EUg(1,1,i-2),DD(1,1,i-1),EUgD(1,1,i-2))
call matmat2(EUgder(1,1,i-2),DD(1,1,i-1),EUgDder(1,1,i-2))
endif
+#endif
enddo
+#ifdef FOURBODY
C Matrices dependent on two consecutive virtual-bond dihedrals.
C The order of matrices is from left to right.
if (wcorr4.gt.0.0d0 .or. wcorr5.gt.0.0d0 .or.wcorr6.gt.0.0d0)
call matmat2(auxmat(1,1),EUg(1,1,i),Ug2DtEUgder(1,1,1,i))
enddo
endif
+#endif
#if defined(MPI) && defined(PARMAT)
#ifdef DEBUG
c if (fg_rank.eq.0) then
call MPI_Allgatherv(sintab2(ivec_start),ivec_count(fg_rank1),
& MPI_DOUBLE_PRECISION,sintab2(1),ivec_count(0),ivec_displ(0),
& MPI_DOUBLE_PRECISION,FG_COMM1,IERR)
+#ifdef FOURBODY
if (wcorr4.gt.0.0d0 .or. wcorr5.gt.0.0d0 .or. wcorr6.gt.0.0d0)
& then
call MPI_Allgatherv(Ctobr(1,ivec_start),ivec_count(fg_rank1),
& MPI_MAT2,Ug2DtEUgder(1,1,1,1),ivec_count(0),ivec_displ(0),
& MPI_MAT2,FG_COMM1,IERR)
endif
+#endif
#else
c Passes matrix info through the ring
isend=fg_rank1
& iprev,6600+irecv,FG_COMM,status,IERR)
c write (iout,*) "Gather PRECOMP12"
c call flush(iout)
+#ifdef FOURBODY
if (wcorr4.gt.0.0d0 .or. wcorr5.gt.0.0d0 .or. wcorr6.gt.0.0d0)
& then
call MPI_SENDRECV(ug2db1t(1,ivec_displ(isend)+1),1,
& Ug2DtEUgder(1,1,1,ivec_displ(irecv)+1),1,
& MPI_PRECOMP23(lenrecv),
& iprev,9900+irecv,FG_COMM,status,IERR)
+#endif
c write (iout,*) "Gather PRECOMP23"
c call flush(iout)
endif
cd enddo
return
end
-C--------------------------------------------------------------------------
+C-----------------------------------------------------------------------------
subroutine eelec(ees,evdw1,eel_loc,eello_turn3,eello_turn4)
C
C This subroutine calculates the average interaction energy and its gradient
include 'COMMON.CHAIN'
include 'COMMON.DERIV'
include 'COMMON.INTERACT'
+#ifdef FOURBODY
include 'COMMON.CONTACTS'
+ include 'COMMON.CONTMAT'
+#endif
+ include 'COMMON.CORRMAT'
include 'COMMON.TORSION'
include 'COMMON.VECTORS'
include 'COMMON.FFIELD'
eello_turn3=0.0d0
eello_turn4=0.0d0
ind=0
+#ifdef FOURBODY
do i=1,nres
num_cont_hb(i)=0
enddo
+#endif
cd print '(a)','Enter EELEC'
cd write (iout,*) 'iatel_s=',iatel_s,' iatel_e=',iatel_e
do i=1,nres
num_conti=0
call eelecij(i,i+2,ees,evdw1,eel_loc)
if (wturn3.gt.0.0d0) call eturn3(i,eello_turn3)
+#ifdef FOURBODY
num_cont_hb(i)=num_conti
+#endif
enddo
do i=iturn4_start,iturn4_end
if (i.lt.1) cycle
zmedi=mod(zmedi,boxzsize)
if (zmedi.lt.0) zmedi=zmedi+boxzsize
+#ifdef FOURBODY
num_conti=num_cont_hb(i)
+#endif
c write(iout,*) "JESTEM W PETLI"
call eelecij(i,i+3,ees,evdw1,eel_loc)
if (wturn4.gt.0.0d0 .and. itype(i+2).ne.ntyp1)
& call eturn4(i,eello_turn4)
+#ifdef FOURBODY
num_cont_hb(i)=num_conti
+#endif
enddo ! i
C Loop over all neighbouring boxes
C do xshift=-1,1
c endif
c write (iout,*) 'i',i,' ielstart',ielstart(i),' ielend',ielend(i)
+#ifdef FOURBODY
num_conti=num_cont_hb(i)
+#endif
C I TU KURWA
do j=ielstart(i),ielend(i)
C do j=16,17
&) cycle
call eelecij(i,j,ees,evdw1,eel_loc)
enddo ! j
+#ifdef FOURBODY
num_cont_hb(i)=num_conti
+#endif
enddo ! i
C enddo ! zshift
C enddo ! yshift
end
C-------------------------------------------------------------------------------
subroutine eelecij(i,j,ees,evdw1,eel_loc)
- implicit real*8 (a-h,o-z)
+ implicit none
include 'DIMENSIONS'
#ifdef MPI
include "mpif.h"
include 'COMMON.CHAIN'
include 'COMMON.DERIV'
include 'COMMON.INTERACT'
+#ifdef FOURBODY
include 'COMMON.CONTACTS'
+ include 'COMMON.CONTMAT'
+#endif
+ include 'COMMON.CORRMAT'
include 'COMMON.TORSION'
include 'COMMON.VECTORS'
include 'COMMON.FFIELD'
include 'COMMON.TIME1'
include 'COMMON.SPLITELE'
include 'COMMON.SHIELD'
- dimension ggg(3),gggp(3),gggm(3),erij(3),dcosb(3),dcosg(3),
+ double precision ggg(3),gggp(3),gggm(3),erij(3),dcosb(3),dcosg(3),
& erder(3,3),uryg(3,3),urzg(3,3),vryg(3,3),vrzg(3,3)
double precision acipa(2,2),agg(3,4),aggi(3,4),aggi1(3,4),
& aggj(3,4),aggj1(3,4),a_temp(2,2),muij(4),gmuij1(4),gmuji1(4),
& gmuij2(4),gmuji2(4)
+ double precision dxi,dyi,dzi
+ double precision dx_normi,dy_normi,dz_normi,aux
+ integer j1,j2,lll,num_conti
common /locel/ a_temp,agg,aggi,aggi1,aggj,aggj1,a22,a23,a32,a33,
& dxi,dyi,dzi,dx_normi,dy_normi,dz_normi,xmedi,ymedi,zmedi,
& num_conti,j1,j2
+ integer k,i,j,iteli,itelj,kkk,l,kkll,m,isubchap,ilist,iresshield
+ double precision ael6i,rrmij,rmij,r0ij,fcont,fprimcont,ees0tmp
+ double precision ees,evdw1,eel_loc,aaa,bbb,ael3i
+ double precision dxj,dyj,dzj,dx_normj,dy_normj,dz_normj,xj,yj,zj,
+ & rij,r3ij,r6ij,cosa,cosb,cosg,fac,ev1,ev2,fac3,fac4,
+ & evdwij,el1,el2,eesij,ees0ij,facvdw,facel,fac1,ecosa,
+ & ecosb,ecosg,ury,urz,vry,vrz,facr,a22der,a23der,a32der,
+ & a33der,eel_loc_ij,cosa4,wij,cosbg1,cosbg2,ees0pij,
+ & ees0pij1,ees0mij,ees0mij1,fac3p,ees0mijp,ees0pijp,
+ & ecosa1,ecosb1,ecosg1,ecosa2,ecosb2,ecosg2,ecosap,ecosbp,
+ & ecosgp,ecosam,ecosbm,ecosgm,ghalf,rlocshield
+ double precision a22,a23,a32,a33,geel_loc_ij,geel_loc_ji
+ double precision dist_init,xj_safe,yj_safe,zj_safe,
+ & xj_temp,yj_temp,zj_temp,dist_temp,xmedi,ymedi,zmedi
+ double precision sscale,sscagrad,scalar
+
c 4/26/02 - AL scaling factor for 1,4 repulsive VDW interactions
#ifdef MOMENT
double precision scal_el /1.0d0/
C zj=zj-zmedi
rij=xj*xj+yj*yj+zj*zj
- sss=sscale(sqrt(rij))
- sssgrad=sscagrad(sqrt(rij))
+ sss=sscale(dsqrt(rij),r_cut_int)
+ if (sss.eq.0.0d0) return
+ sssgrad=sscagrad(dsqrt(rij),r_cut_int)
c if (sss.gt.0.0d0) then
rrmij=1.0D0/rij
rij=dsqrt(rij)
fac_shield(i)=1.0
fac_shield(j)=1.0
eesij=(el1+el2)
- ees=ees+eesij
+ ees=ees+eesij*sss
endif
evdw1=evdw1+evdwij*sss
cd write(iout,'(2(2i3,2x),7(1pd12.4)/2(3(1pd12.4),5x)/)')
cd & xmedi,ymedi,zmedi,xj,yj,zj
if (energy_dec) then
- write (iout,'(a6,2i5,0pf7.3,2i5,3e11.3)')
- &'evdw1',i,j,evdwij
- &,iteli,itelj,aaa,evdw1,sss
- write (iout,'(a6,2i5,0pf7.3,2f8.3)') 'ees',i,j,eesij,
- &fac_shield(i),fac_shield(j)
+ write (iout,'(a6,2i5,0pf7.3,2i5,e11.3,3f10.5)')
+ & 'evdw1',i,j,evdwij,iteli,itelj,aaa,evdw1,sss,rij
+ write (iout,'(a6,2i5,0pf7.3,2f8.3)') 'ees',i,j,eesij,
+ & fac_shield(i),fac_shield(j)
endif
C
*
* Radial derivatives. First process both termini of the fragment (i,j)
*
- ggg(1)=facel*xj
- ggg(2)=facel*yj
- ggg(3)=facel*zj
+ aux=facel*sss+rmij*sssgrad*eesij
+ ggg(1)=aux*xj
+ ggg(2)=aux*yj
+ ggg(3)=aux*zj
if ((fac_shield(i).gt.0).and.(fac_shield(j).gt.0).and.
& (shield_mode.gt.0)) then
C print *,i,j
iresshield=shield_list(ilist,j)
do k=1,3
rlocshield=grad_shield_side(k,ilist,j)*eesij/fac_shield(j)
- & *2.0
+ & *2.0*sss
gshieldx(k,iresshield)=gshieldx(k,iresshield)+
& rlocshield
- & +grad_shield_loc(k,ilist,j)*eesij/fac_shield(j)*2.0
+ & +grad_shield_loc(k,ilist,j)*eesij/fac_shield(j)*2.0*sss
gshieldc(k,iresshield-1)=gshieldc(k,iresshield-1)+rlocshield
C & +grad_shield_loc(k,ilist,j)*eesij/fac_shield(j)
do k=1,3
gshieldc(k,i)=gshieldc(k,i)+
- & grad_shield(k,i)*eesij/fac_shield(i)*2.0
+ & grad_shield(k,i)*eesij/fac_shield(i)*2.0*sss
gshieldc(k,j)=gshieldc(k,j)+
- & grad_shield(k,j)*eesij/fac_shield(j)*2.0
+ & grad_shield(k,j)*eesij/fac_shield(j)*2.0*sss
gshieldc(k,i-1)=gshieldc(k,i-1)+
- & grad_shield(k,i)*eesij/fac_shield(i)*2.0
+ & grad_shield(k,i)*eesij/fac_shield(i)*2.0*sss
gshieldc(k,j-1)=gshieldc(k,j-1)+
- & grad_shield(k,j)*eesij/fac_shield(j)*2.0
+ & grad_shield(k,j)*eesij/fac_shield(j)*2.0*sss
enddo
endif
cgrad gelc(l,k)=gelc(l,k)+ggg(l)
cgrad enddo
cgrad enddo
- if (sss.gt.0.0) then
- ggg(1)=facvdw*xj+sssgrad*rmij*evdwij*xj
- ggg(2)=facvdw*yj+sssgrad*rmij*evdwij*yj
- ggg(3)=facvdw*zj+sssgrad*rmij*evdwij*zj
- else
- ggg(1)=0.0
- ggg(2)=0.0
- ggg(3)=0.0
- endif
+ facvdw=facvdw+sssgrad*rmij*evdwij
+ ggg(1)=facvdw*xj
+ ggg(2)=facvdw*yj
+ ggg(3)=facvdw*zj
c do k=1,3
c ghalf=0.5D0*ggg(k)
c gvdwpp(k,i)=gvdwpp(k,i)+ghalf
cgrad enddo
#else
C MARYSIA
- facvdw=(ev1+evdwij)*sss
+ facvdw=(ev1+evdwij)
facel=(el1+eesij)
fac1=fac
- fac=-3*rrmij*(facvdw+facvdw+facel)
+ fac=-3*rrmij*(facvdw+facvdw+facel)*sss
+ & +(evdwij+eesij)*sssgrad*rrmij
erij(1)=xj*rmij
erij(2)=yj*rmij
erij(3)=zj*rmij
cd & (dcosg(k),k=1,3)
do k=1,3
ggg(k)=(ecosb*dcosb(k)+ecosg*dcosg(k))*
- & fac_shield(i)**2*fac_shield(j)**2
+ & fac_shield(i)**2*fac_shield(j)**2*sss
enddo
c do k=1,3
c ghalf=0.5D0*ggg(k)
do k=1,3
gelc(k,i)=gelc(k,i)
& +((ecosa*(dc_norm(k,j)-cosa*dc_norm(k,i))
- & + ecosb*(erij(k)-cosb*dc_norm(k,i)))*vbld_inv(i+1))
+ & + ecosb*(erij(k)-cosb*dc_norm(k,i)))*vbld_inv(i+1))*sss
& *fac_shield(i)**2*fac_shield(j)**2
gelc(k,j)=gelc(k,j)
& +((ecosa*(dc_norm(k,i)-cosa*dc_norm(k,j))
- & + ecosg*(erij(k)-cosg*dc_norm(k,j)))*vbld_inv(j+1))
+ & + ecosg*(erij(k)-cosg*dc_norm(k,j)))*vbld_inv(j+1))*sss
& *fac_shield(i)**2*fac_shield(j)**2
gelc_long(k,j)=gelc_long(k,j)+ggg(k)
gelc_long(k,i)=gelc_long(k,i)-ggg(k)
C fac_shield(j)=0.6
endif
eel_loc_ij=eel_loc_ij
- & *fac_shield(i)*fac_shield(j)
+ & *fac_shield(i)*fac_shield(j)*sss
c if (energy_dec) write (iout,'(a6,2i5,0pf7.3)')
c & 'eelloc',i,j,eel_loc_ij
C Now derivative over eel_loc
& +a23*gmuij1(2)
& +a32*gmuij1(3)
& +a33*gmuij1(4))
- & *fac_shield(i)*fac_shield(j)
+ & *fac_shield(i)*fac_shield(j)*sss
c write(iout,*) "derivative over thatai"
c write(iout,*) a22*gmuij1(1), a23*gmuij1(2) ,a32*gmuij1(3),
c & a33*gmuij1(4)
& +a33*gmuij2(4)
gloc(nphi+i-1,icg)=gloc(nphi+i-1,icg)+
& geel_loc_ij*wel_loc
- & *fac_shield(i)*fac_shield(j)
+ & *fac_shield(i)*fac_shield(j)*sss
c Derivative over j residue
geel_loc_ji=a22*gmuji1(1)
gloc(nphi+j,icg)=gloc(nphi+j,icg)+
& geel_loc_ji*wel_loc
- & *fac_shield(i)*fac_shield(j)
+ & *fac_shield(i)*fac_shield(j)*sss
geel_loc_ji=
& +a22*gmuji2(1)
c & a33*gmuji2(4)
gloc(nphi+j-1,icg)=gloc(nphi+j-1,icg)+
& geel_loc_ji*wel_loc
- & *fac_shield(i)*fac_shield(j)
+ & *fac_shield(i)*fac_shield(j)*sss
#endif
cd write (iout,*) 'i',i,' j',j,' eel_loc_ij',eel_loc_ij
& gel_loc_loc(i-1)=gel_loc_loc(i-1)+
& (a22*muder(1,i)*mu(1,j)+a23*muder(1,i)*mu(2,j)
& +a32*muder(2,i)*mu(1,j)+a33*muder(2,i)*mu(2,j))
- & *fac_shield(i)*fac_shield(j)
+ & *fac_shield(i)*fac_shield(j)*sss
gel_loc_loc(j-1)=gel_loc_loc(j-1)+
& (a22*mu(1,i)*muder(1,j)+a23*mu(1,i)*muder(2,j)
& +a32*mu(2,i)*muder(1,j)+a33*mu(2,i)*muder(2,j))
- & *fac_shield(i)*fac_shield(j)
+ & *fac_shield(i)*fac_shield(j)*sss
C Derivatives of eello in DC(i+1) thru DC(j-1) or DC(nres-2)
+ aux=eel_loc_ij/sss*sssgrad*rmij
+ ggg(1)=aux*xj
+ ggg(2)=aux*yj
+ ggg(3)=aux*zj
do l=1,3
- ggg(l)=(agg(l,1)*muij(1)+
+ ggg(l)=ggg(l)+(agg(l,1)*muij(1)+
& agg(l,2)*muij(2)+agg(l,3)*muij(3)+agg(l,4)*muij(4))
- & *fac_shield(i)*fac_shield(j)
+ & *fac_shield(i)*fac_shield(j)*sss
gel_loc_long(l,j)=gel_loc_long(l,j)+ggg(l)
gel_loc_long(l,i)=gel_loc_long(l,i)-ggg(l)
cgrad ghalf=0.5d0*ggg(l)
do l=1,3
gel_loc(l,i)=gel_loc(l,i)+(aggi(l,1)*muij(1)+
& aggi(l,2)*muij(2)+aggi(l,3)*muij(3)+aggi(l,4)*muij(4))
- & *fac_shield(i)*fac_shield(j)
+ & *fac_shield(i)*fac_shield(j)*sss
gel_loc(l,i+1)=gel_loc(l,i+1)+(aggi1(l,1)*muij(1)+
& aggi1(l,2)*muij(2)+aggi1(l,3)*muij(3)+aggi1(l,4)*muij(4))
- & *fac_shield(i)*fac_shield(j)
+ & *fac_shield(i)*fac_shield(j)*sss
gel_loc(l,j)=gel_loc(l,j)+(aggj(l,1)*muij(1)+
& aggj(l,2)*muij(2)+aggj(l,3)*muij(3)+aggj(l,4)*muij(4))
- & *fac_shield(i)*fac_shield(j)
+ & *fac_shield(i)*fac_shield(j)*sss
gel_loc(l,j1)=gel_loc(l,j1)+(aggj1(l,1)*muij(1)+
& aggj1(l,2)*muij(2)+aggj1(l,3)*muij(3)+aggj1(l,4)*muij(4))
- & *fac_shield(i)*fac_shield(j)
+ & *fac_shield(i)*fac_shield(j)*sss
enddo
ENDIF
C Change 12/26/95 to calculate four-body contributions to H-bonding energy
c if (j.gt.i+1 .and. num_conti.le.maxconts) then
+#ifdef FOURBODY
if (wcorr+wcorr4+wcorr5+wcorr6.gt.0.0d0
& .and. num_conti.le.maxconts) then
c write (iout,*) i,j," entered corr"
C fac_shield(j)=0.6d0
endif
ees0p(num_conti,i)=0.5D0*fac3*(ees0pij+ees0mij)
- & *fac_shield(i)*fac_shield(j)
+ & *fac_shield(i)*fac_shield(j)*sss
ees0m(num_conti,i)=0.5D0*fac3*(ees0pij-ees0mij)
- & *fac_shield(i)*fac_shield(j)
+ & *fac_shield(i)*fac_shield(j)*sss
C Diagnostics. Comment out or remove after debugging!
c ees0p(num_conti,i)=0.5D0*fac3*ees0pij
c ees0m(num_conti,i)=0.5D0*fac3*ees0mij
gggm(k)=ecosbm*dcosb(k)+ecosgm*dcosg(k)
enddo
gggp(1)=gggp(1)+ees0pijp*xj
+ & +ees0p(num_conti,i)/sss*rmij*xj*sssgrad
gggp(2)=gggp(2)+ees0pijp*yj
+ & +ees0p(num_conti,i)/sss*rmij*yj*sssgrad
gggp(3)=gggp(3)+ees0pijp*zj
+ & +ees0p(num_conti,i)/sss*rmij*zj*sssgrad
gggm(1)=gggm(1)+ees0mijp*xj
+ & +ees0m(num_conti,i)/sss*rmij*xj*sssgrad
gggm(2)=gggm(2)+ees0mijp*yj
+ & +ees0m(num_conti,i)/sss*rmij*yj*sssgrad
gggm(3)=gggm(3)+ees0mijp*zj
+ & +ees0m(num_conti,i)/sss*rmij*zj*sssgrad
C Derivatives due to the contact function
gacont_hbr(1,num_conti,i)=fprimcont*xj
gacont_hbr(2,num_conti,i)=fprimcont*yj
gacontp_hb1(k,num_conti,i)=!ghalfp
& +(ecosap*(dc_norm(k,j)-cosa*dc_norm(k,i))
& + ecosbp*(erij(k)-cosb*dc_norm(k,i)))*vbld_inv(i+1)
- & *fac_shield(i)*fac_shield(j)
+ & *fac_shield(i)*fac_shield(j)*sss
gacontp_hb2(k,num_conti,i)=!ghalfp
& +(ecosap*(dc_norm(k,i)-cosa*dc_norm(k,j))
& + ecosgp*(erij(k)-cosg*dc_norm(k,j)))*vbld_inv(j+1)
- & *fac_shield(i)*fac_shield(j)
+ & *fac_shield(i)*fac_shield(j)*sss
gacontp_hb3(k,num_conti,i)=gggp(k)
- & *fac_shield(i)*fac_shield(j)
+ & *fac_shield(i)*fac_shield(j)*sss
gacontm_hb1(k,num_conti,i)=!ghalfm
& +(ecosam*(dc_norm(k,j)-cosa*dc_norm(k,i))
& + ecosbm*(erij(k)-cosb*dc_norm(k,i)))*vbld_inv(i+1)
- & *fac_shield(i)*fac_shield(j)
+ & *fac_shield(i)*fac_shield(j)*sss
gacontm_hb2(k,num_conti,i)=!ghalfm
& +(ecosam*(dc_norm(k,i)-cosa*dc_norm(k,j))
& + ecosgm*(erij(k)-cosg*dc_norm(k,j)))*vbld_inv(j+1)
- & *fac_shield(i)*fac_shield(j)
+ & *fac_shield(i)*fac_shield(j)*sss
gacontm_hb3(k,num_conti,i)=gggm(k)
- & *fac_shield(i)*fac_shield(j)
+ & *fac_shield(i)*fac_shield(j)*sss
enddo
C Diagnostics. Comment out or remove after debugging!
endif ! num_conti.le.maxconts
endif ! fcont.gt.0
endif ! j.gt.i+1
+#endif
if (wturn3.gt.0.0d0 .or. wturn4.gt.0.0d0) then
do k=1,4
do l=1,3
include 'COMMON.CHAIN'
include 'COMMON.DERIV'
include 'COMMON.INTERACT'
- include 'COMMON.CONTACTS'
+ include 'COMMON.CORRMAT'
include 'COMMON.TORSION'
include 'COMMON.VECTORS'
include 'COMMON.FFIELD'
include 'COMMON.CHAIN'
include 'COMMON.DERIV'
include 'COMMON.INTERACT'
- include 'COMMON.CONTACTS'
+ include 'COMMON.CORRMAT'
include 'COMMON.TORSION'
include 'COMMON.VECTORS'
include 'COMMON.FFIELD'
C peptide-group centers and side chains and its gradient in virtual-bond and
C side-chain vectors.
C
- implicit real*8 (a-h,o-z)
+ implicit none
include 'DIMENSIONS'
include 'COMMON.GEO'
include 'COMMON.VAR'
include 'COMMON.CONTROL'
include 'COMMON.SPLITELE'
integer xshift,yshift,zshift
- dimension ggg(3)
+ double precision ggg(3)
+ integer i,iint,j,k,iteli,itypj,subchap
+ double precision xi,yi,zi,xj,yj,zj,rrij,sss1,sssgrad1,
+ & fac,e1,e2,rij
+ double precision evdw2,evdw2_14,evdwij
+ double precision xj_safe,yj_safe,zj_safe,xj_temp,yj_temp,zj_temp,
+ & dist_temp, dist_init
+ double precision sscale,sscagrad
evdw2=0.0D0
evdw2_14=0.0d0
c print *,boxxsize,boxysize,boxzsize,'wymiary pudla'
C do xshift=-1,1
C do yshift=-1,1
C do zshift=-1,1
- if (energy_dec) write (iout,*) "escp:",r_cut,rlamb
+ if (energy_dec) write (iout,*) "escp:",r_cut_int,rlamb
do i=iatscp_s,iatscp_e
if (itype(i).eq.ntyp1 .or. itype(i+1).eq.ntyp1) cycle
iteli=itel(i)
c print *,xj,yj,zj,'polozenie j'
rrij=1.0D0/(xj*xj+yj*yj+zj*zj)
c print *,rrij
- sss=sscale(1.0d0/(dsqrt(rrij)))
+ sss=sscale(1.0d0/(dsqrt(rrij)),r_cut_int)
c print *,r_cut,1.0d0/dsqrt(rrij),sss,'tu patrz'
c if (sss.eq.0) print *,'czasem jest OK'
if (sss.le.0.0d0) cycle
- sssgrad=sscagrad(1.0d0/(dsqrt(rrij)))
+ sssgrad=sscagrad(1.0d0/(dsqrt(rrij)),r_cut_int)
fac=rrij**expon2
e1=fac*fac*aad(itypj,iteli)
e2=fac*bad(itypj,iteli)
endif
evdwij=e1+e2
evdw2=evdw2+evdwij*sss
- if (energy_dec) write (iout,'(a6,2i5,0pf7.3,2i3,3e11.3)')
- & 'evdw2',i,j,evdwij,iteli,itypj,fac,aad(itypj,iteli),
+ if (energy_dec) write (iout,'(a6,2i5,3f7.3,2i3,3e11.3)')
+ & 'evdw2',i,j,1.0d0/dsqrt(rrij),sss,
+ & evdwij,iteli,itypj,fac,aad(itypj,iteli),
& bad(itypj,iteli)
C
C Calculate contributions to the gradient in the virtual-bond and SC vectors.
estr=0.0d0
estr1=0.0d0
do i=ibondp_start,ibondp_end
+c 3/4/2020 Adam: removed dummy bond graient if Calpha and SC coords are
+c used
+#ifdef FIVEDIAG
+ if (itype(i-1).eq.ntyp1 .or. itype(i).eq.ntyp1) cycle
+ diff = vbld(i)-vbldp0
+#else
if (itype(i-1).eq.ntyp1 .and. itype(i).eq.ntyp1) cycle
c estr1=estr1+gnmr1(vbld(i),-1.0d0,distchainmax)
c do j=1,3
c & "estr1",i,gnmr1(vbld(i),-1.0d0,distchainmax)
c else
C Checking if it involves dummy (NH3+ or COO-) group
- if (itype(i-1).eq.ntyp1 .or. itype(i).eq.ntyp1) then
+ if (itype(i-1).eq.ntyp1 .or. itype(i).eq.ntyp1) then
C YES vbldpDUM is the equlibrium length of spring for Dummy atom
- diff = vbld(i)-vbldpDUM
- if (energy_dec) write(iout,*) "dum_bond",i,diff
- else
-C NO vbldp0 is the equlibrium lenght of spring for peptide group
- diff = vbld(i)-vbldp0
- endif
- if (energy_dec) write (iout,'(a7,i5,4f7.3)')
+ diff = vbld(i)-vbldpDUM
+ if (energy_dec) write(iout,*) "dum_bond",i,diff
+ else
+C NO vbldp0 is the equlibrium length of spring for peptide group
+ diff = vbld(i)-vbldp0
+ endif
+#endif
+ if (energy_dec) write (iout,'(a7,i5,4f7.3)')
& "estr bb",i,vbld(i),vbldp0,diff,AKP*diff*diff
estr=estr+diff*diff
do j=1,3
& phii1*rad2deg,ethetai
c lprn1=.false.
etheta=etheta+ethetai
- if (energy_dec) write (iout,'(a6,i5,0pf7.3)')
- & 'ebend',i,ethetai
if (i.gt.3) gloc(i-3,icg)=gloc(i-3,icg)+wang*dephii
if (i.lt.nres) gloc(i-2,icg)=gloc(i-2,icg)+wang*dephii1
gloc(nphi+i-2,icg)=gloc(nphi+i-2,icg)+wang*dethetai
c & dscp1,dscp2,sumene
c sumene = enesc(x,xx,yy,zz,cost2tab(i+1),sint2tab(i+1))
escloc = escloc + sumene
- if (energy_dec) write (iout,'(a6,i5,0pf7.3)')
- & 'escloc',i,sumene
-c write (2,*) "i",i," escloc",sumene,escloc,it,itype(i)
+ if (energy_dec) write (2,*) "i",i," itype",itype(i)," it",it,
+ & " escloc",sumene,escloc,it,itype(i)
c & ,zz,xx,yy
c#define DEBUG
#ifdef DEBUG
include 'COMMON.IOUNITS'
include 'COMMON.FFIELD'
include 'COMMON.TORCNSTR'
- include 'COMMON.CONTROL'
logical lprn
C Set lprn=.true. for debugging
lprn=.false.
& (itype(i-1).eq.ntyp1).or.(itype(i).eq.ntyp1).or.
& (itype(i+1).eq.ntyp1)) cycle
C In current verion the ALL DUMMY ATOM POTENTIALS ARE OFF
- etors_d_ii=0.0D0
itori=itortyp(itype(i-2))
itori1=itortyp(itype(i-1))
itori2=itortyp(itype(i))
sinphi2=dsin(j*phii1)
etors_d=etors_d+v1cij*cosphi1+v1sij*sinphi1+
& v2cij*cosphi2+v2sij*sinphi2
- if (energy_dec) etors_d_ii=etors_d_ii+
- & v1cij*cosphi1+v1sij*sinphi1+v2cij*cosphi2+v2sij*sinphi2
gloci1=gloci1+j*(v1sij*cosphi1-v1cij*sinphi1)
gloci2=gloci2+j*(v2sij*cosphi2-v2cij*sinphi2)
enddo
sinphi1m2=dsin(l*phii-(k-l)*phii1)
etors_d=etors_d+v1cdij*cosphi1p2+v2cdij*cosphi1m2+
& v1sdij*sinphi1p2+v2sdij*sinphi1m2
- if (energy_dec) etors_d_ii=etors_d_ii+
- & v1cdij*cosphi1p2+v2cdij*cosphi1m2+
- & v1sdij*sinphi1p2+v2sdij*sinphi1m2
gloci1=gloci1+l*(v1sdij*cosphi1p2+v2sdij*cosphi1m2
& -v1cdij*sinphi1p2-v2cdij*sinphi1m2)
gloci2=gloci2+(k-l)*(v1sdij*cosphi1p2-v2sdij*cosphi1m2
& -v1cdij*sinphi1p2+v2cdij*sinphi1m2)
enddo
enddo
- if (energy_dec) write (iout,'(a6,i5,0pf7.3)')
- & 'etor_d',i,etors_d_ii
gloc(i-3,icg)=gloc(i-3,icg)+wtor_d*gloci1
gloc(i-2,icg)=gloc(i-2,icg)+wtor_d*gloci2
enddo
c----------------------------------------------------------------------------
c MODELLER restraint function
subroutine e_modeller(ehomology_constr)
- implicit real*8 (a-h,o-z)
+ implicit none
include 'DIMENSIONS'
- integer nnn, i, j, k, ki, irec, l
+ double precision ehomology_constr
+ integer nnn,i,ii,j,k,ijk,jik,ki,kk,nexl,irec,l
integer katy, odleglosci, test7
real*8 odleg, odleg2, odleg3, kat, kat2, kat3, gdih(max_template)
real*8 Eval,Erot
double precision, dimension (max_template) ::
& gtheta,dscdiff,uscdiffk,guscdiff2,guscdiff3,
& theta_diff
+ double precision sum_godl,sgodl,grad_odl3,ggodl,sum_gdih,
+ & sum_guscdiff,sum_sgdih,sgdih,grad_dih3,usc_diff_i,dxx,dyy,dzz,
+ & betai,sum_sgodl,dij
+ double precision dist,pinorm
c
-
include 'COMMON.SBRIDGE'
include 'COMMON.CHAIN'
include 'COMMON.GEO'
include 'COMMON.INTERACT'
include 'COMMON.VAR'
include 'COMMON.IOUNITS'
- include 'COMMON.MD'
+c include 'COMMON.MD'
include 'COMMON.CONTROL'
+ include 'COMMON.HOMOLOGY'
+ include 'COMMON.QRESTR'
c
c From subroutine Econstr_back
c
esccor=0.0D0
do i=itau_start,itau_end
if ((itype(i-2).eq.ntyp1).or.(itype(i-1).eq.ntyp1)) cycle
+ esccor_ii=0.0D0
isccori=isccortyp(itype(i-2))
isccori1=isccortyp(itype(i-1))
c write (iout,*) "EBACK_SC_COR",i,nterm_sccor(isccori,isccori1)
phii=phi(i)
do intertyp=1,3 !intertyp
- esccor_ii=0.0D0
cc Added 09 May 2012 (Adasko)
cc Intertyp means interaction type of backbone mainchain correlation:
c 1 = SC...Ca...Ca...Ca
v2ij=v2sccor(j,intertyp,isccori,isccori1)
cosphi=dcos(j*tauangle(intertyp,i))
sinphi=dsin(j*tauangle(intertyp,i))
- if (energy_dec) esccor_ii=esccor_ii+v1ij*cosphi+v2ij*sinphi
esccor=esccor+v1ij*cosphi+v2ij*sinphi
gloci=gloci+j*(v2ij*cosphi-v1ij*sinphi)
enddo
- if (energy_dec) write (iout,'(a6,i5,i2,0pf7.3)')
- & 'esccor',i,intertyp,esccor_ii
c write (iout,*) "EBACK_SC_COR",i,v1ij*cosphi+v2ij*sinphi,intertyp
gloc_sc(intertyp,i-3,icg)=gloc_sc(intertyp,i-3,icg)+wsccor*gloci
if (lprn)
return
end
+#ifdef FOURBODY
c----------------------------------------------------------------------------
subroutine multibody(ecorr)
C This subroutine calculates multi-body contributions to energy following
include 'COMMON.DERIV'
include 'COMMON.INTERACT'
include 'COMMON.CONTACTS'
+ include 'COMMON.CONTMAT'
+ include 'COMMON.CORRMAT'
double precision gx(3),gx1(3)
logical lprn
include 'COMMON.DERIV'
include 'COMMON.INTERACT'
include 'COMMON.CONTACTS'
+ include 'COMMON.CONTMAT'
+ include 'COMMON.CORRMAT'
include 'COMMON.SHIELD'
double precision gx(3),gx1(3)
logical lprn
include 'COMMON.DERIV'
include 'COMMON.INTERACT'
include 'COMMON.CONTACTS'
+ include 'COMMON.CONTMAT'
+ include 'COMMON.CORRMAT'
include 'COMMON.CONTROL'
include 'COMMON.LOCAL'
double precision gx(3),gx1(3),time00
parameter (max_cont=maxconts)
parameter (max_dim=26)
include "COMMON.CONTACTS"
+ include 'COMMON.CONTMAT'
+ include 'COMMON.CORRMAT'
double precision zapas(max_dim,maxconts,max_fg_procs),
& zapas_recv(max_dim,maxconts,max_fg_procs)
common /przechowalnia/ zapas
include 'COMMON.LOCAL'
include 'COMMON.INTERACT'
include 'COMMON.CONTACTS'
+ include 'COMMON.CONTMAT'
+ include 'COMMON.CORRMAT'
include 'COMMON.CHAIN'
include 'COMMON.CONTROL'
include 'COMMON.SHIELD'
parameter (max_cont=maxconts)
parameter (max_dim=70)
include "COMMON.CONTACTS"
+ include 'COMMON.CONTMAT'
+ include 'COMMON.CORRMAT'
double precision zapas(max_dim,maxconts,max_fg_procs),
& zapas_recv(max_dim,maxconts,max_fg_procs)
common /przechowalnia/ zapas
include 'COMMON.DERIV'
include 'COMMON.INTERACT'
include 'COMMON.CONTACTS'
+ include 'COMMON.CONTMAT'
+ include 'COMMON.CORRMAT'
include 'COMMON.SHIELD'
include 'COMMON.CONTROL'
double precision gx(3),gx1(3)
include 'COMMON.DERIV'
include 'COMMON.INTERACT'
include 'COMMON.CONTACTS'
+ include 'COMMON.CONTMAT'
+ include 'COMMON.CORRMAT'
include 'COMMON.TORSION'
include 'COMMON.VAR'
include 'COMMON.GEO'
include 'COMMON.DERIV'
include 'COMMON.INTERACT'
include 'COMMON.CONTACTS'
+ include 'COMMON.CONTMAT'
+ include 'COMMON.CORRMAT'
include 'COMMON.TORSION'
include 'COMMON.VAR'
include 'COMMON.GEO'
include 'COMMON.DERIV'
include 'COMMON.INTERACT'
include 'COMMON.CONTACTS'
+ include 'COMMON.CONTMAT'
+ include 'COMMON.CORRMAT'
include 'COMMON.TORSION'
include 'COMMON.VAR'
include 'COMMON.GEO'
include 'COMMON.DERIV'
include 'COMMON.INTERACT'
include 'COMMON.CONTACTS'
+ include 'COMMON.CONTMAT'
+ include 'COMMON.CORRMAT'
include 'COMMON.TORSION'
include 'COMMON.VAR'
include 'COMMON.GEO'
include 'COMMON.DERIV'
include 'COMMON.INTERACT'
include 'COMMON.CONTACTS'
+ include 'COMMON.CONTMAT'
+ include 'COMMON.CORRMAT'
include 'COMMON.TORSION'
include 'COMMON.VAR'
include 'COMMON.GEO'
include 'COMMON.DERIV'
include 'COMMON.INTERACT'
include 'COMMON.CONTACTS'
+ include 'COMMON.CONTMAT'
+ include 'COMMON.CORRMAT'
include 'COMMON.TORSION'
include 'COMMON.VAR'
include 'COMMON.GEO'
include 'COMMON.DERIV'
include 'COMMON.INTERACT'
include 'COMMON.CONTACTS'
+ include 'COMMON.CONTMAT'
+ include 'COMMON.CORRMAT'
include 'COMMON.TORSION'
include 'COMMON.VAR'
include 'COMMON.GEO'
include 'COMMON.DERIV'
include 'COMMON.INTERACT'
include 'COMMON.CONTACTS'
+ include 'COMMON.CONTMAT'
+ include 'COMMON.CORRMAT'
include 'COMMON.TORSION'
include 'COMMON.VAR'
include 'COMMON.GEO'
include 'COMMON.DERIV'
include 'COMMON.INTERACT'
include 'COMMON.CONTACTS'
+ include 'COMMON.CONTMAT'
+ include 'COMMON.CORRMAT'
include 'COMMON.TORSION'
include 'COMMON.VAR'
include 'COMMON.GEO'
include 'COMMON.DERIV'
include 'COMMON.INTERACT'
include 'COMMON.CONTACTS'
+ include 'COMMON.CONTMAT'
+ include 'COMMON.CORRMAT'
include 'COMMON.TORSION'
include 'COMMON.VAR'
include 'COMMON.GEO'
cd write (2,*) 'eel_turn6',ekont*eel_turn6
return
end
-
C-----------------------------------------------------------------------------
+#endif
double precision function scalar(u,v)
!DIR$ INLINEALWAYS scalar
#ifndef OSF
include 'COMMON.INTERACT'
include 'COMMON.VAR'
include 'COMMON.IOUNITS'
- include 'COMMON.MD'
+c include 'COMMON.MD'
+#ifdef LANG0
+#ifdef FIVEDIAG
+ include 'COMMON.LANGEVIN.lang0.5diag'
+#else
+ include 'COMMON.LANGEVIN.lang0'
+#endif
+#else
+ include 'COMMON.LANGEVIN'
+#endif
include 'COMMON.CONTROL'
+ include 'COMMON.SAXS'
include 'COMMON.NAMES'
include 'COMMON.TIME1'
include 'COMMON.FFIELD'
include 'COMMON.INTERACT'
include 'COMMON.VAR'
include 'COMMON.IOUNITS'
- include 'COMMON.MD'
+c include 'COMMON.MD'
+#ifdef LANG0
+#ifdef FIVEDIAG
+ include 'COMMON.LANGEVIN.lang0.5diag'
+#else
+ include 'COMMON.LANGEVIN.lang0'
+#endif
+#else
+ include 'COMMON.LANGEVIN'
+#endif
include 'COMMON.CONTROL'
+ include 'COMMON.SAXS'
include 'COMMON.NAMES'
include 'COMMON.TIME1'
include 'COMMON.FFIELD'