!----------------------------------------------------------------------- ! SUBROUTINE tree_gen ! ! !----------------------------------------------------------------------- ! ! ! Subroutine to insert the bodies into the tree. ! Active bodies are those which are not yet in ! place in the tree, as leaves. The local variables pm1 and nindex ! are used to convert back and forth between physical coordinates ! and subcell coordinates. !----------------------------------------------------------------------- USE fly_h implicit none INCLUDE 'mpif.h' ! Declaration of local variables. ! ------------------------------- REAL(KIND=8), DIMENSION (nsubcell,ndim) :: pm1 REAL(KIND=8), DIMENSION (:) , ALLOCATABLE:: pl,pcl,pcl_par REAL(KIND=8) :: cl_size INTEGER(KIND=4),SAVE :: ns_treegen,icicl,jj INTEGER(KIND=4) :: nbodlist,nclist,indcell,p,nbodlist_old,plocal INTEGER(KIND=4) :: bodycella, counter INTEGER(KIND=4), DIMENSION (nb_loc) :: sub_temp1,par_temp1 INTEGER(KIND=4), DIMENSION (nb_loc) :: temp1 INTEGER(KIND=4), DIMENSION (ndim) :: nindex INTEGER(KIND=4), DIMENSION (:) , ALLOCATABLE ::isub,sub_app INTEGER :: i,j,k,ir,ic,ind_l1, il_sh, pbefore INTEGER(KIND=4), DIMENSION (:) , ALLOCATABLE ::pcg,globalcounter INTEGER :: TID, NTID,nt INTEGER:: status DATA pm1/4*-1.,4*1.,2*-1.,2*1.,2*-1.,2*1.,-1.,1.,-1.,1., & -1.,1.,-1.,1./,nindex/4,2,1/ !----------------------------------------------------------------------- NTID = 1 NTID = OMP_GET_MAX_THREADS() ALLOCATE(pcg(NTID), STAT=status) ALLOCATE(globalcounter(NTID), STAT=status) root=nbodsmax+1 indcell=0 mark_gr_cell=0 !---------------------------------------------------------------------- ! Deallocate old tree, compute coordinates of center of root cell. !---------------------------------------------------------------------- incells=1 subp=0 !zeroing all the tree DO k=1,ndim pos_cell(k,1)=rmin(k)+0.5*rsize ENDDO !----------------------------------------------------------------------- ! Place all bodies on active body list, having root as parent; place ! root on active cell list. !----------------------------------------------------------------------- temp1=0 par_temp1=0 par_temp1(1:nb_res_loc(me+1))=root DO i=1,nb_res_loc(me+1) temp1(i)=i ENDDO nbodlist=nb_res_loc(me+1) nclist=1 icicl=0 lmax=1 cell_ss(1,1)=root cell_ss(1,2)=root size_level(1)=rsize !----------------------------------------------------------------------- ! Loop until no bodies are left active. !----------------------------------------------------------------------- DO WHILE(nclist.GT.0) !Build a new tree levels icicl=icicl+1 IF(me.EQ.0) THEN write(uterm,*)'cycle=',icicl,' incells=',incells,' nclist=',nclist,' ',& ' size=',size_level(lmax) call flush(uterm) ENDIF !----------------------------------------------------------------------- ! Compute subindices for all active bodies. !----------------------------------------------------------------------- !$OMP PARALLEL PRIVATE(i,ic,il_sh, pcl, pl,k,status) ALLOCATE(pl(ndim), STAT=status) ALLOCATE(pcl(ndim), STAT=status) !$OMP DO DO i=1,nbodlist sub_temp1(i)=1 il_sh=par_temp1(i)-nbodsmax pcl(1:ndim)=pos_cell(1:ndim,il_sh) pl(1:3)= pos(1:3,temp1(i)) DO k=1,ndim IF(pl(k).GE.pcl(k)) sub_temp1(i)=sub_temp1(i)+nindex(k) ENDDO !----------------------------------------------------------------------- ! Compute number of bodies in each subcell. !----------------------------------------------------------------------- ic=sub_temp1(i) ! write(uterm,*) "TP0.1 PE=",me," TID=",TID," i=",i," pcl=",pcl," pl=",pl," sub_temp1-ic=",sub_temp1(i)," il_sh=",il_sh," subp(1:8,il_sh)=",subp(1:8,il_sh) IF(subp(ic,il_sh).gt.2) CYCLE !$OMP ATOMIC subp(ic,il_sh)=subp(ic,il_sh)+1 !$OMP END ATOMIC ! write(uterm,*) "TP0.2 PE=",me," TID=",TID," i=",i," pcl=",pcl," pl=",pl," sub_temp1-ic=",sub_temp1(i)," il_sh=",il_sh," subp(1:8,il_sh)=",subp(1:8,il_sh) ENDDO ! do i=1,nbodlist !$OMP END DO DEALLOCATE(pl) DEALLOCATE(pcl) !$OMP END PARALLEL !----------------------------------------------------------------------- ! At the end of bodies analysis flush all the cache !----------------------------------------------------------------------- ! write(uterm,*) "TP1.0 PE=",me," lmax=",lmax," subp=", subp(1:8,cell_ss(lmax,1)-nbodsmax:cell_ss(lmax,2)-nbodsmax) IF(me.eq.0) write(6,*)"PE=",me," cycle=",icicl," nbodlist=",nbodlist indcell=incells+nbodsmax !----------------------------------------------------------------------- ! Prepare arrays for find_group !----------------------------------------------------------------------- !$OMP PARALLEL PRIVATE(il_sh,i,j,bodycella,isub,status) ALLOCATE(isub(nsubcell), STAT=status) IF(lmax.GT.ncrit) THEN !loop only on the local residing cells: !check my starting address on the cell list !$OMP DO DO i=cell_ss(lmax,1),cell_ss(lmax,2) il_sh=i-nbodsmax IF (mark_gr_cell(il_sh).ne.0) CYCLE isub(1:8)=subp(1:8,il_sh) bodycella=0 DO j=1,2**ndim bodycella=bodycella+isub(j) END DO IF(bodycella.le.nbodcrit .and. bodycella.gt.2) mark_gr_cell(il_sh)=bodycella ENDDO !$OMP END DO END IF !if(lmax....) DEALLOCATE(isub) !$OMP END PARALLEL !----------------------------------------------------------------------- ! Create the sub_cells: subp array pointers !----------------------------------------------------------------------- cl_size=size_level(lmax)*.5 il_sh=cell_ss(lmax,1)- nbodsmax counter=0 p=indcell pcg=0 !$OMP PARALLEL PRIVATE(il_sh,TID,nt,pbefore,bodycella,isub,sub_app,pcl_par,pcl,counter,ind_l1,status,i,j,plocal) ALLOCATE(isub(nsubcell), STAT=status) ALLOCATE(sub_app(nsubcell), STAT=status) ALLOCATE(pcl_par(ndim), STAT=status) ALLOCATE(pcl(ndim), STAT=status) !----------------------------------------------------------------------- ! Ogni thread conta quanti indici di sottocella genera in pcg !----------------------------------------------------------------------- TID = 0 TID = OMP_GET_THREAD_NUM() pbefore=0 !$OMP DO DO i = cell_ss(lmax, 1), cell_ss(lmax, 2) il_sh=i-nbodsmax isub(1:nsubcell)=subp(1:nsubcell,il_sh) DO j=1,nsubcell IF (isub(j) .GT. 1) THEN pcg(TID+1)=pcg(TID+1)+1 ENDIF ENDDO ENDDO !$OMP END DO DO nt=1,TID pbefore=pbefore+pcg(nt) ENDDO plocal=p+pbefore ! write(uterm,*) "TP1 CK PE=",me," TID=",TID," of ", NTID," pbefore=",pbefore," plocal=",p ! call flush(uterm) !$OMP DO DO i = cell_ss(lmax, 1), cell_ss(lmax, 2) il_sh=i-nbodsmax isub(1:nsubcell)=subp(1:nsubcell,il_sh) pcl_par(1:ndim)=pos_cell(1:ndim,il_sh) sub_app=isub DO j=1,nsubcell IF (isub(j) .GT. 1) THEN plocal=plocal+1 sub_app(j)=plocal IF(plocal.gt.nbodsmax+ncells) THEN call error('tree_gen error 1: increase ncells in fly_h.F90') ENDIF ind_l1=plocal-nbodsmax pcl=pcl_par pcl(1)=pcl(1)+pm1(j,1)*0.5*cl_size pcl(2)=pcl(2)+pm1(j,2)*0.5*cl_size pcl(3)=pcl(3)+pm1(j,3)*0.5*cl_size pos_cell(1:3,ind_l1)=pcl(1:3) IF(mark_gr_cell(il_sh).ne.0) mark_gr_cell(ind_l1)=-1 !used in find_group ENDIF !if(isub...) ENDDO !DO j=1,nsubcell subp(1:nsubcell,il_sh)=sub_app(1:nsubcell) ! if(me.eq.0) write(uterm,*) "TP1.2 PE=",me," TID=",TID," i=",i," lmax=",lmax," il_sh+nbodsmax=",il_sh+nbodsmax," subp(1:8,il_sh)=",subp(1:8,il_sh)," pos_cell(1:3,il_sh)=",pos_cell(1:3,il_sh) ENDDO !DO i= !$OMP END DO DEALLOCATE(isub) DEALLOCATE(sub_app) DEALLOCATE(pcl_par) DEALLOCATE(pcl) !$OMP END PARALLEL DO i=1,NTID p=p+pcg(i) ENDDO ns_treegen=p-indcell !----------------------------------------------------------------------- ! Find all subcells with one body; add bodies to tree. !----------------------------------------------------------------------- nbodlist_old=nbodlist globalcounter=0 !$OMP PARALLEL PRIVATE(ir,ic,TID,isub,counter,status,i) TID = 0 TID = OMP_GET_THREAD_NUM() ALLOCATE(isub(nsubcell), STAT=status) !$OMP DO DO i=1,nbodlist call flush(uterm) ir=par_temp1(i)-nbodsmax ic=sub_temp1(i) isub(1:nsubcell)=subp(1:nsubcell,ir) IF(isub(ic).GT.nbodsmax) THEN par_temp1(i)=isub(ic) ELSE globalcounter(TID+1)=globalcounter(TID+1)+1 ! write(uterm,*)"TP2 PE=",me," TID=",TID," founded body: gbc(TID+1)=",globalcounter(TID+1) call flush(uterm) subp(ic,ir)=temp1(i) subp(ic,ir)=temp1(i) temp1(i)=0 par_temp1(i)=0 ENDIF ENDDO !$OMP END DO DEALLOCATE(isub) !$OMP END PARALLEL counter=0 DO i=1,NTID counter=counter+globalcounter(i) ! write(uterm,*)"TP2.1 PE=",me," i=",i," counter=",counter," ",globalcounter(i) ENDDO nbodlist=nbodlist-counter !----------------------------------------------------------------------- ! Prepare the arrays to start new cycle. !----------------------------------------------------------------------- IF(ns_treegen.GT.0) THEN lmax=lmax+1 IF(lmax.GT.nmax_level) THEN write(uterm,*)'tree_gen error 1: Max level reached:', lmax,' greater then nmax_level' STOP ENDIF cell_ss(lmax,1)=cell_ss(lmax-1,2)+1 cell_ss(lmax,2)=cell_ss(lmax,1)+ns_treegen-1 size_level(lmax)=cl_size ENDIF incells=incells+ns_treegen nclist=ns_treegen IF(incells.GT.ncells) CALL error('tree_gen error 2: overflow') ns_treegen=0 !----------------------------------------------------------------------- ! Pack of arrays temp1 and par_temp1 !----------------------------------------------------------------------- sub_temp1=0 j=0 DO i=1,nbodlist_old IF(temp1(i).gt.0) THEN j=j+1 sub_temp1(j)=temp1(i) ENDIF ENDDO temp1=sub_temp1 sub_temp1=0 j=0 DO i=1,nbodlist_old IF(par_temp1(i).gt.0) THEN j=j+1 sub_temp1(j)=par_temp1(i) ENDIF ENDDO par_temp1=sub_temp1 sub_temp1=0 ENDDO ! do while (nclist.gt.0) DEALLOCATE(pcg) DEALLOCATE(globalcounter) CALL MPI_BARRIER(MPI_COMM_WORLD,ierror) RETURN END