Factorizations: LU, QR, Cholesky, singular values, eigenvalues, Hessenberg

Functions
array	lu (const array &in)
	LU factorization (packed). More...
void	lu (array &lower, array &upper, const array &in)
	LU factorization. More...
void	lu (array &lower, array &upper, array &pivot, const array &in)
	LU factorization (with pivoting). More...
array	qr (const array &in)
	QR factorization (packed). More...
void	qr (array &q, array &r, const array &in)
	QR factorization. More...
void	qr (array &q, array &r, array &tau, const array &in)
	QR factorization with tau. More...
array	cholesky (const array &X, bool is_upper=true)
	Cholesky decomposition ("Y^T * Y == X"). More...
array	cholesky (unsigned &info, const array &X, bool is_upper=true)
	Cholesky decomposition ("Y^T * Y == X"). More...
array	hessenberg (const array &in)
	Hessenberg matrix form. More...
void	hessenberg (array &h, array &q, const array &in)
	Hessenberg matrix h with unitary permutation matrix q. More...
array	eigen (const array &in, bool is_diag=false)
	Compute eigenvalues. More...
void	eigen (array &values, array &vectors, const array &in)
	Eigenvalues and eigenvectors. More...
array	svd (const array &in, bool is_diag=false)
	Singular values. More...
void	svd (array &s, array &u, array &v, const array &in)
	Singular values with unitary bases: in = u * s * v. More...

Detailed Description

Function Documentation

array af::lu

(

const array &

in

)

LU factorization (packed).

Double-precision or complex values require ArrayFire Pro.

Parameters

[in]

in

Returns

lower and upper matrix combined (getrf)

Examples:

examples/getting_started/lin_algebra.cpp.

void af::lu	(	array &	lower,
		array &	upper,
		const array &	in
	)

LU factorization.

Double-precision or complex values require ArrayFire Pro.

array l, u;
array in = randu(5, 5);
lu(l, u, in); // Unpacked output
print(l);
print(u);

Parameters

[out]	lower	triangular matrix
[out]	upper	triangular matrix
[in]	in

void af::lu	(	array &	lower,
		array &	upper,
		array &	pivot,
		const array &	in
	)

LU factorization (with pivoting).

This function decomposes a matrix A as a product of a lower-triangular matrix l and an upper-triangular matrix u. Note that in general, the product of l and u gives us a matrix which is a row-permutation of the original matrix A, obtained by pre-multiplying A with a permutation matrix P.

The following relationship is obtained as a result:

lu = PA

The permutation matrix P is simply a row-permutation of the identity matrix. The location of non-zero entries (1's) in the P matrix can be obtained via a permutation vector p returned by this function, as explained in the code example below.

Double-precision or complex values require ArrayFire Pro.

#include<stdio.h>
#include <arrayfire.h>
using namespace af;

int main()
{
    // create a lower triangular matrix L from a host buffer
    float h_buff[] = {1, 1, 2, 0, 4, 5, 0, 0, 8};

    array L(3, 3, h_buff);        // L = 1  0  0
    print(L);                     //     1  4  0
                                  //     2  5  8

    // create an upper triangular matrix U from a host buffer
    float g_buff[] = {7, 0, 0, 3, 4, 0, 6, 7, 8};

    array U(3, 3, g_buff);        // U = 7  3  6
    print(U);                     //     0  4  7
                                  //     0  0  8

    // matrix multiplication
    array C = matmul(L, U);       // C = LU =  7    3    6
    print(C);                     //           7   19   34
                                  //          14   26  111

    // reverse engineer l and u from C
    array l, u, p;
    lu(l, u, p, C);

    print(l);   // lower triangular matrix with diagonals  = 1
    print(u);   // upper triangular matrix
    print(p);   // row location of 1's for each column in the permutation matrix

    // l =   1.0000        0.0000        0.0000
    //       0.5000        1.0000        0.0000
    //       0.5000       -0.6000        1.0000

    // u =  14.0000       26.0000      111.0000
    //       0.0000      -10.0000      -49.5000
    //       0.0000        0.0000      -51.2000

    // p =  2   0   1

    // create permutation matrix based on the indices in p
    //
    //                             row  row  row
    //                              2    0    1
    //
    //                              |    |    |
    //                              |    |    |
    //
    float permutation_matrix[] = {  0,   1,   0,
                                    0,   0,   1,
                                    1,   0,   0  };

    array P(3, 3, permutation_matrix);
    print(P);

    // check if lu = PC

    array lu = matmul(l, u);
    print(lu);

    // lu =      14.0000       26.0000      111.0000
    //            7.0000        3.0000        6.0000
    //            7.0000       19.0000       34.0000

    array PC = matmul(P, C);
    print(PC);

    // PC =      14.0000       26.0000      111.0000
    //            7.0000        3.0000        6.0000
    //            7.0000       19.0000       34.0000

    return 0;
}

Parameters

[out]	lower	triangular matrix
[out]	upper	triangular matrix
[out]	pivot	indices
[in]	in

array af::qr

(

const array &

in

)

QR factorization (packed).

Double-precision or complex values require ArrayFire Pro.

array in = randu(5, 5);
array out = qr(in);
print(out);

Parameters

[in]

in

Returns

Q,R factors (geqrf)

void af::qr	(	array &	q,
		array &	r,
		const array &	in
	)

QR factorization.

Double-precision or complex values require ArrayFire Pro.

array q, r;
array in = randu(5, 5);
qr(q, r, in); // Unpacked output
print(q);
print(r);

Parameters

[out]	q	unitary rotation matrix
[out]	r	upper triangular
[in]	in

void af::qr	(	array &	q,
		array &	r,
		array &	tau,
		const array &	in
	)

QR factorization with tau.

Double-precision or complex values require ArrayFire Pro.

Parameters

[out]	q	orthogonal matrix
[out]	r	upper triangular matrix
[out]	tau	Additional information about q
[in]	in

array af::cholesky	(	const array &	X,
		bool	is_upper = true
	)

Cholesky decomposition ("Y^T * Y == X").

Double-precision or complex values require ArrayFire Pro.

float h_x[]={ 123.1029, -15.4097, -15.4097, 62.9234 };
array x = array(2,2,h_x);
array y = cholesky(x,true);
print(x);
print(y);
// x =
//  123.1029      -15.4097
//  -15.4097       62.9234
// y =
//  11.0952       -1.3889
//  0.0000        7.8099

Parameters

[in]	X
[in]	is_upper	true for upper triangular solution, false for lower triangular solution

Returns

upper or lower triangular factor

Note

Exception is throw in case input is not positive definite.

array af::cholesky	(	unsigned &	info,
		const array &	X,
		bool	is_upper = true
	)

Cholesky decomposition ("Y^T * Y == X").

Double-precision or complex values require ArrayFire Pro.

Parameters

[out]	info	details on result of decomposition. If successful info is 0 indicating input is positive definite. Else a positive value representing the row at which decomposition failed.
[in]	X
[in]	is_upper	true for upper triangular solution, false for lower triangular solution

Returns

upper or lower triangular factor

Note

No exception is thrown for this entry point.

array af::hessenberg

(

const array &

in

)

Hessenberg matrix form.

Double-precision or complex values require ArrayFire Pro.

Parameters

[in]

in

Returns

Hessenberg form

void af::hessenberg	(	array &	h,
		array &	q,
		const array &	in
	)

Hessenberg matrix h with unitary permutation matrix q.

such that in = q * h * q^T. Double-precision or complex values require ArrayFire Pro.

Parameters

[out]	h	Hessenberg matrix
[out]	q	unitary matrix
[in]	in

array af::eigen	(	const array &	in,
		bool	is_diag = false
	)

Compute eigenvalues.

Requires ArrayFire Pro.

Parameters

[in]	in
[in]	is_diag	true then produce diagonal matrix of eigenvalues, false then vector of eigenvalues

Returns

eigenvalues

Examples:

examples/getting_started/lin_algebra.cpp.

void af::eigen	(	array &	values,
		array &	vectors,
		const array &	in
	)

Eigenvalues and eigenvectors.

Requires ArrayFire Pro.

This function returns eigenvalues and eigenvectors for the input matrix

#include<stdio.h>
#include <arrayfire.h>
using namespace af;

int main(void)
{
    try {

        // create a coefficient matrix A
        float h_buffer[] = {1, -4, 6, -4, 2, -2, 6, -2, 1};
        array A(3, 3, h_buffer);
        print(A);

        // A =  1  -4   6
        //     -4   2  -2
        //      6  -2   1

        // calculate eigenvalues and eigenvectors
        array values, vectors;
        eigen(values, vectors, A);

        // In this case, the three (scalar) Eigenvalues
        // are located on the diagonal

        print(values);  // values  = -5.3401  0.0000  0.0000
                        //            0.0000 -0.1188  0.0000
                        //            0.0000  0.0000  9.4588

        // The Eigenvectors are returned as column vectors
        // Each vector corresponds to the Eigenvalue in that
        // particular column

        print(vectors); // vectors =  0.7420  0.1805 -0.6457
                        //            0.2330  0.8336  0.5008
                        //           -0.6286  0.5221 -0.5764

    } catch(af::exception &ae) {
        std::cout << ae.what() << std::endl;
    }

    return 0;
}

In the above code example, we obtain three distinct eigenvalues and three distinct eigenvectors, as illustrated in the figure below. As a check, if we pre-multiply each eigenvector with the input matrix, the effect is identical to scaling the eigenvector by a non-zero scalar. This scalar is the associated eigenvalue.

Parameters

[out]	values
[out]	vectors
[in]	in

array af::svd	(	const array &	in,
		bool	is_diag = false
	)

Singular values.

Double-precision or complex values require ArrayFire Pro.

array in = randu(5, 5);
array s = svd(in);
print(s);
// s =
//   2.6075
//   0.9519
//   0.4444
//   0.1295
//   0.1070

Parameters

[in]	in
[in]	is_diag	true then produce diagonal matrix of eigenvalues, false then vector of eigenvalues

Returns

singular values

void af::svd	(	array &	s,
		array &	u,
		array &	v,
		const array &	in
	)

Singular values with unitary bases: in = u * s * v.

Double-precision or complex values require ArrayFire Pro.

array s, u, v;
array in = randu(5, 5);
svd(s, u, v, in); // Unpacked output
print(s); // 5x5 diag matrix
print(u);
print(v);

Parameters

[out]	s	singular values
[out]	u	left unitary matrix
[out]	v	right unitary matrix
[in]	in	The input matrix.