single-complex precision
[Sparse auxiliary]

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.

Function Documentation

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.

Parameters:
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) ) .

Parameters:
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] /

Parameters:
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] /

Parameters:
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

Parameters:
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

Parameters:
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

Parameters:
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

Parameters:
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.

Parameters:
[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

  • = 0: successful exit.
  • < 0: if INFO = -i, the i-th argument had an illegal value
  • = 1: an entry of the matrix A is greater than the COMPLEX overflow threshold, in this case, the content of SA in exit is unspecified.
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.

Parameters:
[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

  • = 0: successful exit.
  • < 0: if INFO = -i, the i-th argument had an illegal value
  • = 1: an entry of the matrix A is greater than the SINGLE PRECISION overflow threshold, in this case, the content of SA in exit is unspecified.

Generated on 17 Sep 2014 for MAGMA by  doxygen 1.6.1