12 integer,
allocatable :: gsq_sfc(:^D&),seq_sfc(:),seq_ig^D(:)
13 integer :: ig^D, ngsq^D, isq, total_number
15 logical,
allocatable :: in_domain(:)
18 ngsq^d=2**ceiling(log(real(
ng^d(1)))/log(2.0));
20 {ngsq^d=max(ngsq^dd) \}
22 total_number={ngsq^d|*}
24 allocate(gsq_sfc(ngsq^d))
26 allocate(seq_sfc(total_number))
28 {
allocate(seq_ig^d(total_number))\}
29 allocate(in_domain(total_number))
34 seq_sfc(gsq_sfc(ig^d))=gsq_sfc(ig^d)
35 {seq_ig^d(gsq_sfc(ig^dd))=ig^d \}
39 if (seq_ig^d(isq)>ng^d(1)|.or.)
then
40 seq_sfc(isq:total_number)=seq_sfc(isq:total_number)-1
41 in_domain(isq)=.false.
45 if(.not.
allocated(iglevel1_sfc))
allocate(iglevel1_sfc(ng^d(1)))
46 if(.not.
allocated(sfc_iglevel1))
allocate(sfc_iglevel1(ndim,nglev1))
48 if(in_domain(isq))
then
49 iglevel1_sfc(seq_ig^d(isq))=seq_sfc(isq)
50 {sfc_iglevel1(^d,seq_sfc(isq))=seq_ig^d(isq) \}
54 deallocate(gsq_sfc,seq_sfc,seq_ig^d,in_domain)
59 use iso_fortran_env,
only : int64
61 integer(kind=4),
intent(in) :: ig^d,ndim
63 integer(kind=8) :: answer, lg^d
72 answer=ior(answer,ior(ishft(iand(lg1,ishft(1_int64,i)),i),&
73 ishft(iand(lg2,ishft(1_int64,i)),i+1)))
76 answer=ior(answer,ior(ishft(iand(lg1,ishft(1_int64,i)),2*i),ior(ishft(&
77 iand(lg2,ishft(1_int64,i)),2*i+1),ishft(iand(lg3,ishft(1_int64,i)),2*i+2))))
89 integer :: ig^D, Morton_no, isfc
98 if (morton_no/=
nleafs)
then
99 call mpistop(
"error in amr_Morton_order: Morton_no/=nleafs")
110 if (tree%node%leaf)
then
111 morton_no=morton_no+1
112 sfc(1,morton_no)=tree%node%igrid
113 sfc(2,morton_no)=tree%node%ipe
114 if (tree%node%active)
then
135 integer :: ipe, blocks_left, procs_left, num_blocks
138 {
#IFDEF EVOLVINGBOUNDARY
152 procs_left =
npe - ipe
153 num_blocks = ceiling(blocks_left / dble(procs_left))
155 blocks_left = blocks_left - num_blocks
159 {
#IFDEF EVOLVINGBOUNDARY
174 integer,
parameter :: wa=3, wp=1
189 integer :: ipe, Morton_no
190 integer :: Mtot, Mstop, Mcurr
192 integer :: nactive(0:npe-1),npassive(0:npe-1)
196 {
#IFDEF EVOLVINGBOUNDARY
210 mstop = (ipe+1)*int(mtot/npe)+min(ipe+1,mod(mtot,npe))
212 mcurr = mcurr + (wa*
sfc(3,morton_no)+wp*(1-
sfc(3,morton_no)))
214 if (
sfc(3,morton_no)==1)
then
215 nactive(ipe) = nactive(ipe) +1
217 npassive(ipe) = npassive(ipe) +1
220 if (mcurr >= mstop)
then
228 xmemory=dble(maxval(npassive+nactive))/&
229 dble(minval(npassive+nactive))
230 xload=dble(maxval(nactive))/&
231 dble(minval(nactive))
239 {
#IFDEF EVOLVINGBOUNDARY
subroutine mpistop(message)
Exit MPI-AMRVAC with an error message.
recursive subroutine get_morton_number(tree)
Module with basic grid data structures.
integer, dimension(:), allocatable, save sfc_to_igrid
Go from a Morton number to an igrid index (for a single processor)
integer, dimension(:), allocatable, save sfc_phybound
Space filling curve used for physical boundary blocks.
integer, dimension(:), allocatable, save morton_start
First Morton number per processor.
integer, save nleafs
Number of leaf block.
integer, dimension(:), allocatable, save igrid_to_sfc
Go from a grid index to Morton number (for a single processor)
integer, dimension(:,:), allocatable, save sfc_iglevel1
Space filling curve for level 1 grid. sfc_iglevel1(^D, MN) gives ig^D (the spatial index of the grid)
integer, dimension(:), allocatable, save morton_stop
Last Morton number per processor.
type(tree_node_ptr), dimension(:^d &), allocatable, save tree_root
Pointers to the coarse grid.
integer, dimension(:,:), allocatable, save sfc
Array to go from a Morton number to an igrid and processor index. Sfc(1:3, MN) contains [igrid,...
This module contains definitions of global parameters and variables and some generic functions/subrou...
double precision xload
Stores the memory and load imbalance, used in printlog.
integer, parameter ndim
Number of spatial dimensions for grid variables.
integer, dimension(:), allocatable ng
number of grid blocks in domain per dimension, in array over levels
integer mype
The rank of the current MPI task.
integer npe
The number of MPI tasks.
subroutine get_morton_range_active
subroutine get_morton_range
Set the Morton range for each processor.
integer(kind=8) function mortonencode(igD, ndim)
subroutine amr_morton_order
Construct Morton-order as a global recursive lexicographic ordering.
subroutine level1_morton_order
build Morton space filling curve for level 1 grid