SHTns  3.4.6
Functions
Initialization functions.

Functions

long nlm_calc (long lmax, long mmax, long mres)
 compute number of spherical harmonics modes (l,m) to represent a real-valued spatial field for given size parameters. Does not require any previous setup. More...
 
long nlm_cplx_calc (long lmax, long mmax, long mres)
 compute number of spherical harmonics modes (l,m) to represent a complex-valued spatial field, for given size parameters. Does not require any previous setup. More...
 
shtns_cfg shtns_create (int lmax, int mmax, int mres, enum shtns_norm norm)
 Defines the sizes of the spectral description. Use for advanced initialization. More...
 
shtns_cfg shtns_create_with_grid (shtns_cfg base, int mmax, int nofft)
 Copy a given config but allow a different (smaller) mmax and the possibility to enable/disable fft (beta). More...
 
void shtns_unset_grid (shtns_cfg shtns)
 release all resources allocated by a grid. More...
 
void shtns_destroy (shtns_cfg shtns)
 release all resources allocated by a given shtns_cfg. More...
 
void shtns_reset ()
 clear all allocated memory (hopefully) and go back to 0 state. More...
 
void * shtns_malloc (size_t size)
 alloc appropriate memory (pinned for gpu, aligned for avx, ...). Use shtns_free to free it.
 
void shtns_free (void *p)
 free memory allocated with shtns_malloc
 
int shtns_set_grid_auto (shtns_cfg shtns, enum shtns_type flags, double eps, int nl_order, int *nlat, int *nphi)
 Precompute everything and choose the optimal nlat and nphi for a given non-linear order. More...
 
int shtns_set_grid (shtns_cfg shtns, enum shtns_type flags, double eps, int nlat, int nphi)
 Precompute everything for a given spatial grid. Use for advanced initialization, after shtns_create. More...
 
shtns_cfg shtns_init (enum shtns_type flags, int lmax, int mmax, int mres, int nlat, int nphi)
 
void shtns_robert_form (shtns_cfg shtns, int robert)
 set the use of Robert form. If robert != 0, the vector synthesis returns a field multiplied by sin(theta), while the analysis divides by sin(theta) before the transform.
 
int shtns_use_threads (int num_threads)
 Enables OpenMP parallel transforms, if available (see Compiling and installing SHTns). More...
 
int shtns_use_gpu (int device_id)
 Enables parallel transforms on selected GPU device, if available (see Compiling and installing SHTns). More...
 
int shtns_gauss_wts (shtns_cfg shtns, double *wts)
 fill the given array with Gauss weights. More...
 

Detailed Description

Function Documentation

◆ nlm_calc()

long nlm_calc ( long  lmax,
long  mmax,
long  mres 
)

compute number of spherical harmonics modes (l,m) to represent a real-valued spatial field for given size parameters. Does not require any previous setup.

return (mmax+1)*(lmax+1) - mres*(mmax*(mmax+1))/2;

*‍/

◆ nlm_cplx_calc()

long nlm_cplx_calc ( long  lmax,
long  mmax,
long  mres 
)

compute number of spherical harmonics modes (l,m) to represent a complex-valued spatial field, for given size parameters. Does not require any previous setup.

return mmax*mmax*mres + (2*mmax+1)*(lmax-mmax*mres+1);

*‍/

◆ shtns_create()

shtns_cfg shtns_create ( int  lmax,
int  mmax,
int  mres,
enum shtns_norm  norm 
)

Defines the sizes of the spectral description. Use for advanced initialization.

This sets the description of spherical harmonic coefficients. It tells SHTns how to interpret spherical harmonic coefficient arrays, and it sets usefull arrays. Returns the configuration to be passed to subsequent transform functions, which is basicaly a pointer to a shtns_info struct.

Parameters
lmax: maximum SH degree that we want to describe.
mmax: number of azimutal wave numbers.
mres: 2.pi/mres is the azimutal periodicity. mmax*mres is the maximum SH order.
norm: define the normalization of the spherical harmonics (shtns_norm)

Condon-Shortley phase (-1)^m is used by default.

renormalization of m>0.

◆ shtns_create_with_grid()

shtns_cfg shtns_create_with_grid ( shtns_cfg  base,
int  mmax,
int  nofft 
)

Copy a given config but allow a different (smaller) mmax and the possibility to enable/disable fft (beta).

Copy a given config but allow a different (smaller) mmax and the possibility to enable/disable fft.

◆ shtns_destroy()

void shtns_destroy ( shtns_cfg  shtns)

release all resources allocated by a given shtns_cfg.

free memory of given config, which cannot be used afterwards.

◆ shtns_gauss_wts()

int shtns_gauss_wts ( shtns_cfg  shtns,
double *  wts 
)

fill the given array with Gauss weights.

fill the given array with Gauss weights. returns the number of weights written (0 if not a Gauss grid).

returns the number of weights written, which may be zero if the grid is not a Gauss grid.

◆ shtns_init()

shtns_cfg shtns_init ( enum shtns_type  flags,
int  lmax,
int  mmax,
int  mres,
int  nlat,
int  nphi 
)

Simple initialization of Spherical Harmonic transforms (backward and forward, vector and scalar, ...) of given size. This function sets all global variables by calling shtns_create followed by shtns_set_grid, with the default normalization and the default polar optimization (see sht_config.h). Returns the configuration to be passed to subsequent transform functions, which is basicaly a pointer to a shtns_info struct.

Parameters
lmax: maximum SH degree that we want to describe.
mmax: number of azimutal wave numbers.
mres: 2.pi/mres is the azimutal periodicity. mmax*mres is the maximum SH order.
nlat,nphi: respectively the number of latitudinal and longitudinal grid points.
flagsallows to choose the type of transform (see shtns_type) and the spatial data layout (see Spatial data layouts and grids used by SHTns)

◆ shtns_reset()

void shtns_reset ( void  )

clear all allocated memory (hopefully) and go back to 0 state.

destroy all configs, free memory, and go back to initial state.

◆ shtns_set_grid()

int shtns_set_grid ( shtns_cfg  shtns,
enum shtns_type  flags,
double  eps,
int  nlat,
int  nphi 
)

Precompute everything for a given spatial grid. Use for advanced initialization, after shtns_create.

Initialization of Spherical Harmonic transforms (backward and forward, vector and scalar, ...) of given size. This function must be called after shtns_create and before any SH transform. and sets all global variables. returns the required number of doubles to be allocated for a spatial field.

Parameters
shtnsis the config created by shtns_create for which the grid will be set.
flagsallows to choose the type of transform (see shtns_type) and the spatial data layout (see Spatial data layouts and grids used by SHTns)
epspolar optimization threshold : polar values of Legendre Polynomials below that threshold are neglected (for high m), leading to increased performance (a few percents) 0 = no polar optimization; 1.e-14 = VERY safe; 1.e-10 = safe; 1.e-6 = aggresive, but still good accuracy.
nlat,nphirespectively the number of latitudinal and longitudinal grid points.

◆ shtns_set_grid_auto()

int shtns_set_grid_auto ( shtns_cfg  shtns,
enum shtns_type  flags,
double  eps,
int  nl_order,
int *  nlat,
int *  nphi 
)

Precompute everything and choose the optimal nlat and nphi for a given non-linear order.

Initialization of Spherical Harmonic transforms (backward and forward, vector and scalar, ...) of given size. This function must be called after shtns_create and before any SH transform. and sets all global variables and internal data. returns the required number of doubles to be allocated for a spatial field.

Parameters
shtnsis the config created by shtns_create for which the grid will be set.
nlat,nphipointers to the number of latitudinal and longitudinal grid points respectively. If 0, they are set to optimal values.
nl_orderdefines the maximum SH degree to be resolved by analysis : lmax_analysis = lmax*nl_order. It is used to set an optimal and anti-aliasing nlat. If 0, the default SHT_DEFAULT_NL_ORDER is used.
flagsallows to choose the type of transform (see shtns_type) and the spatial data layout (see Spatial data layouts and grids used by SHTns)
epspolar optimization threshold : polar values of Legendre Polynomials below that threshold are neglected (for high m), leading to increased performance (a few percents) 0 = no polar optimization; 1.e-14 = VERY safe; 1.e-10 = safe; 1.e-6 = aggresive, but still good accuracy.

◆ shtns_unset_grid()

void shtns_unset_grid ( shtns_cfg  shtns)

release all resources allocated by a grid.

unset the grid.

◆ shtns_use_gpu()

int shtns_use_gpu ( int  device_id)

Enables parallel transforms on selected GPU device, if available (see Compiling and installing SHTns).

Selects the gpu device (device_id % Num_devices). Must be called BEFORE any initialization. Internally calls cudaSetDevice(). Returns the actual device or -1 when no device found.

Call BEFORE any initialization of shtns to select a GPU device. Returns the actual device id used, or -1 if no device found.

  • If device_id >= 0, try to use device with number device_id % device_count.
  • If device_d < 0, do not try to use GPU. WARNING: Calls cudaSetDevice() internally, so it changes the current device for the entire host thread.

◆ shtns_use_threads()

int shtns_use_threads ( int  num_threads)

Enables OpenMP parallel transforms, if available (see Compiling and installing SHTns).

Enables multi-thread transform using OpenMP with num_threads (if available). Returns number of threads that will be used.

Call before any initialization of shtns to use multiple threads. Returns the actual number of threads.

  • If num_threads > 0, specifies the maximum number of threads that should be used.
  • If num_threads <= 0, maximum number of threads is automatically set to the number of processors.
  • If num_threads == 1, openmp will be disabled.