Module maxvol: compute the maximal-volume submatrix

Package teneva, module maxvol: maxvol-like algorithms.

This module contains the implementation of the maxvol algorithm (function “maxvol”) and rect_maxvol algorithm (function “maxvol_rect”) for matrices. The corresponding functions find in a given matrix square and rectangular maximal-volume submatrix, respectively (for the case of square submatrix, it has approximately the maximum value of the modulus of the determinant).



teneva_jax.maxvol.maxvol(A, e=1.05, k=100)[source]

Compute the maximal-volume submatrix for given tall matrix.

Parameters:

  • A (jnp.ndarray) – tall matrix of the shape [n, r] (n > r).

  • e (float) – accuracy parameter (should be >= 1). If the parameter is equal to 1, then the maximum number of iterations will be performed until true convergence is achieved. If the value is greater than one, the algorithm will complete its work faster, but the accuracy will be slightly lower (in most cases, the optimal value is within the range of 1.01 - 1.1).

  • k (int) – maximum number of iterations (should be >= 1).

Returns:

the row numbers I containing the maximal volume submatrix in the form of 1D array of length r and coefficient matrix B in the form of 2D array of shape [n, r], such that A = B A[I, :] and A (A[I, :])^{-1} = B.

Return type:

(jnp.ndarray, jnp.ndarray)

Note

The description of the basic implementation of this algorithm is presented in the work: Goreinov S., Oseledets, I., Savostyanov, D., Tyrtyshnikov, E., Zamarashkin, N. “How to find a good submatrix”. Matrix Methods: Theory, Algorithms And Applications: Dedicated to the Memory of Gene Golub (2010): 247-256.

Examples:

n = 5000                           # Number of rows
r = 50                             # Number of columns
rng, key = jax.random.split(rng)
A = jax.random.normal(key, (n, r)) # Random tall matrix

e = 1.01  # Accuracy parameter
k = 500   # Maximum number of iterations

# Compute row numbers and coefficient matrix:
I, B = teneva.maxvol(A, e, k)

# Maximal-volume square submatrix:
C = A[I, :]

print(f'|Det C|        : {jnp.abs(jnp.linalg.det(C)):-10.2e}')
print(f'Max |B|        : {jnp.max(jnp.abs(B)):-10.2e}')
print(f'Max |A - B C|  : {jnp.max(jnp.abs(A - B @ C)):-10.2e}')
print(f'Selected rows  : {I.size:-10d} > ', jnp.sort(I))

# >>> ----------------------------------------
# >>> Output:

# |Det C|        :   1.29e+40
# Max |B|        :   1.00e+00
# Max |A - B C|  :   9.10e-15
# Selected rows  :         50 >  [ 120  315  571  798 1037 1049 1098 1250 1286 1304 1309 1419 1444 1604
#  1610 1766 1835 1887 1956 2085 2324 2327 2458 2602 2817 2926 3119 3242
#  3322 3497 3508 3705 3715 3722 3743 3771 3811 3904 3973 4068 4101 4165
#  4310 4321 4399 4439 4544 4771 4871 4938]
#


teneva_jax.maxvol.maxvol_rect(A, e=1.1, dr_min=0, dr_max=None, e0=1.05, k0=10)[source]

Compute the maximal-volume rectangular submatrix for given tall matrix.

DRAFT (works with error now) !!!

Within the framework of this function, the original maxvol algorithm is first called (see function “maxvol”). Then additional rows of the matrix are greedily added until the maximum allowed number is reached or until convergence.

Parameters:

  • A (jnp.ndarray) – tall matrix of the shape [n, r] (n > r).

  • e (float) – accuracy parameter.

  • dr_min (int) – minimum number of added rows (should be >= 0 and <= n-r).

  • dr_max (int) – maximum number of added rows (should be >= 0). If the value is not specified, then the number of added rows will be determined by the precision parameter e, while the resulting submatrix can even has the same size as the original matrix A. If r + dr_max is greater than n, then dr_max will be set such that r + dr_max = n.

  • e0 (float) – accuracy parameter for the original maxvol algorithm (should be >= 1). See function “maxvol” for details.

  • k0 (int) – maximum number of iterations for the original maxvol algorithm (should be >= 1). See function “maxvol” for details.

Returns:

the row numbers I containing the rectangular maximal-volume submatrix in the form of 1D array of length r + dr, where dr is a number of additional selected rows (dr >= dr_min and dr <= dr_max) and coefficient matrix B in the form of 2D array of shape [n, r+dr], such that A = B A[I, :].

Return type:

(jnp.ndarray, jnp.ndarray)

Note

The description of the basic implementation of this algorithm is presented in the work: Mikhalev A, Oseledets I., “Rectangular maximum-volume submatrices and their applications.” Linear Algebra and its Applications 538 (2018): 187-211.

Examples:

n = 5000                           # Number of rows
r = 50                             # Number of columns
rng, key = jax.random.split(rng)
A = jax.random.normal(key, (n, r)) # Random tall matrix

e = 1.01    # Accuracy parameter
dr_min = 2  # Minimum number of added rows
dr_max = 8  # Maximum number of added rows
e0 = 1.05   # Accuracy parameter for the original maxvol algorithm
k0 = 50     # Maximum number of iterations for the original maxvol algorithm

THIS IS DRAFT !!!

# Row numbers and coefficient matrix:
I, B = teneva.maxvol_rect(A, e,
    dr_min, dr_max, e0, k0)

# Maximal-volume rectangular submatrix:
C = A[I, :]

# >>> ----------------------------------------
# >>> Output:

# /Users/andrei/opt/anaconda3/envs/teneva_jax/lib/python3.8/site-packages/jax-0.4.8-py3.8.egg/jax/_src/ops/scatter.py:89: FutureWarning: scatter inputs have incompatible types: cannot safely cast value from dtype=float64 to dtype=int64. In future JAX releases this will result in an error.
#   warnings.warn("scatter inputs have incompatible types: cannot safely cast "
#