API

Warning

Page under construction.

ecTrans API

General notes

Note

ecTrans is a legacy code with an accumulated 30 years of history. Over this time certain features enabled through optional arguments will have fallen out of use. We are currently reviewing all options to identify those that can be safely deleted, but this takes time. In the mean time, we have tagged below all options we deem to be "potentially deprecatable".

Variable names

ecTrans in principle follows the coding standard and conventions outlined in the IFS Documentation - Part VI: Technical and Computational Procedures section 1.5. Following these standards, all variable names must begin with a one- or two-character prefix denoting their scope (module level, dummy argument, local variables, loop index, or parameter) and type. These are outlined in Table 1.2. Dummy variables have the following prefixes:

  • K - integer
  • P - real (single or double precision)
  • LD - logical
  • CD - character
  • YD - derived type

KIND parameters

As with the IFS, integer and real variables in ecTrans always have an explicit KIND specification. These are defined in the PARKIND1 module which is part of the Fiat library (a dependency of ecTrans). To understand the subroutines described here, only two must be considered:

  • INTEGER, PARAMETER :: JPIM = SELECTED_INT_KIND(9) (i.e. 4-byte integer)
  • INTEGER, PARAMETER :: JPRD = SELECTED_REAL_KIND(13,300) (i.e. 8-byte float)

SETUP_TRANS0

Signature

SUBROUTINE SETUP_TRANS0(KOUT, KERR, KPRINTLEV, KMAX_RESOL, KPROMATR, KPRGPNS, KPRGPEW, &
  &                     KPRTRW, KCOMBFLEN, LDMPOFF, LDSYNC_TRANS, KTRANS_SYNC_LEVEL, &
  &                     LDEQ_REGIONS, K_REGIONS_NS, K_REGIONS_EW, K_REGIONS, PRAD, LDALLOPERM, &
  &                     KOPT_MEMORY_TR)

Purpose

This subroutine initialises all resolution-agnostic aspects of the ecTrans library. It must be called before any other ecTrans subroutines.

OPTIONAL, INTENT(IN) arguments

  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(IN) :: KOUT
    Fortran unit number for listing output.
    Default: 6 (i.e. STDOUT)
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(IN) :: KERR
    Unit number for error messages.
    Default: 0 (i.e. STDERR)
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(IN) :: KPRINTLEV
    Verbosity of standard output (0: output disabled, 1: normal, 2: debug).
    Default: 0
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(IN) :: KMAX_RESOL
    Maximum number of different resolutions for which spectral transforms will be performed. This is required in advance so ecTrans knows how to correctly size the work arrays.
    Default: 1
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(IN) :: KPROMATR
    Batch size for splitting of vertical field set. This allows one to transform one batch of vertical
    levels at a time rather than transforming all in one go. A value of 0 disables this feature.
    Default: 0
    Potentially deprecatable.
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(IN) :: KPRGPNS
    Splitting level in North-South direction in grid-point space.
    Default: 1
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(IN) :: KPRGPEW
    Splitting level in East-West direction in grid-point space.
    Default: 1
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(IN) :: KPRTRW
    Splitting level in wave direction in spectral space.
    Default: 1
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(IN) :: KCOMBFLEN
    Size of communication buffer [1800000 (8bytes) ].
    This argument is deprecated and will be removed in the next release.
  • LOGICAL, OPTIONAL, INTENT(IN) :: LDMPOFF
    Switch off message passing.
    Default: .FALSE.
    Potentially deprecatable.
  • LOGICAL, OPTIONAL, INTENT(IN) :: LDSYNC_TRANS
    Switch to activate barriers in M->L and L->M transposition routines.
    Default: .FALSE.
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(IN) :: KTRANS_SYNC_LEVEL
    Control parameter for synchronization and blocking in communication routines.
    Default: 0
  • LOGICAL, OPTIONAL, INTENT(IN) :: LDEQ_REGIONS
    Enable EQ_REGIONS partitioning mode.
    Default: .FALSE.
  • LOGICAL, OPTIONAL, INTENT(IN) :: LDALLOPERM
    Allocate certain arrays permanently.
    Default: .FALSE.
  • REAL(KIND=JPRD), OPTIONAL, INTENT(IN) :: PRAD
    Radius of the planet (Earth) in metres.
    Default: 6371229.0
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(IN) :: KOPT_MEMORY_TR
    Memory strategy for gridpoint transpositions (0: heap, 1: stack).
    Default: 0

OPTIONAL, INTENT(OUT) arguments

  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(OUT) :: K_REGIONS(:)
    Number of regions returned by EQ_REGIONS algorithm.
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(OUT) :: K_REGIONS_NS
    Maximum number of North-South partitions, as determined by EQ_REGIONS algorithm.
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(OUT) :: K_REGIONS_EW
    Maximum number of East-West partitions, as determined by EQ_REGIONS algorithm.

SETUP_TRANS

Signature

SUBROUTINE SETUP_TRANS(KSMAX, KDGL, KDLON, KLOEN, LDSPLIT, PSTRET, KTMAX, KRESOL, &
  &                    PWEIGHT, LDGRIDONLY, LDUSERPNM, LDKEEPRPNM, LDUSEFLT, &
  &                    LDSPSETUPONLY, LDPNMONLY, LDUSEFFTW, LDLL, LDSHIFTLL, &
  &                    CDIO_LEGPOL, CDLEGPOLFNAME, KLEGPOLPTR, KLEGPOLPTR_LEN)

Purpose

This subroutine establishes the environment for performing a spectral transform. Each call to this subroutine creates work arrays for a certain resolution. One can do this for NMAX_RESOL different resolutions (internal parameter corresponding to KMAX_RESOL argument of SETUP_TRANS0). One must call SETUP_TRANS0 before calling this subroutine for the first time.

INTENT(IN) arguments

  • INTEGER(KIND=JPIM), INTENT(IN) :: KSMAX
    Spectral truncation to perform the transform up to.
  • INTEGER(KIND=JPIM), INTENT(IN) :: KDGL
    The number of Gaussian latitudes in grid point space.

OPTIONAL, INTENT(IN) arguments

  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(IN) :: KDLON
    Maximum number of points on any latitude (usually the latitude nearest to the Equator).
    If not provided, it will take the value of 2 * KDGL, unless LDLL is .TRUE. in which case
    it will take the value of 2 * KDGL + 2, unless KLOEN is provided in which case it will take
    the value of MAXVAL(KLOEN).
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(IN) :: KLOEN(:)
    An array giving the number of points on each Gaussian latitude (dimension)
    If this is not provided it will be assumed that a full Gaussian grid is used in grid point space,
    for which every latitude will have KDLON points.
  • LOGICAL, OPTIONAL, INTENT(IN) :: LDSPLIT
    Split latitude in grid point space.
    Default: .FALSE.
  • REAL(KIND=JPRD), OPTIONAL, INTENT(IN) :: PSTRET
    Stretching factor for when the Legendre polynomials are computed on the stretched sphere.
    Default: 1.0
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(IN) :: KTMAX
    Spectral truncation to be applied for tendencies.
    Default: KSMAX
  • REAL(KIND=JPRD), OPTIONAL, INTENT(IN) :: PWEIGHT(:)
    The weight per grid-point (for a weighted distribution). If this argument is not provided, a weighted distribution will not be used.
  • LOGICAL, OPTIONAL, INTENT(IN) :: LDGRIDONLY
    Only provide grid point space results.
    Default: .FALSE.
    Potentially deprecatable.
  • LOGICAL, OPTIONAL, INTENT(IN) :: LDUSEFLT
    Use the fast Legendre transform algorithm.
    Default: .FALSE.
  • LOGICAL, OPTIONAL, INTENT(IN) :: LDUSERPNM
    Use the Belusov algorithm to compute the Legendre polynomials.
    Default: .TRUE.
  • LOGICAL, OPTIONAL, INTENT(IN) :: LDKEEPRPNM
    Store Legendre Polynomials instead of recomputing (only applicable when using the fast Legendre transform, otherwise these are always stored).
    Default: .FALSE.
  • LOGICAL, OPTIONAL, INTENT(IN) :: LDSPSETUPONLY
    Only create distributed spectral setup.
    Default: .FALSE.
  • LOGICAL, OPTIONAL, INTENT(IN) :: LDPNMONLY
    Compute the Legendre polynomials only, not the FFTs.
    Default: .FALSE.
  • LOGICAL, OPTIONAL, INTENT(IN) :: LDLL
    Setup second set of input/output latitudes.
    Default: .FALSE.
  • LOGICAL, OPTIONAL, INTENT(IN) :: LDSHIFTLL
    Shift output lon/lat data by 0.5*dx and 0.5*dy.
    Default: .FALSE.
  • CHARACTER(LEN=*), OPTIONAL, INTENT(IN) :: CDIO_LEGPOL
    String config argument to determine several options relating to I/O operations relevant to the Legendre polynomials. If 'readf' is passed, Legendre polynomials will be read from the file given to CDLEGPOLFNAME. If 'writef' is passed, Legendre polynomials will be written to the file given to CDLEGPOLFNAME. If 'membuf' is passed, Legendre polynomials will be read from shared memory.
  • CHARACTER(LEN=*), OPTIONAL, INTENT(IN) :: CDLEGPOLFNAME
    Filename from which to read or to which to write the Legendre polynomials. See CDIO_LEGPOL.
  • TYPE(C_PTR), OPTIONAL, INTENT(IN) :: KLEGPOLPTR
    Pointer to memory segment containing Legendre polynomials. See CDIO_LEGPOL.
  • INTEGER(C_SIZE_T), OPTIONAL, INTENT(IN) :: KLEGPOLPTR_LEN
    Length of memory segment containing Legendre polynomials. See CDIO_LEGPOL.

OPTIONAL, INTENT(OUT) arguments

  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(OUT) :: KRESOL
    A unique integer identifier to act as a handle for the work arrays corresponding to the provided transform resolution.

DIR_TRANS

Signature

SUBROUTINE DIR_TRANS(PSPVOR, PSPDIV, PSPSCALAR, PSPSC3A, PSPSC3B, PSPSC2, LDLATLON, KPROMA, &
  &                  KVSETUV, KVSETSC, KRESOL, KVSETSC3A, KVSETSC3B, KVSETSC2, PGP, PGPUV, &
  &                  PGP3A, PGP3B, PGP2)

Purpose

This subroutine performs a direct spectral transform. DIR_TRANS supports two modes of passing arrays (which we call "call modes" below):

  1. Pass in PGP, receive PSPVOR, PSPDIV, PSPSCALAR.
  2. Pass in PGPUV, PGP3A, PGP3B, PGP2, receive PSPVOR, PSPDIV, PSPSC3A, PSPSC3B, PSPSC2.

Although call mode 1 is simpler, it can often entail having to create array copies in order to prepare the different elements of PGP. Hence, with call mode 2, the input and output arrays are categorized according to whether they are wind related, 3D scalar, or 2D scalar.

Note that PSPSC3A/PGP3A and PSPSC3B/PGP3B are essentially identical in type and size. There are cases where it is useful to have the ability to pass two sets of scalar fields to be transformed independently.

OPTIONAL, INTENT(IN) arguments

  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(IN) :: KPROMA
    Required blocking factor for grid point output.
    Default: D%NGPTOT
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(IN) :: KVSETUV(:)
    Array which maps each vertical level of the output vorticity/divergence to its corresponding
    member of the V-set.
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(IN) :: KVSETSC(:)
    Array which maps each of the output scalar fields to its corresponding member of the V-set
    (call mode 1).
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(IN) :: KVSETSC3A(:)
    Array which maps each vertical level of the output 3D scalar fields to its corresponding
    member of the V-set (call mode 2).
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(IN) :: KVSETSC3B(:)
    Array which maps each vertical level of the output 3D scalar fields to its corresponding
    member of the V-set (call mode 2).
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(IN) :: KVSETSC2(:)
    Array which maps each output 2D scalar field to its corresponding member of the V-set
    (call mode 2).
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(IN) :: KRESOL
    Resolution handle returned by original call to SETUP_TRANS.
  • LOGICAL, OPTIONAL, INTENT(IN) :: LDLATLON
    TODO: what is this?
  • REAL(KIND=JPRB), OPTIONAL, INTENT(IN) :: PGP(:,:,:)
    Array containing all grid point fields (call mode 1).
    Dimensions: (NPROMA, number of grid point fields, number of NPROMA blocks).
    The ordering of grid point fields in this array is as follows:
    1. U wind
    2. V wind
    3. scalar fields
  • REAL(KIND=JPRB), OPTIONAL, INTENT(IN) :: PGPUV(:,:,:,:)
    Array containing wind-related grid point fields (call mode 2).
    Dimensions: (NPROMA, vertical levels, number of wind fields, number of NPROMA blocks).
    DIR_TRANS only operates on U/V winds, not vorticity/divergence. Hence, the third dimension
    of PGPUV must be equal to 2 (one for U and one for V).
  • REAL(KIND=JPRB), OPTIONAL, INTENT(IN) :: PGP3A(:,:,:,:)
    Array containing 3D scalar grid point fields, corresponding to PSPSC3A.
    Dimensions: (NPROMA, vertical levels, number of scalar fields, number of NPROMA blocks).
  • REAL(KIND=JPRB), OPTIONAL, INTENT(IN) :: PGP3B(:,:,:,:)
    Array containing 3D scalar grid point fields, corresponding to PSPSC3B.
    Dimensions: (NPROMA, vertical levels, number of scalar fields, number of NPROMA blocks).
  • REAL(KIND=JPRB), OPTIONAL, INTENT(IN) :: PGP2(:,:,:)
    Array containing 2D scalar grid point fields, corresponding to PSPSC2.
    Dimensions: (NPROMA, number of scalar fields, number of NPROMA blocks).

OPTIONAL, INTENT(OUT) arguments

  • REAL(KIND=JPRB), OPTIONAL, INTENT(OUT) :: PSPVOR(:,:)
    Spectral space vorticity. Dimensions: (vertical levels, spectral coefficients).
  • REAL(KIND=JPRB), OPTIONAL, INTENT(OUT) :: PSPDIV(:,:)
    Spectral space divergence. Dimensions: (vertical levels, spectral coefficients).
  • REAL(KIND=JPRB), OPTIONAL, INTENT(OUT) :: PSPSCALAR(:,:)
    Spectral space scalar fields (call mode 1).
    Dimensions: (total number of 2D scalar fields, spectral coefficients).
  • REAL(KIND=JPRB), OPTIONAL, INTENT(OUT) :: PSPSC3A(:,:,:)
    Spectral space 3D scalar fields (call mode 2) corresponding to PGP3A.
    Dimensions: (vertical levels, spectral coefficients, number of 3D scalar fields).
  • REAL(KIND=JPRB), OPTIONAL, INTENT(OUT) :: PSPSC3B(:,:,:)
    Spectral space 3D scalar fields (call mode 2) corresponding to PGP3B.
    Dimensions: (vertical levels, spectral coefficients, number of 3D scalar fields).
  • REAL(KIND=JPRB), OPTIONAL, INTENT(OUT) :: PSPSC2(:,:)
    Spectral space 2D scalar fields (call mode 2) corresponding to PGP2.

INV_TRANS

Signature

SUBROUTINE INV_TRANS(PSPVOR, PSPDIV, PSPSCALAR, PSPSC3A, PSPSC3B, PSPSC2, &
  &                  FSPGL_PROC, LDSCDERS, LDVORGP, LDDIVGP, LDUVDER, LDLATLON, KPROMA, KVSETUV, &
  &                  KVSETSC, KRESOL, KVSETSC3A, KVSETSC3B, KVSETSC2, PGP, PGPUV, PGP3A, PGP3B, &
  &                  PGP2)

Purpose

This subroutine performs an inverse spectral transform. INV_TRANS supports two modes of passing arrays (which we call "call modes" below):

  1. Pass in PSPVOR, PSPDIV, PSPSCALAR, receive PGP.
  2. Pass in PSPVOR, PSPDIV, PSPSC3A, PSPSC3B, PSPSC2, receive PGPUV, PGP3A, PGP3B, PGP2

Although call mode 1 is simpler, it can often entail having to create array copies in order to subsequently operate on the different elements of PGP. Hence, with call mode 2, the output arrays are already categorized according to whether they are wind related, 3D scalar, or 2D scalar.

Note that PSPSC3A/PGP3A and PSPSC3B/PGP3B are essentially identical in type and size. There are cases where it is useful to have the ability to pass two sets of scalar fields to be transformed independently.

OPTIONAL, INTENT(IN) arguments

  • REAL(KIND=JPRB), OPTIONAL, INTENT(IN) :: PSPVOR(:,:)
    Spectral space vorticity.
    Dimensions: (vertical levels, spectral coefficients).
  • REAL(KIND=JPRB), OPTIONAL, INTENT(IN) :: PSPDIV(:,:)
    Spectral space divergence.
    Dimensions: (vertical levels, spectral coefficients).
  • REAL(KIND=JPRB), OPTIONAL, INTENT(IN) :: PSPSCALAR(:,:)
    Spectral space scalar fields (call mode 1).
    Dimensions: (total number of 2D scalar fields, spectral coefficients).
  • REAL(KIND=JPRB), OPTIONAL, INTENT(IN) :: PSPSC3A(:,:,:)
    Spectral space 3D scalar fields (call mode 2) corresponding to PGP3A.
    Dimensions: (vertical levels, spectral coefficients, number of 3D scalar fields).
  • REAL(KIND=JPRB), OPTIONAL, INTENT(IN) :: PSPSC3B(:,:,:)
    Spectral space 3D scalar fields (call mode 2) corresponding to PGP3B.
    Dimensions: (vertical levels, spectral coefficients, number of 3D scalar fields).
  • REAL(KIND=JPRB), OPTIONAL, INTENT(IN) :: PSPSC2(:,:)
    Spectral space 2D scalar fields (call mode 2) corresponding to PGP2.
    Dimensions: (number of 2D scalar fields, spectral coefficients).
  • LOGICAL, OPTIONAL, INTENT(IN) :: LDSCDERS
    Indicates if derivatives of scalar variables are required.
    Default: .FALSE.
  • LOGICAL, OPTIONAL, INTENT(IN) :: LDVORGP
    Indicates if grid point space vorticity is required.
    Default: .FALSE.
  • LOGICAL, OPTIONAL, INTENT(IN) :: LDDIVGP
    Indicates if grid point space divergence is required.
    Default: .FALSE.
  • LOGICAL, OPTIONAL, INTENT(IN) :: LDUVDER
    Indicates of East-West derivatives of U and V are required.
    Default: .FALSE.
  • LOGICAL, OPTIONAL, INTENT(IN) :: LDLATLON
    Indicates if regular latitude-longitude output is required.
    Default: .FALSE.
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(IN) :: KPROMA
    Required blocking factor for grid point output.
    Default: D%NGPTOT
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(IN) :: KVSETUV(:)
    Array which maps each vertical level of the input vorticity/divergence to its corresponding
    member of the V-set.
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(IN) :: KVSETSC(:)
    Array which maps each of the input scalar fields to its corresponding member of the V-set
    (call mode 1).
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(IN) :: KVSETSC3A(:)
    Array which maps each vertical level of the input 3D scalar fields to its corresponding
    member of the V-set (call mode 2).
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(IN) :: KVSETSC3B(:)
    Array which maps each vertical level of the input 3D scalar fields to its corresponding
    member of the V-set (call mode 2).
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(IN) :: KVSETSC2(:)
    Array which maps each 2D scalar field to its corresponding member of the V-set (call mode 2).
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(IN) :: KRESOL
    Resolution handle returned by original call to SETUP_TRANS.

OPTIONAL, EXTERNAL arguments

  • OPTIONAL, EXTERNAL :: FSPGL_PROC
    External procedure to be executed in Fourier space before transposition.

OPTIONAL, INTENT(OUT) arguments

  • REAL(KIND=JPRB), OPTIONAL, INTENT(OUT) :: PGP(:,:,:)
    Array containing all grid point fields (call mode 1).
    Dimensions: (NPROMA, number of grid point fields, number of NPROMA blocks).
    The ordering of grid point fields in this array is as follows:
    1. vorticity (if PSPVOR is present and LDVORGP)
    2. divergence (if PSPDIV is present and LDDIVGP)
    3. U wind (if PSPVOR and PSPDIV are present)
    4. V wind (if PSPVOR and PSPDIV are present)
    5. scalar fields (if PSPSCALAR is present)
    6. North-South derivative of scalar fields (if PSPSCALAR is present and LDSCDERS)
    7. East-West derivative of U wind (if PSPVOR and PSPDIV are present and LDUVDER)
    8. East-West derivative of V wind (if PSPVOR and PSPDIV are present and LDUVDER)
    9. East-West derivative of scalar fields (if PSPSCALAR is present and LDSCDERS)
  • REAL(KIND=JPRB), OPTIONAL, INTENT(OUT) :: PGPUV(:,:,:,:)
    Array containing grid point fields relating to wind (call mode 2).
    Dimensions: (NPROMA, vertical levels, number of wind-related fields, number of NPROMA blocks).
    The ordering of grid point fields in this array is as follows:
    1. vorticity (if PSPVOR is present and LDVORGP)
    2. divergence (if PSPDIV is present and LDDIVGP)
    3. U wind (if PSPVOR and PSPDIV are present)
    4. V wind (if PSPVOR and PSPDIV are present)
    5. East-West derivative of U wind (if PSPVOR and PSPDIV are present and LDUVDER)
    6. East-West derivative of V wind (if PSPVOR and PSPDIV are present and LDUVDER)
  • REAL(KIND=JPRB), OPTIONAL, INTENT(OUT) :: PGP3A(:,:,:,:)
    Array containing 3D scalar grid point fields corresponding to PSPSC3A.
    Dimensions: (NPROMA, vertical levels, number of 3D scalar fields, number of NPROMA blocks).
    The ordering of grid point fields in this array is as follows:
    1. scalar fields
    2. North-South derivatives of scalar fields (if LDSCDERS)
    3. East-West derivatives of scalar fields (if LDSCDERS)
  • REAL(KIND=JPRB), OPTIONAL, INTENT(OUT) :: PGP3B(:,:,:,:)
    Array containing 3D scalar grid point fields corresponding to PSPSC3B.
    Dimensions: (NPROMA, vertical levels, number of 3D scalar fields, number of NPROMA blocks).
    The ordering of grid point fields in this array is as follows:
    1. scalar fields
    2. North-South derivatives of scalar fields (if LDSCDERS)
    3. East-West derivatives of scalar fields (if LDSCDERS)
  • REAL(KIND=JPRB), OPTIONAL, INTENT(OUT) :: PGP2(:,:,:)
    Array containing 2D scalar grid point fields corresponding to PSPSC2.
    Dimensions: (NPROMA, number of 2D scalar fields, number of NPROMA blocks).
    The ordering of grid point fields in this array is as follows:
    1. scalar fields
    2. North-South derivatives of scalar fields (if LDSCDERS)
    3. East-West derivatives of scalar fields (if LDSCDERS)

TRANS_RELEASE

Signature

SUBROUTINE TRANS_RELEASE(KRESOL)

Purpose

This subroutine releases (i.e. deallocates and nullifies) all arrays related to the given resolution tag KRESOL.

INTENT(IN) arguments

  • INTEGER(KIND=JPIM), INTENT(IN) :: KRESOL
    The handle for the resolution you want to release.

TRANS_END

Signature

SUBROUTINE TRANS_END(CDMODE)

Purpose

This subroutine terminates the transform package by releasing (i.e. deallocating and nullifying) all allocated arrays.

OPTIONAL, INTENT(IN) arguments

  • CHARACTER*5, OPTIONAL, INTENT(IN) :: CDMODE
    A string parameter indicating which mode to finalise with. Either 'FINAL' to finalise
    everything or 'INTER' just to deallocate N_REGIONS and NPRCIDS.
    Default: 'FINAL'

TRANS_INQ

Signature

SUBROUTINE TRANS_INQ(KRESOL, KSPEC, KSPEC2, KSPEC2G, KSPEC2MX, KNUMP, KGPTOT, KGPTOTG, KGPTOTMX, &
  &                  KGPTOTL, KMYMS, KASM0, KUMPP, KPOSSP, KPTRMS, KALLMS, KDIM0G, KFRSTLAT, &
  &                  KLSTLAT, KFRSTLOFF, KPTRLAT, KPTRFRSTLAT, KPTRLSTLAT, KPTRFLOFF, KSTA, &
  &                  KONL, KULTPP, KPTRLS, KNMENG, KPRTRW, KMYSETW, KMYSETV, KMY_REGION_NS, &
  &                  KMY_REGION_EW, LDSPLITLAT, KSMAX, PLAPIN, KNVALUE, KDEF_RESOL, LDLAM, PMU, &
  &                  PGW, PRPNM, KLEI3, KSPOLEGL, KPMS, KDGLU)

Purpose

This subroutine allows one to obtain information on the transform configuration for a given KRESOL handle (or the default handle if this argument is omitted).

OPTIONAL, INTENT(IN) arguments

  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(IN) :: KRESOL
    Resolution handle for which information is required. Default: The first defined resolution handle.

OPTIONAL, INTENT(OUT) arguments

  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(OUT) :: KSPEC
    Number of complex spectral coefficients on this MPI task.
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(OUT) :: KSPEC2
    Number of real and imaginary spectral coefficients on this MPI task. Equal to 2 * KSPEC.
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(OUT) :: KSPEC2G
    Number of real and imaginary spectral coefficients across all MPI tasks.
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(OUT) :: KSPEC2MX
    The maximum value of KSPEC2 over all MPI tasks.
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(OUT) :: KNUMP
    Number of zonal wave numbers handled by this MPI task for the Legendre transform.
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(OUT) :: KGPTOT
    Number of grid columns on this MPI task.
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(OUT) :: KGPTOTG
    Total number of grid columns across all MPI tasks.
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(OUT) :: KGPTOTMX
    Maximum number of grid columns on any of the MPI tasks.
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(OUT) :: KGPTOTL(:,:)
    TODO: check. Number of grid columns one each PE (dimension N_REGIONS_NS:N_REGIONS_EW)
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(OUT) :: KMYMS(:)
    The zonal wavenumbers handled by this MPI task for the Legendre transform.
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(OUT) :: KASM0(0:)
    Address in a spectral array of (m, n=m).
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(OUT) :: KUMPP(:)
    TODO: check this. Number of wave numbers each wave set is responsible for
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(OUT) :: KPOSSP(:)
    Defines partitioning of global spectral fields among PEs
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(OUT) :: KPTRMS(:)
    Pointer to the first wave number of a given A-set
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(OUT) :: KALLMS(:)
    Wave numbers for all wave-set concatenated together to give all wave numbers in wave-set order
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(OUT) :: KDIM0G(0:)
    Defines partitioning of global spectral fields among PEs
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(OUT) :: KFRSTLAT(:)
    First latitude of each a-set in grid-point space
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(OUT) :: KLSTLAT(:)
    Last latitude of each a-set in grid-point space
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(OUT) :: KFRSTLOFF
    Offset for first lat of own a-set in grid-point space
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(OUT) :: KPTRLAT(:)
    Pointer to the start of each latitude
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(OUT) :: KPTRFRSTLAT(:)
    Pointer to the first latitude of each a-set in NSTA and NONL arrays
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(OUT) :: KPTRLSTLAT(:)
    Pointer to the last latitude of each a-set in NSTA and NONL arrays
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(OUT) :: KPTRFLOFF
    Offset for pointer to the first latitude of own a-set NSTA and NONL arrays, i.e. nptrfrstlat
    (myseta)-1
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(OUT) :: KSTA(:,:)
    Position of first grid column for the latitudes on a processor. The information is available
    for all processors. The b-sets are distinguished by the last dimension of nsta(). The latitude
    band for each a-set is addressed by nptrfrstlat(jaset), nptrlstlat(jaset), and
    nptrfloff=nptrfrstlat(myseta) on this processors a-set. Each split latitude has two entries in
    nsta(,:) which necessitates the rather complex addressing of nsta(,:) and the overdimensioning
    of nsta by N_REGIONS_NS.
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(OUT) :: KONL(:,:)
    Number of grid columns for the latitudes on a processor. Similar to nsta() in data structure.
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(OUT) :: KULTPP(:)
    number of latitudes for which each a-set is calculating the FFTs.
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(OUT) :: KPTRLS(:)
    pointer to first global latitude of each a-set for which it performs the Fourier calculations
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(OUT) :: KNMENG(:)
    associated (with NLOENG) cut-off zonal wavenumber
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(OUT) :: KPRTRW
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(OUT) :: KMYSETW
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(OUT) :: KMYSETV
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(OUT) :: KMY_REGION_NS
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(OUT) :: KMY_REGION_EW
  • LOGICAL, OPTIONAL, INTENT(OUT) :: LDSPLITLAT(:)
    TRUE if latitude is split in grid point space over two a-sets
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(OUT) :: KSMAX
    spectral truncation
  • REAL(KIND=JPRB), OPTIONAL, INTENT(OUT) :: PLAPIN(-1:)
    Eigenvalues of the inverse Laplace operator
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(OUT) :: KNVALUE(:)
    n value for each KSPEC2 spectral coeffient
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(OUT) :: KDEF_RESOL
    number or resolutions defined
  • LOGICAL, OPTIONAL, INTENT(OUT) :: LDLAM
    .T. if the corresponding resolution is LAM, .F. if it is global
  • REAL(KIND=JPRD), OPTIONAL, INTENT(OUT) :: PMU(:)
    sin(Gaussian latitudes)
  • REAL(KIND=JPRB), OPTIONAL, INTENT(OUT) :: PGW(:)
    Gaussian weights
  • REAL(KIND=JPRB), OPTIONAL, INTENT(OUT) :: PRPNM(:,:)
    Legendre polynomials
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(OUT) :: KLEI3
    First dimension of Legendre polynomials
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(OUT) :: KSPOLEGL
    Second dimension of Legendre polynomials
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(OUT) :: KPMS(0:)
    Address for legendre polynomial for given M (NSMAX)
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(OUT) :: KDGLU(0:)
    Number of active points in an hemisphere for a given wavenumber "m"

DIR_TRANSAD

Signature

SUBROUTINE DIR_TRANSAD(PSPVOR, PSPDIV, PSPSCALAR, PSPSC3A, PSPSC3B, PSPSC2, KPROMA, KVSETUV, &
  &                    KVSETSC, KRESOL, KVSETSC3A, KVSETSC3B, KVSETSC2, PGP, PGPUV, PGP3A, &
  &                    PGP3B, PGP2)

Purpose

This is the adjoint version of the subroutine for performing a direct spectral transform. Like DIR_TRANS, DIR_TRANSAD supports two modes of passing arrays (which we call "call modes" below):

  1. Pass in PGP, receive PSPVOR, PSPDIV, PSPSCALAR.
  2. Pass in PGPUV, PGP3A, PGP3B, PGP2, receive PSPVOR, PSPDIV, PSPSC3A, PSPSC3B, PSPSC2.

Although call mode 1 is simpler, it can often entail having to create array copies in order to prepare the different elements of PGP. Hence, with call mode 2, the input and output arrays are categorized according to whether they are wind related, 3D scalar, or 2D scalar.

Note that PSPSC3A/PGP3A and PSPSC3B/PGP3B are essentially identical in type and size. There are cases where it is useful to have the ability to pass two sets of scalar fields to be transformed independently.

OPTIONAL, INTENT(IN) arguments

  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(IN) :: KPROMA
    Required blocking factor for grid point output.
    Default: D%NGPTOT
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(IN) :: KVSETUV(:)
    Array which maps each vertical level of the output vorticity/divergence to its corresponding
    member of the V-set.
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(IN) :: KVSETSC(:)
    Array which maps each of the output scalar fields to its corresponding member of the V-set
    (call mode 1).
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(IN) :: KVSETSC3A(:)
    Array which maps each vertical level of the output 3D scalar fields to its corresponding
    member of the V-set (call mode 2).
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(IN) :: KVSETSC3B(:)
    Array which maps each vertical level of the output 3D scalar fields to its corresponding
    member of the V-set (call mode 2).
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(IN) :: KVSETSC2(:)
    Array which maps each output 2D scalar field to its corresponding member of the V-set
    (call mode 2).
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(IN) :: KRESOL
    Resolution handle returned by original call to SETUP_TRANS.

OPTIONAL, INTENT(OUT) arguments

  • REAL(KIND=JPRB), OPTIONAL, INTENT(OUT) :: PGP(:,:,:)
    Array containing all grid point fields (call mode 1).
    Dimensions: (NPROMA, number of grid point fields, number of NPROMA blocks).
    The ordering of grid point fields in this array is as follows:
    1. U wind
    2. V wind
    3. scalar fields
  • REAL(KIND=JPRB), OPTIONAL, INTENT(OUT) :: PGPUV(:,:,:,:)
    Array containing wind-related grid point fields (call mode 2).
    Dimensions: (NPROMA, vertical levels, number of wind fields, number of NPROMA blocks).
    DIR_TRANS only operates on U/V winds, not vorticity/divergence. Hence, the third dimension
    of PGPUV must be equal to 2 (one for U and one for V).
  • REAL(KIND=JPRB), OPTIONAL, INTENT(OUT) :: PGP3A(:,:,:,:)
    Array containing 3D scalar grid point fields, corresponding to PSPSC3A.
    Dimensions: (NPROMA, vertical levels, number of scalar fields, number of NPROMA blocks).
  • REAL(KIND=JPRB), OPTIONAL, INTENT(OUT) :: PGP3B(:,:,:,:)
    Array containing 3D scalar grid point fields, corresponding to PSPSC3B.
    Dimensions: (NPROMA, vertical levels, number of scalar fields, number of NPROMA blocks).
  • REAL(KIND=JPRB), OPTIONAL, INTENT(OUT) :: PGP2(:,:,:)
    Array containing 2D scalar grid point fields, corresponding to PSPSC2.
    Dimensions: (NPROMA, number of scalar fields, number of NPROMA blocks).

OPTIONAL, INTENT(INOUT) arguments

  • REAL(KIND=JPRB), OPTIONAL, INTENT(INOUT) :: PSPVOR(:,:)
    Spectral space vorticity. Dimensions: (vertical levels, spectral coefficients).
  • REAL(KIND=JPRB), OPTIONAL, INTENT(INOUT) :: PSPDIV(:,:)
    Spectral space divergence. Dimensions: (vertical levels, spectral coefficients).
  • REAL(KIND=JPRB), OPTIONAL, INTENT(INOUT) :: PSPSCALAR(:,:)
    Spectral space scalar fields (call mode 1).
    Dimensions: (total number of 2D scalar fields, spectral coefficients).
  • REAL(KIND=JPRB), OPTIONAL, INTENT(INOUT) :: PSPSC3A(:,:,:)
    Spectral space 3D scalar fields (call mode 2) corresponding to PGP3A.
    Dimensions: (vertical levels, spectral coefficients, number of 3D scalar fields).
  • REAL(KIND=JPRB), OPTIONAL, INTENT(INOUT) :: PSPSC3B(:,:,:)
    Spectral space 3D scalar fields (call mode 2) corresponding to PGP3B.
    Dimensions: (vertical levels, spectral coefficients, number of 3D scalar fields).
  • REAL(KIND=JPRB), OPTIONAL, INTENT(INOUT) :: PSPSC2(:,:)
    Spectral space 2D scalar fields (call mode 2) corresponding to PGP2.

INV_TRANSAD

Signature

SUBROUTINE INV_TRANSAD(PSPVOR, PSPDIV, PSPSCALAR, PSPSC3A, PSPSC3B, PSPSC2, FSPGL_PROC, &
  &                    LDSCDERS, LDVORGP, LDDIVGP, LDUVDER, KPROMA, KVSETUV, KVSETSC, KRESOL, &
  &                    KVSETSC3A, KVSETSC3B, KVSETSC2, PGP, PGPUV, PGP3A, PGP3B, &
  &                    PGP2)

Purpose

This is the adjoint version of the subroutine for performing an inverse spectral transform. Like INV_TRANS, INV_TRANSAD supports two modes of passing arrays (which we call "call modes" below):

  1. Pass in PSPVOR, PSPDIV, PSPSCALAR, receive PGP.
  2. Pass in PSPVOR, PSPDIV, PSPSC3A, PSPSC3B, PSPSC2, receive PGPUV, PGP3A, PGP3B, PGP2

Although call mode 1 is simpler, it can often entail having to create array copies in order to subsequently operate on the different elements of PGP. Hence, with call mode 2, the output arrays are already categorized according to whether they are wind related, 3D scalar, or 2D scalar.

Note that PSPSC3A/PGP3A and PSPSC3B/PGP3B are essentially identical in type and size. There are cases where it is useful to have the ability to pass two sets of scalar fields to be transformed independently.

OPTIONAL, INTENT(IN) arguments

  • LOGICAL, OPTIONAL, INTENT(IN) :: LDSCDERS
    Indicates if derivatives of scalar variables are required.
    Default: .FALSE.
  • LOGICAL, OPTIONAL, INTENT(IN) :: LDVORGP
    Indicates if grid point space vorticity is required.
    Default: .FALSE.
  • LOGICAL, OPTIONAL, INTENT(IN) :: LDDIVGP
    Indicates if grid point space divergence is required.
    Default: .FALSE.
  • LOGICAL, OPTIONAL, INTENT(IN) :: LDUVDER
    Indicates of East-West derivatives of U and V are required.
    Default: .FALSE.
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(IN) :: KPROMA
    Required blocking factor for grid point output.
    Default: D%NGPTOT
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(IN) :: KVSETUV(:)
    Array which maps each vertical level of the input vorticity/divergence to its corresponding
    member of the V-set.
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(IN) :: KVSETSC(:)
    Array which maps each of the input scalar fields to its corresponding member of the V-set
    (call mode 1).
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(IN) :: KVSETSC3A(:)
    Array which maps each vertical level of the input 3D scalar fields to its corresponding
    member of the V-set (call mode 2).
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(IN) :: KVSETSC3B(:)
    Array which maps each vertical level of the input 3D scalar fields to its corresponding
    member of the V-set (call mode 2).
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(IN) :: KVSETSC2(:)
    Array which maps each 2D scalar field to its corresponding member of the V-set (call mode 2).
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(IN) :: KRESOL
    Resolution handle returned by original call to SETUP_TRANS.
  • REAL(KIND=JPRB), OPTIONAL, INTENT(IN) :: PGP(:,:,:)
    Array containing all grid point fields (call mode 1).
    Dimensions: (NPROMA, number of grid point fields, number of NPROMA blocks).
    The ordering of grid point fields in this array is as follows:
    1. vorticity (if PSPVOR is present and LDVORGP)
    2. divergence (if PSPDIV is present and LDDIVGP)
    3. U wind (if PSPVOR and PSPDIV are present)
    4. V wind (if PSPVOR and PSPDIV are present)
    5. scalar fields (if PSPSCALAR is present)
    6. North-South derivative of scalar fields (if PSPSCALAR is present and LDSCDERS)
    7. East-West derivative of U wind (if PSPVOR and PSPDIV are present and LDUVDER)
    8. East-West derivative of V wind (if PSPVOR and PSPDIV are present and LDUVDER)
    9. East-West derivative of scalar fields (if PSPSCALAR is present and LDSCDERS)
  • REAL(KIND=JPRB), OPTIONAL, INTENT(IN) :: PGPUV(:,:,:,:)
    Array containing grid point fields relating to wind (call mode 2).
    Dimensions: (NPROMA, vertical levels, number of wind-related fields, number of NPROMA blocks).
    The ordering of grid point fields in this array is as follows:
    1. vorticity (if PSPVOR is present and LDVORGP)
    2. divergence (if PSPDIV is present and LDDIVGP)
    3. U wind (if PSPVOR and PSPDIV are present)
    4. V wind (if PSPVOR and PSPDIV are present)
    5. East-West derivative of U wind (if PSPVOR and PSPDIV are present and LDUVDER)
    6. East-West derivative of V wind (if PSPVOR and PSPDIV are present and LDUVDER)
  • REAL(KIND=JPRB), OPTIONAL, INTENT(IN) :: PGP3A(:,:,:,:)
    Array containing 3D scalar grid point fields corresponding to PSPSC3A.
    Dimensions: (NPROMA, vertical levels, number of 3D scalar fields, number of NPROMA blocks).
    The ordering of grid point fields in this array is as follows:
    1. scalar fields
    2. North-South derivatives of scalar fields (if LDSCDERS)
    3. East-West derivatives of scalar fields (if LDSCDERS)
  • REAL(KIND=JPRB), OPTIONAL, INTENT(IN) :: PGP3B(:,:,:,:)
    Array containing 3D scalar grid point fields corresponding to PSPSC3B.
    Dimensions: (NPROMA, vertical levels, number of 3D scalar fields, number of NPROMA blocks).
    The ordering of grid point fields in this array is as follows:
    1. scalar fields
    2. North-South derivatives of scalar fields (if LDSCDERS)
    3. East-West derivatives of scalar fields (if LDSCDERS)
  • REAL(KIND=JPRB), OPTIONAL, INTENT(IN) :: PGP2(:,:,:)
    Array containing 2D scalar grid point fields corresponding to PSPSC2.
    Dimensions: (NPROMA, number of 2D scalar fields, number of NPROMA blocks).
    The ordering of grid point fields in this array is as follows:
    1. scalar fields
    2. North-South derivatives of scalar fields (if LDSCDERS)
    3. East-West derivatives of scalar fields (if LDSCDERS)

OPTIONAL, EXTERNAL arguments

  • OPTIONAL, EXTERNAL :: FSPGL_PROC
    External procedure to be executed in Fourier space before transposition.

OPTIONAL, INTENT(OUT) arguments

  • REAL(KIND=JPRB), OPTIONAL, INTENT(OUT) :: PSPVOR(:,:)
    Spectral space vorticity.
    Dimensions: (vertical levels, spectral coefficients).
  • REAL(KIND=JPRB), OPTIONAL, INTENT(OUT) :: PSPDIV(:,:)
    Spectral space divergence.
    Dimensions: (vertical levels, spectral coefficients).
  • REAL(KIND=JPRB), OPTIONAL, INTENT(OUT) :: PSPSCALAR(:,:)
    Spectral space scalar fields (call mode 1).
    Dimensions: (total number of 2D scalar fields, spectral coefficients).
  • REAL(KIND=JPRB), OPTIONAL, INTENT(OUT) :: PSPSC3A(:,:,:)
    Spectral space 3D scalar fields (call mode 2) corresponding to PGP3A.
    Dimensions: (vertical levels, spectral coefficients, number of 3D scalar fields).
  • REAL(KIND=JPRB), OPTIONAL, INTENT(OUT) :: PSPSC3B(:,:,:)
    Spectral space 3D scalar fields (call mode 2) corresponding to PGP3B.
    Dimensions: (vertical levels, spectral coefficients, number of 3D scalar fields).
  • REAL(KIND=JPRB), OPTIONAL, INTENT(OUT) :: PSPSC2(:,:)
    Spectral space 2D scalar fields (call mode 2) corresponding to PGP2.
    Dimensions: (number of 2D scalar fields, spectral coefficients).

DIST_GRID

Signature

SUBROUTINE DIST_GRID(PGPG, KPROMA, KFDISTG, KFROM, KRESOL, PGP, KSORT)

Purpose

The subroutine distributes one or more fields resident entirely on a single MPI task among all other tasks according to the grid point decomposition used in ecTrans.

INTENT(IN) arguments

  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(IN) :: KFDISTG
    The number of fields to be distributed.
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(IN) :: KFROM
    Array which, for each field to be distributed, which MPI task will be sending the field.

OPTIONAL, INTENT(IN) arguments

  • REAL(KIND=JPRB), OPTIONAL, INTENT(IN) :: PGPG
    Array containing fields to be distributed.
    Dimensions: (number of global grid points, number of fields on this MPI task). Note that this is optional, because some tasks may only receive fields (and so they wouldn't have any fields to offer through PGPG).
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(IN) :: KPROMA
    "Blocking factor" used for grid point output. Default: D%NGPTOT
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(IN) :: KRESOL
    Resolution handle returned by original call to SETUP_TRANS.
    Default: 1 (i.e. first resolution handle)
  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(IN) :: KSORT
    Array allowing to rearrange fields in the output array. For each element, specify which
    element you want the field to end up in.

INTENT(OUT) arguments

  • REAL(KIND=JPRB), OPTIONAL, INTENT(OUT) :: PGP(:,:)
    Array containing grid point output.
    Dimensions: (NPROMA, number of fields, number of NPROMA blocks).

DIST_GRID_32

Note

This subroutine is deprecated and will be removed in a future release.

GATH_GRID_32

Note

This subroutine is deprecated and will be removed in a future release.

GET_CURRENT

Signature

SUBROUTINE GET_CURRENT(KRESOL,LDLAM)

Purpose

This subroutine allows one to retrieve the current resolution handle, NCUR_RESOL, and the LAM parameter, LDLAM, which indicates if the local area version of ecTrans is being used or not.

Note

LDLAM will always be .FALSE. for the moment. ecTrans currently only supports global spectral transforms, though local (bifourier) transforms may be supported in future.

OPTIONAL, INTENT(OUT) arguments

  • INTEGER(KIND=JPIM), OPTIONAL, INTENT(OUT) :: KRESOL
    The current default resolution handle.
  • LOGICAL, OPTIONAL, INTENT(OUT) :: LDLAM
    Whether the local area version of ecTrans is being used or not.