Functions | |
void | magmablas_clag2z_sparse (magma_int_t M, magma_int_t N, const magmaFloatComplex *SA, magma_int_t ldsa, magmaDoubleComplex *A, magma_int_t lda, magma_int_t *info) |
CLAG2Z converts a COMPLEX matrix SA to a COMPLEX_16 matrix A. | |
magma_int_t | magma_clobpcg_res (magma_int_t num_rows, magma_int_t num_vecs, float *evalues, magmaFloatComplex *X, magmaFloatComplex *R, float *res) |
This routine computes for Block-LOBPCG, the set of residuals. | |
magma_int_t | magma_clobpcg_shift (magma_int_t num_rows, magma_int_t num_vecs, magma_int_t shift, magmaFloatComplex *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. | |
magma_int_t | magma_ccopyscale (int n, int k, magmaFloatComplex *r, magmaFloatComplex *v, magmaFloatComplex *skp) |
Computes the correction term of the pipelined GMRES according to P. | |
magma_int_t | magma_c_spmv_shift (magmaFloatComplex alpha, magma_c_sparse_matrix A, magmaFloatComplex lambda, magma_c_vector x, magmaFloatComplex beta, magma_int_t offset, magma_int_t blocksize, magma_index_t *add_rows, magma_c_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. | |
magma_int_t | magma_vector_clag2z (magma_c_vector x, magma_z_vector *y) |
convertes magma_c_vector from C to Z | |
magma_int_t | magma_sparse_matrix_clag2z (magma_c_sparse_matrix A, magma_z_sparse_matrix *B) |
convertes magma_c_sparse_matrix from C to Z | |
magma_int_t | magma_vector_slag2d (magma_s_vector x, magma_d_vector *y) |
convertes magma_s_vector from C to Z | |
magma_int_t | magma_sparse_matrix_slag2d (magma_s_sparse_matrix A, magma_d_sparse_matrix *B) |
convertes magma_s_sparse_matrix from C to Z | |
void | magmablas_slag2d_sparse (magma_int_t M, magma_int_t N, const float *SA, magma_int_t ldsa, double *A, magma_int_t lda, magma_int_t *info) |
SLAG2D converts a SINGLE PRECISION matrix SA to a DOUBLE PRECISION matrix A. |
magma_int_t magma_c_spmv_shift | ( | magmaFloatComplex | alpha, | |
magma_c_sparse_matrix | A, | |||
magmaFloatComplex | lambda, | |||
magma_c_vector | x, | |||
magmaFloatComplex | beta, | |||
magma_int_t | offset, | |||
magma_int_t | blocksize, | |||
magma_index_t * | add_rows, | |||
magma_c_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 | magmaFloatComplex scalar alpha | |
A | magma_c_sparse_matrix sparse matrix A | |
lambda | magmaFloatComplex scalar lambda | |
x | magma_c_vector input vector x | |
beta | magmaFloatComplex 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_c_vector output vector y |
magma_int_t magma_ccopyscale | ( | int | n, | |
int | k, | |||
magmaFloatComplex * | r, | |||
magmaFloatComplex * | v, | |||
magmaFloatComplex * | 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 | magmaFloatComplex* vector of length n | |
v | magmaFloatComplex* vector of length n | |
skp | magmaFloatComplex* array of parameters |
magma_int_t magma_clobpcg_res | ( | magma_int_t | num_rows, | |
magma_int_t | num_vecs, | |||
float * | evalues, | |||
magmaFloatComplex * | X, | |||
magmaFloatComplex * | 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_caxpy(m, MAGMA_C_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 | magmaFloatComplex* block of eigenvector approximations | |
R | magmaFloatComplex* block of residuals | |
res | float* array of residuals |
magma_int_t magma_clobpcg_shift | ( | magma_int_t | num_rows, | |
magma_int_t | num_vecs, | |||
magma_int_t | shift, | |||
magmaFloatComplex * | 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 | magmaFloatComplex* input/output vector x |
magma_int_t magma_sparse_matrix_clag2z | ( | magma_c_sparse_matrix | A, | |
magma_z_sparse_matrix * | B | |||
) |
convertes magma_c_sparse_matrix from C to Z
A | magma_c_sparse_matrix input matrix descriptor | |
B | magma_z_sparse_matrix* output matrix descriptor |
magma_int_t magma_sparse_matrix_slag2d | ( | magma_s_sparse_matrix | A, | |
magma_d_sparse_matrix * | B | |||
) |
convertes magma_s_sparse_matrix from C to Z
A | magma_s_sparse_matrix input matrix descriptor | |
B | magma_d_sparse_matrix* output matrix descriptor |
magma_int_t magma_vector_clag2z | ( | magma_c_vector | x, | |
magma_z_vector * | y | |||
) |
convertes magma_c_vector from C to Z
x | magma_c_vector input vector descriptor | |
y | magma_z_vector* output vector descriptor |
magma_int_t magma_vector_slag2d | ( | magma_s_vector | x, | |
magma_d_vector * | y | |||
) |
convertes magma_s_vector from C to Z
x | magma_s_vector input vector descriptor | |
y | magma_d_vector* output vector descriptor |
void magmablas_clag2z_sparse | ( | magma_int_t | M, | |
magma_int_t | N, | |||
const magmaFloatComplex * | SA, | |||
magma_int_t | ldsa, | |||
magmaDoubleComplex * | A, | |||
magma_int_t | lda, | |||
magma_int_t * | info | |||
) |
CLAG2Z converts a COMPLEX matrix SA to a COMPLEX_16 matrix A.
RMAX is the overflow for the COMPLEX arithmetic. CLAG2Z checks that all the entries of A are between -RMAX and RMAX. If not the convertion is aborted and a flag is raised.
[in] | M | INTEGER The number of lines of the matrix A. M >= 0. |
[in] | N | INTEGER The number of columns of the matrix A. N >= 0. |
[in] | SA | COMPLEX array, dimension (LDSA,N) On entry, the M-by-N coefficient matrix SA. |
[in] | ldsa | INTEGER The leading dimension of the array A. LDA >= max(1,M). |
[out] | A | COMPLEX_16 array, dimension (LDA,N) On exit, if INFO=0, the M-by-N coefficient matrix A; if INFO>0, the content of A is unspecified. |
[in] | lda | INTEGER The leading dimension of the array SA. LDSA >= max(1,M). |
[out] | info | INTEGER
|
void magmablas_slag2d_sparse | ( | magma_int_t | M, | |
magma_int_t | N, | |||
const float * | SA, | |||
magma_int_t | ldsa, | |||
double * | A, | |||
magma_int_t | lda, | |||
magma_int_t * | info | |||
) |
SLAG2D converts a SINGLE PRECISION matrix SA to a DOUBLE PRECISION matrix A.
RMAX is the overflow for the SINGLE PRECISION arithmetic. SLAG2D checks that all the entries of A are between -RMAX and RMAX. If not the convertion is aborted and a flag is raised.
[in] | M | INTEGER The number of lines of the matrix A. M >= 0. |
[in] | N | INTEGER The number of columns of the matrix A. N >= 0. |
[in] | SA | SINGLE PRECISION array, dimension (LDSA,N) On entry, the M-by-N coefficient matrix SA. |
[in] | ldsa | INTEGER The leading dimension of the array A. LDA >= max(1,M). |
[out] | A | DOUBLE PRECISION array, dimension (LDA,N) On exit, if INFO=0, the M-by-N coefficient matrix A; if INFO>0, the content of A is unspecified. |
[in] | lda | INTEGER The leading dimension of the array SA. LDSA >= max(1,M). |
[out] | info | INTEGER
|