Octopus
centroids_oct_m Module Reference

Data Types

type  centroids_t
 Encapsulate centroid points and their indexing across a domain-decomposed mesh. More...
 

Functions/Subroutines

subroutine centroids_init (this, mesh, positions)
 Initialise a centroids_t instance. More...
 
pure integer function centroids_get_n_centroids (this)
 Number of centroids in local domain. More...
 
pure integer function centroids_get_n_centroids_global (this)
 Total number of centroids in system. More...
 
pure integer(int64) function centroids_global_centroid_to_global_mesh_index (this, ic)
 Given centroid index ic, return the global mesh index ipg. More...
 
pure integer(int64) function, dimension(:), allocatable centroids_global_centroid_to_global_mesh_indices (this)
 Return global mesh indices for all centroids in this domain. More...
 
pure integer function centroids_get_local_mesh_index (this, ic)
 Get the mesh index ip of centroid ic in this domain. More...
 
pure integer function, dimension(:), allocatable centroids_get_local_mesh_indices (this)
 Return mesh indices ip for all centroids in this domain. More...
 
pure integer function centroids_local_index_to_global_index (this, ic)
 Get the global centroid index, given centroid ic this domain. More...
 
pure integer function, dimension(:), allocatable centroids_local_index_to_global_indices (this)
 Return global centroid indices for all centroids in this domain. More...
 
pure real(real64) function, dimension(this%ndim) centroids_get_local_position (this, mesh, ic)
 Get the position of the centroid ic. More...
 
pure real(real64) function, dimension(:, :), allocatable centroids_get_all_positions (this)
 Return a copy of all centroid positions. More...
 
subroutine centroids_finalize (this)
 Finalize a centroids instance. More...
 
subroutine exact_centroids_init_index_maps (this, mesh)
 Initialise index maps with centroid points that are commensurate with the grid. More...
 

Function/Subroutine Documentation

◆ centroids_init()

subroutine centroids_oct_m::centroids_init ( class(centroids_t), intent(out)  this,
class(mesh_t), intent(in)  mesh,
real(real64), dimension(:, :), intent(in)  positions 
)
private

Initialise a centroids_t instance.

Centroid positions must correspond to grid points.

Parameters
[out]thisCentroid instance
[in]meshMesh or grid instance
[in]positionsAll centroid cartesian positions

Definition at line 187 of file centroids.F90.

◆ centroids_get_n_centroids()

pure integer function centroids_oct_m::centroids_get_n_centroids ( class(centroids_t), intent(in)  this)
private

Number of centroids in local domain.

Definition at line 205 of file centroids.F90.

◆ centroids_get_n_centroids_global()

pure integer function centroids_oct_m::centroids_get_n_centroids_global ( class(centroids_t), intent(in)  this)
private

Total number of centroids in system.

Definition at line 215 of file centroids.F90.

◆ centroids_global_centroid_to_global_mesh_index()

pure integer(int64) function centroids_oct_m::centroids_global_centroid_to_global_mesh_index ( class(centroids_t), intent(in)  this,
integer, intent(in)  ic 
)
private

Given centroid index ic, return the global mesh index ipg.

Note
One could alternatively implement this with:
ip = this%ic_to_ip(ic)
ipg = mesh_local2global(mesh, ip)
which would remove the need to store icg_to_ipg, but this would require mesh to be an argument of the function.

Definition at line 234 of file centroids.F90.

◆ centroids_global_centroid_to_global_mesh_indices()

pure integer(int64) function, dimension(:), allocatable centroids_oct_m::centroids_global_centroid_to_global_mesh_indices ( class(centroids_t), intent(in)  this)
private

Return global mesh indices for all centroids in this domain.

Definition at line 247 of file centroids.F90.

◆ centroids_get_local_mesh_index()

pure integer function centroids_oct_m::centroids_get_local_mesh_index ( class(centroids_t), intent(in)  this,
integer, intent(in)  ic 
)
private

Get the mesh index ip of centroid ic in this domain.

Definition at line 263 of file centroids.F90.

◆ centroids_get_local_mesh_indices()

pure integer function, dimension(:), allocatable centroids_oct_m::centroids_get_local_mesh_indices ( class(centroids_t), intent(in)  this)
private

Return mesh indices ip for all centroids in this domain.

Definition at line 274 of file centroids.F90.

◆ centroids_local_index_to_global_index()

pure integer function centroids_oct_m::centroids_local_index_to_global_index ( class(centroids_t), intent(in)  this,
integer, intent(in)  ic 
)
private

Get the global centroid index, given centroid ic this domain.

Definition at line 284 of file centroids.F90.

◆ centroids_local_index_to_global_indices()

pure integer function, dimension(:), allocatable centroids_oct_m::centroids_local_index_to_global_indices ( class(centroids_t), intent(in)  this)
private

Return global centroid indices for all centroids in this domain.

Definition at line 295 of file centroids.F90.

◆ centroids_get_local_position()

pure real(real64) function, dimension(this%ndim) centroids_oct_m::centroids_get_local_position ( class(centroids_t), intent(in)  this,
class(mesh_t), intent(in)  mesh,
integer, intent(in)  ic 
)
private

Get the position of the centroid ic.

Definition at line 305 of file centroids.F90.

◆ centroids_get_all_positions()

pure real(real64) function, dimension(:, :), allocatable centroids_oct_m::centroids_get_all_positions ( class(centroids_t), intent(in)  this)
private

Return a copy of all centroid positions.

Definition at line 318 of file centroids.F90.

◆ centroids_finalize()

subroutine centroids_oct_m::centroids_finalize ( type(centroids_t), intent(inout)  this)
private

Finalize a centroids instance.

Parameters
[in,out]thisCentroid instance

Definition at line 328 of file centroids.F90.

◆ exact_centroids_init_index_maps()

subroutine centroids_oct_m::exact_centroids_init_index_maps ( type(centroids_t), intent(inout)  this,
class(mesh_t), intent(in)  mesh 
)
private

Initialise index maps with centroid points that are commensurate with the grid.

This routine breaks single responsibility, but it makes sense to reuse the the various indices as they are computed. This routine computes:

  • Convert centroid positions to global mesh index, ipg
  • Set number of centroids in the domain of process
  • Map global centroid index icg to centroid index in domain process
  • Map centroid index ic in domain process, to mesh index ip in the same domain
    Parameters
    [in]meshMesh or grid instance

Definition at line 349 of file centroids.F90.