MAGMA  1.5.0
Matrix Algebra for GPU and Multicore Architectures
 All Functions Groups
double-complex precision

Functions

magma_int_t magma_znan_inf (magma_uplo_t uplo, magma_int_t m, magma_int_t n, const magmaDoubleComplex *A, magma_int_t lda, magma_int_t *cnt_nan, magma_int_t *cnt_inf)
 magma_znan_inf checks a matrix that is located on the CPU host for NAN (not-a-number) and INF (infinity) values. More...
 
magma_int_t magma_znan_inf_gpu (magma_uplo_t uplo, magma_int_t m, magma_int_t n, const magmaDoubleComplex *dA, magma_int_t ldda, magma_int_t *cnt_nan, magma_int_t *cnt_inf)
 magma_znan_inf checks a matrix that is located on the CPU host for NAN (not-a-number) and INF (infinity) values. More...
 
void magma_zprint (magma_int_t m, magma_int_t n, const magmaDoubleComplex *A, magma_int_t lda)
 magma_zprint prints a matrix that is located on the CPU host. More...
 
void magma_zprint_gpu (magma_int_t m, magma_int_t n, const magmaDoubleComplex *dA, magma_int_t ldda)
 magma_zprint_gpu prints a matrix that is located on the GPU device. More...
 
void magmablas_zclaswp (magma_int_t n, magmaDoubleComplex *a, magma_int_t lda, magmaFloatComplex *sa, magma_int_t m, const magma_int_t *ipiv, magma_int_t incx)
 Row i of A is cast to single precision in row ipiv[i] of SA, for 0 <= i < M. More...
 
void magmablas_zgeadd (magma_int_t m, magma_int_t n, magmaDoubleComplex alpha, const magmaDoubleComplex *dA, magma_int_t ldda, magmaDoubleComplex *dB, magma_int_t lddb)
 ZGEADD adds two matrices, dB = alpha*dA + dB. More...
 
void magmablas_zgeadd_batched (magma_int_t m, magma_int_t n, magmaDoubleComplex alpha, const magmaDoubleComplex *const *dAarray, magma_int_t ldda, magmaDoubleComplex **dBarray, magma_int_t lddb, magma_int_t batchCount)
 ZGEADD adds two sets of matrices, dAarray[i] = alpha*dAarray[i] + dBarray[i], for i = 0, ..., batchCount-1. More...
 
void magmablas_zlacpy (magma_uplo_t uplo, magma_int_t m, magma_int_t n, const magmaDoubleComplex *dA, magma_int_t ldda, magmaDoubleComplex *dB, magma_int_t lddb)
 

Note

More...
 
void magmablas_zlacpy_batched (magma_uplo_t uplo, magma_int_t m, magma_int_t n, const magmaDoubleComplex *const *dAarray, magma_int_t ldda, magmaDoubleComplex **dBarray, magma_int_t lddb, magma_int_t batchCount)
 

Note

More...
 
void magmablas_zlag2c (magma_int_t m, magma_int_t n, const magmaDoubleComplex *A, magma_int_t lda, magmaFloatComplex *SA, magma_int_t ldsa, magma_int_t *info)
 ZLAG2C converts a double-complex matrix, A, to a single-complex matrix, SA. More...
 
double magmablas_zlange (magma_norm_t norm, magma_int_t m, magma_int_t n, const magmaDoubleComplex *A, magma_int_t lda, double *dwork)
 ZLANGE returns the value of the one norm, or the Frobenius norm, or the infinity norm, or the element of largest absolute value of a real matrix A. More...
 
double magmablas_zlanhe (magma_norm_t norm, magma_uplo_t uplo, magma_int_t n, const magmaDoubleComplex *A, magma_int_t lda, double *dwork)
 ZLANHE returns the value of the one norm, or the Frobenius norm, or the infinity norm, or the element of largest absolute value of a complex Hermitian matrix A. More...
 
void magmablas_zlaset_stream (magma_uplo_t uplo, magma_int_t m, magma_int_t n, magmaDoubleComplex offdiag, magmaDoubleComplex diag, magmaDoubleComplex *dA, magma_int_t ldda, magma_queue_t stream)
 ZLASET_STREAM initializes a 2-D array A to DIAG on the diagonal and OFFDIAG on the off-diagonals. More...
 
void magmablas_zlaset (magma_uplo_t uplo, magma_int_t m, magma_int_t n, magmaDoubleComplex offdiag, magmaDoubleComplex diag, magmaDoubleComplex *dA, magma_int_t ldda)
 
void magmablas_zlaset_band_stream (magma_uplo_t uplo, magma_int_t m, magma_int_t n, magma_int_t k, magmaDoubleComplex offdiag, magmaDoubleComplex diag, magmaDoubleComplex *dA, magma_int_t ldda, magma_queue_t stream)
 ZLASET_BAND_STREAM initializes the main diagonal of dA to DIAG, and the K-1 sub- or super-diagonals to OFFDIAG. More...
 
void magmablas_zlaset_band (magma_uplo_t uplo, magma_int_t m, magma_int_t n, magma_int_t k, magmaDoubleComplex offdiag, magmaDoubleComplex diag, magmaDoubleComplex *dA, magma_int_t ldda)
 
void magmablas_zswapdblk (magma_int_t n, magma_int_t nb, magmaDoubleComplex *dA1, magma_int_t ldda1, magma_int_t inca1, magmaDoubleComplex *dA2, magma_int_t ldda2, magma_int_t inca2)
 This is an auxiliary MAGMA routine. More...
 
void magmablas_zsymmetrize (magma_uplo_t uplo, magma_int_t m, magmaDoubleComplex *dA, magma_int_t ldda)
 ZSYMMETRIZE copies lower triangle to upper triangle, or vice-versa, to make dA a general representation of a symmetric matrix. More...
 
void magmablas_zsymmetrize_tiles (magma_uplo_t uplo, magma_int_t m, magmaDoubleComplex *dA, magma_int_t ldda, magma_int_t ntile, magma_int_t mstride, magma_int_t nstride)
 ZSYMMETRIZE_TILES copies lower triangle to upper triangle, or vice-versa, to make some blocks of dA into general representations of a symmetric block. More...
 
void magmablas_ztranspose_stream (magma_int_t m, magma_int_t n, const magmaDoubleComplex *dA, magma_int_t ldda, magmaDoubleComplex *dAT, magma_int_t lddat, magma_queue_t stream)
 ztranspose_stream copies and transposes a matrix dA to matrix dAT. More...
 
void magmablas_ztranspose (magma_int_t m, magma_int_t n, const magmaDoubleComplex *dA, magma_int_t ldda, magmaDoubleComplex *dAT, magma_int_t lddat)
 
void magmablas_ztranspose_inplace_stream (magma_int_t n, magmaDoubleComplex *dA, magma_int_t ldda, magma_queue_t stream)
 ztranspose_inplace_stream transposes a square N-by-N matrix in-place. More...
 
void magmablas_ztranspose_inplace (magma_int_t n, magmaDoubleComplex *dA, magma_int_t ldda)
 

Detailed Description

Function Documentation

magma_int_t magma_znan_inf ( magma_uplo_t  uplo,
magma_int_t  m,
magma_int_t  n,
const magmaDoubleComplex *  A,
magma_int_t  lda,
magma_int_t *  cnt_nan,
magma_int_t *  cnt_inf 
)

magma_znan_inf checks a matrix that is located on the CPU host for NAN (not-a-number) and INF (infinity) values.

NAN is created by 0/0 and similar. INF is created by x/0 and similar, where x != 0.

Parameters
[in]uplomagma_uplo_t Specifies what part of the matrix A to check.
  • = MagmaUpper: Upper triangular part of A
  • = MagmaLower: Lower triangular part of A
  • = MagmaFull: All of A
[in]mINTEGER The number of rows of the matrix A. M >= 0.
[in]nINTEGER The number of columns of the matrix A. N >= 0.
[in]ACOMPLEX_16 array, dimension (LDA,N), on the CPU host. The M-by-N matrix to be printed.
[in]ldaINTEGER The leading dimension of the array A. LDA >= max(1,M).
[out]cnt_nanINTEGER* If non-NULL, on exit contains the number of NAN values in A.
[out]cnt_infINTEGER* If non-NULL, on exit contains the number of INF values in A.
Returns
  • >= 0: Returns number of NAN + number of INF values.
  • < 0: If it returns -i, the i-th argument had an illegal value, or another error occured, such as memory allocation failed.
magma_int_t magma_znan_inf_gpu ( magma_uplo_t  uplo,
magma_int_t  m,
magma_int_t  n,
const magmaDoubleComplex *  dA,
magma_int_t  ldda,
magma_int_t *  cnt_nan,
magma_int_t *  cnt_inf 
)

magma_znan_inf checks a matrix that is located on the CPU host for NAN (not-a-number) and INF (infinity) values.

NAN is created by 0/0 and similar. INF is created by x/0 and similar, where x != 0.

Parameters
[in]uplomagma_uplo_t Specifies what part of the matrix A to check.
  • = MagmaUpper: Upper triangular part of A
  • = MagmaLower: Lower triangular part of A
  • = MagmaFull: All of A
[in]mINTEGER The number of rows of the matrix A. M >= 0.
[in]nINTEGER The number of columns of the matrix A. N >= 0.
[in,out]dACOMPLEX_16 array, dimension (LDDA,N), on the GPU device. The M-by-N matrix to be printed.
[in]lddaINTEGER The leading dimension of the array A. LDDA >= max(1,M).
[out]cnt_nanINTEGER* If non-NULL, on exit contains the number of NAN values in A.
[out]cnt_infINTEGER* If non-NULL, on exit contains the number of INF values in A.
Returns
  • >= 0: Returns number of NAN + number of INF values.
  • < 0: If it returns -i, the i-th argument had an illegal value, or another error occured, such as memory allocation failed.
void magma_zprint ( magma_int_t  m,
magma_int_t  n,
const magmaDoubleComplex *  A,
magma_int_t  lda 
)

magma_zprint prints a matrix that is located on the CPU host.

The output is intended to be Matlab compatible, to be useful in debugging.

Parameters
[in]mINTEGER The number of rows of the matrix A. M >= 0.
[in]nINTEGER The number of columns of the matrix A. N >= 0.
[in]ACOMPLEX_16 array, dimension (LDA,N), on the CPU host. The M-by-N matrix to be printed.
[in]ldaINTEGER The leading dimension of the array A. LDA >= max(1,M).
void magma_zprint_gpu ( magma_int_t  m,
magma_int_t  n,
const magmaDoubleComplex *  dA,
magma_int_t  ldda 
)

magma_zprint_gpu prints a matrix that is located on the GPU device.

Internally, it allocates CPU memory and copies the matrix to the CPU. The output is intended to be Matlab compatible, to be useful in debugging.

Parameters
[in]mINTEGER The number of rows of the matrix A. M >= 0.
[in]nINTEGER The number of columns of the matrix A. N >= 0.
[in]dACOMPLEX_16 array, dimension (LDDA,N), on the GPU device. The M-by-N matrix to be printed.
[in]lddaINTEGER The leading dimension of the array A. LDDA >= max(1,M).
void magmablas_zclaswp ( magma_int_t  n,
magmaDoubleComplex *  a,
magma_int_t  lda,
magmaFloatComplex *  sa,
magma_int_t  m,
const magma_int_t *  ipiv,
magma_int_t  incx 
)

Row i of A is cast to single precision in row ipiv[i] of SA, for 0 <= i < M.

N - (input) INTEGER. On entry, N specifies the number of columns of the matrix A.

A - (input) DOUBLE PRECISION array on the GPU, dimension (LDA,N) On entry, the M-by-N matrix to which the row interchanges will be applied.

LDA - (input) INTEGER. LDA specifies the leading dimension of A.

SA - (output) REAL array on the GPU, dimension (LDA,N) On exit, the single precision, permuted matrix.

M - (input) The number of rows to be interchanged.

IPIV - (input) INTEGER array on the GPU, dimension (M) The vector of pivot indices. Row i of A is cast to single precision in row ipiv[i] of SA, for 0 <= i < m.

INCX - (input) INTEGER If IPIV is negative, the pivots are applied in reverse order, otherwise in straight-forward order.

void magmablas_zgeadd ( magma_int_t  m,
magma_int_t  n,
magmaDoubleComplex  alpha,
const magmaDoubleComplex *  dA,
magma_int_t  ldda,
magmaDoubleComplex *  dB,
magma_int_t  lddb 
)

ZGEADD adds two matrices, dB = alpha*dA + dB.

Parameters
[in]mINTEGER The number of rows of the matrix dA. M >= 0.
[in]nINTEGER The number of columns of the matrix dA. N >= 0.
[in]alphaCOMPLEX_16 The scalar alpha.
[in]dACOMPLEX_16 array, dimension (LDDA,N) The m by n matrix dA.
[in]lddaINTEGER The leading dimension of the array dA. LDDA >= max(1,M).
[in,out]dBCOMPLEX_16 array, dimension (LDDB,N) The m by n matrix dB.
[in]lddbINTEGER The leading dimension of the array dB. LDDB >= max(1,M).
void magmablas_zgeadd_batched ( magma_int_t  m,
magma_int_t  n,
magmaDoubleComplex  alpha,
const magmaDoubleComplex *const *  dAarray,
magma_int_t  ldda,
magmaDoubleComplex **  dBarray,
magma_int_t  lddb,
magma_int_t  batchCount 
)

ZGEADD adds two sets of matrices, dAarray[i] = alpha*dAarray[i] + dBarray[i], for i = 0, ..., batchCount-1.

Parameters
[in]mINTEGER The number of rows of each matrix dAarray[i]. M >= 0.
[in]nINTEGER The number of columns of each matrix dAarray[i]. N >= 0.
[in]alphaCOMPLEX_16 The scalar alpha.
[in]dAarrayarray on GPU, dimension(batchCount), of pointers to arrays, with each array a COMPLEX_16 array, dimension (LDDA,N) The m by n matrices dAarray[i].
[in]lddaINTEGER The leading dimension of each array dAarray[i]. LDDA >= max(1,M).
[in,out]dBarrayarray on GPU, dimension(batchCount), of pointers to arrays, with each array a COMPLEX_16 array, dimension (LDDB,N) The m by n matrices dBarray[i].
[in]lddbINTEGER The leading dimension of each array dBarray[i]. LDDB >= max(1,M).
[in]batchCountINTEGER The number of matrices to add; length of dAarray and dBarray. batchCount >= 0.
void magmablas_zlacpy ( magma_uplo_t  uplo,
magma_int_t  m,
magma_int_t  n,
const magmaDoubleComplex *  dA,
magma_int_t  ldda,
magmaDoubleComplex *  dB,
magma_int_t  lddb 
)

Note

  • UPLO Parameter is disabled
  • Do we want to provide a generic function to the user with all the options?

ZLACPY copies all or part of a two-dimensional matrix dA to another matrix dB.

Parameters
[in]uplomagma_uplo_t Specifies the part of the matrix dA to be copied to dB.
  • = MagmaUpper: Upper triangular part
  • = MagmaLower: Lower triangular part Otherwise: All of the matrix dA
[in]mINTEGER The number of rows of the matrix dA. M >= 0.
[in]nINTEGER The number of columns of the matrix dA. N >= 0.
[in]dACOMPLEX_16 array, dimension (LDDA,N) The m by n matrix dA. If UPLO = MagmaUpper, only the upper triangle or trapezoid is accessed; if UPLO = MagmaLower, only the lower triangle or trapezoid is accessed.
[in]lddaINTEGER The leading dimension of the array dA. LDDA >= max(1,M).
[out]dBCOMPLEX_16 array, dimension (LDDB,N) The m by n matrix dB. On exit, dB = dA in the locations specified by UPLO.
[in]lddbINTEGER The leading dimension of the array dB. LDDB >= max(1,M).
void magmablas_zlacpy_batched ( magma_uplo_t  uplo,
magma_int_t  m,
magma_int_t  n,
const magmaDoubleComplex *const *  dAarray,
magma_int_t  ldda,
magmaDoubleComplex **  dBarray,
magma_int_t  lddb,
magma_int_t  batchCount 
)

Note

  • UPLO Parameter is disabled
  • Do we want to provide a generic function to the user with all the options?

ZLACPY copies all or part of a set of two-dimensional matrices dAarray[i] to another set of matrices dBarray[i], for i = 0, ..., batchCount-1.

Parameters
[in]uplomagma_uplo_t Specifies the part of each matrix dAarray[i] to be copied to dBarray[i].
  • = MagmaUpper: Upper triangular part
  • = MagmaLower: Lower triangular part Otherwise: All of each matrix dAarray[i]
[in]mINTEGER The number of rows of each matrix dAarray[i]. M >= 0.
[in]nINTEGER The number of columns of each matrix dAarray[i]. N >= 0.
[in]dAarrayarray on GPU, dimension(batchCount), of pointers to arrays, with each array a COMPLEX_16 array, dimension (LDDA,N) The m by n matrices dAarray[i]. If UPLO = MagmaUpper, only the upper triangle or trapezoid is accessed; if UPLO = MagmaLower, only the lower triangle or trapezoid is accessed.
[in]lddaINTEGER The leading dimension of each array dAarray[i]. LDDA >= max(1,M).
[out]dBarrayarray on GPU, dimension(batchCount), of pointers to arrays, with each array a COMPLEX_16 array, dimension (LDDB,N) The m by n matrices dBarray[i]. On exit, matrix dBarray[i] = matrix dAarray[i] in the locations specified by UPLO.
[in]lddbINTEGER The leading dimension of each array dBarray[i]. LDDB >= max(1,M).
[in]batchCountINTEGER The number of matrices to add; length of dAarray and dBarray. batchCount >= 0.
void magmablas_zlag2c ( magma_int_t  m,
magma_int_t  n,
const magmaDoubleComplex *  A,
magma_int_t  lda,
magmaFloatComplex *  SA,
magma_int_t  ldsa,
magma_int_t *  info 
)

ZLAG2C converts a double-complex matrix, A, to a single-complex matrix, SA.

RMAX is the overflow for the single-complex arithmetic. ZLAG2C checks that all the entries of A are between -RMAX and RMAX. If not, the conversion is aborted and a flag is raised.

Parameters
[in]mINTEGER The number of lines of the matrix A. m >= 0.
[in]nINTEGER The number of columns of the matrix A. n >= 0.
[in]ACOMPLEX_16 array, dimension (LDA,n) On entry, the m-by-n coefficient matrix A.
[in]ldaINTEGER The leading dimension of the array A. LDA >= max(1,m).
[out]SACOMPLEX array, dimension (LDSA,n) On exit, if INFO=0, the m-by-n coefficient matrix SA; if INFO>0, the content of SA is unspecified.
[in]ldsaINTEGER The leading dimension of the array SA. LDSA >= max(1,m).
[out]infoINTEGER
  • = 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.
double magmablas_zlange ( magma_norm_t  norm,
magma_int_t  m,
magma_int_t  n,
const magmaDoubleComplex *  A,
magma_int_t  lda,
double *  dwork 
)

ZLANGE returns the value of the one norm, or the Frobenius norm, or the infinity norm, or the element of largest absolute value of a real matrix A.

Description

ZLANGE returns the value

ZLANGE = ( max(abs(A(i,j))), NORM = 'M' or 'm' ** not yet supported ( ( norm1(A), NORM = '1', 'O' or 'o' ** not yet supported ( ( normI(A), NORM = 'I' or 'i' ( ( normF(A), NORM = 'F', 'f', 'E' or 'e' ** not yet supported

where norm1 denotes the one norm of a matrix (maximum column sum), normI denotes the infinity norm of a matrix (maximum row sum) and normF denotes the Frobenius norm of a matrix (square root of sum of squares). Note that max(abs(A(i,j))) is not a consistent matrix norm.

Parameters
[in]normCHARACTER*1 Specifies the value to be returned in ZLANGE as described above.
[in]mINTEGER The number of rows of the matrix A. M >= 0. When M = 0, ZLANGE is set to zero.
[in]nINTEGER The number of columns of the matrix A. N >= 0. When N = 0, ZLANGE is set to zero.
[in]ADOUBLE PRECISION array on the GPU, dimension (LDA,N) The m by n matrix A.
[in]ldaINTEGER The leading dimension of the array A. LDA >= max(M,1).
dwork(workspace) DOUBLE PRECISION array on the GPU, dimension (MAX(1,LWORK)), where LWORK >= M when NORM = 'I'; otherwise, WORK is not referenced.
double magmablas_zlanhe ( magma_norm_t  norm,
magma_uplo_t  uplo,
magma_int_t  n,
const magmaDoubleComplex *  A,
magma_int_t  lda,
double *  dwork 
)

ZLANHE returns the value of the one norm, or the Frobenius norm, or the infinity norm, or the element of largest absolute value of a complex Hermitian matrix A.

ZLANHE = ( max(abs(A(i,j))), NORM = 'M' or 'm' ( ( norm1(A), NORM = '1', 'O' or 'o' ** supported only for (PRECISION_s || PRECISION_d || PRECISION_c || CUDA_ARCH >= 200) ( ( normI(A), NORM = 'I' or 'i' ** supported only for (PRECISION_s || PRECISION_d || PRECISION_c || CUDA_ARCH >= 200) ( ( normF(A), NORM = 'F', 'f', 'E' or 'e' ** not yet supported

where norm1 denotes the one norm of a matrix (maximum column sum), normI denotes the infinity norm of a matrix (maximum row sum) and normF denotes the Frobenius norm of a matrix (square root of sum of squares). Note that max(abs(A(i,j))) is not a consistent matrix norm.

Returns ZLANHE < 0: if ZLANHE = -i, the i-th argument had an illegal value.

Arguments:

Parameters
[in]normCHARACTER*1 Specifies the value to be returned in ZLANHE as described above.
[in]uplomagma_uplo_t Specifies whether the upper or lower triangular part of the Hermitian matrix A is to be referenced.
  • = MagmaUpper: Upper triangular part of A is referenced
  • = MagmaLower: Lower triangular part of A is referenced
[in]nINTEGER The order of the matrix A. N >= 0. When N = 0, ZLANHE is set to zero.
[in]ACOMPLEX*16 array on the GPU, dimension (LDA,N) The Hermitian matrix A. If UPLO = MagmaUpper, the leading n by n upper triangular part of A contains the upper triangular part of the matrix A, and the strictly lower triangular part of A is not referenced. If UPLO = MagmaLower, the leading n by n lower triangular part of A contains the lower triangular part of the matrix A, and the strictly upper triangular part of A is not referenced. Note that the imaginary parts of the diagonal elements need not be set and are assumed to be zero.
[in]ldaINTEGER The leading dimension of the array A. LDA >= max(N,1).
dwork(workspace) DOUBLE PRECISION array on the GPU, dimension (MAX(1,LWORK)), where LWORK >= N. NOTE: this is different than LAPACK, where WORK is required only for norm1 and normI. Here max-norm also requires work.
void magmablas_zlaset ( magma_uplo_t  uplo,
magma_int_t  m,
magma_int_t  n,
magmaDoubleComplex  offdiag,
magmaDoubleComplex  diag,
magmaDoubleComplex *  dA,
magma_int_t  ldda 
)
void magmablas_zlaset_band ( magma_uplo_t  uplo,
magma_int_t  m,
magma_int_t  n,
magma_int_t  k,
magmaDoubleComplex  offdiag,
magmaDoubleComplex  diag,
magmaDoubleComplex *  dA,
magma_int_t  ldda 
)
void magmablas_zlaset_band_stream ( magma_uplo_t  uplo,
magma_int_t  m,
magma_int_t  n,
magma_int_t  k,
magmaDoubleComplex  offdiag,
magmaDoubleComplex  diag,
magmaDoubleComplex *  dA,
magma_int_t  ldda,
magma_queue_t  stream 
)

ZLASET_BAND_STREAM initializes the main diagonal of dA to DIAG, and the K-1 sub- or super-diagonals to OFFDIAG.

This is the same as ZLASET_BAND, but adds stream argument.

Parameters
[in]uplomagma_uplo_t Specifies the part of the matrix dA to be set.
  • = MagmaUpper: Upper triangular part
  • = MagmaLower: Lower triangular part
[in]mINTEGER The number of rows of the matrix dA. M >= 0.
[in]nINTEGER The number of columns of the matrix dA. N >= 0.
[in]kINTEGER The number of diagonals to set, including the main diagonal. K >= 0. Currently, K <= 1024 due to CUDA restrictions (max. number of threads per block).
[in]offdiagCOMPLEX_16 Off-diagonal elements in the band are set to OFFDIAG.
[in]diagCOMPLEX_16 All the main diagonal elements are set to DIAG.
[in]dACOMPLEX_16 array, dimension (LDDA,N) The M-by-N matrix dA. If UPLO = MagmaUpper, only the upper triangle or trapezoid is accessed; if UPLO = MagmaLower, only the lower triangle or trapezoid is accessed. On exit, A(i,j) = ALPHA, 1 <= i <= m, 1 <= j <= n where i != j, abs(i-j) < k; A(i,i) = BETA , 1 <= i <= min(m,n)
[in]lddaINTEGER The leading dimension of the array dA. LDDA >= max(1,M).
[in]streammagma_queue_t Stream to execute ZLASET in.
Author
Raffaele Solca
Mark Gates
void magmablas_zlaset_stream ( magma_uplo_t  uplo,
magma_int_t  m,
magma_int_t  n,
magmaDoubleComplex  offdiag,
magmaDoubleComplex  diag,
magmaDoubleComplex *  dA,
magma_int_t  ldda,
magma_queue_t  stream 
)

ZLASET_STREAM initializes a 2-D array A to DIAG on the diagonal and OFFDIAG on the off-diagonals.

This is the same as ZLASET, but adds stream argument.

Parameters
[in]uplomagma_uplo_t Specifies the part of the matrix dA to be set.
  • = MagmaUpper: Upper triangular part
  • = MagmaLower: Lower triangular part Otherwise: All of the matrix dA is set.
[in]mINTEGER The number of rows of the matrix dA. M >= 0.
[in]nINTEGER The number of columns of the matrix dA. N >= 0.
[in]offdiagCOMPLEX_16 The scalar OFFDIAG. (In LAPACK this is called ALPHA.)
[in]diagCOMPLEX_16 The scalar DIAG. (In LAPACK this is called BETA.)
[in]dACOMPLEX_16 array, dimension (LDDA,N) The M-by-N matrix dA. If UPLO = MagmaUpper, only the upper triangle or trapezoid is accessed; if UPLO = MagmaLower, only the lower triangle or trapezoid is accessed. On exit, A(i,j) = OFFDIAG, 1 <= i <= m, 1 <= j <= n, i != j; A(i,i) = DIAG, 1 <= i <= min(m,n)
[in]lddaINTEGER The leading dimension of the array dA. LDDA >= max(1,M).
[in]streammagma_queue_t Stream to execute in.
void magmablas_zswapdblk ( magma_int_t  n,
magma_int_t  nb,
magmaDoubleComplex *  dA1,
magma_int_t  ldda1,
magma_int_t  inca1,
magmaDoubleComplex *  dA2,
magma_int_t  ldda2,
magma_int_t  inca2 
)

This is an auxiliary MAGMA routine.

It swaps diagonal blocks of size nb x nb between matrices dA1 and dA2 on the GPU.

The number of blocks swapped is (n-1)/nb. For i = 1 .. (n-1)/nb matrices dA1 + i * nb * (ldda1 + inca1) and dA2 + i * nb * (ldda2 + inca2) are swapped.

void magmablas_zsymmetrize ( magma_uplo_t  uplo,
magma_int_t  m,
magmaDoubleComplex *  dA,
magma_int_t  ldda 
)

ZSYMMETRIZE copies lower triangle to upper triangle, or vice-versa, to make dA a general representation of a symmetric matrix.

Parameters
[in]uplomagma_uplo_t Specifies the part of the matrix dA that is valid on input.
  • = MagmaUpper: Upper triangular part
  • = MagmaLower: Lower triangular part
[in]mINTEGER The number of rows of the matrix dA. M >= 0.
[in,out]dACOMPLEX_16 array, dimension (LDDA,N) The m by m matrix dA.
[in]lddaINTEGER The leading dimension of the array dA. LDDA >= max(1,M).
void magmablas_zsymmetrize_tiles ( magma_uplo_t  uplo,
magma_int_t  m,
magmaDoubleComplex *  dA,
magma_int_t  ldda,
magma_int_t  ntile,
magma_int_t  mstride,
magma_int_t  nstride 
)

ZSYMMETRIZE_TILES copies lower triangle to upper triangle, or vice-versa, to make some blocks of dA into general representations of a symmetric block.

This processes NTILE blocks, typically the diagonal blocks. Each block is offset by mstride rows and nstride columns from the previous block.

Parameters
[in]uplomagma_uplo_t Specifies the part of the matrix dA that is valid on input.
  • = MagmaUpper: Upper triangular part
  • = MagmaLower: Lower triangular part
[in]mINTEGER The number of rows & columns of each square block of dA. M >= 0.
[in,out]dACOMPLEX_16 array, dimension (LDDA,N) The matrix dA. N = m + nstride*(ntile-1).
[in]lddaINTEGER The leading dimension of the array dA. LDDA >= max(1, m + mstride*(ntile-1)).
[in]ntileINTEGER Number of blocks to symmetrize.
[in]mstrideINTEGER Row offset from start of one block to start of next block.
[in]nstrideINTEGER Column offset from start of one block to start of next block.
void magmablas_ztranspose ( magma_int_t  m,
magma_int_t  n,
const magmaDoubleComplex *  dA,
magma_int_t  ldda,
magmaDoubleComplex *  dAT,
magma_int_t  lddat 
)
void magmablas_ztranspose_inplace ( magma_int_t  n,
magmaDoubleComplex *  dA,
magma_int_t  ldda 
)
void magmablas_ztranspose_inplace_stream ( magma_int_t  n,
magmaDoubleComplex *  dA,
magma_int_t  ldda,
magma_queue_t  stream 
)

ztranspose_inplace_stream transposes a square N-by-N matrix in-place.

Same as ztranspose_inplace, but adds stream argument.

Parameters
[in]nINTEGER The number of rows & columns of the matrix dA. N >= 0.
[in]dACOMPLEX_16 array, dimension (LDDA,N) The N-by-N matrix dA. On exit, dA(j,i) = dA_original(i,j), for 0 <= i,j < N.
[in]lddaINTEGER The leading dimension of the array dA. LDDA >= N.
[in]streammagma_queue_t Stream to execute in.
void magmablas_ztranspose_stream ( magma_int_t  m,
magma_int_t  n,
const magmaDoubleComplex *  dA,
magma_int_t  ldda,
magmaDoubleComplex *  dAT,
magma_int_t  lddat,
magma_queue_t  stream 
)

ztranspose_stream copies and transposes a matrix dA to matrix dAT.

Same as ztranspose, but adds stream argument.

Parameters
[in]mINTEGER The number of rows of the matrix dA. M >= 0.
[in]nINTEGER The number of columns of the matrix dA. N >= 0.
[in]dACOMPLEX_16 array, dimension (LDDA,N) The M-by-N matrix dA.
[in]lddaINTEGER The leading dimension of the array dA. LDDA >= M.
[in]dATCOMPLEX_16 array, dimension (LDDA,N) The N-by-M matrix dAT.
[in]lddatINTEGER The leading dimension of the array dAT. LDDAT >= N.
[in]streammagma_queue_t Stream to execute in.