![]() |
MAGMA
1.5.0
Matrix Algebra for GPU and Multicore Architectures
|
Functions | |
magma_int_t | magma_s_precond (magma_s_sparse_matrix A, magma_s_vector b, magma_s_vector *x, magma_s_preconditioner precond) |
For a given input matrix A and vectors x, y and the preconditioner parameters, the respective preconditioner is chosen. More... | |
magma_int_t | magma_s_precondsetup (magma_s_sparse_matrix A, magma_s_vector b, magma_s_preconditioner *precond) |
For a given input matrix A and vectors x, y and the preconditioner parameters, the respective preconditioner is preprocessed. More... | |
magma_int_t | magma_s_applyprecond (magma_s_sparse_matrix A, magma_s_vector b, magma_s_vector *x, magma_s_preconditioner *precond) |
For a given input matrix A and vectors x, y and the preconditioner parameters, the respective preconditioner is applied. More... | |
magma_int_t | magma_s_applyprecond_left (magma_s_sparse_matrix A, magma_s_vector b, magma_s_vector *x, magma_s_preconditioner *precond) |
For a given input matrix A and vectors x, y and the preconditioner parameters, the respective left preconditioner is applied. More... | |
magma_int_t | magma_s_applyprecond_right (magma_s_sparse_matrix A, magma_s_vector b, magma_s_vector *x, magma_s_preconditioner *precond) |
For a given input matrix A and vectors x, y and the preconditioner parameters, the respective right-preconditioner is applied. More... | |
magma_int_t | magma_sresidual (magma_s_sparse_matrix A, magma_s_vector b, magma_s_vector x, float *res) |
Computes the residual ||b-Ax|| for a solution approximation x. More... | |
magma_int_t | magma_ccsrsplit (magma_int_t bsize, magma_c_sparse_matrix A, magma_c_sparse_matrix *D, magma_c_sparse_matrix *R) |
– MAGMA (version 1.5.0-beta3) – Univ. More... | |
magma_int_t | magma_dcsrsplit (magma_int_t bsize, magma_d_sparse_matrix A, magma_d_sparse_matrix *D, magma_d_sparse_matrix *R) |
– MAGMA (version 1.5.0-beta3) – Univ. More... | |
magma_int_t | magma_snonlinres (magma_s_sparse_matrix A, magma_s_sparse_matrix L, magma_s_sparse_matrix U, magma_s_sparse_matrix *LU, real_Double_t *res) |
Computes the nonlinear residual A- LU and returns the difference as well es the Frobenius norm of the difference. More... | |
magma_int_t | magma_silures (magma_s_sparse_matrix A, magma_s_sparse_matrix L, magma_s_sparse_matrix U, magma_s_sparse_matrix *LU, real_Double_t *res) |
Computes the ILU residual A- LU and returns the difference as well es the Frobenius norm of the difference. More... | |
magma_int_t | magma_s_vfree (magma_s_vector *x) |
Free the memory of a magma_s_vector. More... | |
magma_int_t | magma_s_mfree (magma_s_sparse_matrix *A) |
Free the memory of a magma_s_sparse_matrix. More... | |
magma_int_t | magma_s_vinit (magma_s_vector *x, magma_location_t mem_loc, magma_int_t num_rows, float values) |
Initialize a magma_s_vector. More... | |
magma_int_t | magma_srowentries (magma_s_sparse_matrix *A) |
Checks the maximal number of nonzeros in a row of matrix A. More... | |
magma_int_t | magma_sdiameter (magma_s_sparse_matrix *A) |
Computes the diameter of a sparse matrix and stores the value in diameter. More... | |
magma_int_t | magma_s_csr_compressor (float **val, magma_index_t **row, magma_index_t **col, float **valn, magma_index_t **rown, magma_index_t **coln, magma_int_t *n) |
Helper function to compress CSR containing zero-entries. More... | |
magma_int_t | magma_s_mconvert (magma_s_sparse_matrix A, magma_s_sparse_matrix *B, magma_storage_t old_format, magma_storage_t new_format) |
Converter between different sparse storage formats. More... | |
magma_int_t | magma_s_LUmergein (magma_s_sparse_matrix L, magma_s_sparse_matrix U, magma_s_sparse_matrix *B) |
Merges an ILU factorization into one matrix. More... | |
magma_int_t | magma_s_mtransfer (magma_s_sparse_matrix A, magma_s_sparse_matrix *B, magma_location_t src, magma_location_t dst) |
Copies a matrix from memory location src to memory location dst. More... | |
magma_int_t | magma_s_vtransfer (magma_s_vector x, magma_s_vector *y, magma_location_t src, magma_location_t dst) |
Copies a vector from memory location src to memory location dst. More... | |
magma_int_t | magma_s_vvisu (magma_s_vector x, magma_int_t offset, magma_int_t visulen) |
Visualizes part of a vector of type magma_s_vector. More... | |
magma_int_t | magma_s_vread (magma_s_vector *x, magma_int_t length, char *filename) |
Reads in a float vector of length "length". More... | |
magma_int_t | magma_scsrsplit (magma_int_t bsize, magma_s_sparse_matrix A, magma_s_sparse_matrix *D, magma_s_sparse_matrix *R) |
– MAGMA (version 1.5.0-beta3) – Univ. More... | |
magma_int_t | magma_smgenerator (magma_int_t n, magma_int_t offdiags, magma_index_t *diag_offset, float *diag_vals, magma_s_sparse_matrix *A) |
Generate a symmetric n x n CSR matrix for a stencil. More... | |
magma_int_t | magma_silustruct (magma_s_sparse_matrix *A, magma_int_t levels) |
This routine computes the fill-in structure of an ILU(levels) factorization based on the successive multiplication of upper and lower triangular factors using the CUSPARSE library. More... | |
magma_int_t | magma_smhom (magma_s_sparse_matrix A, magma_int_t b, magma_index_t *p) |
– MAGMA (version 1.5.0-beta3) – Univ. More... | |
magma_int_t | magma_smscale (magma_s_sparse_matrix *A, magma_scale_t scaling) |
– MAGMA (version 1.5.0-beta3) – Univ. More... | |
magma_int_t | magma_smdiagadd (magma_s_sparse_matrix *A, float add) |
– MAGMA (version 1.5.0-beta3) – Univ. More... | |
magma_int_t | magma_smsort (magma_s_sparse_matrix *A) |
Sorts columns and rows for a matrix in COO or CSRCOO format. More... | |
magma_int_t | magma_s_initP2P (magma_int_t *bw_bmark, magma_int_t *num_gpus) |
Initializes P2P communication between GPUs. More... | |
magma_int_t | magma_ssolverinfo (magma_s_solver_par *solver_par, magma_s_preconditioner *precond_par) |
– MAGMA (version 1.5.0-beta3) – Univ. More... | |
magma_int_t | magma_ssolverinfo_free (magma_s_solver_par *solver_par, magma_s_preconditioner *precond_par) |
– MAGMA (version 1.5.0-beta3) – Univ. More... | |
magma_int_t | magma_ssolverinfo_init (magma_s_solver_par *solver_par, magma_s_preconditioner *precond_par) |
– MAGMA (version 1.5.0-beta3) – Univ. More... | |
magma_int_t | s_transpose_csr (magma_int_t n_rows, magma_int_t n_cols, magma_int_t nnz, float *val, magma_index_t *row, magma_index_t *col, magma_int_t *new_n_rows, magma_int_t *new_n_cols, magma_int_t *new_nnz, float **new_val, magma_index_t **new_row, magma_index_t **new_col) |
Transposes a matrix stored in CSR format. More... | |
magma_int_t | magma_s_csrtranspose (magma_s_sparse_matrix A, magma_s_sparse_matrix *B) |
Helper function to transpose CSR matrix. More... | |
magma_int_t | magma_s_cucsrtranspose (magma_s_sparse_matrix A, magma_s_sparse_matrix *B) |
Helper function to transpose CSR matrix. More... | |
magma_int_t | magma_zcsrsplit (magma_int_t bsize, magma_z_sparse_matrix A, magma_z_sparse_matrix *D, magma_z_sparse_matrix *R) |
– MAGMA (version 1.5.0-beta3) – Univ. More... | |
magma_int_t | read_s_csr_from_binary (magma_int_t *n_row, magma_int_t *n_col, magma_int_t *nnz, float **val, magma_index_t **row, magma_index_t **col, const char *filename) |
Reads in a matrix stored in coo format from a binary and converts it into CSR format. More... | |
magma_int_t | read_s_csr_from_mtx (magma_storage_t *type, magma_location_t *location, magma_int_t *n_row, magma_int_t *n_col, magma_int_t *nnz, float **val, magma_index_t **row, magma_index_t **col, const char *filename) |
Reads in a matrix stored in coo format from a Matrix Market (.mtx) file and converts it into CSR format. More... | |
magma_int_t | write_s_csr_mtx (magma_int_t n_row, magma_int_t n_col, magma_int_t nnz, float **val, magma_index_t **row, magma_index_t **col, magma_order_t MajorType, const char *filename) |
Writes a CSR matrix to a file using Matrix Market format. More... | |
magma_int_t | print_s_csr_mtx (magma_int_t n_row, magma_int_t n_col, magma_int_t nnz, float **val, magma_index_t **row, magma_index_t **col, magma_order_t MajorType) |
Prints a CSR matrix in Matrix Market format. More... | |
magma_int_t | print_s_csr (magma_int_t n_row, magma_int_t n_col, magma_int_t nnz, float **val, magma_index_t **row, magma_index_t **col) |
Prints a CSR matrix in CSR format. More... | |
magma_int_t | magma_s_mvisu (magma_s_sparse_matrix A) |
Prints a sparse matrix in CSR format. More... | |
magma_int_t | magma_s_csr_mtx (magma_s_sparse_matrix *A, const char *filename) |
Reads in a matrix stored in coo format from a Matrix Market (.mtx) file and converts it into CSR format. More... | |
magma_int_t | magma_s_csr_mtxsymm (magma_s_sparse_matrix *A, const char *filename) |
– MAGMA (version 1.5.0-beta3) – Univ. More... | |
magma_int_t | magma_s_spmv_shift (float alpha, magma_s_sparse_matrix A, float lambda, magma_s_vector x, float beta, magma_int_t offset, magma_int_t blocksize, magma_index_t *add_rows, magma_s_vector y) |
For a given input matrix A and vectors x, y and scalars alpha, beta the wrapper determines the suitable SpMV computing y = alpha * ( A - lambda I ) * x + beta * y. More... | |
magma_int_t | magma_slobpcg_res (magma_int_t num_rows, magma_int_t num_vecs, float *evalues, float *X, float *R, float *res) |
This routine computes for Block-LOBPCG, the set of residuals. More... | |
magma_int_t | magma_slobpcg_shift (magma_int_t num_rows, magma_int_t num_vecs, magma_int_t shift, float *x) |
For a Block-LOBPCG, the set of residuals (entries consecutive in memory) shrinks and the vectors are shifted in case shift residuals drop below threshold. More... | |
magma_int_t | magma_scopyscale (int n, int k, float *r, float *v, float *skp) |
Computes the correction term of the pipelined GMRES according to P. More... | |
magma_int_t magma_ccsrsplit | ( | magma_int_t | bsize, |
magma_c_sparse_matrix | A, | ||
magma_c_sparse_matrix * | D, | ||
magma_c_sparse_matrix * | R | ||
) |
– MAGMA (version 1.5.0-beta3) – Univ.
of Tennessee, Knoxville Univ. of California, Berkeley Univ. of Colorado, Denver November 2011
Splits a CSR matrix into two matrices, one containing the diagonal blocks with the diagonal element stored first, one containing the rest of the original matrix.
bsize | magma_int_t size of the diagonal blocks |
A | magma_c_sparse_matrix CSR input matrix |
D | magma_c_sparse_matrix* CSR matrix containing diagonal blocks |
R | magma_c_sparse_matrix* CSR matrix containing rest |
magma_int_t magma_dcsrsplit | ( | magma_int_t | bsize, |
magma_d_sparse_matrix | A, | ||
magma_d_sparse_matrix * | D, | ||
magma_d_sparse_matrix * | R | ||
) |
– MAGMA (version 1.5.0-beta3) – Univ.
of Tennessee, Knoxville Univ. of California, Berkeley Univ. of Colorado, Denver November 2011
Splits a CSR matrix into two matrices, one containing the diagonal blocks with the diagonal element stored first, one containing the rest of the original matrix.
bsize | magma_int_t size of the diagonal blocks |
A | magma_d_sparse_matrix CSR input matrix |
D | magma_d_sparse_matrix* CSR matrix containing diagonal blocks |
R | magma_d_sparse_matrix* CSR matrix containing rest |
magma_int_t magma_s_applyprecond | ( | magma_s_sparse_matrix | A, |
magma_s_vector | b, | ||
magma_s_vector * | x, | ||
magma_s_preconditioner * | precond | ||
) |
For a given input matrix A and vectors x, y and the preconditioner parameters, the respective preconditioner is applied.
E.g. for Jacobi: the scaling-vetor, for ILU the triangular solves.
A | magma_s_sparse_matrix sparse matrix A |
x | magma_s_vector input vector x |
y | magma_s_vector input vector y |
precond | magma_s_preconditioner preconditioner |
magma_int_t magma_s_applyprecond_left | ( | magma_s_sparse_matrix | A, |
magma_s_vector | b, | ||
magma_s_vector * | x, | ||
magma_s_preconditioner * | precond | ||
) |
For a given input matrix A and vectors x, y and the preconditioner parameters, the respective left preconditioner is applied.
E.g. for Jacobi: the scaling-vetor, for ILU the left triangular solve.
A | magma_s_sparse_matrix sparse matrix A |
x | magma_s_vector input vector x |
y | magma_s_vector input vector y |
precond | magma_s_preconditioner preconditioner |
magma_int_t magma_s_applyprecond_right | ( | magma_s_sparse_matrix | A, |
magma_s_vector | b, | ||
magma_s_vector * | x, | ||
magma_s_preconditioner * | precond | ||
) |
For a given input matrix A and vectors x, y and the preconditioner parameters, the respective right-preconditioner is applied.
E.g. for Jacobi: the scaling-vetor, for ILU the right triangular solve.
A | magma_s_sparse_matrix sparse matrix A |
x | magma_s_vector input vector x |
y | magma_s_vector input vector y |
precond | magma_s_preconditioner preconditioner |
magma_int_t magma_s_csr_compressor | ( | float ** | val, |
magma_index_t ** | row, | ||
magma_index_t ** | col, | ||
float ** | valn, | ||
magma_index_t ** | rown, | ||
magma_index_t ** | coln, | ||
magma_int_t * | n | ||
) |
Helper function to compress CSR containing zero-entries.
val | float** input val pointer to compress |
row | magma_int_t** input row pointer to modify |
col | magma_int_t** input col pointer to compress |
valn | float** output val pointer |
rown | magma_int_t** output row pointer |
coln | magma_int_t** output col pointer |
n | magma_int_t* number of rows in matrix |
magma_int_t magma_s_csr_mtx | ( | magma_s_sparse_matrix * | A, |
const char * | filename | ||
) |
Reads in a matrix stored in coo format from a Matrix Market (.mtx) file and converts it into CSR format.
It duplicates the off-diagonal entries in the symmetric case.
A | magma_s_sparse_matrix* matrix in magma sparse matrix format |
filename | const char* filname of the mtx matrix |
magma_int_t magma_s_csr_mtxsymm | ( | magma_s_sparse_matrix * | A, |
const char * | filename | ||
) |
– MAGMA (version 1.5.0-beta3) – Univ.
of Tennessee, Knoxville Univ. of California, Berkeley Univ. of Colorado, Denver November 2011
Reads in a SYMMETRIC matrix stored in coo format from a Matrix Market (.mtx) file and converts it into CSR format. It does not duplicate the off-diagonal entries!
A | magma_s_sparse_matrix* matrix in magma sparse matrix format |
filename | const char* filname of the mtx matrix |
magma_int_t magma_s_csrtranspose | ( | magma_s_sparse_matrix | A, |
magma_s_sparse_matrix * | B | ||
) |
Helper function to transpose CSR matrix.
A | magma_s_sparse_matrix input matrix (CSR) |
B | magma_s_sparse_matrix* output matrix (CSR) |
magma_int_t magma_s_cucsrtranspose | ( | magma_s_sparse_matrix | A, |
magma_s_sparse_matrix * | B | ||
) |
Helper function to transpose CSR matrix.
Using the CUSPARSE CSR2CSC function.
A | magma_s_sparse_matrix input matrix (CSR) |
B | magma_s_sparse_matrix* output matrix (CSR) |
magma_int_t magma_s_initP2P | ( | magma_int_t * | bw_bmark, |
magma_int_t * | num_gpus | ||
) |
Initializes P2P communication between GPUs.
bw_bmark | magma_int_t* input: run the benchmark (1/0) |
num_gpus | magma_int_t* output: number of GPUs |
magma_int_t magma_s_LUmergein | ( | magma_s_sparse_matrix | L, |
magma_s_sparse_matrix | U, | ||
magma_s_sparse_matrix * | B | ||
) |
Merges an ILU factorization into one matrix.
works only for the symmetric case!!!
L | magma_s_sparse_matrix sparse matrix L |
U | magma_s_sparse_matrix sparse matrix U |
B | magma_s_sparse_matrix* output sparse matrix B |
magma_int_t magma_s_mconvert | ( | magma_s_sparse_matrix | A, |
magma_s_sparse_matrix * | B, | ||
magma_storage_t | old_format, | ||
magma_storage_t | new_format | ||
) |
Converter between different sparse storage formats.
A | magma_s_sparse_matrix sparse matrix A |
B | magma_s_sparse_matrix* copy of A in new format |
old_format | magma_storage_t original storage format |
new_format | magma_storage_t new storage format |
magma_int_t magma_s_mfree | ( | magma_s_sparse_matrix * | A | ) |
Free the memory of a magma_s_sparse_matrix.
A | magma_s_sparse_matrix* matrix to free |
magma_int_t magma_s_mtransfer | ( | magma_s_sparse_matrix | A, |
magma_s_sparse_matrix * | B, | ||
magma_location_t | src, | ||
magma_location_t | dst | ||
) |
Copies a matrix from memory location src to memory location dst.
A | magma_s_sparse_matrix sparse matrix A |
B | magma_s_sparse_matrix* copy of A |
src | magma_location_t original location A |
dst | magma_location_t location of the copy of A |
magma_int_t magma_s_mvisu | ( | magma_s_sparse_matrix | A | ) |
Prints a sparse matrix in CSR format.
A | magma_s_sparse_matrix sparse matrix in Magma_CSR format |
magma_int_t magma_s_precond | ( | magma_s_sparse_matrix | A, |
magma_s_vector | b, | ||
magma_s_vector * | x, | ||
magma_s_preconditioner | precond | ||
) |
For a given input matrix A and vectors x, y and the preconditioner parameters, the respective preconditioner is chosen.
It approximates x for A x = y.
A | magma_s_sparse_matrix sparse matrix A |
x | magma_s_vector input vector x |
y | magma_s_vector input vector y |
precond | magma_s_preconditioner preconditioner |
magma_int_t magma_s_precondsetup | ( | magma_s_sparse_matrix | A, |
magma_s_vector | b, | ||
magma_s_preconditioner * | precond | ||
) |
For a given input matrix A and vectors x, y and the preconditioner parameters, the respective preconditioner is preprocessed.
E.g. for Jacobi: the scaling-vetor, for ILU the factorization.
A | magma_s_sparse_matrix sparse matrix A |
x | magma_s_vector input vector x |
y | magma_s_vector input vector y |
precond | magma_s_preconditioner preconditioner |
magma_int_t magma_s_spmv_shift | ( | float | alpha, |
magma_s_sparse_matrix | A, | ||
float | lambda, | ||
magma_s_vector | x, | ||
float | beta, | ||
magma_int_t | offset, | ||
magma_int_t | blocksize, | ||
magma_index_t * | add_rows, | ||
magma_s_vector | y | ||
) |
For a given input matrix A and vectors x, y and scalars alpha, beta the wrapper determines the suitable SpMV computing y = alpha * ( A - lambda I ) * x + beta * y.
alpha | float scalar alpha |
A | magma_s_sparse_matrix sparse matrix A |
lambda | float scalar lambda |
x | magma_s_vector input vector x |
beta | float scalar beta |
offset | magma_int_t in case not the main diagonal is scaled |
blocksize | magma_int_t in case of processing multiple vectors |
add_rows | magma_int_t* in case the matrixpowerskernel is used |
y | magma_s_vector output vector y |
magma_int_t magma_s_vfree | ( | magma_s_vector * | x | ) |
Free the memory of a magma_s_vector.
x | magma_s_vector* vector to free |
magma_int_t magma_s_vinit | ( | magma_s_vector * | x, |
magma_location_t | mem_loc, | ||
magma_int_t | num_rows, | ||
float | values | ||
) |
Initialize a magma_s_vector.
x | magma_s_vector vector to initialize |
mem_loc | magma_location_t memory for vector |
num_rows | magma_int_t desired length of vector |
values | float entries in vector |
magma_int_t magma_s_vread | ( | magma_s_vector * | x, |
magma_int_t | length, | ||
char * | filename | ||
) |
Reads in a float vector of length "length".
x | magma_s_vector vector to read in |
length | magma_int_t length of vector |
filename | char* file where vector is stored |
magma_int_t magma_s_vtransfer | ( | magma_s_vector | x, |
magma_s_vector * | y, | ||
magma_location_t | src, | ||
magma_location_t | dst | ||
) |
Copies a vector from memory location src to memory location dst.
x | magma_s_vector vector x |
y | magma_s_vector* copy of x |
src | magma_location_t original location x |
dst | magma_location_t location of the copy of x |
magma_int_t magma_s_vvisu | ( | magma_s_vector | x, |
magma_int_t | offset, | ||
magma_int_t | visulen | ||
) |
Visualizes part of a vector of type magma_s_vector.
With input vector x , offset, visulen, the entries offset - (offset + visulen) of x are visualized.
x | magma_s_vector vector to visualize |
offset | magma_int_t start inex of visualization |
visulen | magma_int_t number of entries to visualize |
magma_int_t magma_scopyscale | ( | int | n, |
int | k, | ||
float * | r, | ||
float * | v, | ||
float * | skp | ||
) |
Computes the correction term of the pipelined GMRES according to P.
Ghysels and scales and copies the new search direction
Returns the vector v = r/ ( skp[k] - (sum_i=1^k skp[i]^2) ) .
n | int length of v_i |
k | int skp entries v_i^T * r ( without r ) |
r | float* vector of length n |
v | float* vector of length n |
skp | float* array of parameters |
magma_int_t magma_scsrsplit | ( | magma_int_t | bsize, |
magma_s_sparse_matrix | A, | ||
magma_s_sparse_matrix * | D, | ||
magma_s_sparse_matrix * | R | ||
) |
– MAGMA (version 1.5.0-beta3) – Univ.
of Tennessee, Knoxville Univ. of California, Berkeley Univ. of Colorado, Denver November 2011
Splits a CSR matrix into two matrices, one containing the diagonal blocks with the diagonal element stored first, one containing the rest of the original matrix.
bsize | magma_int_t size of the diagonal blocks |
A | magma_s_sparse_matrix CSR input matrix |
D | magma_s_sparse_matrix* CSR matrix containing diagonal blocks |
R | magma_s_sparse_matrix* CSR matrix containing rest |
magma_int_t magma_sdiameter | ( | magma_s_sparse_matrix * | A | ) |
Computes the diameter of a sparse matrix and stores the value in diameter.
A | magma_s_sparse_matrix* sparse matrix |
magma_int_t magma_silures | ( | magma_s_sparse_matrix | A, |
magma_s_sparse_matrix | L, | ||
magma_s_sparse_matrix | U, | ||
magma_s_sparse_matrix * | LU, | ||
real_Double_t * | res | ||
) |
Computes the ILU residual A- LU and returns the difference as well es the Frobenius norm of the difference.
A | magma_s_sparse_matrix input sparse matrix in CSR |
L | magma_s_sparse_matrix input sparse matrix in CSR |
U | magma_s_sparse_matrix input sparse matrix in CSR |
LU | magma_s_sparse_matrix* output sparse matrix in A-LU in CSR |
res | real_Double_t* residual |
magma_int_t magma_silustruct | ( | magma_s_sparse_matrix * | A, |
magma_int_t | levels | ||
) |
This routine computes the fill-in structure of an ILU(levels) factorization based on the successive multiplication of upper and lower triangular factors using the CUSPARSE library.
A | magma_s_sparse_matrix* matrix in magma sparse matrix format |
levels | magma_int_t fill in level |
magma_int_t magma_slobpcg_res | ( | magma_int_t | num_rows, |
magma_int_t | num_vecs, | ||
float * | evalues, | ||
float * | X, | ||
float * | R, | ||
float * | res | ||
) |
This routine computes for Block-LOBPCG, the set of residuals.
R = Ax - x evalues It replaces: for(int i=0; i < n; i++){ magma_saxpy(m, MAGMA_S_MAKE(-evalues[i],0),blockX+i*m,1,blockR+i*m,1); } The memory layout of x is:
/ x1[0] x2[0] x3[0] \ | x1[1] x2[1] x3[1] |
x = | x1[2] x2[2] x3[2] | = x1[0] x1[1] x1[2] x1[3] x1[4] x2[0] x2[1] . | x1[3] x2[3] x3[3] | \ x1[4] x2[4] x3[4] /
num_rows | magma_int_t number of rows |
num_vecs | magma_int_t number of vectors |
evalues | float* array of eigenvalues/approximations |
X | float* block of eigenvector approximations |
R | float* block of residuals |
res | float* array of residuals |
magma_int_t magma_slobpcg_shift | ( | magma_int_t | num_rows, |
magma_int_t | num_vecs, | ||
magma_int_t | shift, | ||
float * | x | ||
) |
For a Block-LOBPCG, the set of residuals (entries consecutive in memory) shrinks and the vectors are shifted in case shift residuals drop below threshold.
The memory layout of x is:
/ x1[0] x2[0] x3[0] \ | x1[1] x2[1] x3[1] | x = | x1[2] x2[2] x3[2] | = x1[0] x2[0] x3[0] x1[1] x2[1] x3[1] x1[2] . | x1[3] x2[3] x3[3] | \ x1[4] x2[4] x3[4] /
num_rows | magma_int_t number of rows |
num_vecs | magma_int_t number of vectors |
shift | magma_int_t shift number |
x | float* input/output vector x |
magma_int_t magma_smdiagadd | ( | magma_s_sparse_matrix * | A, |
float | add | ||
) |
– MAGMA (version 1.5.0-beta3) – Univ.
of Tennessee, Knoxville Univ. of California, Berkeley Univ. of Colorado, Denver November 2011
Adds a multiple of the Identity matrix to a matrix: A = A+add * I
A | magma_s_sparse_matrix* input/output matrix |
add | float scaling for the identity matrix |
magma_int_t magma_smgenerator | ( | magma_int_t | n, |
magma_int_t | offdiags, | ||
magma_index_t * | diag_offset, | ||
float * | diag_vals, | ||
magma_s_sparse_matrix * | A | ||
) |
Generate a symmetric n x n CSR matrix for a stencil.
n | magma_int_t number of rows |
offdiags | magma_int_t number of offdiagonals |
diag_offset | magma_int_t* array containing the offsets (length offsets+1) |
diag_vals | float* array containing the values (length offsets+1) |
A | magma_s_sparse_matrix* matrix to generate |
magma_int_t magma_smhom | ( | magma_s_sparse_matrix | A, |
magma_int_t | b, | ||
magma_index_t * | p | ||
) |
– MAGMA (version 1.5.0-beta3) – Univ.
of Tennessee, Knoxville Univ. of California, Berkeley Univ. of Colorado, Denver November 2011
Takes a matrix and a blocksize b to generate a homomorphism that orders the matrix entries according to the subdomains of size b x b. Returns p on the device
example:
/ a 0 0 b 0 \ | 0 c 0 d 0 |
A= | 0 e f g 0 | b = 2 | h 0 0 0 0 | \ i j 0 0 0 /
will generate the projection:
0 2 1 3 4 7 8 9 10 11
according to
a c b d e h f g i j
A | magma_s_sparse_matrix input/output matrix |
b | magma_int_t blocksize |
p | magma_index_t* homomorphism vector containing the indices |
magma_int_t magma_smscale | ( | magma_s_sparse_matrix * | A, |
magma_scale_t | scaling | ||
) |
– MAGMA (version 1.5.0-beta3) – Univ.
of Tennessee, Knoxville Univ. of California, Berkeley Univ. of Colorado, Denver November 2011
Scales a matrix.
A | magma_s_sparse_matrix* input/output matrix |
scaling | magma_scale_t scaling type (unit rownorm / unit diagonal) |
magma_int_t magma_smsort | ( | magma_s_sparse_matrix * | A | ) |
Sorts columns and rows for a matrix in COO or CSRCOO format.
A | magma_s_sparse_matrix* matrix in magma sparse matrix format |
magma_int_t magma_snonlinres | ( | magma_s_sparse_matrix | A, |
magma_s_sparse_matrix | L, | ||
magma_s_sparse_matrix | U, | ||
magma_s_sparse_matrix * | LU, | ||
real_Double_t * | res | ||
) |
Computes the nonlinear residual A- LU and returns the difference as well es the Frobenius norm of the difference.
A | magma_s_sparse_matrix input sparse matrix in CSR |
L | magma_s_sparse_matrix input sparse matrix in CSR |
U | magma_s_sparse_matrix input sparse matrix in CSR |
LU | magma_s_sparse_matrix* output sparse matrix in A-LU in CSR |
res | real_Double_t* residual |
magma_int_t magma_sresidual | ( | magma_s_sparse_matrix | A, |
magma_s_vector | b, | ||
magma_s_vector | x, | ||
float * | res | ||
) |
Computes the residual ||b-Ax|| for a solution approximation x.
A | magma_s_sparse_matrix input matrix A |
b | magma_s_vector RHS b |
x | magma_s_vector solution approximation |
res | float* return residual |
magma_int_t magma_srowentries | ( | magma_s_sparse_matrix * | A | ) |
Checks the maximal number of nonzeros in a row of matrix A.
Inserts the data into max_nnz_row.
A | magma_s_sparse_matrix* sparse matrix |
magma_int_t magma_ssolverinfo | ( | magma_s_solver_par * | solver_par, |
magma_s_preconditioner * | precond_par | ||
) |
– MAGMA (version 1.5.0-beta3) – Univ.
of Tennessee, Knoxville Univ. of California, Berkeley Univ. of Colorado, Denver November 2011
Prints information about a previously called solver.
solver_par | magma_s_solver_par* structure containing all solver information |
precond_par | magma_s_preconditioner* structure containing all preconditioner information |
magma_int_t magma_ssolverinfo_free | ( | magma_s_solver_par * | solver_par, |
magma_s_preconditioner * | precond_par | ||
) |
– MAGMA (version 1.5.0-beta3) – Univ.
of Tennessee, Knoxville Univ. of California, Berkeley Univ. of Colorado, Denver November 2011
Frees any memory assocoiated with the verbose mode of solver_par. The other values are set to default.
solver_par | magma_s_solver_par* structure containing all solver information |
precond_par | magma_s_preconditioner* structure containing all preconditioner information |
magma_int_t magma_ssolverinfo_init | ( | magma_s_solver_par * | solver_par, |
magma_s_preconditioner * | precond_par | ||
) |
– MAGMA (version 1.5.0-beta3) – Univ.
of Tennessee, Knoxville Univ. of California, Berkeley Univ. of Colorado, Denver November 2011
Initializes all solver and preconditioner parameters.
solver_par | magma_s_solver_par* structure containing all solver information |
precond_par | magma_s_preconditioner* structure containing all preconditioner information |
magma_int_t magma_zcsrsplit | ( | magma_int_t | bsize, |
magma_z_sparse_matrix | A, | ||
magma_z_sparse_matrix * | D, | ||
magma_z_sparse_matrix * | R | ||
) |
– MAGMA (version 1.5.0-beta3) – Univ.
of Tennessee, Knoxville Univ. of California, Berkeley Univ. of Colorado, Denver November 2011
Splits a CSR matrix into two matrices, one containing the diagonal blocks with the diagonal element stored first, one containing the rest of the original matrix.
bsize | magma_int_t size of the diagonal blocks |
A | magma_z_sparse_matrix CSR input matrix |
D | magma_z_sparse_matrix* CSR matrix containing diagonal blocks |
R | magma_z_sparse_matrix* CSR matrix containing rest |
magma_int_t print_s_csr | ( | magma_int_t | n_row, |
magma_int_t | n_col, | ||
magma_int_t | nnz, | ||
float ** | val, | ||
magma_index_t ** | row, | ||
magma_index_t ** | col | ||
) |
Prints a CSR matrix in CSR format.
n_row | magma_int_t* number of rows in matrix |
n_col | magma_int_t* number of columns in matrix |
nnz | magma_int_t* number of nonzeros in matrix |
val | float** value array of CSR |
row | magma_index_t** row pointer of CSR |
col | magma_index_t** column indices of CSR |
magma_int_t print_s_csr_mtx | ( | magma_int_t | n_row, |
magma_int_t | n_col, | ||
magma_int_t | nnz, | ||
float ** | val, | ||
magma_index_t ** | row, | ||
magma_index_t ** | col, | ||
magma_order_t | MajorType | ||
) |
Prints a CSR matrix in Matrix Market format.
n_row | magma_int_t* number of rows in matrix |
n_col | magma_int_t* number of columns in matrix |
nnz | magma_int_t* number of nonzeros in matrix |
val | float** value array of CSR |
row | magma_index_t** row pointer of CSR |
col | magma_index_t** column indices of CSR |
MajorType | magma_index_t Row or Column sort default: 0 = RowMajor, 1 = ColMajor |
magma_int_t read_s_csr_from_binary | ( | magma_int_t * | n_row, |
magma_int_t * | n_col, | ||
magma_int_t * | nnz, | ||
float ** | val, | ||
magma_index_t ** | row, | ||
magma_index_t ** | col, | ||
const char * | filename | ||
) |
Reads in a matrix stored in coo format from a binary and converts it into CSR format.
It duplicates the off-diagonal entries in the symmetric case.
n_row | magma_int_t* number of rows in matrix |
n_col | magma_int_t* number of columns in matrix |
nnz | magma_int_t* number of nonzeros in matrix |
val | float** value array of CSR output |
row | magma_index_t** row pointer of CSR output |
col | magma_index_t** column indices of CSR output |
filename | const char* filname of the mtx matrix |
magma_int_t read_s_csr_from_mtx | ( | magma_storage_t * | type, |
magma_location_t * | location, | ||
magma_int_t * | n_row, | ||
magma_int_t * | n_col, | ||
magma_int_t * | nnz, | ||
float ** | val, | ||
magma_index_t ** | row, | ||
magma_index_t ** | col, | ||
const char * | filename | ||
) |
Reads in a matrix stored in coo format from a Matrix Market (.mtx) file and converts it into CSR format.
It duplicates the off-diagonal entries in the symmetric case.
type | magma_storage_t* storage type of matrix |
location | magma_location_t* location of matrix |
n_row | magma_int_t* number of rows in matrix |
n_col | magma_int_t* number of columns in matrix |
nnz | magma_int_t* number of nonzeros in matrix |
val | float** value array of CSR output |
row | magma_index_t** row pointer of CSR output |
col | magma_index_t** column indices of CSR output |
filename | const char* filname of the mtx matrix |
magma_int_t s_transpose_csr | ( | magma_int_t | n_rows, |
magma_int_t | n_cols, | ||
magma_int_t | nnz, | ||
float * | val, | ||
magma_index_t * | row, | ||
magma_index_t * | col, | ||
magma_int_t * | new_n_rows, | ||
magma_int_t * | new_n_cols, | ||
magma_int_t * | new_nnz, | ||
float ** | new_val, | ||
magma_index_t ** | new_row, | ||
magma_index_t ** | new_col | ||
) |
Transposes a matrix stored in CSR format.
n_rows | magma_int_t number of rows in input matrix |
n_cols | magma_int_t number of columns in input matrix |
nnz | magma_int_t number of nonzeros in input matrix |
val | float* value array of input matrix |
row | magma_index_t* row pointer of input matrix |
col | magma_index_t* column indices of input matrix |
new_n_rows | magma_index_t* number of rows in transposed matrix |
new_n_cols | magma_index_t* number of columns in transposed matrix |
new_nnz | magma_index_t* number of nonzeros in transposed matrix |
new_val | float** value array of transposed matrix |
new_row | magma_index_t** row pointer of transposed matrix |
new_col | magma_index_t** column indices of transposed matrix |
magma_int_t write_s_csr_mtx | ( | magma_int_t | n_row, |
magma_int_t | n_col, | ||
magma_int_t | nnz, | ||
float ** | val, | ||
magma_index_t ** | row, | ||
magma_index_t ** | col, | ||
magma_order_t | MajorType, | ||
const char * | filename | ||
) |
Writes a CSR matrix to a file using Matrix Market format.
n_row | magma_int_t* number of rows in matrix |
n_col | magma_int_t* number of columns in matrix |
nnz | magma_int_t* number of nonzeros in matrix |
val | float** value array of CSR |
row | magma_index_t** row pointer of CSR |
col | magma_index_t** column indices of CSR |
MajorType | magma_index_t Row or Column sort default: 0 = RowMajor, 1 = ColMajor |
filename | const char* output-filname of the mtx matrix |