Taylor models bounding

#include <audi/taylor_model_bounding.hpp>

std::vector<std::vector<int>> audi::generate_combinations(const std::vector<int> &limits, bool cap_sum_indices = false)

Function for the combinations of given maximum indices.

Function gives either all possible combinations given a list of maximum limits, or (if cap_sum_indices = true), the total sum of the indices of a given combination are limited by the maximum index given by the input maximum indices.

The cap_sum_indices defaults to off, and is generally only set to true when retrieving the polynomial from a gdual (to prevent a monomial from exceeding the degree of the gdual).

Parameters:

limits – std::vector<int> the maximum degree limits per dimension
cap_sum_indices – bool whether or not to limit the combinations by the maximum of the maximum degrees

Returns:

a nested vector of combinations MxN where M is the number of combinations, and N is the number of dimensions that then have a given degree.

template<typename T> std::pair<std::vector<T>, std::vector<std::vector<int>>> audi::get_poly(const audi::gdual<T> &tpol)

Extract the polynomial representation of a gdual.

Retrieves the coefficients and corresponding exponents of the nonzero monomials contained in a gdual.

For univariate polynomials, the exponents are generated from 0 up to the degree of the gdual. For multivariate polynomials, all valid exponent combinations are generated such that the sum of indices does not exceed the degree of the gdual.

Template Parameters:

T – numeric type of the coefficients (e.g., double, long double)

Parameters:

tpol – the input truncated polynomial (audi::gdual<T>) from which the polynomial representation is extracted

Returns:

a pair (coeffs, exponents) where:

coeffs is a vector of coefficients of type T
exponents is a vector of integer vectors, each containing the exponents of a corresponding monomial term in coeffs

unsigned int audi::get_ndim(const std::vector<double> &coeffs, const std::vector<std::vector<int>> &exps)

Get the number of dimensions in a polynomial representation.

Determines the dimensionality of a polynomial given its coefficients and exponent matrix.

Parameters:

coeffs – vector of coefficients
exps – vector of exponent rows, each row representing the exponents of a monomial term

Throws:

std::invalid_argument – if the number of coefficients does not match the number of exponent rows, or if exponent rows are inconsistent in length

Returns:

the number of dimensions (variables) of the polynomial

std::vector<int> audi::get_max_degrees(const std::vector<std::vector<int>> &exponents, unsigned int ndim)

Get the maximum degree in each dimension.

Computes the maximum degree reached in each variable of a polynomial, given its exponent representation. The result is a vector of per-dimension maximum exponents.

For multivariate polynomials (ndim > 1), returns one maximum per variable.
For univariate polynomials (ndim == 1), returns a single-element vector containing the maximum degree.
For zero-dimensional input (ndim == 0), returns {0}.

Parameters:

exponents – vector of exponent rows, each row representing the exponents of a monomial term
ndim – number of dimensions (variables) of the polynomial

Throws:

std::invalid_argument – if exponent rows are inconsistent in length or if ndim < 0

Returns:

vector of maximum degrees, one per dimension

template<typename T> T audi::get_coefficient(const std::vector<int> &comb, const std::vector<T> &coeffs, const std::vector<std::vector<int>> &exps)

Retrieve a coefficient by exponent combination.

Finds the coefficient of a monomial term in a polynomial representation given its exponent combination. If the exponent vector is not found in the exponent list, the coefficient is assumed to be zero.

Template Parameters:

T – numeric type of the coefficients (e.g., double, long double)

Parameters:

comb – exponent combination to look up (multi-index)
coeffs – vector of coefficients
exps – vector of exponent rows corresponding to the coefficients

Returns:

the coefficient corresponding to the given exponent combination, or 0 if not found

template<typename T> std::vector<std::vector<T>> audi::get_a_matrix_vec(const std::vector<T> &coeffs, const std::vector<std::vector<int>> &exps)

Construct the coefficient matrix representation of a polynomial.

Builds a 2D matrix A from a polynomial given in coefficient–exponent form.

The rows of A correspond to the degree of the first variable.
The columns of A represent the flattened index of all remaining variables.

More precisely, the matrix has size:

(max_degrees[0] + 1) rows
∏_{i=1..ndim-1} (max_degrees[i] + 1) columns

Each entry A[i][j] contains the coefficient of the monomial with multi-index (i, *rest) mapped into column index j.

Template Parameters:

T – numeric type of the coefficients (e.g., double, long double)

Parameters:

coeffs – vector of coefficients
exps – vector of exponent rows corresponding to the coefficients

Returns:

a 2D vector A holding the coefficients arranged as a matrix

unsigned long long audi::binomial(int n, int k)

Compute a binomial coefficient.

Computes the number of ways to choose k elements from n without repetition:

\[ \binom{n}{k} = \frac{n!}{k!(n-k)!} \]

Uses an iterative approach that avoids intermediate factorial computation and exploits the symmetry relation \(\binom{n}{k} = \binom{n}{n-k}\).

Parameters:

n – nonnegative integer
k – integer satisfying \(0 \leq k \leq n\)

Returns:

the binomial coefficient as an unsigned long long

unsigned long long audi::binom_product(const std::vector<int> &n, const std::vector<int> &k)

Compute the product of binomial coefficients.

Given two vectors n and k of the same length, computes the product

\[ \prod_{i=1}^m \binom{n_i}{k_i} \]

where \(m = n.size() = k.size()\).

Parameters:

n – vector of nonnegative integers
k – vector of nonnegative integers of the same size as n

Throws:

std::runtime_error – if n and k do not have the same size

Returns:

the product of binomial coefficients as an unsigned long long

std::vector<std::vector<double>> audi::get_d_acc_matrix(const std::vector<int> &max_degrees, int dim)

Construct the diagonal accumulation matrix for a given dimension.

Builds a diagonal matrix of size (n+1) x (n+1) where n is the maximum degree in the specified dimension. The diagonal entries are given by

\[ (D_{\text{acc}})_{kk} = \frac{1}{\binom{n}{k}}, \quad k = 0, \dots, n. \]

Parameters:

max_degrees – vector of maximum degrees per dimension
dim – the dimension index to construct the matrix for

Throws:

std::out_of_range – if dim is negative or not less than max_degrees.size()

Returns:

a (n+1) x (n+1) diagonal matrix of scaling factors

std::vector<std::vector<double>> audi::get_d_matrix(const std::vector<int> &max_degrees, int dim, double t)

Construct the diagonal evaluation matrix for a given dimension.

Builds a diagonal matrix of size (n+1) x (n+1) where n is the maximum degree in the specified dimension. The diagonal entries are

\[\begin{split} (D)_{kk} = \begin{cases} 1 & \text{if } k = 0, \\ t^k & \text{if } k \geq 1. \end{cases} \end{split}\]

This corresponds to evaluating powers of t up to the maximum degree.

Parameters:

max_degrees – vector of maximum degrees per dimension
dim – the dimension index to construct the matrix for
t – the evaluation point

Returns:

a (n+1) x (n+1) diagonal matrix of powers of t

std::vector<std::vector<double>> audi::get_lower_pascal_matrix(const std::vector<int> &max_degrees, int dim)

Construct the lower Pascal matrix for a given dimension.

Builds a lower-triangular Pascal matrix of size (n+1) x (n+1) where n is the maximum degree in the specified dimension. The entries are defined by

\[\begin{split} P_{ij} = \begin{cases} \binom{i}{j} & \text{if } j \leq i, \\ 0 & \text{if } j > i. \end{cases} \end{split}\]

This matrix encodes binomial coefficients in triangular form.

Parameters:

max_degrees – vector of maximum degrees per dimension
dim – the dimension index to construct the matrix for

Returns:

a (n+1) x (n+1) lower-triangular Pascal matrix

std::unordered_map<int, std::pair<std::string, int_d>> audi::build_dim_to_var_map(const var_map_i &domain)

Build a dimension-to-variable mapping.

Constructs a mapping from dimension index to (variable name, interval) pairs based on the given domain.

Parameters:: domain – a mapping from variable name to its interval
Returns:: unordered map from dimension index (0-based) to (var_name, interval)

std::vector<std::vector<double>> audi::get_q_matrix(const std::vector<int> &max_degrees, const var_map_i &domain, int dim)

Construct the Q-matrix for a given dimension.

Builds the Q-matrix (from Titi (2018)) associated with a given variable domain and maximum degree. This function basically maps the domain to [0, 1] for the subsequent matrix compuations.

If the lower bound a of the variable interval is zero, the matrix reduces to a diagonal matrix \(D(b)\) where diagonal entries are \(b^k\).
Otherwise, the matrix is constructed as

\[ Q = D\!\left(\frac{b-a}{a}\right) \cdot P^T \cdot D(a), \]

where:

P is the lower Pascal matrix of size (n+1) x (n+1)
D(x) is the diagonal evaluation matrix with entries \(x^k\)
n = max_degrees[dim]

Parameters:

max_degrees – vector of maximum degrees per dimension
domain – variable domain mapping (variable name → interval)
dim – the dimension index

Throws:

std::out_of_range – if dim is invalid
std::runtime_error – if dim is not found in the domain map

Returns:

a (n+1) x (n+1) Q-matrix for the specified dimension

std::vector<int> audi::shift_indices(const std::vector<int> &indices)

Shift indices cyclically to the left.

Rotates a vector of indices one position to the left.

Example:

Input: {2, 3, 4}
Output: {3, 4, 2}

Parameters:: indices – vector of indices
Returns:: shifted vector with elements rotated left by one

std::vector<std::vector<double>> audi::get_cycled_2d_array_vec(const std::vector<std::vector<double>> &arr, const std::vector<int> &max_degrees, int dim)

Cycle a 2D array representation along dimensions.

Rearranges a 2D array representation of polynomial coefficients by cycling the role of dimensions. This operation effectively rotates the “leading” dimension used in the row index.

If the number of dimensions is one, the input array is returned unchanged.

Parameters:

arr – 2D array of shape (rows x cols) representing coefficients
max_degrees – vector of maximum degrees per dimension
dim – the current dimension index to cycle

Returns:

a 2D array with cycled dimension ordering; new shape is (max_degrees[next_dim]+1, ∏_{m≠next_dim} (max_degrees[m]+1))

std::vector<std::vector<double>> audi::get_titi_base_lambda_generalbox(const std::vector<double> &coeffs, const std::vector<std::vector<int>> &exps, const var_map_i &domain)

Construct the base Lambda matrix for the general box domain.

Computes the base \(\Lambda\) matrix representation for a polynomial expressed in the Bernstein basis over a general (non-unit) domain.

Iteratively constructs \(\Lambda\) by applying, for each dimension:

diagonal accumulation matrix \(D_{\text{acc}}\)
Q-matrix (interval scaling/translation)
cycling of dimensions

Intermediate Lambda matrices are stored in a dictionary by dimension.

Parameters:

coeffs – vector of polynomial coefficients
exps – exponent matrix (each row a monomial’s exponents)
domain – variable domain mapping (variable name → interval)

Returns:

the final \(\Lambda\) matrix of shape ∏ (deg_i+1) x ∏ (deg_i+1) for all dimensions

template<typename T> std::vector<std::vector<T>> audi::get_titi_bernstein_patch_ndim_generalbox(const std::vector<T> &coeffs, const std::vector<std::vector<int>> &exps, const var_map_i &domain)

Construct the Bernstein patch for an n-dimensional polynomial (with a general box domain)

Computes the Bernstein patch (control points) of a multivariate polynomial defined on a general hyper-rectangular domain. Uses audi::get_titi_base_lambda_generalbox(const std::vector<double>&coeffs, const std::vector<std::vector<int>> &exps, const var_map_i &domain) followed by repeated multiplication with lower Pascal matrices and cycling of dimensions.

Template Parameters:

T – numeric type of coefficients

Parameters:

coeffs – vector of polynomial coefficients
exps – exponent matrix (each row a monomial’s exponents)
domain – variable domain mapping (variable name → interval)

Returns:

a 2D array of control points in the Bernstein basis, reshaped according to (deg_1+1, deg_2+1, ..., deg_n+1)