!----------------------------------------------------------------------- ! ! SUBROUTINE ilist_group(p,nterms,nterms_gr,bcount_ele,iterms,iterms_gr, & pmass, pmass_gr, pquad, pquad_gr, drdotdr_gr,dx_gr,dy_gr,dz_gr) ! ! !----------------------------------------------------------------------- ! ! ! Subroutine to form the interaction list of a grouping cell. ! ! P= number of cell group ! nterms= number of elements INside Sphere ! nterms_gr= number of elements OUTside Sphere ! bcount_ele= number of elements INside C_group ! iterms_gr=list of global elements outside the spere ! !----------------------------------------------------------------------- USE fly_h implicit none ! Declaration of local variables. ! ------------------------------- INTEGER(KIND=4) :: bcount_ele INTEGER(KIND=4) :: j,ind_in,ind_in3,indx,lmax1,ii INTEGER(KIND=4) ::p,i,nnodes,nsubdiv,nterms,nterms_gr,nl INTEGER(KIND=4), DIMENSION (nsubcell) :: isub INTEGER(KIND=4), DIMENSION (ndim) ::dire INTEGER(KIND=4), DIMENSION (maxilf)::nodelist,twasubp,twisubset INTEGER :: il_sh LOGICAL ::tolcrit REAL(KIND=8) :: dp_x,dp_y,dp_z,cc_test REAL(KIND=8) :: ddx,ddy,ddz,dd2,sphere,cl2,cl_s,mass_s REAL(KIND=8) :: ddr,dd_min REAL(KIND=8), DIMENSION (ndim) ::pos_l, pos_cl,d_c,dr2,dp REAL(KIND=8), DIMENSION (ndim) ::sh REAL(KIND=8), DIMENSION (ndim,ndim) ::dd_r REAL(KIND=8), DIMENSION (2*ndim-1) ::quad_s INTEGER(KIND=4), DIMENSION(maxnterm) :: iterms ,iterms_gr REAL(KIND=8), DIMENSION(maxnterm) :: pmass ,pmass_gr REAL(KIND=8), DIMENSION(2*ndim-1,maxnterm) :: pquad ,pquad_gr REAL(KIND=8), DIMENSION(maxnterm)::drdotdr_gr ,dx_gr ,dy_gr ,dz_gr !======================================================================= ! DO i=nb_loc+1,cell_ss(lmax,2) ! write(140+me,*)"il1 PE=",me," cella=",i," pos=",pos_cell(1:3,i-nb_loc) ! ENDDO ! STOP !----------------------------------------------------------------------- ! Group-Cell element properties !----------------------------------------------------------------------- il_sh = p - nb_loc do i=1,lmax if ((p .ge. cell_ss(i,1)) .and. (p .le. cell_ss(i,2))) then cl2=size_level(i) exit endif end do cl2=cl2*cl2 sphere=(27./4.)*cl2 bcount_ele=0 !----------------------------------------------------------------------- ! Initialize list of cells to examine. !----------------------------------------------------------------------- nterms=0 nterms_gr=0 nnodes=1 nodelist(1)=root lmax1=0 ! CACHE pos_l(1:ndim)=pos_cell(1:ndim,il_sh) !----------------------------------------------------------------------- ! Loop until no cells are left to examine. !-----------------------------------------------------------------------! DO WHILE (nnodes.GT.0) lmax1 = lmax1 + 1 cl_s = size_level(lmax1) nsubdiv=0 !----------------------------------------------------------------------- ! Apply tolerance criterion to list of cells. !----------------------------------------------------------------------- DO i=1,nnodes nl=nodelist(i) IF( nl .le. nb_loc) THEN pos_cl(1:ndim)=pos(1:ndim,nl) ddx=pos_l(1)-pos_cl(1) ddy=pos_l(2)-pos_cl(2) ddz=pos_l(3)-pos_cl(3) dd2=ddx**2+ddy**2+ddz**2 IF(dd2.GT.sphere) THEN nterms_gr=nterms_gr+1 iterms_gr(nterms_gr)=nl pmass_gr(nterms_gr)=mass_read dx_gr(nterms_gr)=ddx dy_gr(nterms_gr)=ddy dz_gr(nterms_gr)=ddz drdotdr_gr(nterms_gr)=dd2 ELSE nterms=nterms+1 iterms(nterms)=nl pmass(nterms)=mass_read ENDIF CYCLE ! goto next nodelist(i) element ELSE !IF( nl .le. nb_loc) il_sh=nl-nb_loc pos_cl(1:ndim)=pos_cell(1:ndim,il_sh) cc_test = cl_s**2*tol2inv ddx=pos_l(1)-pos_cl(1) ddy=pos_l(2)-pos_cl(2) ddz=pos_l(3)-pos_cl(3) dd2=ddx**2+ddy**2+ddz**2 tolcrit= dd2 .GE. cc_test !----------------------------------------------------------------------- ! Ewald boundary periodic conditions section !----------------------------------------------------------------------- IF(tolcrit) THEN !----------------------------------------------------------------------- ! Ewald: Search for "mirror" cells. ! ! For each simmetry plane checks whether the ! distance body[pos(p)]--cell[pos(pc)] is such ! that the reflexion w.r.t. direction i ! will make a nearest cell to body pos(p). !----------------------------------------------------------------------- dp(1)=0. dp(2)=0. dp(3)=0. ind_in=0 DO ii=1, ndim sh(ii)=pos_l(ii)-pos_cl(ii) ddr=abs(sh(ii)) if(ddr.gt.L2)then ind_in=ind_in+1 dire(ind_in)=ii ! ind_in3=INT(sh(ii)/ddr) ind_in3=1 IF(sh(ii).lt.0) ind_in3=-1 d_c(ind_in)=ind_in3*Lbox endif ENDDO ! If no reflection w.r.t. any plane can produce nearest cell, go to end IF(ind_in.eq.0) THEN mass_s=mass_cell(il_sh) quad_s(1:2*ndim-1)=quad(1:2*ndim-1,il_sh) IF(dd2.GT.sphere) THEN nterms_gr=nterms_gr+1 iterms_gr(nterms_gr)=nl pmass_gr(nterms_gr)=mass_s pquad_gr(1:2*ndim-1,nterms_gr)=quad_s(1:2*ndim-1) dx_gr(nterms_gr)=ddx dy_gr(nterms_gr)=ddy dz_gr(nterms_gr)=ddz drdotdr_gr(nterms_gr)=dd2 ELSE IF(nl.NE.p) THEN nterms=nterms+1 iterms(nterms)=nl pmass(nterms)=mass_s pquad(1:2*ndim-1,nterms)=quad_s(1:2*ndim-1) ENDIF ENDIF CYCLE ! goto next nodelist(i) element ENDIF ! i=index of the mirror cell j=direction (x,y or z) do ii=1, ind_in do j=1, ndim dd_r(ii,j)=pos_cl(j) if(j.eq.dire(ii))dd_r(ii,j)=dd_r(ii,j)+d_c(ii) enddo enddo ! Search the nearest cell among the 3 reflexions if(ind_in.eq.1)then indx=ind_in go to 200 endif DO ii=1, ind_in dr2(ii)=(pos_l(1)-dd_r(ii,1))* & (pos_l(1)-dd_r(ii,1)) + & (pos_l(2)-dd_r(ii,2))*(pos_l(2)-dd_r(ii,2))+ & (pos_l(3)-dd_r(ii,3))*(pos_l(3)-dd_r(ii,3)) ENDDO indx=1 dd_min=dr2(1) if(ind_in.eq.2)then if(dr2(2) .le. dd_min) indx=2 endif if(ind_in.eq.3)then if(dr2(2).le.dd_min)indx=2 if(dr2(3).le.dd_min.and.dr2(3).le.dr2(2))indx=3 endif ! Assigns to dp(i=1,ndim) the coordinate shift to reach ! from cell "pc" the nearest (to particle "p") cell outside the domain. ! 200 continue DO ii=1, ndim dp(ii)=dd_r(indx,ii)-pos_cl(ii) ENDDO dp_x=pos_l(1)-(pos_cl(1)+dp(1)) dp_y=pos_l(2)-(pos_cl(2)+dp(2)) dp_z=pos_l(3)-(pos_cl(3)+dp(3)) tolcrit=(dp_x*dp_x+dp_y*dp_y+dp_z*dp_z).GE.cc_test ENDIF !----------------------------------------------------------------------- ! End of Ewald section !----------------------------------------------------------------------- IF(tolcrit) THEN mass_s=mass_cell(il_sh) quad_s(1:2*ndim-1)=quad(1:2*ndim-1,il_sh) IF(dd2.GT.sphere) THEN nterms_gr=nterms_gr+1 iterms_gr(nterms_gr)=nl pmass_gr(nterms_gr)=mass_s pquad_gr(1:2*ndim-1,nterms_gr)=quad_s(1:2*ndim-1) dx_gr(nterms_gr)=ddx dy_gr(nterms_gr)=ddy dz_gr(nterms_gr)=ddz drdotdr_gr(nterms_gr)=dd2 ELSE IF(nl.NE.p) THEN nterms=nterms+1 iterms(nterms)=nl pmass(nterms)=mass_s pquad(1:2*ndim-1,nterms)=quad_s(1:2*ndim-1) ENDIF ENDIF ELSE !IF(tolcrit) IF(nl.ne.p) THEN nsubdiv=nsubdiv+1 twisubset(nsubdiv)=i !nodelist index: cell to be opened ENDIF ENDIF !IF(tolcrit) ENDIF !IF( nl .le. nb_loc) END DO ! DO i=1,nnodes !----------------------------------------------------------------------- ! Add cells which satisfy criterion to interaction list. Note that, ! depending on theta, self-interaction term will be included. !----------------------------------------------------------------------- IF(nterms_gr.GT.maxnterm) THEN write(uterm,*)'PE=',me,' nterms_gr=',nterms_gr, & ' p=',p,' nnodes=',nnodes call flush(uterm) CALL error('ilist_group error 1: overflow') ENDIF IF(nterms.GT.maxnterm) THEN write(uterm,*)'PE=',me,' nterms=',nterms, & ' p=',p,' nnodes=',nnodes call flush(uterm) CALL error('ilist_group error 2: overflow') ENDIF !----------------------------------------------------------------------- ! Add subcells of cells which fail tolerance criterion to list of ! cells to examine. !----------------------------------------------------------------------- IF(8*nsubdiv.GT.maxilf.OR.8*nsubdiv.GT.ncells) THEN write(uterm,*)'Errore PE=',me,' nterms=',nterms, & ' p=',p CALL error('ilist_group error 3: overflow') ENDIF DO 40 i=1,nsubdiv il_sh=nodelist(twisubset(i))-nb_loc isub(1:nsubcell)=subp(1:nsubcell,il_sh) twasubp(i)=isub(1) twasubp(i+nsubdiv)=isub(2) twasubp(i+2*nsubdiv)=isub(3) twasubp(i+3*nsubdiv)=isub(4) twasubp(i+4*nsubdiv)=isub(5) twasubp(i+5*nsubdiv)=isub(6) twasubp(i+6*nsubdiv)=isub(7) twasubp(i+7*nsubdiv)=isub(8) 40 CONTINUE nnodes=0 DO i=1,8*nsubdiv IF(twasubp(i).NE.0) THEN nnodes=nnodes+1 twisubset(nnodes)=i ENDIF ENDDO DO 60 i=1,nnodes nodelist(i)=twasubp(twisubset(i)) !list of the new cell to examin 60 CONTINUE ENDDO ! DO WHILE (nnodes.GT.0) ! start new cycle with new nodelist(nnodes) !----------------------------------------------------------------------- ! Examin the group cell 'p' and put particles in the il_near list !----------------------------------------------------------------------- nsubdiv=1 nodelist(1)=p 100 CONTINUE DO WHILE (nsubdiv.GT.0) DO i = 1, nsubdiv il_sh=nodelist(i)-nb_loc isub(1:nsubcell)=subp(1:nsubcell,il_sh) twasubp(i)=isub(1) twasubp(i+nsubdiv)=isub(2) twasubp(i+2*nsubdiv)=isub(3) twasubp(i+3*nsubdiv)=isub(4) twasubp(i+4*nsubdiv)=isub(5) twasubp(i+5*nsubdiv)=isub(6) twasubp(i+6*nsubdiv)=isub(7) twasubp(i+7*nsubdiv)=isub(8) END DO nnodes=0 DO i=1,8*nsubdiv IF(twasubp(i).GT.nb_loc) THEN nnodes=nnodes+1 nodelist(nnodes)=twasubp(i) ELSE IF(twasubp(i).GT.0) THEN bcount_ele=bcount_ele+1 nterms=nterms+1 iterms(nterms)=twasubp(i) pmass(nterms)=mass_read ENDIF END DO nsubdiv=nnodes IF(nterms.GT.maxnterm) THEN write(uterm,*)'PE=',me,' nterms=',nterms, & ' p=',p,' nnodes=',nnodes call flush(uterm) CALL error('ilist_group error 4: overflow') ENDIF IF(8*nsubdiv.GT.maxilf .or. 8*nsubdiv.GT.ncells) THEN write(uterm,*)'PE=',me,' nsubiv=',nsubdiv, & ' p=',p call flush(uterm) CALL error('ilist_group error 5: overflow') ENDIF ENDDO !dowhile RETURN END