charm_crd

charm_crd#

Module to work with the coordinates of evaluation points/cells:

defines charm_point and charm_cell structures,
allocates, initializes and frees charm_point and charm_cell,
computes charm_point for a few pre-defined grid types.

References:

Sneeuw, N. (1994) Global spherical harmonic analysis by least-squares and numerical quadrature methods in historical perspective. Geophysical Journal International 118:707-716
Driscoll, J. R., Healy, D. M. (1994) Computing Fourier transforms and convolutions on the 2-sphere. Advances in Applied Mathematics 15:202-250
Wieczorek, M. A., Meschede, M. (2018) SHTools: Tools for Working with Spherical Harmonics. Geochemistry, Geophysics, Geosystems 19:2574-2592

Note

This documentation is written for double precision version of CHarm.

Allocate, initialize and free the charm_point structure

These functions allocate, initialize and free the charm_point structure, which is designed to store evaluation points.

charm_point *charm_crd_point_malloc(int type, size_t nlat, size_t nlon)#

Allocates evaluation points of a given type with nlat latitudes, nlon longitudes and nlat spherical radii. If type is CHARM_CRD_POINT_GRID_GL, CHARM_CRD_POINT_GRID_DH1 or CHARM_CRD_POINT_GRID_DH2, nlat integration weights are allocated in addition. The memory of all array elements is uninitialized, so their values are undefined.

Valid values of type are constants CHARM_CRD_POINT_SCATTERED, CHARM_CRD_POINT_GRID, CHARM_CRD_POINT_GRID_GL, CHARM_CRD_POINT_GRID_DH1 and CHARM_CRD_POINT_GRID_DH2.

For charm_point *pnt returned by this function with type being CHARM_CRD_POINT_GRID, CHARM_CRD_POINT_GRID_GL, CHARM_CRD_POINT_GRID_DH1 or CHARM_CRD_POINT_GRID_DH2, the coordinates of the grid point on the i-th latitude parallel and the j-th meridian have to be accessed as
- pnt->lat[i] for the grid latitude,
- pnt->lon[j] for the grid longitude, and
- pnt->r[i] for the spherical radii
with i = 0, 1, ..., nlat - 1 and j = 0, 1, ..., nlon - 1.
For type set to CHARM_CRD_POINT_SCATTERED, the coordinates of the i-th point have to be accessed as
- pnt->lat[i] for the grid latitude,
- pnt->lon[i] for the grid longitude, and
- pnt->r[i] for the spherical radii
with i = 0, 1, ..., nlat - 1, nlat == nlon.

Warning

The structure returned must be deallocated by calling charm_crd_point_free(). The usual deallocation with free will lead to memory leaks.

Returns:: On success, returned is a pointer to the charm_point structure. On error, NULL is returned. Errors include memory allocation failures and incorrect input arguments (e.g., unsupported value of type).

charm_point *charm_crd_point_calloc(int type, size_t nlat, size_t nlon)#: The same as charm_crd_point_malloc() but all array elements are initialized to zero.

charm_point *charm_crd_point_init(int type, size_t nlat, size_t nlon, double *lat, double *lon, double *r)#

For a given point type, takes nlat spherical latitudes from the array pointed to by lat, nlon spherical longitudes from lon and nlat spherical radii from r (shallow copy).

Valid values of type are constants CHARM_CRD_POINT_SCATTERED, CHARM_CRD_POINT_GRID, CHARM_CRD_POINT_GRID_GL, CHARM_CRD_POINT_GRID_DH1 and CHARM_CRD_POINT_GRID_DH2.

If type is CHARM_CRD_POINT_SCATTERED,
- nlat must be equal to nlon and
- each of lat, lon and r must have access to nlat elements.
If type is CHARM_CRD_POINT_GRID, CHARM_CRD_POINT_GRID_GL, CHARM_CRD_POINT_GRID_DH1 or CHARM_CRD_POINT_GRID_DH2,
- lat must have access to nlat elements,
- lon must have access to nlon elements, and
- r must have access to nlat elements.

Note

See charm_crd_point_malloc() to find out how to access the coordinates stored in the returned charm_point structure.

Note

See charm_shc_init() for the rationale behind the shallow copies of input arrays lat, lon and r. The memory of these arrays must be treated in the same fashion as shown therein with the c and s arrays.

Warning

The structure returned must be deallocated using charm_crd_point_free(). The usual deallocation with free will lead to memory leaks.

Returns:: On success, returned is a pointer to the charm_point structure. On error, NULL is returned. Errors include memory allocation failures and incorrect input arguments (e.g., unsupported value of type).

void charm_crd_point_free(charm_point *pnt)#

Frees the memory associated with pnt. No operation is performed if pnt is NULL.

If pnt->owner is 1, the function releases all the memory that is associated with pnt, including that of its array members. If pnt->owner is 0, the arrays are not released from the memory, because they were not allocated by CHarm.

Quadrature point grids

These functions compute the Gauss-Legendre and Driscoll-Healy quadrature grids and the associated integration weights.

charm_point *charm_crd_point_gl(unsigned long nmax, double r)#

Computes the Gauss-Legendre grid associated with the harmonic degree nmax (Sneeuw, 1994) on the sphere with the radius r. The integration weights are computed on the unit sphere.

Loops are parallelized using OpenMP.

Note

The Gauss-Legendre grid is non-equiangular in latitude. The longitudinal step is constant.

Warning

The structure returned must be deallocated by calling charm_crd_point_free(). The usual deallocation with free will lead to memory leaks.

Returns:: On success, returned is a pointer to the charm_point structure. On error, NULL is returned. In addition to the memory allocation failure, the error may be due to the exceeded maximum number of iterations to compute the latitudes or due to the overflow problem. The latter may happen in single precision if nmax is larger than 7200 or so. The overflow check is not performed if the output from charm_misc_buildopt_isfinite() is zero.

void charm_crd_point_gl_shape(unsigned long nmax, size_t *nlat, size_t *nlon)#: Computes the number of latitudes nlat and the number of longitudes nlon of the Gauss-Legendre grid that corresponds to the maximum harmonic degree nmax.

charm_point *charm_crd_point_dh1(unsigned long nmax, double r)#

Computes the non-equiangular Driscoll-Healy grid (as defined by Driscoll and Healy, 1994) associated with the maximum spherical harmonic degree nmax on the sphere with the radius r. The integration weights are computed on the unit sphere.

Loops are parallelized using OpenMP.

Warning

The structure returned must be deallocated by calling charm_crd_point_free(). The usual deallocation with free will lead to memory leaks.

Returns:: On success, returned is a pointer to the charm_point structure. On error, NULL is returned.

void charm_crd_point_dh1_shape(unsigned long nmax, size_t *nlat, size_t *nlon)#: Computes the number of latitudes nlat and the number of longitudes nlon of the non-equiangular Driscoll-Healy grid that corresponds to the maximum harmonic degree nmax.

charm_point *charm_crd_point_dh2(unsigned long nmax, double r)#: The same as charm_crd_point_dh1() but for the equiangular modification of the Driscoll-Healy grid after Wieczorek and Meschede (2018).

void charm_crd_point_dh2_shape(unsigned long nmax, size_t *nlat, size_t *nlon)#: Computes the number of latitudes nlat and the number of longitudes nlon of the equiangular Driscoll-Healy grid that corresponds to the maximum harmonic degree nmax.

Allocate, initialize and free the charm_cell structure

These functions allocate, initialize and free the charm_cell structure, which is designed to store evaluation cells.

charm_cell *charm_crd_cell_malloc(int type, size_t nlat, size_t nlon)#

Allocates evaluation cells of a given type with nlat minimum and maximum latitudes, nlon minimum and maximum longitudes and nlat spherical radii. The memory of all array elements is uninitialized, so their values are undefined.

Valid values of type are constants CHARM_CRD_CELL_GRID and CHARM_CRD_CELL_SCATTERED.

For charm_cell *cell returned by this function with type being CHARM_CRD_CELL_GRID, the coordinates of the grid cell on the i-th latitude parallel and the j-th meridian have to be accessed as
- cell->latmin[i] for the minimum cell latitude,
- cell->latmax[i] for the maximum cell latitude,
- cell->lonmin[j] for the minimum cell longitude,
- cell->lonmax[j] for the maximum cell longitude,
- cell->r[i] for the spherical radii
with i = 0, 1, ..., nlat - 1 and j = 0, 1, ..., nlon - 1.
For type set to CHARM_CRD_CELL_SCATTERED, the coordinates of the i-th cell have to be accessed as
- cell->latmin[i] for the minimum cell latitude,
- cell->latmax[i] for the maximum cell latitude,
- cell->lonmin[i] for the minimum cell longitude,
- cell->lonmax[i] for the maximum cell longitude,
- cell->r[i] for the spherical radii
with i = 0, 1, ..., nlat - 1, nlat == nlon.

Warning

The structure returned must be deallocated by calling charm_crd_cell_free(). The usual deallocation with free will lead to memory leaks.

Returns:: On success, returned is a pointer to the charm_cell structure. On error, NULL is returned. Errors include memory allocation failures and incorrect input arguments (e.g., unsupported value of type).

charm_cell *charm_crd_cell_calloc(int type, size_t nlat, size_t nlon)#: The same as charm_crd_cell_malloc() but the memory of all array elements is initialized to zero.

charm_cell *charm_crd_cell_init(int type, size_t nlat, size_t nlon, double *latmin, double *latmax, double *lonmin, double *lonmax, double *r)#

For a given type, takes nlat minimum and maximum spherical latitudes from the arrays pointed to by latmin and latmax, nlon minimum and maximum spherical longitudes from lonmin and lonmax and nlat spherical radii from r, respectively (shallow copy).

Accepted values of type are CHARM_CRD_CELL_GRID and CHARM_CRD_CELL_SCATTERED.

If type is CHARM_CRD_CELL_SCATTERED,
- nlat must be equal to nlon and
- each of latmin, latmax, lonmin, lonmax and r must have access to nlat elements.
If type is CHARM_CRD_CELL_GRID,
- latmin and latmax must each have access to nlat elements,
- lonmin and lonmax must each have access to nlon elements,
- r must have access to nlat elements.

Note

See charm_crd_cell_malloc() to find out how to access the coordinates stored in the returned charm_cell structure.

Note

See charm_shc_init() for the rationale behind the shallow copies of input arrays latmin, latmax, lonmin, lonmax and r. The memory of these arrays must be treated in the same fashion as shown therein with the c and s arrays.

Warning

The structure created by this function must be deallocated using charm_crd_cell_free(). The usual deallocation with free will lead to memory leaks.

Returns:: On success, returned is a pointer to the charm_cell structure. On error, NULL is returned. Errors include memory allocation failures and incorrect input arguments (e.g., unsupported value of type).

void charm_crd_cell_free(charm_cell *cell)#

Frees the memory associated with cell. No operation is performed if cell is NULL.

If cell->owner is 1, the function releases all the memory that is associated with cell, including that of its array members. If cell->owner is 0, the arrays are not released from the memory, because they were not allocated by CHarm.

Enums

Defines supported values of charm_point.type to CHARM_CRD_POINT_SCATTERED, CHARM_CRD_POINT_GRID, CHARM_CRD_POINT_GRID_GL, CHARM_CRD_POINT_GRID_DH1, CHARM_CRD_POINT_GRID_DH2. For charm_cell.type, valid constants are CHARM_CRD_CELL_SCATTERED and CHARM_CRD_CELL_GRID.

Values:

enumerator CHARM_CRD_CELL_GRID = -2#: Custom grid of cells.

enumerator CHARM_CRD_CELL_SCATTERED#: Scattered cells.

enumerator CHARM_CRD_POINT_SCATTERED = 1#: Scattered points.

enumerator CHARM_CRD_POINT_GRID#: Custom point grid.

enumerator CHARM_CRD_POINT_GRID_GL#: Gauss-Legendre point grid.

enumerator CHARM_CRD_POINT_GRID_DH1#: Driscoll-Healy point grid as defined by Driscoll and Healy (1994) (non-equiangular).

enumerator CHARM_CRD_POINT_GRID_DH2#: The equiangular modification of the Driscoll-Healy point grid after Wieczorek and Meschede (2018).

struct charm_point#

Structure to store evaluation points.

If CHarm is compiled with the MPI support, the data of the structure can be distributed across MPI processes. This is useful if your dataset is so large that it would not fit into the RAM of a single computing node. The members of the structure are therefore divided into two groups.

Members common to both non-distributed and distributed structures. These are always accessible.
Members specific to distributed structures only. These are accessible only when CHarm was compiled with the MPI support.

Members common to non-distributed and distributed structures

int type#

Organization of points in the structure.

The following options are available:

CHARM_CRD_POINT_SCATTERED for user-defined scattered points,
CHARM_CRD_POINT_GRID for a user-defined grid,
CHARM_CRD_POINT_GRID_GL for the Gauss-Legendre grid (Sneeuw, 1994) with nmax + 1 latitudes and 2 * nmax + 2 longitudes (nmax is maximum harmonic degree),
CHARM_CRD_POINT_GRID_DH1 for the non-equiangular Driscoll-Healy grid as defined by Driscoll and Healy (1994) with 2 * nmax + 2 latitudes and 2 * nmax + 2 longitudes, and
CHARM_CRD_POINT_GRID_DH2 for an equiangular modification of the Driscoll-Healy grid (Wieczorek and Meschede, 2018) with 2 * nmax + 2 latitudes and 4 * nmax + 4 longitudes.

size_t nlat#

Total number of points in the latitudinal direction. If 0, no points are stored in the structure.

If the structure is distributed, this member represents the sum of charm_point.local_nlat across all MPI processes within charm_point.comm.

size_t nlon#

Total number of points in the longitudinal direction. If 0, no points are stored in the structure.

If the structure is distributed and charm_point.type == CHARM_CRD_POINT_SCATTERED, then this member represents the sum of charm_point.local_nlon across all MPI processes within charm_point.comm. This means that each MPI process may hold different number of longitudes, hence scattered points.

If the structure is distributed and charm_point.type is either CHARM_CRD_POINT_GRID, CHARM_CRD_POINT_GRID_GL, CHARM_CRD_POINT_GRID_DH1 or CHARM_CRD_POINT_GRID_DH2, then this member has always the same value as charm_point.local_nlon. This means that all point grids must always be distributed across MPI processes in latitudinal chunks, that is, each sub-grid must consist of the same number of longitudes, so that charm_point.local_nlon and charm_point.nlon are equal.

size_t npoint#

Total number of points.

For scattered points ( charm_point.type == CHARM_CRD_POINT_SCATTERED), equals to charm_point.nlat and to charm_point.nlon.
For grids (charm_point.type is any of CHARM_CRD_POINT_GRID*), equals to charm_point.nlat * charm_point.nlon.

If the structure is distributed, this member represents the sum of charm_point.local_npoint across all MPI processes within charm_point.comm.

double *lat#

Pointer to an array of latitudes in radians.

For non-distributed structures, charm_point.nlat array elements are associated with the pointer. For distributed structures, charm_point.local_nlat elements are associated.

double *lon#

Pointer to an array of longitudes in radians.

For non-distributed structures, charm_point.nlon array elements are associated with the pointer. For distributed structures, charm_point.local_nlon elements are associated.

double *r#

Pointer to an array of spherical radii in metres.

The same number of array elements is associated with this pointer as in the case of charm_point.lat.

double *w#

Pointer to an array of integration weights on the unit sphere. The pointer is used (i.e. is not NULL) only for charm_point structures returned by the charm_crd_point_gl(), charm_crd_point_dh1() and charm_crd_point_dh2() functions, that is, if charm_point.type is CHARM_CRD_POINT_GRID_GL, CHARM_CRD_POINT_GRID_DH1, or CHARM_CRD_POINT_GRID_DH2. Otherwise, the pointer is NULL and is never used.

The same number of array elements is associated with this pointer as in the case of charm_point.lat.

_Bool owner#

If 1, the memory associated with charm_point.lat, charm_point.lon, charm_point.r and charm_point.w was allocated by CHarm, so charm_crd_point_free() deallocates it.
If 0, the memory associated with these arrays was allocated outside CHarm, so charm_crd_point_free() does not deallocate it (the user allocated this memory outside CHarm, so the user should free the memory outside CHarm as well).

_Bool distributed#

If 0, the structure is not distributed, meaning that each process (be it an OpenMP or MPI process) accesses the same data through charm_point.lat, charm_point.lon, charm_point.r and charm_point.w.
If 1, the structure is distributed across MPI processes, that is, each MPI process accesses different data through charm_point.lat, charm_point.lon, charm_point.r and charm_point.w.

Note

If CHarm was compiled without the MPI support (default), charm_point.distributed is always 0.

Members available to distributed structures only

Note

The members that follow are available only when CHarm is compiled with the MPI support (--enable-mpi, refer to charm_mpi for further details).

size_t local_nlat#: Number of points in the latitudinal direction locally accessible to an MPI process. If 0, no points are stored locally in the structure.

size_t local_nlon#: Number of points in the longitudinal direction locally accessible to an MPI process. If 0, no points are stored locally in the structure.

size_t local_npoint#

Total number of points locally accessible to an MPI process.

It is computed in the same way as charm_point.npoint, but using charm_point.local_nlat and charm_point.local_nlon instead of charm_point.nlat and charm_point.nlon.

size_t local_0_start#

Index, at which the local portions of charm_point.lat, charm_point.r and charm_point.w start in the charm_point structure. Any value is accepted if charm_point.local_nlat is 0.

Assume that charm_point *pnt with pnt->distributed == 0 holds all your data points. Now assume that you distribute the same data in chunks across MPI processes, obtaining some charm_point *pnt_dist with pnt_dist->distributed == 1. For an MPI process, pnt_dist.local_0_start is the index of pnt_dist->lat[0] in the pnt->lat array. This means that pnt->lat[pnt_dist->local_0_start + idx] == pnt_dist->[idx], where idx < pnt_dist->local_nlat. If the grid is symmetric with respect to the equator, the index in pnt_dist.local_0_start is always divided by 2.

Note

Sounds complicated? It isn’t! Have a look at the charm_mpi module and the cookbook.

MPI_Comm comm#: MPI communicator defining the group of MPI processes to distribute the structure across.

struct charm_cell#

Structure to store evaluation cells.

This structure cannot be distributed across MPI processes.

Public Members

int type#

Organization of cells in the structure.

The following options are available:

CHARM_CRD_CELL_GRID for a user-defined grid of cells,
CHARM_CRD_CELL_SCATTERED for a user-defined set of scattered cells.

size_t nlat#: Total number of cells in the latitudinal direction.

size_t nlon#: Total number of cells in the longitudinal direction.

size_t ncell#

Total number of cells.

For scattered cells (charm_cell.type set to CHARM_CRD_CELL_SCATTERED), charm_cell.ncell equals to charm_cell.nlat and to charm_cell.nlon.
For grids (charm_cell.type set to CHARM_CRD_CELL_GRID), charm_cell.ncell is given as charm_cell.nlat * charm_cell.nlon.

double *latmin#: Pointer to a charm_cell.nlat array elements representing minimum cell latitudes in radians.

double *latmax#: Pointer to a charm_cell.nlat array elements representing maximum cell latitudes in radians.

double *lonmin#: Pointer to a charm_cell.nlon array elements representing minimum cell longitudes in radians.

double *lonmax#: Pointer to a charm_cell.nlon array elements representing maximum cell longitudes in radians.

double *r#: Pointer to a charm_cell.nlat array elements representing spherical radii in metres.

_Bool owner#: If 1, the memory associated with charm_cell.latmin, charm_cell.latmax, charm_cell.lonmin, charm_cell.lonmax and charm_cell.r was allocated by CHarm, so charm_crd_cell_free() deallocates it. If charm_cell.owner is 0, the memory associated with these arrays was allocated outside CHarm, so charm_crd_cell_free() does not deallocate it (the user allocated this memory outside CHarm, so the user should free the memory outside CHarm as well).

charm_crd

Contents

charm_crd#