integer ishield_listbuf(-1:nres), &
shield_listbuf(maxcontsshi,-1:nres),k,j,i,iii,impishi,mojint,jjj
! print *,"I START ENERGY"
- imatupdate=1
+ imatupdate=100
! if (mod(itime_mat,imatupdate).eq.0) call make_SCSC_inter_list
! real(kind=8), dimension(:),allocatable:: fac_shieldbuf
! real(kind=8), dimension(:,:,:),allocatable:: &
time_Bcastw=time_Bcastw+MPI_Wtime()-time00
! call chainbuild_cart
endif
-! if (mod(itime_mat,imatupdate).eq.0) call make_SCSC_inter_list
+! print *,"itime_mat",itime_mat,imatupdate
+ if (nfgtasks.gt.1) then
+ call MPI_Bcast(itime_mat,1,MPI_INT,king,FG_COMM,IERROR)
+ endif
+ if (mod(itime_mat,imatupdate).eq.0) call make_SCp_inter_list
+ if (mod(itime_mat,imatupdate).eq.0) call make_SCSC_inter_list
+ if (mod(itime_mat,imatupdate).eq.0) call make_pp_inter_list
! print *,'Processor',myrank,' calling etotal ipot=',ipot
! print *,'Processor',myrank,' nnt=',nnt,' nct=',nct
etors_nucl=0.0d0
estr_nucl=0.0d0
ecorr3_nucl=0.0d0
+ ecorr_nucl=0.0d0
ebe_nucl=0.0d0
evdwsb=0.0d0
eelsb=0.0d0
eelpsb=0.0d0
evdwpp=0.0d0
eespp=0.0d0
+ etors_d_nucl=0.0d0
endif
! write(iout,*) ecorr_nucl,"ecorr_nucl",nres_molec(2)
- print *,"before ecatcat",wcatcat
+! print *,"before ecatcat",wcatcat
+ if (nres_molec(5).gt.0) then
if (nfgtasks.gt.1) then
if (fg_rank.eq.0) then
call ecatcat(ecationcation)
else
call ecats_prot_amber(ecation_prot)
endif
- if (nres_molec(2).gt.0) then
+ else
+ ecationcation=0.0d0
+ ecation_prot=0.0d0
+ endif
+ if ((nres_molec(2).gt.0).and.(nres_molec(1).gt.0)) then
call eprot_sc_base(escbase)
call epep_sc_base(epepbase)
call eprot_sc_phosphate(escpho)
! include 'COMMON.SBRIDGE'
logical :: lprn
!el local variables
- integer :: iint,itypi,itypi1,itypj,subchap
+ integer :: iint,itypi,itypi1,itypj,subchap,icont
real(kind=8) :: rrij,xi,yi,zi,sig,rij_shift,fac,e1,e2,sigm,epsi
real(kind=8) :: evdw,sig0ij
real(kind=8) :: xj_safe,yj_safe,zj_safe,xj_temp,yj_temp,zj_temp,&
dPOLdOM1=0.0d0
- do i=iatsc_s,iatsc_e
+ do icont=g_listscsc_start,g_listscsc_end
+ i=newcontlisti(icont)
+ j=newcontlistj(icont)
+
+! do i=iatsc_s,iatsc_e
!C print *,"I am in EVDW",i
itypi=iabs(itype(i,1))
! if (i.ne.47) cycle
!
! Calculate SC interaction energy.
!
- do iint=1,nint_gr(i)
- do j=istart(i,iint),iend(i,iint)
+! do iint=1,nint_gr(i)
+! do j=istart(i,iint),iend(i,iint)
IF (dyn_ss_mask(i).and.dyn_ss_mask(j)) THEN
call dyn_ssbond_ene(i,j,evdwij)
evdw=evdw+evdwij
! Calculate angular part of the gradient.
call sc_grad
ENDIF ! dyn_ss
- enddo ! j
- enddo ! iint
+! enddo ! j
+! enddo ! iint
enddo ! i
! print *,"ZALAMKA", evdw
! write (iout,*) "Number of loop steps in EGB:",ind
#endif
#else
if (i.gt. nnt+2 .and. i.lt.nct+2) then
-! write(iout,*) "i,",molnum(i)
+! write(iout,*) "i,",molnum(i),nloctyp
! print *, "i,",molnum(i),i,itype(i-2,1)
if (molnum(i).eq.1) then
+ if (itype(i-2,1).eq.ntyp1) then
+ iti=nloctyp
+ else
iti = itype2loc(itype(i-2,1))
+ endif
else
iti=nloctyp
endif
0.0d0,1.0d0,0.0d0,&
0.0d0,0.0d0,1.0d0/),shape(unmat))
!el local variables
- integer :: i,k,j
+ integer :: i,k,j,icont
real(kind=8) :: ees,evdw1,eel_loc,eello_turn3,eello_turn4
real(kind=8) :: fac,t_eelecij,fracinbuf
! 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
+! do i=iatel_s,iatel_e
+! JPRDLC
+ do icont=g_listpp_start,g_listpp_end
+ i=newcontlistppi(icont)
+ j=newcontlistppj(icont)
if (itype(i,1).eq.ntyp1 .or. itype(i+1,1).eq.ntyp1) cycle
dxi=dc(1,i)
dyi=dc(2,i)
! write (iout,*) 'i',i,' ielstart',ielstart(i),' ielend',ielend(i)
num_conti=num_cont_hb(i)
- do j=ielstart(i),ielend(i)
+! do j=ielstart(i),ielend(i)
! write (iout,*) i,j,itype(i,1),itype(j,1)
if (itype(j,1).eq.ntyp1.or. itype(j+1,1).eq.ntyp1) cycle
call eelecij(i,j,ees,evdw1,eel_loc)
- enddo ! j
+! enddo ! j
num_cont_hb(i)=num_conti
enddo ! i
! write (iout,*) "Number of loop steps in EELEC:",ind
! sss_ele_grad=0.0d0
! print *,sss_ele_cut,sss_ele_grad,&
! (rij),r_cut_ele,rlamb_ele
-! if (sss_ele_cut.le.0.0) go to 128
+ if (sss_ele_cut.le.0.0) go to 128
rmij=1.0D0/rij
r3ij=rrmij*rmij
!grad enddo
!grad enddo
! 9/28/08 AL Gradient compotents will be summed only at the end
- ggg(1)=facvdw*xj &
+ ggg(1)=facvdw*xj+sss_ele_grad*rmij*evdwij*xj &
*((sslipi+sslipj)/2.0d0*lipscale**2+1.0d0)
- ggg(2)=facvdw*yj &
+ ggg(2)=facvdw*yj+sss_ele_grad*rmij*evdwij*yj &
*((sslipi+sslipj)/2.0d0*lipscale**2+1.0d0)
- ggg(3)=facvdw*zj &
+ ggg(3)=facvdw*zj+sss_ele_grad*rmij*evdwij*zj &
*((sslipi+sslipj)/2.0d0*lipscale**2+1.0d0)
do k=1,3
+a32*gmuij1(3)&
+a33*gmuij1(4))&
*fac_shield(i)*fac_shield(j)&
- *sss_ele_cut
+ *sss_ele_cut &
+ *((sslipi+sslipj)/2.0d0*lipscale+1.0d0)
+
!c write(iout,*) "derivative over thatai"
!c write(iout,*) a22*gmuij1(1), a23*gmuij1(2) ,a32*gmuij1(3),
gloc(nphi+i-1,icg)=gloc(nphi+i-1,icg)+&
geel_loc_ij*wel_loc&
*fac_shield(i)*fac_shield(j)&
- *sss_ele_cut
+ *sss_ele_cut &
+ *((sslipi+sslipj)/2.0d0*lipscale+1.0d0)
!c Derivative over j residue
gloc(nphi+j,icg)=gloc(nphi+j,icg)+&
geel_loc_ji*wel_loc&
*fac_shield(i)*fac_shield(j)&
- *sss_ele_cut
+ *sss_ele_cut &
+ *((sslipi+sslipj)/2.0d0*lipscale+1.0d0)
geel_loc_ji=&
gloc(nphi+j-1,icg)=gloc(nphi+j-1,icg)+&
geel_loc_ji*wel_loc&
*fac_shield(i)*fac_shield(j)&
- *sss_ele_cut
+ *sss_ele_cut &
+ *((sslipi+sslipj)/2.0d0*lipscale+1.0d0)
+
#endif
! write (iout,*) 'i',i,' j',j,' eel_loc_ij',eel_loc_ij
ees0p(num_conti,i)=0.5D0*fac3*(ees0pij+ees0mij) &
*sss_ele_cut &
*fac_shield(i)*fac_shield(j)
+! *((sslipi+sslipj)/2.0d0*lipscale+1.0d0)
ees0m(num_conti,i)=0.5D0*fac3*(ees0pij-ees0mij) &
*sss_ele_cut &
*fac_shield(i)*fac_shield(j)
+! *((sslipi+sslipj)/2.0d0*lipscale+1.0d0)
! Diagnostics. Comment out or remove after debugging!
! ees0p(num_conti,i)=0.5D0*fac3*ees0pij
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) &
- *sss_ele_cut*fac_shield(i)*fac_shield(j)
+ *sss_ele_cut*fac_shield(i)*fac_shield(j) ! &
+! *((sslipi+sslipj)/2.0d0*lipscale+1.0d0)
+
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)&
- *sss_ele_cut*fac_shield(i)*fac_shield(j)
+ *sss_ele_cut*fac_shield(i)*fac_shield(j)! &
+! *((sslipi+sslipj)/2.0d0*lipscale+1.0d0)
+
gacontp_hb3(k,num_conti,i)=gggp(k) &
*sss_ele_cut*fac_shield(i)*fac_shield(j)
+! *((sslipi+sslipj)/2.0d0*lipscale+1.0d0)
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) &
*sss_ele_cut*fac_shield(i)*fac_shield(j)
+! *((sslipi+sslipj)/2.0d0*lipscale+1.0d0)
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) &
*sss_ele_cut*fac_shield(i)*fac_shield(j)
+! *((sslipi+sslipj)/2.0d0*lipscale+1.0d0)
gacontm_hb3(k,num_conti,i)=gggm(k) &
*sss_ele_cut*fac_shield(i)*fac_shield(j)
+! *((sslipi+sslipj)/2.0d0*lipscale+1.0d0)
enddo
! Diagnostics. Comment out or remove after debugging!
!C Derivatives in theta
gloc(nphi+i,icg)=gloc(nphi+i,icg) &
+0.5d0*(gpizda1(1,1)+gpizda1(2,2))*wturn3&
- *fac_shield(i)*fac_shield(j)
+ *fac_shield(i)*fac_shield(j) &
+ *((sslipi+sslipj)/2.0d0*lipscale+1.0d0)
+
gloc(nphi+i+1,icg)=gloc(nphi+i+1,icg)&
+0.5d0*(gpizda2(1,1)+gpizda2(2,2))*wturn3&
- *fac_shield(i)*fac_shield(j)
+ *fac_shield(i)*fac_shield(j) &
+ *((sslipi+sslipj)/2.0d0*lipscale+1.0d0)
+
+
!C#endif
! include 'COMMON.CONTROL'
real(kind=8),dimension(3) :: ggg
!el local variables
- integer :: i,iint,j,k,iteli,itypj,subchap
+ integer :: i,iint,j,k,iteli,itypj,subchap,icont
real(kind=8) :: evdw2,evdw2_14,xi,yi,zi,xj,yj,zj,rrij,fac,&
e1,e2,evdwij,rij
real(kind=8) :: xj_safe,yj_safe,zj_safe,xj_temp,yj_temp,zj_temp,&
evdw2_14=0.0d0
!d print '(a)','Enter ESCP'
!d write (iout,*) 'iatscp_s=',iatscp_s,' iatscp_e=',iatscp_e
- do i=iatscp_s,iatscp_e
+! do i=iatscp_s,iatscp_e
+ do icont=g_listscp_start,g_listscp_end
+ i=newcontlistscpi(icont)
+ j=newcontlistscpj(icont)
if (itype(i,1).eq.ntyp1 .or. itype(i+1,1).eq.ntyp1) cycle
iteli=itel(i)
xi=0.5D0*(c(1,i)+c(1,i+1))
zi=mod(zi,boxzsize)
if (zi.lt.0) zi=zi+boxzsize
- do iint=1,nscp_gr(i)
+! do iint=1,nscp_gr(i)
- do j=iscpstart(i,iint),iscpend(i,iint)
+! do j=iscpstart(i,iint),iscpend(i,iint)
itypj=iabs(itype(j,1))
if (itypj.eq.ntyp1) cycle
! Uncomment following three lines for SC-p interactions
gvdwc_scpp(k,i)=gvdwc_scpp(k,i)-ggg(k)
gvdwc_scp(k,j)=gvdwc_scp(k,j)+ggg(k)
enddo
- enddo
+! enddo
- enddo ! iint
+! enddo ! iint
enddo ! i
do i=1,nct
do j=1,3
! & dscp1,dscp2,sumene
! sumene = enesc(x,xx,yy,zz,cost2tab(i+1),sint2tab(i+1))
escloc = escloc + sumene
+ if (energy_dec) write (2,*) "i",i," itype",itype(i,1)," it",it, &
+ " escloc",sumene,escloc,it,itype(i,1)
! write (2,*) "i",i," escloc",sumene,escloc,it,itype(i,1)
! & ,zz,xx,yy
!#define DEBUG
allocate(uygrad(3,3,2,nres))
allocate(uzgrad(3,3,2,nres))
!(3,3,2,maxres)
+! allocateion of lists JPRDLA
+ allocate(newcontlistppi(200*nres))
+ allocate(newcontlistscpi(200*nres))
+ allocate(newcontlisti(200*nres))
+ allocate(newcontlistppj(200*nres))
+ allocate(newcontlistscpj(200*nres))
+ allocate(newcontlistj(200*nres))
return
end subroutine alloc_ener_arrays
! r06 = rcat0**6
! r012 = r06**2
! k0 = 332*(2*2)/80
- Evan1cat=epscalc*(r012/rcal**6)
- Evan2cat=epscalc*2*(r06/rcal**3)
+ Evan1cat=epscalc*(r012/(rcal**6))
+ Evan2cat=epscalc*2*(r06/(rcal**3))
Eeleccat=k0/ract
r7 = rcal**7
r4 = rcal**4
gradcatcat(k,i)=gradcatcat(k,i)-gg(k)
gradcatcat(k,j)=gradcatcat(k,j)+gg(k)
enddo
-
+ if (energy_dec) write (iout,*) i,j,Evan1cat,Evan2cat,Eeleccat,&
+ r012,rcal**6,ichargecat(itypi)*ichargecat(itypj)
! write(iout,*) "ecatcat",i,j, ecationcation,xj,yj,zj
ecationcation=ecationcation+Evan1cat+Evan2cat+Eeleccat
enddo
!c! No charges - do nothing
eheadtail = 0.0d0
- ELSE IF (isel.eq.1 .and. iabs(Qj).eq.1) THEN
+ ELSE IF (isel.eq.1) THEN
!c! Nonpolar-charge interactions
if ((itype(i,1).eq.27).or.(itype(i,1).eq.26).or.(itype(i,1).eq.25)) then
Qi=Qi*2
eheadtail = epol
! eheadtail = 0.0d0
- ELSE IF (isel.eq.3 .and. icharge(itypj).eq.2) THEN
+ ELSE IF (isel.eq.3) THEN
!c! Dipole-charge interactions
if ((itype(i,1).eq.27).or.(itype(i,1).eq.26).or.(itype(i,1).eq.25)) then
Qi=Qi*2
Qj=Qj*2
Qij=Qij*2
endif
+ write(iout,*) "KURWA0",d1
+
CALL edq_cat(ecl, elj, epol)
eheadtail = ECL + elj + epol
! eheadtail = 0.0d0
- ELSE IF ((isel.eq.2.and. &
- iabs(Qi).eq.1).and. &
- nstatecat(itypi,itypj).eq.1) THEN
+ ELSE IF ((isel.eq.2)) THEN
!c! Same charge-charge interaction ( +/+ or -/- )
if ((itype(i,1).eq.27).or.(itype(i,1).eq.26).or.(itype(i,1).eq.25)) then
real(kind=8) :: facd4, adler, Fgb, facd3
integer troll,jj,istate
real (kind=8) :: dcosom1(3),dcosom2(3)
+ evdw=0.0d0
eps_out=80.0d0
sss_ele_cut=1.0d0
! print *,"EVDW KURW",evdw,nres
!c!-------------------------------------------------------------------
!c! ecl
+ write(iout,*) "KURWA2",Rhead
sparrow = w1 * Qj * om1
hawk = w2 * Qj * Qj * (1.0d0 - sqom2)
ECL = sparrow / Rhead**2.0d0 &
include 'mpif.h'
real*8 :: xi,yi,zi,xj,yj,zj,xj_safe,yj_safe,zj_safe,xj_temp,yj_temp,zj_temp
real*8 :: dist_init, dist_temp,r_buff_list
- integer:: contlisti(50*nres),contlistj(50*nres)
- integer :: newcontlisti(50*nres),newcontlistj(50*nres)
+ integer:: contlisti(200*nres),contlistj(200*nres)
+! integer :: newcontlisti(200*nres),newcontlistj(200*nres)
integer i,j,itypi,itypj,subchap,xshift,yshift,zshift,iint,ilist_sc,g_ilist_sc
integer displ(0:nprocs),i_ilist_sc(0:nprocs),ierr
! print *,"START make_SC"
+ r_buff_list=5.0
ilist_sc=0
do i=iatsc_s,iatsc_e
itypi=iabs(itype(i,1))
write (iout,*) i,newcontlisti(i),newcontlistj(i)
enddo
#endif
-
+ call int_bounds(g_ilist_sc,g_listscsc_start,g_listscsc_end)
return
end subroutine make_SCSC_inter_list
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
subroutine make_SCp_inter_list
+ use MD_data, only: itime_mat
+
include 'mpif.h'
real*8 :: xi,yi,zi,xj,yj,zj,xj_safe,yj_safe,zj_safe,xj_temp,yj_temp,zj_temp
real*8 :: dist_init, dist_temp,r_buff_list
- integer:: contlistscpi(50*nres),contlistscpj(50*nres)
- integer :: newcontlistscpi(50*nres),newcontlistscpj(50*nres)
+ integer:: contlistscpi(200*nres),contlistscpj(200*nres)
+! integer :: newcontlistscpi(200*nres),newcontlistscpj(200*nres)
integer i,j,itypi,itypj,subchap,xshift,yshift,zshift,iint,ilist_scp,g_ilist_scp
integer displ(0:nprocs),i_ilist_scp(0:nprocs),ierr
- print *,"START make_SC"
+! print *,"START make_SC"
+ r_buff_list=5.0
ilist_scp=0
do i=iatscp_s,iatscp_e
if (itype(i,1).eq.ntyp1 .or. itype(i+1,1).eq.ntyp1) cycle
yj=yj_safe-yi
zj=zj_safe-zi
endif
+#ifdef DEBUG
+ ! r_buff_list is a read value for a buffer
+ if ((sqrt(dist_init).le.(r_cut_ele)).and.(ifirstrun.eq.0)) then
+! Here the list is created
+ ilist_scp_first=ilist_scp_first+1
+! this can be substituted by cantor and anti-cantor
+ contlistscpi_f(ilist_scp_first)=i
+ contlistscpj_f(ilist_scp_first)=j
+ endif
+#endif
! r_buff_list is a read value for a buffer
if (sqrt(dist_init).le.(r_cut_ele+r_buff_list)) then
! Here the list is created
do i=1,g_ilist_scp
write (iout,*) i,newcontlistscpi(i),newcontlistscpj(i)
enddo
+
+! if (ifirstrun.eq.0) ifirstrun=1
+! do i=1,ilist_scp_first
+! do j=1,g_ilist_scp
+! if ((newcontlistscpi(j).eq.contlistscpi_f(i)).and.&
+! (newcontlistscpj(j).eq.contlistscpj_f(i))) go to 126
+! enddo
+! print *,itime_mat,"ERROR matrix needs updating"
+! print *,contlistscpi_f(i),contlistscpj_f(i)
+! 126 continue
+! enddo
#endif
+ call int_bounds(g_ilist_scp,g_listscp_start,g_listscp_end)
+
return
end subroutine make_SCp_inter_list
!-----------------------------------------------------------------------------
+ subroutine make_pp_inter_list
+ include 'mpif.h'
+ real*8 :: xi,yi,zi,xj,yj,zj,xj_safe,yj_safe,zj_safe,xj_temp,yj_temp,zj_temp
+ real*8 :: xmedj,ymedj,zmedj
+ real*8 :: dist_init, dist_temp,r_buff_list,dxi,dyi,dzi,xmedi,ymedi,zmedi
+ real*8 :: dx_normi,dy_normi,dz_normi,dxj,dyj,dzj,dx_normj,dy_normj,dz_normj
+ integer:: contlistppi(200*nres),contlistppj(200*nres)
+! integer :: newcontlistppi(200*nres),newcontlistppj(200*nres)
+ integer i,j,itypi,itypj,subchap,xshift,yshift,zshift,iint,ilist_pp,g_ilist_pp
+ integer displ(0:nprocs),i_ilist_pp(0:nprocs),ierr
+! print *,"START make_SC"
+ ilist_pp=0
+ r_buff_list=5.0
+ do i=iatel_s,iatel_e
+ if (itype(i,1).eq.ntyp1 .or. itype(i+1,1).eq.ntyp1) cycle
+ dxi=dc(1,i)
+ dyi=dc(2,i)
+ dzi=dc(3,i)
+ dx_normi=dc_norm(1,i)
+ dy_normi=dc_norm(2,i)
+ dz_normi=dc_norm(3,i)
+ xmedi=c(1,i)+0.5d0*dxi
+ ymedi=c(2,i)+0.5d0*dyi
+ zmedi=c(3,i)+0.5d0*dzi
+ xmedi=dmod(xmedi,boxxsize)
+ if (xmedi.lt.0) xmedi=xmedi+boxxsize
+ ymedi=dmod(ymedi,boxysize)
+ if (ymedi.lt.0) ymedi=ymedi+boxysize
+ zmedi=dmod(zmedi,boxzsize)
+ if (zmedi.lt.0) zmedi=zmedi+boxzsize
+ do j=ielstart(i),ielend(i)
+! write (iout,*) i,j,itype(i,1),itype(j,1)
+ if (itype(j,1).eq.ntyp1.or. itype(j+1,1).eq.ntyp1) cycle
+
+! 1,j)
+ dxj=dc(1,j)
+ dyj=dc(2,j)
+ dzj=dc(3,j)
+ dx_normj=dc_norm(1,j)
+ dy_normj=dc_norm(2,j)
+ dz_normj=dc_norm(3,j)
+! xj=c(1,j)+0.5D0*dxj-xmedi
+! yj=c(2,j)+0.5D0*dyj-ymedi
+! zj=c(3,j)+0.5D0*dzj-zmedi
+ xj=c(1,j)+0.5D0*dxj
+ yj=c(2,j)+0.5D0*dyj
+ zj=c(3,j)+0.5D0*dzj
+ xj=mod(xj,boxxsize)
+ if (xj.lt.0) xj=xj+boxxsize
+ yj=mod(yj,boxysize)
+ if (yj.lt.0) yj=yj+boxysize
+ zj=mod(zj,boxzsize)
+ if (zj.lt.0) zj=zj+boxzsize
+
+ dist_init=(xj-xmedi)**2+(yj-ymedi)**2+(zj-zmedi)**2
+ xj_safe=xj
+ yj_safe=yj
+ zj_safe=zj
+ do xshift=-1,1
+ do yshift=-1,1
+ do zshift=-1,1
+ xj=xj_safe+xshift*boxxsize
+ yj=yj_safe+yshift*boxysize
+ zj=zj_safe+zshift*boxzsize
+ dist_temp=(xj-xmedi)**2+(yj-ymedi)**2+(zj-zmedi)**2
+ if(dist_temp.lt.dist_init) then
+ dist_init=dist_temp
+ xj_temp=xj
+ yj_temp=yj
+ zj_temp=zj
+ endif
+ enddo
+ enddo
+ enddo
+
+ if (sqrt(dist_init).le.(r_cut_ele+r_buff_list)) then
+! Here the list is created
+ ilist_pp=ilist_pp+1
+! this can be substituted by cantor and anti-cantor
+ contlistppi(ilist_pp)=i
+ contlistppj(ilist_pp)=j
+ endif
+ enddo
+ enddo
+! enddo
+#ifdef DEBUG
+ write (iout,*) "before MPIREDUCE",ilist_pp
+ do i=1,ilist_pp
+ write (iout,*) i,contlistppi(i),contlistppj(i)
+ enddo
+#endif
+ if (nfgtasks.gt.1)then
+
+ call MPI_Reduce(ilist_pp,g_ilist_pp,1,&
+ MPI_INTEGER,MPI_SUM,king,FG_COMM,IERR)
+! write(iout,*) "before bcast",g_ilist_sc
+ call MPI_Gather(ilist_pp,1,MPI_INTEGER,&
+ i_ilist_pp,1,MPI_INTEGER,king,FG_COMM,IERR)
+ displ(0)=0
+ do i=1,nfgtasks-1,1
+ displ(i)=i_ilist_pp(i-1)+displ(i-1)
+ enddo
+! write(iout,*) "before gather",displ(0),displ(1)
+ call MPI_Gatherv(contlistppi,ilist_pp,MPI_INTEGER,&
+ newcontlistppi,i_ilist_pp,displ,MPI_INTEGER,&
+ king,FG_COMM,IERR)
+ call MPI_Gatherv(contlistppj,ilist_pp,MPI_INTEGER,&
+ newcontlistppj,i_ilist_pp,displ,MPI_INTEGER,&
+ king,FG_COMM,IERR)
+ call MPI_Bcast(g_ilist_pp,1,MPI_INT,king,FG_COMM,IERR)
+! write(iout,*) "before bcast",g_ilist_sc
+! call MPI_Bcast(g_ilist_sc,1,MPI_INT,king,FG_COMM)
+ call MPI_Bcast(newcontlistppi,g_ilist_pp,MPI_INT,king,FG_COMM,IERR)
+ call MPI_Bcast(newcontlistppj,g_ilist_pp,MPI_INT,king,FG_COMM,IERR)
+
+! call MPI_Bcast(g_ilist_sc,1,MPI_INT,king,FG_COMM)
+
+ else
+ g_ilist_pp=ilist_pp
+
+ do i=1,ilist_pp
+ newcontlistppi(i)=contlistppi(i)
+ newcontlistppj(i)=contlistppj(i)
+ enddo
+ endif
+ call int_bounds(g_ilist_pp,g_listpp_start,g_listpp_end)
+#ifdef DEBUG
+ write (iout,*) "after MPIREDUCE",g_ilist_pp
+ do i=1,g_ilist_pp
+ write (iout,*) i,newcontlistppi(i),newcontlistppj(i)
+ enddo
+#endif
+ return
+ end subroutine make_pp_inter_list
+
+!-----------------------------------------------------------------------------
+!-----------------------------------------------------------------------------
+
+
+
end module energy