arrayfire.lapack module¶
Dense Linear Algebra functions (solve, inverse, etc).
- arrayfire.lapack.cholesky(A, is_upper=True)[source]¶
Cholesky decomposition
- Parameters:
- A: af.Array
A 2 dimensional, symmetric, positive definite matrix.
- is_upper: optional: bool. default: True
Specifies if output R is upper triangular (if True) or lower triangular (if False).
- Returns:
- (R,info): tuple of af.Array, int.
R - triangular matrix.
info - 0 if decomposition sucessful.
- Note
- The original matrix A can be reconstructed using the outputs in the following manner.
>>> A = af.matmulNT(R, R) #if R is upper triangular ..
- arrayfire.lapack.cholesky_inplace(A, is_upper=True)[source]¶
In place Cholesky decomposition.
- Parameters:
- A: af.Array
a 2 dimensional, symmetric, positive definite matrix.
Trinangular matrix on exit.
- is_upper: optional: bool. default: True.
Specifies if output R is upper triangular (if True) or lower triangular (if False).
- Returns:
- infoint.
0 if decomposition sucessful.
- arrayfire.lapack.det(A)[source]¶
Determinant of a matrix.
- Parameters:
- A: af.Array
A 2 dimensional arrayfire array
- Returns:
- res: scalar
Determinant of the matrix.
- arrayfire.lapack.inverse(A, options=MATPROP.NONE)[source]¶
Invert a matrix.
- Parameters:
- A: af.Array
A 2 dimensional arrayfire array
- options: optional: af.MATPROP. default: af.MATPROP.NONE.
Additional options to speed up computations.
Currently needs to be one of af.MATPROP.NONE.
- Returns:
- AI: af.Array
A 2 dimensional array that is the inverse of A
- arrayfire.lapack.is_lapack_available()[source]¶
Function to check if the arrayfire library was built with lapack support.
- arrayfire.lapack.lu(A)[source]¶
LU decomposition.
- Parameters:
- A: af.Array
A 2 dimensional arrayfire array.
- Returns:
- (L,U,P): tuple of af.Arrays
L - Lower triangular matrix.
U - Upper triangular matrix.
P - Permutation array.
- arrayfire.lapack.lu_inplace(A, pivot='lapack')[source]¶
In place LU decomposition.
- Parameters:
- A: af.Array
a 2 dimensional arrayfire array on entry.
Contains L in the lower triangle on exit.
Contains U in the upper triangle on exit.
- Returns:
- P: af.Array
Permutation array.
- arrayfire.lapack.norm(A, norm_type=NORM.VECTOR_2, p=1.0, q=1.0)[source]¶
Norm of an array or a matrix.
- Parameters:
- A: af.Array
A 1 or 2 dimensional arrayfire array
- norm_type: optional: af.NORM. default: af.NORM.EUCLID.
Type of norm to be calculated.
- p: scalar. default 1.0.
Used only if norm_type is one of af.NORM.VECTOR_P, af.NORM_MATRIX_L_PQ
- q: scalar. default 1.0.
Used only if norm_type is af.NORM_MATRIX_L_PQ
- Returns:
- res: scalar
norm of the input
- arrayfire.lapack.qr(A)[source]¶
QR decomposition.
- Parameters:
- A: af.Array
A 2 dimensional arrayfire array.
- Returns:
- (Q,R,T): tuple of af.Arrays
Q - Orthogonal matrix.
R - Upper triangular matrix.
T - Vector containing additional information to solve a least squares problem.
- arrayfire.lapack.qr_inplace(A)[source]¶
In place QR decomposition.
- Parameters:
- A: af.Array
a 2 dimensional arrayfire array on entry.
Packed Q and R matrices on exit.
- Returns:
- T: af.Array
Vector containing additional information to solve a least squares problem.
- arrayfire.lapack.rank(A, tol=1e-05)[source]¶
Rank of a matrix.
- Parameters:
- A: af.Array
A 2 dimensional arrayfire array
- tol: optional: scalar. default: 1E-5.
Tolerance for calculating rank
- Returns:
- r: int
Rank of A within the given tolerance
- arrayfire.lapack.solve(A, B, options=MATPROP.NONE)[source]¶
Solve a system of linear equations.
- Parameters:
- A: af.Array
A 2 dimensional arrayfire array representing the coefficients of the system.
- B: af.Array
A 1 or 2 dimensional arrayfire array representing the constants of the system.
- options: optional: af.MATPROP. default: af.MATPROP.NONE.
Additional options to speed up computations.
Currently needs to be one of af.MATPROP.NONE, af.MATPROP.LOWER, af.MATPROP.UPPER.
- Returns:
- X: af.Array
A 1 or 2 dimensional arrayfire array representing the unknowns in the system.
- arrayfire.lapack.solve_lu(A, P, B, options=MATPROP.NONE)[source]¶
Solve a system of linear equations, using LU decomposition.
- Parameters:
- A: af.Array
A 2 dimensional arrayfire array representing the coefficients of the system.
This matrix should be decomposed previously using lu_inplace(A).
- P: af.Array
Permutation array.
This array is the output of an earlier call to lu_inplace(A)
- B: af.Array
A 1 or 2 dimensional arrayfire array representing the constants of the system.
- Returns:
- X: af.Array
A 1 or 2 dimensional arrayfire array representing the unknowns in the system.