double-complex precision
[GPU kernels for sparse LA]

Functions

magma_int_t magma_zbajac_csr (magma_int_t localiters, magma_z_sparse_matrix D, magma_z_sparse_matrix R, magma_z_vector b, magma_z_vector *x)
 This routine is a block-asynchronous Jacobi iteration performing s local Jacobi-updates within the block.
magma_int_t magma_zbcsrblockinfo5 (magma_int_t lustep, magma_int_t num_blocks, magma_int_t c_blocks, magma_int_t size_b, magma_index_t *blockinfo, magmaDoubleComplex *val, magmaDoubleComplex **AII)
 For a Block-CSR ILU factorization, this routine copies the filled blocks from the original matrix A and initializes the blocks that will later be filled in the factorization process with zeros.
magma_int_t magma_zbcsrvalcpy (magma_int_t size_b, magma_int_t num_blocks, magma_int_t num_zblocks, magmaDoubleComplex **Aval, magmaDoubleComplex **Bval, magmaDoubleComplex **Bval2)
 For a Block-CSR ILU factorization, this routine copies the filled blocks from the original matrix A and initializes the blocks that will later be filled in the factorization process with zeros.
magma_int_t magma_zbcsrluegemm (magma_int_t size_b, magma_int_t num_brows, magma_int_t kblocks, magmaDoubleComplex **dA, magmaDoubleComplex **dB, magmaDoubleComplex **dC)
 For a Block-CSR ILU factorization, this routine updates all blocks in the trailing matrix.
magma_int_t magma_zbcsrlupivloc (magma_int_t size_b, magma_int_t kblocks, magmaDoubleComplex **dA, magma_int_t *ipiv)
 For a Block-CSR ILU factorization, this routine updates all blocks in the trailing matrix.
magma_int_t magma_zbcsrswp (magma_int_t r_blocks, magma_int_t size_b, magma_int_t *ipiv, magmaDoubleComplex *x)
 For a Block-CSR ILU factorization, this routine swaps rows in the vector *x according to the pivoting in *ipiv.
magma_int_t magma_zbcsrtrsv (magma_uplo_t uplo, magma_int_t r_blocks, magma_int_t c_blocks, magma_int_t size_b, magmaDoubleComplex *A, magma_index_t *blockinfo, magmaDoubleComplex *x)
 For a Block-CSR ILU factorization, this routine performs the triangular solves.
void magma_zcompact (magma_int_t m, magma_int_t n, magmaDoubleComplex *dA, magma_int_t ldda, double *dnorms, double tol, magma_index_t *active, magma_index_t *cBlock)
 ZCOMPACT takes a set of n vectors of size m (in dA) and their norms and compacts them into the cBlock size<=n vectors that have norms > tol.
magma_int_t magma_zjacobisetup_vector_gpu (int num_rows, magmaDoubleComplex *b, magmaDoubleComplex *d, magmaDoubleComplex *c, magmaDoubleComplex *x)
 Prepares the Jacobi Iteration according to x^(k+1) = D^(-1) * b - D^(-1) * (L+U) * x^k x^(k+1) = c - M * x^k.
magma_int_t magma_zlobpcg_maxpy (magma_int_t num_rows, magma_int_t num_vecs, magmaDoubleComplex *X, magmaDoubleComplex *Y)
 This routine computes a axpy for a mxn matrix:.
int magma_zbicgmerge1 (int n, magmaDoubleComplex *skp, magmaDoubleComplex *v, magmaDoubleComplex *r, magmaDoubleComplex *p)
 Mergels multiple operations into one kernel:.
int magma_zbicgmerge2 (int n, magmaDoubleComplex *skp, magmaDoubleComplex *r, magmaDoubleComplex *v, magmaDoubleComplex *s)
 Mergels multiple operations into one kernel:.
int magma_zbicgmerge3 (int n, magmaDoubleComplex *skp, magmaDoubleComplex *p, magmaDoubleComplex *s, magmaDoubleComplex *t, magmaDoubleComplex *x, magmaDoubleComplex *r)
 Mergels multiple operations into one kernel:.
int magma_zbicgmerge4 (int type, magmaDoubleComplex *skp)
 Performs some parameter operations for the BiCGSTAB with scalars on GPU.
magma_int_t magma_zbicgmerge_spmv1 (magma_z_sparse_matrix A, magmaDoubleComplex *d1, magmaDoubleComplex *d2, magmaDoubleComplex *d_p, magmaDoubleComplex *d_r, magmaDoubleComplex *d_v, magmaDoubleComplex *skp)
 Merges the first SpmV using CSR with the dot product and the computation of alpha.
magma_int_t magma_zbicgmerge_spmv2 (magma_z_sparse_matrix A, magmaDoubleComplex *d1, magmaDoubleComplex *d2, magmaDoubleComplex *d_s, magmaDoubleComplex *d_t, magmaDoubleComplex *skp)
 Merges the second SpmV using CSR with the dot product and the computation of omega.
magma_int_t magma_zbicgmerge_xrbeta (int n, magmaDoubleComplex *d1, magmaDoubleComplex *d2, magmaDoubleComplex *rr, magmaDoubleComplex *r, magmaDoubleComplex *p, magmaDoubleComplex *s, magmaDoubleComplex *t, magmaDoubleComplex *x, magmaDoubleComplex *skp)
 Merges the second SpmV using CSR with the dot product and the computation of omega.
magma_int_t magma_zcgmerge_spmv1 (magma_z_sparse_matrix A, magmaDoubleComplex *d1, magmaDoubleComplex *d2, magmaDoubleComplex *d_d, magmaDoubleComplex *d_z, magmaDoubleComplex *skp)
 Merges the first SpmV using different formats with the dot product and the computation of rho.

Function Documentation

magma_int_t magma_zbajac_csr ( magma_int_t  localiters,
magma_z_sparse_matrix  D,
magma_z_sparse_matrix  R,
magma_z_vector  b,
magma_z_vector *  x 
)

This routine is a block-asynchronous Jacobi iteration performing s local Jacobi-updates within the block.

Input format is two CSR matrices, one containing the diagonal blocks, one containing the rest.

Parameters:
localiters magma_int_t number of local Jacobi-like updates
D magma_z_sparse_matrix input matrix with diagonal blocks
R magma_z_sparse_matrix input matrix with non-diagonal parts
b magma_z_vector RHS
x magma_z_vector* iterate/solution
magma_int_t magma_zbcsrblockinfo5 ( magma_int_t  lustep,
magma_int_t  num_blocks,
magma_int_t  c_blocks,
magma_int_t  size_b,
magma_index_t *  blockinfo,
magmaDoubleComplex *  val,
magmaDoubleComplex **  AII 
)

For a Block-CSR ILU factorization, this routine copies the filled blocks from the original matrix A and initializes the blocks that will later be filled in the factorization process with zeros.

Parameters:
lustep magma_int_t lustep
num_blocks magma_int_t number of nonzero blocks
c_blocks magma_int_t number of column-blocks
size_b magma_int_t blocksize
blockinfo magma_int_t* block filled? location?
val magmaDoubleComplex* pointers to the nonzero blocks in A
AII magmaDoubleComplex** pointers to the respective nonzero blocks in B
magma_int_t magma_zbcsrluegemm ( magma_int_t  size_b,
magma_int_t  num_brows,
magma_int_t  kblocks,
magmaDoubleComplex **  dA,
magmaDoubleComplex **  dB,
magmaDoubleComplex **  dC 
)

For a Block-CSR ILU factorization, this routine updates all blocks in the trailing matrix.

Parameters:
size_b magma_int_t blocksize in BCSR
num_brows magma_int_t number of block rows
kblocks magma_int_t number of blocks in row
dA magmaDoubleComplex** input blocks of matrix A
dB magmaDoubleComplex** input blocks of matrix B
dC magmaDoubleComplex** output blocks of matrix C
magma_int_t magma_zbcsrlupivloc ( magma_int_t  size_b,
magma_int_t  kblocks,
magmaDoubleComplex **  dA,
magma_int_t *  ipiv 
)

For a Block-CSR ILU factorization, this routine updates all blocks in the trailing matrix.

Parameters:
size_b magma_int_t blocksize in BCSR
kblocks magma_int_t number of blocks
dA magmaDoubleComplex** matrix in BCSR
ipiv magma_int_t* array containing pivots
magma_int_t magma_zbcsrswp ( magma_int_t  r_blocks,
magma_int_t  size_b,
magma_int_t *  ipiv,
magmaDoubleComplex *  x 
)

For a Block-CSR ILU factorization, this routine swaps rows in the vector *x according to the pivoting in *ipiv.

Parameters:
r_blocks magma_int_t number of blocks
size_b magma_int_t blocksize in BCSR
ipiv magma_int_t* array containing pivots
x magmaDoubleComplex* input/output vector x
magma_int_t magma_zbcsrtrsv ( magma_uplo_t  uplo,
magma_int_t  r_blocks,
magma_int_t  c_blocks,
magma_int_t  size_b,
magmaDoubleComplex *  A,
magma_index_t *  blockinfo,
magmaDoubleComplex *  x 
)

For a Block-CSR ILU factorization, this routine performs the triangular solves.

Parameters:
uplo magma_uplo_t upper/lower fill structure
r_blocks magma_int_t number of blocks in row
c_blocks magma_int_t number of blocks in column
size_b magma_int_t blocksize in BCSR
A magmaDoubleComplex* upper/lower factor
blockinfo magma_int_t* array containing matrix information
x magmaDoubleComplex* input/output vector x
magma_int_t magma_zbcsrvalcpy ( magma_int_t  size_b,
magma_int_t  num_blocks,
magma_int_t  num_zblocks,
magmaDoubleComplex **  Aval,
magmaDoubleComplex **  Bval,
magmaDoubleComplex **  Bval2 
)

For a Block-CSR ILU factorization, this routine copies the filled blocks from the original matrix A and initializes the blocks that will later be filled in the factorization process with zeros.

Parameters:
size_b magma_int_t blocksize in BCSR
num_blocks magma_int_t number of nonzero blocks
num_zblocks magma_int_t number of zero-blocks (will later be filled)
Aval magmaDoubleComplex** pointers to the nonzero blocks in A
Bval magmaDoubleComplex** pointers to the nonzero blocks in B
Bval2 magmaDoubleComplex** pointers to the zero blocks in B
int magma_zbicgmerge1 ( int  n,
magmaDoubleComplex *  skp,
magmaDoubleComplex *  v,
magmaDoubleComplex *  r,
magmaDoubleComplex *  p 
)

Mergels multiple operations into one kernel:.

p = beta*p p = p-omega*beta*v p = p+r

-> p = r + beta * ( p - omega * v )

Parameters:
n int dimension n
skp magmaDoubleComplex* set of scalar parameters
v magmaDoubleComplex* input v
r magmaDoubleComplex* input r
p magmaDoubleComplex* input/output p
int magma_zbicgmerge2 ( int  n,
magmaDoubleComplex *  skp,
magmaDoubleComplex *  r,
magmaDoubleComplex *  v,
magmaDoubleComplex *  s 
)

Mergels multiple operations into one kernel:.

s=r s=s-alpha*v

-> s = r - alpha * v

Parameters:
n int dimension n
skp magmaDoubleComplex* set of scalar parameters
r magmaDoubleComplex* input r
v magmaDoubleComplex* input v
s magmaDoubleComplex* input/output s
int magma_zbicgmerge3 ( int  n,
magmaDoubleComplex *  skp,
magmaDoubleComplex *  p,
magmaDoubleComplex *  s,
magmaDoubleComplex *  t,
magmaDoubleComplex *  x,
magmaDoubleComplex *  r 
)

Mergels multiple operations into one kernel:.

x=x+alpha*p x=x+omega*s r=s r=r-omega*t

-> x = x + alpha * p + omega * s -> r = s - omega * t

Parameters:
n int dimension n
skp magmaDoubleComplex* set of scalar parameters
p magmaDoubleComplex* input p
s magmaDoubleComplex* input s
t magmaDoubleComplex* input t
x magmaDoubleComplex* input/output x
r magmaDoubleComplex* input/output r
int magma_zbicgmerge4 ( int  type,
magmaDoubleComplex *  skp 
)

Performs some parameter operations for the BiCGSTAB with scalars on GPU.

Parameters:
type int kernel type
skp magmaDoubleComplex* vector with parameters
magma_int_t magma_zbicgmerge_spmv1 ( magma_z_sparse_matrix  A,
magmaDoubleComplex *  d1,
magmaDoubleComplex *  d2,
magmaDoubleComplex *  d_p,
magmaDoubleComplex *  d_r,
magmaDoubleComplex *  d_v,
magmaDoubleComplex *  skp 
)

Merges the first SpmV using CSR with the dot product and the computation of alpha.

Parameters:
A magma_z_sparse_matrix system matrix
d1 magmaDoubleComplex* temporary vector
d2 magmaDoubleComplex* temporary vector
d_p magmaDoubleComplex* input vector p
d_r magmaDoubleComplex* input vector r
d_v magmaDoubleComplex* output vector v
skp magmaDoubleComplex* array for parameters ( skp[0]=alpha )
magma_int_t magma_zbicgmerge_spmv2 ( magma_z_sparse_matrix  A,
magmaDoubleComplex *  d1,
magmaDoubleComplex *  d2,
magmaDoubleComplex *  d_s,
magmaDoubleComplex *  d_t,
magmaDoubleComplex *  skp 
)

Merges the second SpmV using CSR with the dot product and the computation of omega.

Parameters:
A magma_z_sparse_matrix input matrix
d1 magmaDoubleComplex* temporary vector
d2 magmaDoubleComplex* temporary vector
d_s magmaDoubleComplex* input vector s
d_t magmaDoubleComplex* output vector t
skp magmaDoubleComplex* array for parameters
magma_int_t magma_zbicgmerge_xrbeta ( int  n,
magmaDoubleComplex *  d1,
magmaDoubleComplex *  d2,
magmaDoubleComplex *  rr,
magmaDoubleComplex *  r,
magmaDoubleComplex *  p,
magmaDoubleComplex *  s,
magmaDoubleComplex *  t,
magmaDoubleComplex *  x,
magmaDoubleComplex *  skp 
)

Merges the second SpmV using CSR with the dot product and the computation of omega.

Parameters:
n int dimension n
d1 magmaDoubleComplex* temporary vector
d2 magmaDoubleComplex* temporary vector
rr magmaDoubleComplex* input vector rr
r magmaDoubleComplex* input/output vector r
p magmaDoubleComplex* input vector p
s magmaDoubleComplex* input vector s
t magmaDoubleComplex* input vector t
x magmaDoubleComplex* output vector x
skp magmaDoubleComplex* array for parameters
magma_int_t magma_zcgmerge_spmv1 ( magma_z_sparse_matrix  A,
magmaDoubleComplex *  d1,
magmaDoubleComplex *  d2,
magmaDoubleComplex *  d_d,
magmaDoubleComplex *  d_z,
magmaDoubleComplex *  skp 
)

Merges the first SpmV using different formats with the dot product and the computation of rho.

Parameters:
A magma_z_sparse_matrix input matrix
d1 magmaDoubleComplex* temporary vector
d2 magmaDoubleComplex* temporary vector
d_d magmaDoubleComplex* input vector d
d_z magmaDoubleComplex* input vector z
skp magmaDoubleComplex* array for parameters ( skp[3]=rho )
void magma_zcompact ( magma_int_t  m,
magma_int_t  n,
magmaDoubleComplex *  dA,
magma_int_t  ldda,
double *  dnorms,
double  tol,
magma_index_t *  active,
magma_index_t *  cBlock 
)

ZCOMPACT takes a set of n vectors of size m (in dA) and their norms and compacts them into the cBlock size<=n vectors that have norms > tol.

The active mask array has 1 or 0, showing if a vector remained or not in the compacted resulting set of vectors.

Parameters:
[in] m INTEGER The number of rows of the matrix dA. M >= 0.
[in] n INTEGER The number of columns of the matrix dA. N >= 0.
[in,out] dA COMPLEX DOUBLE PRECISION array, dimension (LDDA,N) The m by n matrix dA.
[in] ldda INTEGER The leading dimension of the array dA. LDDA >= max(1,M).
[in] dnorms DOUBLE PRECISION array, dimension N The norms of the N vectors in dA
[in] tol DOUBLE PRECISON The tolerance value used in the criteria to compact or not.
[out] active INTEGER array, dimension N A mask of 1s and 0s showing if a vector remains or has been removed
[out] cBlock magma_index_t* The number of vectors that remain in dA (i.e., with norms > tol).
magma_int_t magma_zjacobisetup_vector_gpu ( int  num_rows,
magmaDoubleComplex *  b,
magmaDoubleComplex *  d,
magmaDoubleComplex *  c,
magmaDoubleComplex *  x 
)

Prepares the Jacobi Iteration according to x^(k+1) = D^(-1) * b - D^(-1) * (L+U) * x^k x^(k+1) = c - M * x^k.

Returns the vector c. It calls a GPU kernel

Parameters:
num_rows magma_int_t number of rows
b magma_z_vector RHS b
d magma_z_vector vector with diagonal entries
c magma_z_vector* c = D^(-1) * b
x magma_z_vector* iteration vector
magma_int_t magma_zlobpcg_maxpy ( magma_int_t  num_rows,
magma_int_t  num_vecs,
magmaDoubleComplex *  X,
magmaDoubleComplex *  Y 
)

This routine computes a axpy for a mxn matrix:.

Y = X + Y

It replaces: magma_zaxpy(m*n, c_one, Y, 1, X, 1);

/ 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
X magmaDoubleComplex* input vector X
Y magmaDoubleComplex* input/output vector Y

Generated on 17 Sep 2014 for MAGMA by  doxygen 1.6.1