!-----------------------------NUCLEIC GRADIENT
real(kind=8),dimension(:,:),allocatable ::gradb_nucl,gradbx_nucl, &
gvdwpsb1,gelpp,gvdwpsb,gelsbc,gelsbx,gvdwsbx,gvdwsbc,gsbloc,&
- gsblocx
+ gsblocx,gradcorr_nucl,gradxorr_nucl,gradcorr3_nucl,gradxorr3_nucl
! real(kind=8),dimension(:,:),allocatable :: gloc,gloc_x !(maxvar,2)
real(kind=8),dimension(:,:),allocatable :: gel_loc,gel_loc_long,&
gcorr3_turn,gcorr4_turn,gcorr6_turn,gradb,gradbx !(3,maxres)
real(kind=8),dimension(:,:,:,:),allocatable :: uygrad,uzgrad !(3,3,2,maxres)
!-----------------------------------------------------------------------------
! common /przechowalnia/
- real(kind=8),dimension(:,:,:),allocatable :: zapas !(max_dim,maxconts,max_fg_procs)
+ real(kind=8),dimension(:,:,:),allocatable :: zapas
+ real(kind=8),dimension(:,:,:,:),allocatable ::zapas2 !(max_dim,maxconts,max_fg_procs)
real(kind=8),dimension(:,:,:),allocatable :: fromto !(3,3,maxdim)(maxdim=(maxres-1)*(maxres-2)/2)
!-----------------------------------------------------------------------------
!-----------------------------------------------------------------------------
etube=0.0d0
endif
!--------------------------------------------------------
+ print *,"before",ees,evdw1,ecorr
call ebond_nucl(estr_nucl)
call ebend_nucl(ebe_nucl)
call etor_nucl(etors_nucl)
call esb_gb(evdwsb,eelsb)
-! call multibody_hb(ecorr,ecorr3,n_corr,n_corr1)
call epp_nucl_sub(evdwpp,eespp)
call epsb(evdwpsb,eelpsb)
call esb(esbloc)
+ call multibody_hb_nucl(ecorr_nucl,ecorr3_nucl,n_corr,n_corr1)
print *,"after ebend", ebe_nucl
#ifdef TIMING
!
! Loop over all pairs of interacting peptide groups except i,i+2 and i,i+3
!
+ print *,"iatel_s,iatel_e,",iatel_s,iatel_e
do i=iatel_s,iatel_e
if (itype(i,1).eq.ntyp1 .or. itype(i+1,1).eq.ntyp1) cycle
dxi=dc(1,i)
iteli=itel(i)
itelj=itel(j)
if (j.eq.i+2 .and. itelj.eq.2) iteli=2
- aaa=app_nucl(iteli,itelj)
- bbb=bpp_nucl(iteli,itelj)
- ael6i=ael6_nucl(iteli,itelj)
- ael3i=ael3_nucl(iteli,itelj)
+ aaa=app(iteli,itelj)
+ bbb=bpp(iteli,itelj)
+ ael6i=ael6(iteli,itelj)
+ ael3i=ael3(iteli,itelj)
dxj=dc(1,j)
dyj=dc(2,j)
dzj=dc(3,j)
allocate(gvdwsbc(3,-1:nres))
allocate(gsbloc(3,-1:nres))
allocate(gsblocx(3,-1:nres))
+ allocate(gradcorr_nucl(3,-1:nres))
+ allocate(gradxorr_nucl(3,-1:nres))
+ allocate(gradcorr3_nucl(3,-1:nres))
+ allocate(gradxorr3_nucl(3,-1:nres))
!(3,maxres)
allocate(grad_shield_side(3,50,nres))
subroutine esb_gb(evdwsb,eelsb)
use comm_locel
use calc_data_nucl
- integer :: iint,itypi,itypi1,itypj,subchap
+ integer :: iint,itypi,itypi1,itypj,subchap,num_conti2
real(kind=8) :: xi,yi,zi,sig,rij_shift,fac,e1,e2,sigm,epsi
real(kind=8) :: evdw,sig0iji,evdwsb,eelsb,ecorr,eelij
real(kind=8) :: xj_safe,yj_safe,zj_safe,xj_temp,yj_temp,zj_temp,&
! print *,"iastsc_nucl",iatsc_s_nucl,iatsc_e_nucl
do i=iatsc_s_nucl,iatsc_e_nucl
num_conti=0
+ num_conti2=0
itypi=itype(i,2)
! PRINT *,"I=",i,itypi
if (itypi.eq.ntyp1_molec(2)) cycle
gg(3)=zj*fac
!C Calculate angular part of the gradient.
call sc_grad_nucl
- call eelsbij(eelij)
+ call eelsbij(eelij,num_conti2)
if (energy_dec .and. &
(j.eq.i+1.or.j.eq.nres-i+1.or.j.eq.nres-i.or.j.eq.nres-i+2)) &
write (istat,'(e14.5)') evdwij
eelsb=eelsb+eelij
enddo ! j
enddo ! iint
- num_cont_hb(i)=num_conti
+ num_cont_hb(i)=num_conti2
enddo ! i
!c write (iout,*) "Number of loop steps in EGB:",ind
!cccc energy_dec=.false.
return
end subroutine esb_gb
!-------------------------------------------------------------------------------
- subroutine eelsbij(eesij)
+ subroutine eelsbij(eesij,num_conti2)
use comm_locel
use calc_data_nucl
real(kind=8),dimension(3) :: ggg,gggp,gggm,dcosb,dcosg
real(kind=8),dimension(3,3) :: erder,uryg,urzg,vryg,vrzg
real(kind=8) :: xj_safe,yj_safe,zj_safe,xj_temp,yj_temp,zj_temp,&
dist_temp, dist_init,rlocshield,fracinbuf
- integer xshift,yshift,zshift,ilist,iresshield
+ integer xshift,yshift,zshift,ilist,iresshield,num_conti2
!c 4/26/02 - AL scaling factor for 1,4 repulsive VDW interactions
real(kind=8) scal_el /0.5d0/
gelsbc(k,j)=gelsbc(k,j)+ggg(k)
gelsbc(k,i)=gelsbc(k,i)-ggg(k)
enddo
- IF ( (wcorr_nucl.gt.0.0d0.or.wcorr3_nucl.gt.0.0d0) .and. j.gt.i+1 .and.&
+! IF ( (wcorr_nucl.gt.0.0d0.or.wcorr3_nucl.gt.0.0d0) .and.
+ IF ( j.gt.i+1 .and.&
num_conti.le.maxconts) THEN
!C
!C Calculate the contact function. The ith column of the array JCONT will
!c write (2,*) "fcont",fcont
if (fcont.gt.0.0D0) then
num_conti=num_conti+1
+ num_conti2=num_conti2+1
+
if (num_conti.gt.maxconts) then
write (iout,*) 'WARNING - max. # of contacts exceeded;',&
' will skip next contacts for this conf.'
enesc_nucl=sumene
return
end function enesc_nucl
+!-----------------------------------------------------------------------------
+ subroutine multibody_hb_nucl(ecorr,ecorr3,n_corr,n_corr1)
+#ifdef MPI
+ include 'mpif.h'
+ integer,parameter :: max_cont=2000
+ integer,parameter:: max_dim=2*(8*3+6)
+ integer, parameter :: msglen1=max_cont*max_dim
+ integer,parameter :: msglen2=2*msglen1
+ integer source,CorrelType,CorrelID,Error
+ real(kind=8) :: buffer(max_cont,max_dim)
+ integer status(MPI_STATUS_SIZE)
+ integer :: ierror,nbytes
+#endif
+ real(kind=8),dimension(3):: gx(3),gx1(3)
+ real(kind=8) :: time00
+ logical lprn,ldone
+ integer i,j,i1,j1,jj,kk,num_conti,num_conti1,nn
+ real(kind=8) ecorr,ecorr3
+ integer :: n_corr,n_corr1,mm,msglen
+!C Set lprn=.true. for debugging
+ lprn=.false.
+ n_corr=0
+ n_corr1=0
+#ifdef MPI
+ if(.not.allocated(zapas2)) allocate(zapas2(3,maxconts,nres,8))
+
+ if (nfgtasks.le.1) goto 30
+ if (lprn) then
+ write (iout,'(a)') 'Contact function values:'
+ do i=nnt,nct-1
+ write (iout,'(2i3,50(1x,i2,f5.2))') &
+ i,num_cont_hb(i),(jcont_hb(j,i),facont_hb(j,i), &
+ j=1,num_cont_hb(i))
+ enddo
+ endif
+!C Caution! Following code assumes that electrostatic interactions concerning
+!C a given atom are split among at most two processors!
+ CorrelType=477
+ CorrelID=fg_rank+1
+ ldone=.false.
+ do i=1,max_cont
+ do j=1,max_dim
+ buffer(i,j)=0.0D0
+ enddo
+ enddo
+ mm=mod(fg_rank,2)
+!c write (*,*) 'MyRank',MyRank,' mm',mm
+ if (mm) 20,20,10
+ 10 continue
+!c write (*,*) 'Sending: MyRank',MyRank,' mm',mm,' ldone',ldone
+ if (fg_rank.gt.0) then
+!C Send correlation contributions to the preceding processor
+ msglen=msglen1
+ nn=num_cont_hb(iatel_s_nucl)
+ call pack_buffer(max_cont,max_dim,iatel_s,0,buffer)
+!c write (*,*) 'The BUFFER array:'
+!c do i=1,nn
+!c write (*,'(i2,9(3f8.3,2x))') i,(buffer(i,j),j=1,30)
+!c enddo
+ if (ielstart_nucl(iatel_s_nucl).gt.iatel_s_nucl+ispp) then
+ msglen=msglen2
+ call pack_buffer(max_cont,max_dim,iatel_s+1,30,buffer)
+!C Clear the contacts of the atom passed to the neighboring processor
+ nn=num_cont_hb(iatel_s_nucl+1)
+!c do i=1,nn
+!c write (*,'(i2,9(3f8.3,2x))') i,(buffer(i,j+30),j=1,30)
+!c enddo
+ num_cont_hb(iatel_s_nucl)=0
+ endif
+!cd write (iout,*) 'Processor ',fg_rank,MyRank,
+!cd & ' is sending correlation contribution to processor',fg_rank-1,
+!cd & ' msglen=',msglen
+!c write (*,*) 'Processor ',fg_rank,MyRank,
+!c & ' is sending correlation contribution to processor',fg_rank-1,
+!c & ' msglen=',msglen,' CorrelType=',CorrelType
+ time00=MPI_Wtime()
+ call MPI_Send(buffer,msglen,MPI_DOUBLE_PRECISION,fg_rank-1, &
+ CorrelType,FG_COMM,IERROR)
+ time_sendrecv=time_sendrecv+MPI_Wtime()-time00
+!cd write (iout,*) 'Processor ',fg_rank,
+!cd & ' has sent correlation contribution to processor',fg_rank-1,
+!cd & ' msglen=',msglen,' CorrelID=',CorrelID
+!c write (*,*) 'Processor ',fg_rank,
+!c & ' has sent correlation contribution to processor',fg_rank-1,
+!c & ' msglen=',msglen,' CorrelID=',CorrelID
+!c msglen=msglen1
+ endif ! (fg_rank.gt.0)
+ if (ldone) goto 30
+ ldone=.true.
+ 20 continue
+!c write (*,*) 'Receiving: MyRank',MyRank,' mm',mm,' ldone',ldone
+ if (fg_rank.lt.nfgtasks-1) then
+!C Receive correlation contributions from the next processor
+ msglen=msglen1
+ if (ielend_nucl(iatel_e_nucl).lt.nct_molec(2)-1) msglen=msglen2
+!cd write (iout,*) 'Processor',fg_rank,
+!cd & ' is receiving correlation contribution from processor',fg_rank+1,
+!cd & ' msglen=',msglen,' CorrelType=',CorrelType
+!c write (*,*) 'Processor',fg_rank,
+!c &' is receiving correlation contribution from processor',fg_rank+1,
+!c & ' msglen=',msglen,' CorrelType=',CorrelType
+ time00=MPI_Wtime()
+ nbytes=-1
+ do while (nbytes.le.0)
+ call MPI_Probe(fg_rank+1,CorrelType,FG_COMM,status,IERROR)
+ call MPI_Get_count(status,MPI_DOUBLE_PRECISION,nbytes,IERROR)
+ enddo
+!c print *,'Processor',myrank,' msglen',msglen,' nbytes',nbytes
+ call MPI_Recv(buffer,nbytes,MPI_DOUBLE_PRECISION, &
+ fg_rank+1,CorrelType,FG_COMM,status,IERROR)
+ time_sendrecv=time_sendrecv+MPI_Wtime()-time00
+!c write (*,*) 'Processor',fg_rank,
+!c &' has received correlation contribution from processor',fg_rank+1,
+!c & ' msglen=',msglen,' nbytes=',nbytes
+!c write (*,*) 'The received BUFFER array:'
+!c do i=1,max_cont
+!c write (*,'(i2,9(3f8.3,2x))') i,(buffer(i,j),j=1,60)
+!c enddo
+ if (msglen.eq.msglen1) then
+ call unpack_buffer(max_cont,max_dim,iatel_e_nucl+1,0,buffer)
+ else if (msglen.eq.msglen2) then
+ call unpack_buffer(max_cont,max_dim,iatel_e_nucl,0,buffer)
+ call unpack_buffer(max_cont,max_dim,iatel_e_nucl+1,30,buffer)
+ else
+ write (iout,*) &
+ 'ERROR!!!! message length changed while processing correlations.'
+ write (*,*) &
+ 'ERROR!!!! message length changed while processing correlations.'
+ call MPI_Abort(MPI_COMM_WORLD,Error,IERROR)
+ endif ! msglen.eq.msglen1
+ endif ! fg_rank.lt.nfgtasks-1
+ if (ldone) goto 30
+ ldone=.true.
+ goto 10
+ 30 continue
+#endif
+ if (lprn) then
+ write (iout,'(a)') 'Contact function values:'
+ do i=nnt_molec(2),nct_molec(2)-1
+ write (iout,'(2i3,50(1x,i2,f5.2))') &
+ i,num_cont_hb(i),(jcont_hb(j,i),facont_hb(j,i), &
+ j=1,num_cont_hb(i))
+ enddo
+ endif
+ ecorr=0.0D0
+ ecorr3=0.0d0
+!C Remove the loop below after debugging !!!
+ do i=nnt_molec(2),nct_molec(2)
+ do j=1,3
+ gradcorr_nucl(j,i)=0.0D0
+ gradxorr_nucl(j,i)=0.0D0
+ gradcorr3_nucl(j,i)=0.0D0
+ gradxorr3_nucl(j,i)=0.0D0
+ enddo
+ enddo
+ print *,"iatsc_s_nucl,iatsc_e_nucl",iatsc_s_nucl,iatsc_e_nucl
+!C Calculate the local-electrostatic correlation terms
+ do i=iatsc_s_nucl,iatsc_e_nucl
+ i1=i+1
+ num_conti=num_cont_hb(i)
+ num_conti1=num_cont_hb(i+1)
+ print *,i,num_conti,num_conti1
+ do jj=1,num_conti
+ j=jcont_hb(jj,i)
+ do kk=1,num_conti1
+ j1=jcont_hb(kk,i1)
+!c write (iout,*) 'i=',i,' j=',j,' i1=',i1,' j1=',j1,
+!c & ' jj=',jj,' kk=',kk
+ if (j1.eq.j+1 .or. j1.eq.j-1) then
+!C
+!C Contacts I-J and (I+1)-(J+1) or (I+1)-(J-1) occur simultaneously.
+!C The system gains extra energy.
+!C Tentative expression & coefficients; assumed d(stacking)=4.5 A,
+!C parallel dipoles of stacknig bases and sin(mui)sin(muj)/eps/d^3=0.7
+!C Need to implement full formulas 34 and 35 from Liwo et al., 1998.
+!C
+ ecorr=ecorr+ehbcorr_nucl(i,j,i+1,j1,jj,kk,0.528D0,0.132D0)
+ if (energy_dec) write (iout,'(a6,2i5,0pf7.3)') &
+ 'ecorrh',i,j,ehbcorr_nucl(i,j,i+1,j1,jj,kk,0.528D0,0.132D0)
+ n_corr=n_corr+1
+ else if (j1.eq.j) then
+!C
+!C Contacts I-J and I-(J+1) occur simultaneously.
+!C The system loses extra energy.
+!C Tentative expression & c?oefficients; assumed d(stacking)=4.5 A,
+!C parallel dipoles of stacknig bases and sin(mui)sin(muj)/eps/d^3=0.7
+!C Need to implement full formulas 32 from Liwo et al., 1998.
+!C
+!c write (iout,*) 'ecorr3: i=',i,' j=',j,' i1=',i1,' j1=',j1,
+!c & ' jj=',jj,' kk=',kk
+ ecorr3=ecorr3+ehbcorr3_nucl(i,j,i+1,j,jj,kk,0.310D0,-0.155D0)
+ endif
+ enddo ! kk
+ do kk=1,num_conti
+ j1=jcont_hb(kk,i)
+!c write (iout,*) 'ecorr3: i=',i,' j=',j,' i1=',i1,' j1=',j1,
+!c & ' jj=',jj,' kk=',kk
+ if (j1.eq.j+1) then
+!C Contacts I-J and (I+1)-J occur simultaneously.
+!C The system loses extra energy.
+ ecorr3=ecorr3+ehbcorr3_nucl(i,j,i,j+1,jj,kk,0.310D0,-0.155D0)
+ endif ! j1==j+1
+ enddo ! kk
+ enddo ! jj
+ enddo ! i
+ return
+ end subroutine multibody_hb_nucl
+!-----------------------------------------------------------
+ real(kind=8) function ehbcorr_nucl(i,j,k,l,jj,kk,coeffp,coeffm)
+! implicit real*8 (a-h,o-z)
+! include 'DIMENSIONS'
+! include 'COMMON.IOUNITS'
+! include 'COMMON.DERIV'
+! include 'COMMON.INTERACT'
+! include 'COMMON.CONTACTS'
+ real(kind=8),dimension(3) :: gx,gx1
+ logical :: lprn
+!el local variables
+ integer :: i,j,k,l,jj,kk,ll,ilist,m, iresshield
+ real(kind=8) :: coeffp,coeffm,eij,ekl,ees0pij,ees0pkl,ees0mij,&
+ ees0mkl,ees,coeffpees0pij,coeffmees0mij,&
+ coeffpees0pkl,coeffmees0mkl,gradlongij,gradlongkl, &
+ rlocshield
+
+ lprn=.false.
+ eij=facont_hb(jj,i)
+ ekl=facont_hb(kk,k)
+ ees0pij=ees0p(jj,i)
+ ees0pkl=ees0p(kk,k)
+ ees0mij=ees0m(jj,i)
+ ees0mkl=ees0m(kk,k)
+ ekont=eij*ekl
+ ees=-(coeffp*ees0pij*ees0pkl+coeffm*ees0mij*ees0mkl)
+! print *,"ehbcorr_nucl",ekont,ees
+!cd ees=-(coeffp*ees0pkl+coeffm*ees0mkl)
+!C Following 4 lines for diagnostics.
+!cd ees0pkl=0.0D0
+!cd ees0pij=1.0D0
+!cd ees0mkl=0.0D0
+!cd ees0mij=1.0D0
+!cd write (iout,*)'Contacts have occurred for nucleic bases',
+!cd & i,j,' fcont:',eij,' eij',' eesij',ees0pij,ees0mij,' and ',k,l
+!cd & ,' fcont ',ekl,' eeskl',ees0pkl,ees0mkl,' ees=',ees
+!C Calculate the multi-body contribution to energy.
+! ecorr_nucl=ecorr_nucl+ekont*ees
+!C Calculate multi-body contributions to the gradient.
+ coeffpees0pij=coeffp*ees0pij
+ coeffmees0mij=coeffm*ees0mij
+ coeffpees0pkl=coeffp*ees0pkl
+ coeffmees0mkl=coeffm*ees0mkl
+ do ll=1,3
+ gradxorr_nucl(ll,i)=gradxorr_nucl(ll,i) &
+ -ekont*(coeffpees0pkl*gacontp_hb1(ll,jj,i)+&
+ coeffmees0mkl*gacontm_hb1(ll,jj,i))
+ gradxorr_nucl(ll,j)=gradxorr_nucl(ll,j) &
+ -ekont*(coeffpees0pkl*gacontp_hb2(ll,jj,i)+&
+ coeffmees0mkl*gacontm_hb2(ll,jj,i))
+ gradxorr_nucl(ll,k)=gradxorr_nucl(ll,k) &
+ -ekont*(coeffpees0pij*gacontp_hb1(ll,kk,k)+&
+ coeffmees0mij*gacontm_hb1(ll,kk,k))
+ gradxorr_nucl(ll,l)=gradxorr_nucl(ll,l) &
+ -ekont*(coeffpees0pij*gacontp_hb2(ll,kk,k)+ &
+ coeffmees0mij*gacontm_hb2(ll,kk,k))
+ gradlongij=ees*ekl*gacont_hbr(ll,jj,i)- &
+ ekont*(coeffpees0pkl*gacontp_hb3(ll,jj,i)+ &
+ coeffmees0mkl*gacontm_hb3(ll,jj,i))
+ gradcorr_nucl(ll,j)=gradcorr_nucl(ll,j)+gradlongij
+ gradcorr_nucl(ll,i)=gradcorr_nucl(ll,i)-gradlongij
+ gradlongkl=ees*eij*gacont_hbr(ll,kk,k)- &
+ ekont*(coeffpees0pij*gacontp_hb3(ll,kk,k)+ &
+ coeffmees0mij*gacontm_hb3(ll,kk,k))
+ gradcorr_nucl(ll,l)=gradcorr_nucl(ll,l)+gradlongkl
+ gradcorr_nucl(ll,k)=gradcorr_nucl(ll,k)-gradlongkl
+ gradxorr_nucl(ll,i)=gradxorr_nucl(ll,i)-gradlongij
+ gradxorr_nucl(ll,j)=gradxorr_nucl(ll,j)+gradlongij
+ gradxorr_nucl(ll,k)=gradxorr_nucl(ll,k)-gradlongkl
+ gradxorr_nucl(ll,l)=gradxorr_nucl(ll,l)+gradlongkl
+ enddo
+ ehbcorr_nucl=ekont*ees
+ return
+ end function ehbcorr_nucl
+!-------------------------------------------------------------------------
+
+ real(kind=8) function ehbcorr3_nucl(i,j,k,l,jj,kk,coeffp,coeffm)
+! implicit real*8 (a-h,o-z)
+! include 'DIMENSIONS'
+! include 'COMMON.IOUNITS'
+! include 'COMMON.DERIV'
+! include 'COMMON.INTERACT'
+! include 'COMMON.CONTACTS'
+ real(kind=8),dimension(3) :: gx,gx1
+ logical :: lprn
+!el local variables
+ integer :: i,j,k,l,jj,kk,ll,ilist,m, iresshield
+ real(kind=8) :: coeffp,coeffm,eij,ekl,ees0pij,ees0pkl,ees0mij,&
+ ees0mkl,ees,coeffpees0pij,coeffmees0mij,&
+ coeffpees0pkl,coeffmees0mkl,gradlongij,gradlongkl, &
+ rlocshield
+
+ lprn=.false.
+ eij=facont_hb(jj,i)
+ ekl=facont_hb(kk,k)
+ ees0pij=ees0p(jj,i)
+ ees0pkl=ees0p(kk,k)
+ ees0mij=ees0m(jj,i)
+ ees0mkl=ees0m(kk,k)
+ ekont=eij*ekl
+ ees=-(coeffp*ees0pij*ees0pkl+coeffm*ees0mij*ees0mkl)
+!cd ees=-(coeffp*ees0pkl+coeffm*ees0mkl)
+!C Following 4 lines for diagnostics.
+!cd ees0pkl=0.0D0
+!cd ees0pij=1.0D0
+!cd ees0mkl=0.0D0
+!cd ees0mij=1.0D0
+!cd write (iout,*)'Contacts have occurred for nucleic bases',
+!cd & i,j,' fcont:',eij,' eij',' eesij',ees0pij,ees0mij,' and ',k,l
+!cd & ,' fcont ',ekl,' eeskl',ees0pkl,ees0mkl,' ees=',ees
+!C Calculate the multi-body contribution to energy.
+! ecorr=ecorr+ekont*ees
+!C Calculate multi-body contributions to the gradient.
+ coeffpees0pij=coeffp*ees0pij
+ coeffmees0mij=coeffm*ees0mij
+ coeffpees0pkl=coeffp*ees0pkl
+ coeffmees0mkl=coeffm*ees0mkl
+ do ll=1,3
+ gradxorr3_nucl(ll,i)=gradxorr3_nucl(ll,i) &
+ -ekont*(coeffpees0pkl*gacontp_hb1(ll,jj,i)+&
+ coeffmees0mkl*gacontm_hb1(ll,jj,i))
+ gradxorr3_nucl(ll,j)=gradxorr3_nucl(ll,j) &
+ -ekont*(coeffpees0pkl*gacontp_hb2(ll,jj,i)+ &
+ coeffmees0mkl*gacontm_hb2(ll,jj,i))
+ gradxorr3_nucl(ll,k)=gradxorr3_nucl(ll,k) &
+ -ekont*(coeffpees0pij*gacontp_hb1(ll,kk,k)+ &
+ coeffmees0mij*gacontm_hb1(ll,kk,k))
+ gradxorr3_nucl(ll,l)=gradxorr3_nucl(ll,l) &
+ -ekont*(coeffpees0pij*gacontp_hb2(ll,kk,k)+ &
+ coeffmees0mij*gacontm_hb2(ll,kk,k))
+ gradlongij=ees*ekl*gacont_hbr(ll,jj,i)- &
+ ekont*(coeffpees0pkl*gacontp_hb3(ll,jj,i)+ &
+ coeffmees0mkl*gacontm_hb3(ll,jj,i))
+ gradcorr3_nucl(ll,j)=gradcorr3_nucl(ll,j)+gradlongij
+ gradcorr3_nucl(ll,i)=gradcorr3_nucl(ll,i)-gradlongij
+ gradlongkl=ees*eij*gacont_hbr(ll,kk,k)- &
+ ekont*(coeffpees0pij*gacontp_hb3(ll,kk,k)+ &
+ coeffmees0mij*gacontm_hb3(ll,kk,k))
+ gradcorr3_nucl(ll,l)=gradcorr3_nucl(ll,l)+gradlongkl
+ gradcorr3_nucl(ll,k)=gradcorr3_nucl(ll,k)-gradlongkl
+ gradxorr3_nucl(ll,i)=gradxorr3_nucl(ll,i)-gradlongij
+ gradxorr3_nucl(ll,j)=gradxorr3_nucl(ll,j)+gradlongij
+ gradxorr3_nucl(ll,k)=gradxorr3_nucl(ll,k)-gradlongkl
+ gradxorr3_nucl(ll,l)=gradxorr3_nucl(ll,l)+gradlongkl
+ enddo
+ ehbcorr3_nucl=ekont*ees
+ return
+ end function ehbcorr3_nucl
+#ifdef MPI
+ subroutine pack_buffer(dimen1,dimen2,atom,indx,buffer)
+ integer dimen1,dimen2,atom,indx,numcont,i,ii,k,j,num_kont,num_kont_old
+ real(kind=8):: buffer(dimen1,dimen2)
+ num_kont=num_cont_hb(atom)
+ do i=1,num_kont
+ do k=1,8
+ do j=1,3
+ buffer(i,indx+(k-1)*3+j)=zapas2(j,i,atom,k)
+ enddo ! j
+ enddo ! k
+ buffer(i,indx+25)=facont_hb(i,atom)
+ buffer(i,indx+26)=ees0p(i,atom)
+ buffer(i,indx+27)=ees0m(i,atom)
+ buffer(i,indx+28)=d_cont(i,atom)
+ buffer(i,indx+29)=dfloat(jcont_hb(i,atom))
+ enddo ! i
+ buffer(1,indx+30)=dfloat(num_kont)
+ return
+ end subroutine pack_buffer
+!c------------------------------------------------------------------------------
+ subroutine unpack_buffer(dimen1,dimen2,atom,indx,buffer)
+ integer dimen1,dimen2,atom,indx,numcont,i,ii,k,j,num_kont,num_kont_old
+ real(kind=8):: buffer(dimen1,dimen2)
+! double precision zapas
+! common /contacts_hb/ zapas(3,maxconts,maxres,8),
+! & facont_hb(maxconts,maxres),ees0p(maxconts,maxres),
+! & ees0m(maxconts,maxres),d_cont(maxconts,maxres),
+! & num_cont_hb(maxres),jcont_hb(maxconts,maxres)
+ num_kont=buffer(1,indx+30)
+ num_kont_old=num_cont_hb(atom)
+ num_cont_hb(atom)=num_kont+num_kont_old
+ do i=1,num_kont
+ ii=i+num_kont_old
+ do k=1,8
+ do j=1,3
+ zapas2(j,ii,atom,k)=buffer(i,indx+(k-1)*3+j)
+ enddo ! j
+ enddo ! k
+ facont_hb(ii,atom)=buffer(i,indx+25)
+ ees0p(ii,atom)=buffer(i,indx+26)
+ ees0m(ii,atom)=buffer(i,indx+27)
+ d_cont(i,atom)=buffer(i,indx+28)
+ jcont_hb(ii,atom)=buffer(i,indx+29)
+ enddo ! i
+ return
+ end subroutine unpack_buffer
+!c------------------------------------------------------------------------------
+#endif
!----------------------------------------------------------------------------
!-----------------------------------------------------------------------------
projects / unres4.git / commitdiff