basix::precompute Namespace Reference

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, MDSPAN_IMPL_STANDARD_NAMESPACE::mdspan< const T, MDSPAN_IMPL_STANDARD_NAMESPACE::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, MDSPAN_IMPL_STANDARD_NAMESPACE::mdspan< const T, MDSPAN_IMPL_STANDARD_NAMESPACE::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].

1. 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.
2. 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.
3. 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.
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.
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.
6. 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_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}\).

Note

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

See apply_permutation().

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.

Returns

The three parts of a precomputed representation of the matrix. These are (as described above):

• A permutation (precomputed as in prepare_permutation());
• the vector \(D\);

apply_matrix()

template<typename T , typename E >

void basix::precompute::apply_matrix	(	std::span< const std::size_t >	v_size_t,
		MDSPAN_IMPL_STANDARD_NAMESPACE::mdspan< const T, MDSPAN_IMPL_STANDARD_NAMESPACE::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, this will apply the permutation to every block. The offset is set, this will start applying the permutation at the offsetth 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 vector created by 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,
		MDSPAN_IMPL_STANDARD_NAMESPACE::mdspan< const T, MDSPAN_IMPL_STANDARD_NAMESPACE::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}\)).

Note

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

See apply_matrix().