arb_mat.h

Returns a pointer to the entry in row i and column j of mat; the same as arb_mat_entry, as a function.

Initialises mat as the zero matrix with r rows and c columns.

Frees the memory used by mat.

Swap two matrices. arb_mat_swap swaps the structures (also when the dimensions differ); arb_mat_swap_entrywise swaps the entries one by one and requires equal dimensions, which also works for windows.

Returns the number of bytes allocated on the heap by x, excluding sizeof(arb_mat_struct).

Initialise window as a view of the submatrix of mat formed by rows r_1 \le i < r_2 and columns c_1 \le j < c_2, sharing its entries, and release it again (which frees nothing). Require 0 \le r_1 \le r_2 \le the number of rows and likewise for the columns. A window must not be used after mat has been cleared.

Set dest to src, which must have the same dimensions: exactly, or with the entries rounded to prec bits.

Sets mat to a random matrix of random density whose nonzero entries are generated by arb_randtest or arb_randtest_precise with midpoints of up to prec bits and exponents of up to mag_bits bits.

Set the square matrix mat to a random lower triangular matrix with exact entries and a positive diagonal (a Cholesky factor), or to the symmetric positive definite matrix L L^T for such a random L.

Print mat to standard output or to file, one row per line as [a, b, ...], each entry printed with arb_fprintd to digits decimal digits.

Returns whether the matrices have the same dimensions and identical balls as entries.

Returns whether the matrices have the same dimensions and each entry of mat1 overlaps the corresponding entry of mat2.

Return whether the matrices have the same dimensions and each entry of mat2 is contained in the corresponding entry of mat1.

Return whether mat1 and mat2 certainly represent the same matrix (same dimensions and all entries exact and equal), respectively certainly represent different matrices (different dimensions, or some pair of entries that do not overlap).

Return whether mat has zero rows or zero columns, respectively as many rows as columns.

Return whether all entries of mat have zero radius, are exactly zero, or are finite.

Return whether all entries below, above, respectively off the main diagonal of mat are exactly zero. The matrix need not be square.

Set all entries of mat to zero; the diagonal entries to one and the others to zero; all entries to one; or all entries to indeterminate balls [\operatorname{NaN} \pm \infty].

Sets mat to the Hilbert matrix with entries 1/(j + k + 1), rounded to prec bits. A wrapper of gr_mat_hilbert.

Sets mat to a Pascal matrix of binomial coefficients: for triangular 0 the symmetric matrix with entries \binom{j+k}{j} (the rows of Pascal's triangle as antidiagonals), for 1 the upper triangular matrix with entries \binom{k}{j}, for -1 the lower triangular matrix with entries \binom{j}{k}. A wrapper of gr_mat_pascal; the entries are computed by recurrences and are exact as long as they fit in prec bits.

Sets mat to the matrix of Stirling numbers s(j, k): unsigned of the first kind for kind 0, signed of the first kind for 1, of the second kind for 2. A wrapper of gr_mat_stirling; the entries are computed by recurrences.

Sets mat to the orthogonal DCT-II matrix of order n, the smaller dimension of mat,

which satisfies A^{-1} = A^T (if mat is not square, the formula is used for all entries). type is ignored and should be 0.

Sets dest to the transpose of src, exactly. Aliasing is allowed for square matrices.

Sets b to an upper bound for the infinity norm \max_i \sum_j |a_{ij}| of every matrix in A.

Set res to an enclosure, respectively an upper bound, of the Frobenius norm (\sum_{i,j} a_{ij}^2)^{1/2}.

Sets dest to -\mathit{src}, exactly.

Set res to the sum or difference of matrices of equal dimensions.

Set C (res) to the matrix product AB. The classical version computes each entry with arb_dot; the threaded version does the same with the rows of the output distributed over the threads of flint_get_num_threads(); the block version multiplies the midpoints exactly as integer matrices and bounds the radii separately. The default chooses automatically.

Sets C to the entrywise (Hadamard) product of matrices of equal dimensions.

Set the output to A^2 for a square matrix A. The classical version uses the symmetry of the products a_{ik} a_{kj} to save multiplications.

Sets res to A^e for a square matrix A, by binary exponentiation. A wrapper of gr_mat_pow_ui.

Approximate matrix multiplication: the input radii are ignored and the output is a floating-point approximation of the product of the midpoints. The radii of the output are not necessarily zero (they are whatever the algorithm used leaves), and carry no meaning.

Sets B to 2^c A, exactly.

Set B to B + cA.

Set B to cA.

Set B to A / c.

Set res to the product vA of the row vector v (of length the number of rows of A) by A, respectively the product Av of A by the column vector v, each entry computed by arb_dot. The underscore versions do not allow res to alias v.

Compute an LU decomposition PA = LU of an n \times n matrix A by Gaussian elimination with partial pivoting in ball arithmetic. Return nonzero if n pivots that do not contain zero are found, which proves that A is invertible; then L (unit lower triangular, without its diagonal) and U are stored in LU, and perm[i] is the row of A that became row i. Return zero, leaving perm and LU undefined, otherwise. LU may be the same object as A.

Solve LX = B for a lower triangular L, respectively UX = B for an upper triangular U, by substitution. Only the relevant triangle is read; the diagonal must not contain zero (not checked; otherwise the output is not finite). If unit is nonzero, the diagonal is taken to be all ones and is not read. The recursive versions use block substitution with matrix products; the defaults use them when B has at least 40 rows and 40 columns.

Solves AX = B given an LU decomposition PA = LU computed by arb_mat_lu (assumed valid, not checked). X may be the same object as B but not as LU.

Solve AX = B for a square matrix A and X, B with m columns. Return nonzero if A is proved invertible, in which case X contains the exact solution for every pair of matrices in the balls A and B; return zero (with X undefined) if m > 0 and this could not be proved. The lu version uses Gaussian elimination in ball arithmetic, the precond version preconditions with an approximate inverse [HS1967], and the default uses the lu version if n \le 4 or \mathit{prec} > 10n and the precond version otherwise. See Solving linear systems with ball matrices.

Solves AX = B given an approximate inverse R of A and an approximate solution T, by bounding the error of T (Theorem 10.2 of [Rum2010]). Returns nonzero if \|I - RA\|_\infty < 1 could be proved, which shows that A is invertible; X then contains the exact solution. Returns zero (with X undefined) otherwise, for example if R is not a good enough approximation.

Sets X to A^{-1} by solving AX = I with arb_mat_solve. Returns nonzero if A is proved invertible (and X contains the exact inverse) and zero otherwise, with X undefined.

Set det to an enclosure of the determinant of the square matrix A. The lu version uses Gaussian elimination and bounds the remaining submatrix by Hadamard's inequality when no more pivots can be found; the precond version multiplies A by approximate inverses of the factors of an approximate LU decomposition and encloses the determinant of the result, which is close to the identity [Rum2010]. The default uses direct formulas for n \le 3 and triangular matrices, the lu version if n \le 10 or \mathit{prec} > 10n, and the precond version otherwise; it throws an exception for a non-square matrix.

Approximately solve UX = B or LX = B for a triangular matrix with nonzero diagonal (not checked); with unit set, the diagonal is taken to be ones and not read.

Compute an approximate LU decomposition with partial pivoting, returning zero if an exactly zero pivot column is encountered, and solve with such a decomposition (assumed to come from arb_mat_approx_lu).

Approximately solve AX = B or invert A. Return zero if the approximate LU decomposition fails, nonzero otherwise; a nonzero return value says nothing about the accuracy of the result.

Compute the Cholesky decomposition A = LL^T with L lower triangular, by the Cholesky–Banachiewicz algorithm. Return nonzero if all pivots are certainly positive, which proves that every symmetric matrix represented by the lower triangle of A is positive definite; return zero if A is not positive definite or the precision is insufficient. The underscore version works in place and leaves the strict upper triangle undefined.

Solves AX = B given the Cholesky factor L of A computed by arb_mat_cho (assumed valid, not checked). X may be the same object as B but not as L.

Solves AX = B for a symmetric positive definite matrix A using arb_mat_cho. Returns nonzero if the decomposition succeeds, which proves that A is positive definite and that X contains the exact solution, and zero (if B has at least one column) otherwise, with X undefined.

Set X to A^{-1} for a symmetric positive definite A, given its Cholesky factor L (assumed valid) or computing it. The inverse is obtained by the method of [Kri2013], which computes fewer intermediate results than solving AX = I. arb_mat_spd_inv returns zero (with X undefined) if the Cholesky decomposition fails.

Compute the decomposition A = LDL^T with L unit lower triangular and D diagonal, stored in res: the strict lower triangle holds L and the diagonal holds D. Return nonzero if all entries of D are certainly positive, which proves that A is positive definite, and zero otherwise. The underscore versions work in place and leave the strict upper triangle undefined; the default uses Algorithm 4.1.2 of [GVL1996] (_arb_mat_ldl_golub_and_van_loan).

Solve AX = B, respectively compute A^{-1} (by the method of [Kri2013]), given the LDL^T decomposition of A encoded in L as computed by arb_mat_ldl (assumed valid). X may be the same object as B but not as L.

Set poly to the characteristic polynomial \det(xI - A) of the square n \times n matrix A (an exception is raised otherwise). The underscore version writes n + 1 coefficients. Wrappers of _gr_mat_charpoly_berkowitz, the division-free algorithm of Berkowitz with O(n^4) operations, which avoids divisions by balls that contain zero (Characteristic and minimal polynomials). A non-finite matrix gives indeterminate coefficients.

Set the n \times n matrix mat to the companion matrix of the polynomial of degree n (the underscore version reads n + 1 coefficients). If the leading coefficient cannot be inverted, the result is indeterminate. Wrappers of _gr_mat_companion and gr_mat_companion.

Sets S to the truncated exponential series \sum_{k=0}^{N-1} A^k / k!, using rectangular splitting with O(\sqrt{N}) matrix multiplications. The factorials are generated by their recurrence so that the scalars are small integers, and the divisions are postponed by computing partial factorials of length O(\sqrt{N}) [Joh2014b].

Sets B to the exponential \exp(A) = \sum_{k \ge 0} A^k/k! of the square matrix A.

Sets trace to the sum of the diagonal entries of the square matrix mat.

Set res to the product of the diagonal entries of mat; the underscore version takes the product of the entries with index a \le i < b (which must be in range; the empty product is 1).

Return the number of entries of mat that are exactly zero (arb_is_zero), respectively not exactly zero.

Sets B to the matrix of midpoints of A, with zero radii.

Adds err to the radius of every entry of mat.

Attempts to set B to the symmetric positive definite integer matrix obtained by rounding 2^{\mathit{prec}} A entrywise, using the midpoint of the intersection of a_{jk} and a_{kj} for both entries. Returns 1 on success, and 0 (with B undefined) if some pair a_{jk}, a_{kj} does not overlap or the rounded matrix is not positive definite. As for arf_get_fmpz, huge exponents can exhaust memory.

Given a symmetric positive definite matrix A, sets U to an invertible integer matrix such that U A U^T is close to being LLL-reduced with the parameters \delta, \eta: if arb_mat_spd_get_fmpz_mat succeeds at precision prec, fmpz_lll is applied to the rounded Gram matrix, and otherwise U is the identity. If A is not finite, U is not modified.

Given a symmetric positive definite matrix A (only the lower triangle is read; neither property is checked), returns nonzero if A is certainly LLL-reduced for the parameters \delta, \eta as a Gram matrix: the Gram–Schmidt coefficients satisfy |\mu_{j,k}| \le \eta and \delta \|b_{k-1}^*\|^2 \le \|b_k^*\|^2 + \mu_{k,k-1}^2 \|b_{k-1}^*\|^2.

Swaps rows r and s of mat and, if perm is not NULL, the entries r and s of perm.

Returns the index i with \mathit{start\_row} \le i < \mathit{end\_row} of the entry in column c that does not contain zero and has the largest midpoint in absolute value, or -1 if every entry in that range contains zero.

Adds to the radii of C an upper bound for the product of the nonnegative matrices A (\mathit{ar} \times \mathit{ac}, row-major array) and B (\mathit{ac} \times \mathit{bc}, given as a column-major array). The exponents of all entries must be small (this is not checked), which arb_mat_mul_block ensures. Products are computed in double precision when the exponent range permits and the result is multiplied by 1 + 2.01(\mathit{ac} + 1) 2^{-52} to account for rounding errors.