Matrix and permutation pre-computation. More...

Functions
void	prepare_permutation (std::span< std::size_t > perm)
	Prepare a permutation.

template<typename E >
void	apply_permutation (std::span< const std::size_t > perm, std::span< E > data, std::size_t offset=0, std::size_t n=1)
	Apply a (precomputed) permutation \(v = P u\).

template<typename E >
void	apply_permutation_mapped (std::span< const std::size_t > perm, std::span< E > data, std::span< const int > emap, std::size_t n=1)
	Permutation of mapped data.

template<typename E >
void	apply_inv_permutation_right (std::span< const std::size_t > perm, std::span< E > data, std::size_t offset=0, std::size_t n=1)
	Apply a (precomputed) permutation to some transposed data.

template<std::floating_point T>
std::vector< std::size_t >	prepare_matrix (std::pair< std::vector< T >, std::array< std::size_t, 2 > > &A)
	Prepare a square matrix.

template<typename T , typename E >
void	apply_matrix (std::span< const std::size_t > v_size_t, md::mdspan< const T, md::dextents< std::size_t, 2 > > M, std::span< E > data, std::size_t offset=0, std::size_t n=1)
	Apply a (precomputed) matrix.

template<typename T , typename E >
void	apply_tranpose_matrix_right (std::span< const std::size_t > v_size_t, md::mdspan< const T, md::dextents< std::size_t, 2 > > M, std::span< E > data, std::size_t offset=0, std::size_t n=1)
	Apply a (precomputed) matrix to some transposed data.

Detailed Description

Matrix and permutation pre-computation.

Function Documentation

◆ prepare_permutation()

void basix::precompute::prepare_permutation ( std::span< std::size_t > perm )

Prepare a permutation.

Computes representation of the permutation that allows the permutation to be applied without any temporary memory assignment. In pseudo code, this function does the following:

FOR index, entry IN perm:
    new_index = entry
    WHILE new_index < index:
        new_index = perm[new_index]
    OUT[index] = new_index

Example

Consider the permutation P = [1, 4, 0, 5, 2, 3].

First, we look at the 0th entry. P[0] is 1. This is greater than 0, so the 0th entry of the output is 1.
Next, we look at the 1st entry. P[1] is 4. This is greater than 1, so the 1st entry of the output is 4.
Next, we look at the 2nd entry. P[2] is 0. This is less than 2, so we look at P[0]. P[0] is 1. This is less than 2, so we look at P[1]. P[1] is 4. This is greater than 2, so the 2nd entry of the output is 4.
Next, we look at the 3rd entry. P[3] is 5. This is greater than 3, so the 3rd entry of the output is 5.
Next, we look at the 4th entry. P[4] is 2. This is less than 4, so we look at P[2]. P[2] is 0. This is less than 4, so we look at P[0]. P[0] is 1. This is less than 4, so we look at P[1]. P[1] is 4. This is greater than (or equal to) 4, so the 4th entry of the output is 4.
Next, we look at the 5th entry. P[5] is 3. This is less than 5, so we look at P[3]. P[3] is 5. This is greater than (or equal to) 5, so the 5th entry of the output is 5.

Hence, the output of this function in this case is [1, 4, 4, 5, 4, 5].

For an example of how the permutation in this form is applied, see apply_permutation().

Parameters

[in,out] perm A permutation.

◆ apply_permutation()

template<typename E >

void basix::precompute::apply_permutation	(	std::span< const std::size_t >	perm,
		std::span< E >	data,
		std::size_t	offset = `0`,
		std::size_t	n = `1`
	)

Apply a (precomputed) permutation \(v = P u\).

This uses the representation returned by prepare_permutation() to apply a permutation without needing any temporary memory. In pseudo code, this function does the following:

FOR index, entry IN perm:

SWAP(data[index], data[entry])

If n is set, this will apply the permutation to every block. The offset is set, this will start applying the permutation at the offsetth block.

Example

Consider the permutation P = [1, 4, 0, 5, 2, 3]. In the documentation of prepare_permutation(), we saw that the precomputed representation of this permutation is P2 = [1, 4, 4, 5, 4, 5]. In this example, we look at how this representation can be used to apply this permutation to the array A = [a, b, c, d, e, f].

P2[0] is 1, so we swap A[0] and A[1]. After this, A = [b, a, c, d, e, f].
P2[1] is 4, so we swap A[1] and A[4]. After this, A = [b, e, c, d, a, f].
P2[2] is 4, so we swap A[2] and A[4]. After this, A = [b, e, a, d, c, f].
P2[3] is 5, so we swap A[3] and A[5]. After this, A = [b,e,a,f,c,d].
P2[4] is 4, so we swap A[4] and A[4]. This changes nothing.
P2[5] is 5, so we swap A[5] and A[5]. This changes nothing.

Therefore the result of applying this permutation is [b, e, a, f, c, d] (which is what we get if we apply the permutation directly).

Note: This function is designed to be called at runtime, so its performance is critical.

Parameters

[in]	perm	A permutation in precomputed form (as returned by prepare_permutation()).
[in,out]	data	The data to apply the permutation to. It has shape `(m, n)` (uses row-major storage), where the permutation matrix has shape `(m, m)`.
[in]	offset	The position in the data to start applying the permutation.
[in]	n	The block size of the data.

◆ apply_permutation_mapped()

template<typename E >

void basix::precompute::apply_permutation_mapped	(	std::span< const std::size_t >	perm,
		std::span< E >	data,
		std::span< const int >	emap,
		std::size_t	n = `1`
	)

Permutation of mapped data.

Parameters

[in]	perm	A permutation in precomputed form (as returned by prepare_permutation()).
[in,out]	data	The data to apply the permutation to.
[in]	emap	Map from the permutation's local entity numbering to the position of that entity's block in `data`.
[in]	n	The block size of the data.

◆ apply_inv_permutation_right()

template<typename E >

void basix::precompute::apply_inv_permutation_right	(	std::span< const std::size_t >	perm,
		std::span< E >	data,
		std::size_t	offset = `0`,
		std::size_t	n = `1`
	)

Apply a (precomputed) permutation to some transposed data.

Applies \(v = u P^{T}\).

Unlike apply_permutation(), where the n values of a block are interleaved (block b occupies data[n*(offset+i)+b] for i in the permutation), here each of the n blocks is contiguous: block b occupies data[data_size*b+offset+i] for i in the permutation.

Note: This function is designed to be called at runtime, so its performance is critical.

Parameters

[in]	perm	A permutation in precomputed form (as returned by prepare_permutation()).
[in,out]	data	The data to apply the permutation to.
[in]	offset	The position in each block to start applying the permutation.
[in]	n	The number of blocks in the data.

◆ prepare_matrix()

template<std::floating_point T>

std::vector< std::size_t > basix::precompute::prepare_matrix ( std::pair< std::vector< T >, std::array< std::size_t, 2 > > & A )

Prepare a square matrix.

Computes the LU decomposition of the transpose of the matrix.

This function returns the permutation \(P\)P in the representation \(PA^t=LU\) (precomputed as in prepare_permutation()). The LU decomposition of \(A^t\) is computed in-place

For an example of how the permutation in this form is applied, see apply_matrix().

Parameters

[in,out] A The matrix data. On exit, holds the LU factors of \(A^t\) in place.

Returns: The permutation \(P\) (in prepared form, as produced by prepare_permutation()) satisfying \(PA^t=LU\).

◆ apply_matrix()

template<typename T , typename E >

void basix::precompute::apply_matrix	(	std::span< const std::size_t >	v_size_t,
		md::mdspan< const T, md::dextents< std::size_t, 2 > >	M,
		std::span< E >	data,
		std::size_t	offset = `0`,
		std::size_t	n = `1`
	)

Apply a (precomputed) matrix.

This uses the representation returned by prepare_matrix() to apply a matrix without needing any temporary memory. In pseudo code, this function does the following:

INPUT perm, mat, data
apply_permutation(perm, data)
FOR i IN RANGE(dim):
    FOR j IN RANGE(i+1, dim):
        data[i] += mat[i, j] * data[j]
FOR i IN RANGE(dim - 1, -1, -1):
    data[i] *= M[i, i]
    FOR j in RANGE(i):
        data[i] += mat[i, j] * data[j]

If n is set, the permutation is applied to every block. If offset is set, application starts at the offset-th block.

Note: This function is designed to be called at runtime, so its performance is critical.

Parameters

[in]	v_size_t	A permutation, as computed by prepare_matrix().
[in]	M	The matrix data left in `A` after the call to prepare_matrix().
[in,out]	data	The data to apply the permutation to
[in]	offset	The position in the data to start applying the permutation
[in]	n	The block size of the data

◆ apply_tranpose_matrix_right()

template<typename T , typename E >

void basix::precompute::apply_tranpose_matrix_right	(	std::span< const std::size_t >	v_size_t,
		md::mdspan< const T, md::dextents< std::size_t, 2 > >	M,
		std::span< E >	data,
		std::size_t	offset = `0`,
		std::size_t	n = `1`
	)

Apply a (precomputed) matrix to some transposed data.

Computes \(v^{T} = M u^{T}\) (or equivalently \(v = u M^{T}\)).

Uses the same block-contiguous data layout as apply_inv_permutation_right() (rather than the interleaved layout used by apply_matrix()): block b occupies data[data_size*b+offset+i].

Note: This function is designed to be called at runtime, so its performance is critical.

Parameters

[in]	v_size_t	A permutation, as computed by prepare_matrix().
[in]	M	The matrix data left in `A` after the call to prepare_matrix().
[in,out]	data	The data to apply the permutation to.
[in]	offset	The position in each block to start applying the permutation.
[in]	n	The number of blocks in the data.

Functions

Detailed Description

Function Documentation

◆ prepare_permutation()

Example

◆ apply_permutation()

Example

◆ apply_permutation_mapped()

◆ apply_inv_permutation_right()

◆ prepare_matrix()

◆ apply_matrix()

◆ apply_tranpose_matrix_right()