# ply_polyBaseExc_module Module

Module to change bases functions of a modal representation.

\author{Jens Zudrop}

This module provides routines for a fast basis exchange between Legendre and Chebyshev polynomials. The alogorihtm is fast as it completes this task in O(n log(n)) operations (where n is the number of modal coefficients). The alogrithm is based on the following publication: Alpert, B., & Rokhlin, V. (1991). A fast algorithm for the evaluation of Legendre expansions. SIAM Journal on Scientific and Statistical Computing. There also an alternative implementation with O(n) operations is described.

Some recommendations to achieve a fast transformation. The number of minimal blocks (n/s) should be a power of two plus 1: (n/s) = 2^k + 1 This yields the minimal number of blocks, that need to be computed. If this is not possible, it is probably good to at lease have an odd number of blocks (n/s = 2*k + 1), to reduce the number of smallest blocks. The remainder mod(n, s) + s, should be even, as otherwise there is an additional diagonal that needs to be computed. Similarily also s itself should probably be even.

## Variables

TypeVisibilityAttributesNameInitial
integer, public, parameter:: ply_legToCheb_param =1
integer, public, parameter:: ply_chebToLeg_param =2

## Interfaces

• ### private subroutine Copy_trafo_params(left, right)

#### Arguments

TypeIntentOptionalAttributesName
type(ply_trafo_params_type), intent(out) :: left

fpt to copy to

type(ply_trafo_params_type), intent(in) :: right

fpt to copy from

## Derived Types

### type, public :: ply_trafo_params_type

#### Components

TypeVisibilityAttributesNameInitial
type(ply_sub_vec), private, allocatable:: u(:,:)

Lagrange polynomials evaluated at the Chebyshev points on [0,+1].

real(kind=rk), private, allocatable:: diag(:,:)

The array to store the diagonals of the matrix in.

real(kind=rk), private, allocatable:: adapter(:,:,:)

The array to store the adapters between diagonal and blocks in.

type(ply_submatrix_type), private, allocatable:: sub(:)

Data of all sub matrices (separated from the diagonal). Size is the number of different sub matrix sizes, i.e. h.

integer, private :: nBlocks

Number of blocks in one direction of the matrix.

integer, private :: striplen

Length of stripes to use in the matrix operation.

integer, private :: remainder

Remaining columns close to the diagonal after subdividing the matrix into blocks

integer, private :: nDiagonals

Number of full diagonals that need to be considered close the diagonal of the matrix.

integer, private :: nBlockDiagonals

Number of diagonals in triangle blocks, that remain between blocks and full diagonals.

integer, private :: n

The number of modal coefficients to convert

integer, private :: k

The number of Cheb coefficients to approximate M

integer, private :: s

The size of the smallest subblock of M

integer, private :: h

The number of subblocks (per direction) in M

integer, private :: subblockingWidth

The width of the subblocks used during the unrolled base exchange to ensure a better cache usage.

integer, private :: trafo

The transformation type

type(ply_subvector_type), private, allocatable:: b(:)

Conversion data structure used for fpt.

### type, private :: ply_sub_vec

#### Components

TypeVisibilityAttributesNameInitial
real(kind=rk), private, allocatable:: dat(:)

### type, private :: ply_matrixExpCoeff_type

Expansion coefficients for a certain submatrix.

#### Components

TypeVisibilityAttributesNameInitial
real(kind=rk), private, allocatable:: coeff(:)

Polynomials expansion coefficients.

### type, private :: ply_coldata_type

Information for a set of local rows in the current block

#### Components

TypeVisibilityAttributesNameInitial
type(ply_matrixExpCoeff_type), private, allocatable:: rowDat(:)

The Chebyshev expansion coefficients for a set of block local rows.

### type, private :: ply_rowdata_type

Sparse data for information of a column in a sub matrix

#### Components

TypeVisibilityAttributesNameInitial
type(ply_coldata_type), private, allocatable:: subCol(:)

Column data. Dimension is the different number of column rows in the sub matrix.

### type, private :: ply_submatrix_type

Sparse data for information of a row in a sub matrix.

#### Components

TypeVisibilityAttributesNameInitial
type(ply_rowdata_type), private, allocatable:: subRow(:)

Row data. Dimension is the different number of block rows in the sub matrix.

### type, private :: ply_matrixExpCoeffOddEven_type

Expansion coefficients for a certain submatrix.

#### Components

TypeVisibilityAttributesNameInitial
real(kind=rk), private, allocatable:: coeff(:,:)

Polynomials expansion coefficients.

### type, private :: ply_subvector_type

Sparse data for a subvector

#### Components

TypeVisibilityAttributesNameInitial
type(ply_matrixExpCoeffOddEven_type), private, allocatable:: col(:)

Expansion coefficients for a column

## Functions

### public elemental function ply_lambda(val) result(funcVal)

\todo: as we use a relation of gamma, it might be better to use the gammln function provided by the numerical recipes, and just use the difference in an exponential function.

#### Arguments

TypeIntentOptionalAttributesName
real(kind=rk), intent(in) :: val

### private elemental function ply_m(iReal, jReal) result(mVal)

#### Arguments

TypeIntentOptionalAttributesName
real(kind=rk), intent(in) :: iReal
real(kind=rk), intent(in) :: jReal

### private elemental function ply_m_int(iReal, jReal) result(mVal)

#### Arguments

TypeIntentOptionalAttributesName
integer, intent(in) :: iReal
integer, intent(in) :: jReal

### private elemental function ply_l(iReal, jReal) result(lVal)

#### Arguments

TypeIntentOptionalAttributesName
real(kind=rk), intent(in) :: iReal
real(kind=rk), intent(in) :: jReal

### private elemental function ply_l_int(i, j) result(lVal)

#### Arguments

TypeIntentOptionalAttributesName
integer, intent(in) :: i
integer, intent(in) :: j

## Subroutines

### public subroutine ply_fpt_init(n, params, trafo, blocksize, approx_terms, striplen, subblockingWidth)

#### Arguments

TypeIntentOptionalAttributesName
integer, intent(in) :: n
type(ply_trafo_params_type), intent(inout) :: params
integer, intent(in) :: trafo
integer, intent(in), optional :: blocksize

Smallest block that is approximated by approx_terms coefficients.

Please note, that this has to be larger than 2*approx_terms to result in a reduced number of operations. Default is 64.

integer, intent(in), optional :: approx_terms

Number of approximation terms used to compute off-diagonal products.

Defaults to 18, which is the suggested accuracy for double precision.

integer, intent(in) :: striplen

Length to use in vectorization, this is the number of independent matrix multiplications that are to be done simultaneously.

integer, intent(in), optional :: subblockingWidth

The width of the subblocks used during the unrolled base exchange to ensure a better cache usage.

### public subroutine ply_fpt_exec(alph, gam, params, nIndeps)

Convert strip of coefficients of a modal representation in terms of Legendre polynomials to modal coefficients in terms of Chebyshev polynomials.

#### Arguments

TypeIntentOptionalAttributesName
real(kind=rk), intent(inout) :: alph(:)

Modal coefficients of the Legendre expansion. Size has to be: (1:params%n*nIndeps)

The direction which is to be transformed has to run fastest in the array.

real(kind=rk), intent(out) :: gam(:)

Modal coefficients of the Chebyshev expansion. Size has to be: (1:params%n*nIndeps)

type(ply_trafo_params_type), intent(inout) :: params

The parameters of the fast polynomial transformation.

integer, intent(in) :: nIndeps

Number of values that can be computed independently.

### public subroutine ply_fpt_single(alph, gam, params)

Convert strip of coefficients of a modal representation in terms of Legendre polynomials to modal coefficients in terms of Chebyshev polynomials.

#### Arguments

TypeIntentOptionalAttributesName
real(kind=rk), intent(inout) :: alph(params%n)

Modal coefficients of the Legendre expansion. Size has to be: params%n

The direction which is to be transformed has to run fastest in the array.

real(kind=rk), intent(out) :: gam(params%n)

Modal coefficients of the Chebyshev expansion. Size has to be: params%n

type(ply_trafo_params_type), intent(inout) :: params

The parameters of the fast polynomial transformation.

### public subroutine ply_fpt_exec_striped(nIndeps, alph, gam, params)

Convert coefficients of a modal representation in terms of Legendre polynomials to modal coefficients in terms of Chebyshev polynomials.

#### Arguments

TypeIntentOptionalAttributesName
integer, intent(in) :: nIndeps

Number of values that can be computed independently.

real(kind=rk), intent(in) :: alph(:)

Modal coefficients of the Legendre expansion. Size has to be: (1:params%n*indeps,nVars)

The direction which is to be transformed has to run fastest in the array.

real(kind=rk), intent(out) :: gam(:)

Modal coefficients of the Chebyshev expansion. Size has to be: (1:indeps*params%n,nVars)

Note, that the resulting array will have changed layout, and the transformed direction will run slowest in the array.

type(ply_trafo_params_type), intent(inout) :: params

The parameters of the fast polynomial transformation.

### private subroutine Copy_trafo_params(left, right)

#### Arguments

TypeIntentOptionalAttributesName
type(ply_trafo_params_type), intent(out) :: left

fpt to copy to

type(ply_trafo_params_type), intent(in) :: right

fpt to copy from

### private subroutine ply_calculate_coeff_strip(n, s, gam, matrix, alph, nDiagonals, block_offset, remainder, strip_lb, strip_ub, subblockingWidth)

#### Arguments

TypeIntentOptionalAttributesName
integer, intent(in) :: n
integer, intent(in) :: s
real(kind=rk), intent(inout) :: gam(:)

Modal coefficients of the Chebyshev expansion. Size has to be: (1:indeps*params%n,nVars)

Note, that the resulting array will have changed layout, and the transformed direction will run slowest in the array.

real(kind=rk), intent(in) :: matrix(:,:)

The arraz that holds the coefficients to calculate.

real(kind=rk), intent(in) :: alph(:)

Modal coefficients of the Legendre expansion. Size has to be: (1:params%n*indeps,nVars)

The direction which is to be transformed has to run fastest in the array.

integer, intent(in) :: nDiagonals

The number of diagonals to calculcate

integer, intent(in) :: block_offset

The offset of the block relative to the origin of the whole matrix.

integer, intent(in) :: remainder

The diagonals that are not covered by any block.

integer, intent(in) :: strip_lb

The lower bound of the strip to calculate.

integer, intent(in) :: strip_ub

The upper bound of the strip to calculate.

integer, intent(in) :: subblockingWidth

The subblocking width defines the size of the blocking of the diagonal strides, i.e. in y direction. This subblocking is used to get a better cache locality.