Mesh implementation

The data structures:

Definition of basis_set_abst_t

The mesh descriptor is defined in grid/mesh.F90:

  type, extends(basis_set_abst_t) :: mesh_t
    ! Components are public by default
    class(box_t),               pointer :: box           !< simulation box; polymorphic pointer
    class(coordinate_system_t), pointer :: coord_system  !< the underlying coordinate system (affine, cubic, curvilinear, etc.)
    type(index_t)                       :: idx           !< the indexing scheme.
    logical :: use_curvilinear

    real(real64), allocatable :: spacing(:)         !< the (constant) spacing between the points

    ! When running serially, the local number of points is
    ! equal to the global number of points.
    ! Otherwise, the next two are different on each node.
    integer      :: np                  !< Local number of points in mesh
    integer      :: np_part             !< Local points plus ghost points plus boundary points.
    integer(int64)  :: np_global        !< Global number of points in mesh.
    integer(int64)  :: np_part_global   !< Global number of inner points and boundary points.
    logical         :: parallel_in_domains !< will I run parallel in domains?
    type(mpi_grp_t) :: mpi_grp             !< the mpi group describing parallelization in domains
    type(par_vec_t) :: pv                  !< describes parallel vectors defined on the mesh.
    type(partition_t) :: partition         !< describes how the inner points are assigned to the domains

    real(real64),   allocatable :: x(:,:)            !< The (local) \b points
    real(real64)                :: volume_element    !< The global volume element.
    real(real64),   allocatable :: vol_pp(:)         !< Element of volume for curvilinear coordinates.

    logical            :: masked_periodic_boundaries !< flag whether to use the mask defined below
    character(len=256) :: periodic_boundary_mask     !< mask for which the periodic boundary conditions
    !!                                                  are replaced by zero boundary conditions
  contains
    procedure :: end => mesh_end                !< @copydoc mesh_oct_m::mesh_end
    procedure :: init => mesh_init              !< @copydoc mesh_oct_m::mesh_init
    procedure :: write_info => mesh_write_info  !< @copydoc mesh_oct_m::mesh_write_info
    procedure :: dmesh_allreduce_0, zmesh_allreduce_0, imesh_allreduce_0, lmesh_allreduce_0
    procedure :: dmesh_allreduce_1, zmesh_allreduce_1, imesh_allreduce_1, lmesh_allreduce_1
    procedure :: dmesh_allreduce_2, zmesh_allreduce_2, imesh_allreduce_2, lmesh_allreduce_2
    procedure :: dmesh_allreduce_3, zmesh_allreduce_3, imesh_allreduce_3, lmesh_allreduce_3
    procedure :: dmesh_allreduce_4, zmesh_allreduce_4, imesh_allreduce_4, lmesh_allreduce_4
    procedure :: dmesh_allreduce_5, zmesh_allreduce_5, imesh_allreduce_5, lmesh_allreduce_5
    generic :: allreduce => dmesh_allreduce_0, zmesh_allreduce_0, imesh_allreduce_0, lmesh_allreduce_0
    generic :: allreduce => dmesh_allreduce_1, zmesh_allreduce_1, imesh_allreduce_1, lmesh_allreduce_1
    generic :: allreduce => dmesh_allreduce_2, zmesh_allreduce_2, imesh_allreduce_2, lmesh_allreduce_2
    generic :: allreduce => dmesh_allreduce_3, zmesh_allreduce_3, imesh_allreduce_3, lmesh_allreduce_3
    generic :: allreduce => dmesh_allreduce_4, zmesh_allreduce_4, imesh_allreduce_4, lmesh_allreduce_4
    generic :: allreduce => dmesh_allreduce_5, zmesh_allreduce_5, imesh_allreduce_5, lmesh_allreduce_5
  end type mesh_t

which uses the structure index_t, containing the range and the mapping arrays:

  type index_t
    ! Components are public by default
    integer                 :: dim        !< the dimension
    integer, allocatable    :: nr(:,:)    !< dimensions of the box where the points are contained
    integer, allocatable    :: ll(:)      !< literally nr(2,:) - nr(1,:) + 1 - 2*enlarge(:)
    integer, allocatable    :: enlarge(:) !< number of points to add for boundary conditions
    integer(int64)          :: checksum
    integer(int64), pointer :: grid_to_spatial_global(:) !< map: global grid index -> spatial index
    integer(int64), pointer :: spatial_to_grid_global(:) !< inverse map: spatial index -> global grid index
    integer                 :: type               !< index type: IDX_CUBIC or IDX_HILBERT
    integer                 :: bits               !< bits per dimension for Hilbert index
    integer, allocatable    :: offset(:)          !< offset for getting the indices from the spatial index
    integer, allocatable    :: stride(:)
    type(MPI_Win)           :: window_grid_to_spatial !< handle to shared memory window for map
    type(MPI_Win)           :: window_spatial_to_grid !< handle to shared memory window for inverse map
  end type index_t

About lxyz and lxyz_inv maps:

The direct map:

lxyz(d, i) = R(d)i

where d denotes the dimension (i.e. x, y or z) and i is the index of the point.

The inverse map:

lxyz_inv( Rxi , Ryi , Rzi ) = i

The points defined by lxyz define a rectangular box of d dimensions. The real mesh vectors are related to Ri by multiplication with spacing and possibly distortion for curvilinear coordinates.

Note that this index array is the same in all parallel domains, and relates the integer coordinates to the global index of a given point.

  type mesh_cube_map_t
    ! Components are public by default
    integer              :: nmap      !< The number of maps
    integer, allocatable :: map(:, :)
    type(accel_mem_t)    :: map_buffer
  end type mesh_cube_map_t

Note on packed states:

Note on curvilinear meshes:

The mesh::x(:,:) array always contains a regular mesh, which gets ‘distorded’ to a curvilinear mesh by additional function calls.

Domain decomposition