#include "slu_mt_machines.h"
#include "slu_mt_Cnames.h"
#include "supermatrix.h"
#include "slu_mt_util.h"
#include "pxgstrf_synch.h"
Go to the source code of this file.
Classes | |
struct | GlobalLU_t |
struct | pxgstrf_shared_t |
struct | pdgstrf_threadarg_t |
Typedefs | |
typedef int | int_t |
Functions | |
void | pdgssv (int, SuperMatrix *, int *, int *, SuperMatrix *, SuperMatrix *, SuperMatrix *, int *) |
void | pdgssvx (int, superlumt_options_t *, SuperMatrix *, int *, int *, equed_t *, double *, double *, SuperMatrix *, SuperMatrix *, SuperMatrix *, SuperMatrix *, double *, double *, double *, double *, superlu_memusage_t *, int *) |
void | dgsequ (SuperMatrix *, double *, double *, double *, double *, double *, int *) |
Driver related. | |
void | dlaqgs (SuperMatrix *, double *, double *, double, double, double, equed_t *) |
void | dgscon (char *, SuperMatrix *, SuperMatrix *, double, double *, int *) |
double | dPivotGrowth (int, SuperMatrix *, int *, SuperMatrix *, SuperMatrix *) |
void | dgsrfs (trans_t, SuperMatrix *, SuperMatrix *, SuperMatrix *, int *, int *, equed_t, double *, double *, SuperMatrix *, SuperMatrix *, double *, double *, Gstat_t *, int *) |
int | sp_dtrsv (char *, char *, char *, SuperMatrix *, SuperMatrix *, double *, int *) |
int | sp_dgemv (char *, double, SuperMatrix *, double *, int, double, double *, int) |
Performs one of the matrix-vector operations y := alpha*A*x + beta*y, or y := alpha*A'*x + beta*y,. | |
int | sp_dgemm (char *, int, int, int, double, SuperMatrix *, double *, int, double, double *, int) |
void | pxgstrf_scheduler (const int, const int, const int *, int *, int *, pxgstrf_shared_t *) |
int | dParallelInit (int, pxgstrf_relax_t *, superlumt_options_t *, pxgstrf_shared_t *) |
int | ParallelFinalize () |
void | pdgstrf_StackFree () |
int | queue_init (queue_t *, int) |
int | queue_destroy (queue_t *) |
int | EnqueueRelaxSnode (queue_t *, int, pxgstrf_relax_t *, pxgstrf_shared_t *) |
int | EnqueueDomains (queue_t *, struct Branch *, pxgstrf_shared_t *) |
int | Enqueue (queue_t *, qitem_t) |
int | Dequeue (queue_t *, qitem_t *) |
int | NewNsuper (const int, pxgstrf_shared_t *, int *) |
int | lockon (int *) |
void | PartDomains (const int, const float, SuperMatrix *, int *, int *) |
void | dCreate_CompCol_Matrix (SuperMatrix *, int, int, int, double *, int *, int *, Stype_t, Dtype_t, Mtype_t) |
Supernodal LU factor related. | |
void | dCreate_CompCol_Permuted (SuperMatrix *, int, int, int, double *, int *, int *, int *, Stype_t, Dtype_t, Mtype_t) |
void | dCopy_CompCol_Matrix (SuperMatrix *, SuperMatrix *) |
Copy matrix A into matrix B. | |
void | dCreate_Dense_Matrix (SuperMatrix *, int, int, double *, int, Stype_t, Dtype_t, Mtype_t) |
void | dCreate_SuperNode_Matrix (SuperMatrix *, int, int, int, double *, int *, int *, int *, int *, int *, Stype_t, Dtype_t, Mtype_t) |
void | dCreate_SuperNode_Permuted (SuperMatrix *, int, int, int, double *, int *, int *, int *, int *, int *, int *, int *, int *, Stype_t, Dtype_t, Mtype_t) |
void | dCopy_Dense_Matrix (int, int, double *, int, double *, int) |
void | Destroy_SuperMatrix_Store (SuperMatrix *) |
Deallocate the structure pointing to the actual storage of the matrix. | |
void | Destroy_CompCol_Matrix (SuperMatrix *) |
void | Destroy_CompCol_Permuted (SuperMatrix *) |
A is of type Stype==NCP. | |
void | Destroy_CompCol_NCP (SuperMatrix *) |
void | Destroy_SuperNode_Matrix (SuperMatrix *) |
void | Destroy_SuperNode_SCP (SuperMatrix *) |
void | dallocateA (int, int, double **, int **, int **) |
Allocate storage for original matrix A. | |
void | StatAlloc (const int, const int, const int, const int, Gstat_t *) |
void | StatInit (const int, const int, Gstat_t *) |
void | StatFree (Gstat_t *) |
void | get_perm_c (int, SuperMatrix *, int *) |
void | dsp_colorder (SuperMatrix *, int *, superlumt_options_t *, SuperMatrix *) |
int | sp_coletree (int *, int *, int *, int, int, int *) |
int | dPresetMap (const int, SuperMatrix *, pxgstrf_relax_t *, superlumt_options_t *, GlobalLU_t *) |
int | qrnzcnt (int, int, int *, int *, int *, int *, int *, int *, int *, int *, int *, int *) |
int | DynamicSetMap (const int, const int, const int, pxgstrf_shared_t *) |
void | pdgstrf (superlumt_options_t *, SuperMatrix *, int *, SuperMatrix *, SuperMatrix *, Gstat_t *, int *) |
void | pdgstrf_init (int, fact_t, trans_t, yes_no_t, int, int, double, yes_no_t, double, int *, int *, void *, int, SuperMatrix *, SuperMatrix *, superlumt_options_t *, Gstat_t *) |
pdgstrf_threadarg_t * | pdgstrf_thread_init (SuperMatrix *, SuperMatrix *, SuperMatrix *, superlumt_options_t *, pxgstrf_shared_t *, Gstat_t *, int *) |
void | pdgstrf_thread_finalize (pdgstrf_threadarg_t *, pxgstrf_shared_t *, SuperMatrix *, int *, SuperMatrix *, SuperMatrix *) |
void | pdgstrf_finalize (superlumt_options_t *, SuperMatrix *) |
void | pxgstrf_finalize (superlumt_options_t *, SuperMatrix *) |
void | pdgstrf_relax_snode (const int, superlumt_options_t *, pxgstrf_relax_t *) |
int | pdgstrf_factor_snode (const int, const int, SuperMatrix *, const double, yes_no_t *, int *, int *, int *, int *, int *, int *, double *, double *, pxgstrf_shared_t *, int *) |
void | pxgstrf_mark_busy_descends (int, int, int *, pxgstrf_shared_t *, int *, int *) |
int | pdgstrf_snode_dfs (const int, const int, const int, const int *, const int *, const int *, int *, int *, int *, pxgstrf_shared_t *) |
int | pdgstrf_snode_bmod (const int, const int, const int, const int, double *, double *, GlobalLU_t *, Gstat_t *) |
void | pdgstrf_panel_dfs (const int, const int, const int, const int, SuperMatrix *, int *, int *, int *, int *, int *, int *, int *, int *, int *, int *, int *, int *, int *, double *, GlobalLU_t *) |
void | pdgstrf_panel_bmod (const int, const int, const int, const int, const int, int *, int *, int *, int *, int *, int *, int *, int *, double *, double *, pxgstrf_shared_t *) |
void | pdgstrf_bmod1D (const int, const int, const int, const int, const int, const int, const int, int, int, int *, int *, int *, int *, double *, double *, GlobalLU_t *, Gstat_t *) |
void | pdgstrf_bmod2D (const int, const int, const int, const int, const int, const int, const int, int, int, int *, int *, int *, int *, double *, double *, GlobalLU_t *, Gstat_t *) |
void | pdgstrf_bmod1D_mv2 (const int, const int, const int, const int, const int, const int, const int, int, int, int *, int *, int *, int *, double *, double *, GlobalLU_t *, Gstat_t *) |
void | pdgstrf_bmod2D_mv2 (const int, const int, const int, const int, const int, const int, const int, int, int, int *, int *, int *, int *, double *, double *, GlobalLU_t *, Gstat_t *) |
void | pxgstrf_super_bnd_dfs (const int, const int, const int, const int, const int, SuperMatrix *, int *, int *, int *, int *, int *, int *, int *, pxgstrf_shared_t *) |
int | pdgstrf_column_dfs (const int, const int, const int, const int, int *, int *, int *, int, int *, int *, int *, int *, int *, int *, int *, int *, pxgstrf_shared_t *) |
int | pdgstrf_column_bmod (const int, const int, const int, const int, int *, int *, double *, double *, pxgstrf_shared_t *, Gstat_t *) |
int | pdgstrf_pivotL (const int, const int, const double, yes_no_t *, int *, int *, int *, int *, GlobalLU_t *, Gstat_t *) |
int | pdgstrf_copy_to_ucol (const int, const int, const int, const int *, const int *, const int *, double *, pxgstrf_shared_t *) |
void | pxgstrf_pruneL (const int, const int *, const int, const int, const int *, const int *, int *, int *, GlobalLU_t *) |
void | pxgstrf_resetrep_col (const int, const int *, int *) |
void | countnz (const int, int *, int *, int *, GlobalLU_t *) |
Count the total number of nonzeros in factors L and U, and in the symmetrically reduced L. | |
void | fixupL (const int, const int *, GlobalLU_t *) |
Fix up the data storage lsub for L-subscripts. It removes the subscript sets for structural pruning, and applies permuation to the remaining subscripts. | |
void | compressSUP (const int, GlobalLU_t *) |
int | spcoletree (int *, int *, int *, int, int, int *) |
int * | TreePostorder (int, int *) |
void | dreadmt (int *, int *, int *, double **, int **, int **) |
void | dreadhb (int *, int *, int *, double **, int **, int **) |
void | dGenXtrue (int, int, double *, int) |
void | dFillRHS (trans_t, int, double *, int, SuperMatrix *, SuperMatrix *) |
Let rhs[i] = sum of i-th row of A, so the solution vector is all 1's. | |
void | dgstrs (trans_t, SuperMatrix *, SuperMatrix *, int *, int *, SuperMatrix *, Gstat_t *, int *) |
void | dlsolve (int, int, double *, double *) |
Solves a dense UNIT lower triangular system. | |
void | dusolve (int, int, double *, double *) |
Solves a dense upper triangular system. | |
void | dmatvec (int, int, int, double *, double *, double *) |
Performs a dense matrix-vector multiply: Mxvec = Mxvec + M * vec. | |
int | dgemm_ (char *, char *, int *, int *, int *, double *, double *, int *, double *, int *, double *, double *, int *) |
int | dtrsm_ (char *, char *, char *, char *, int *, int *, double *, double *, int *, double *, int *) |
int | dtrsv_ (char *, char *, char *, int *, double *, int *, double *, int *) |
int | dgemv_ (char *, int *, int *, double *, double *, int *, double *, int *, double *, double *, int *) |
float | pdgstrf_MemInit (int, int, superlumt_options_t *, SuperMatrix *, SuperMatrix *, GlobalLU_t *) |
float | pdgstrf_memory_use (const int, const int, const int) |
int | pdgstrf_WorkInit (int, int, int **, double **) |
void | pxgstrf_SetIWork (int, int, int *, int **, int **, int **, int **, int **, int **, int **) |
void | pdgstrf_SetRWork (int, int, double *, double **, double **) |
void | pdgstrf_WorkFree (int *, double *, GlobalLU_t *) |
int | pdgstrf_MemXpand (int, int, MemType, int *, GlobalLU_t *) |
int * | intMalloc (int) |
int * | intCalloc (int) |
double * | doubleMalloc (int) |
double * | doubleCalloc (int) |
int | memory_usage () |
int | superlu_dQuerySpace (int, SuperMatrix *, SuperMatrix *, int, superlu_memusage_t *) |
int | Glu_alloc (const int, const int, const int, const MemType, int *, pxgstrf_shared_t *) |
double | SuperLU_timer_ () |
Timer function. | |
int | sp_ienv (int) |
double | dlamch_ () |
int | lsame_ (char *, char *) |
int | xerbla_ (char *, int *) |
void | superlu_abort_and_exit (char *) |
Global statistics variale. | |
void | ifill (int *, int, int) |
Fills an integer array with a given value. | |
void | dfill (double *, int, double) |
Fills a double precision array with a given value. | |
void | dinf_norm_error (int, SuperMatrix *, double *) |
Check the inf-norm of the error vector. | |
void | dstat_allocate (int) |
void | snode_profile (int, int *) |
void | super_stats (int, int *, int *) |
void | panel_stats (int, int, int *, Gstat_t *) |
void | PrintSumm (char *, int, int, int) |
Print a summary of the testing results. | |
void | dPrintPerf (SuperMatrix *, SuperMatrix *, superlu_memusage_t *, double, double, double *, double *, char *, Gstat_t *) |
void | dCompRow_to_CompCol (int m, int n, int nnz, double *a, int *colind, int *rowptr, double **at, int **rowind, int **colptr) |
Convert a row compressed storage into a column compressed storage. | |
void | print_lu_col (int, char *, int, int, int, int *, GlobalLU_t *) |
void | print_panel_seg (int, int, int, int, int *, int *) |
Diagnostic print of segment info after panel_dfs(). | |
void | dcheck_zero_vec (int, char *, int, double *) |
void | check_repfnz (int, int, int, int *) |
Check whether repfnz[] == EMPTY after reset. |
typedef int int_t |
void check_repfnz | ( | int | , | |
int | , | |||
int | , | |||
int * | ||||
) |
Check whether repfnz[] == EMPTY after reset.
void compressSUP | ( | const | int, | |
GlobalLU_t * | ||||
) |
void countnz | ( | const | int, | |
int * | , | |||
int * | , | |||
int * | , | |||
GlobalLU_t * | ||||
) |
Count the total number of nonzeros in factors L and U, and in the symmetrically reduced L.
void dallocateA | ( | int | , | |
int | , | |||
double ** | , | |||
int ** | , | |||
int ** | ||||
) |
Allocate storage for original matrix A.
void dcheck_zero_vec | ( | int | , | |
char * | , | |||
int | , | |||
double * | ||||
) |
void dCompRow_to_CompCol | ( | int | m, | |
int | n, | |||
int | nnz, | |||
double * | a, | |||
int * | colind, | |||
int * | rowptr, | |||
double ** | at, | |||
int ** | rowind, | |||
int ** | colptr | |||
) |
Convert a row compressed storage into a column compressed storage.
void dCopy_CompCol_Matrix | ( | SuperMatrix * | , | |
SuperMatrix * | ||||
) |
Copy matrix A into matrix B.
void dCopy_Dense_Matrix | ( | int | , | |
int | , | |||
double * | , | |||
int | , | |||
double * | , | |||
int | ||||
) |
Copies a two-dimensional matrix X to another matrix Y.
void dCreate_CompCol_Matrix | ( | SuperMatrix * | , | |
int | , | |||
int | , | |||
int | , | |||
double * | , | |||
int * | , | |||
int * | , | |||
Stype_t | , | |||
Dtype_t | , | |||
Mtype_t | ||||
) |
Supernodal LU factor related.
void dCreate_CompCol_Permuted | ( | SuperMatrix * | , | |
int | , | |||
int | , | |||
int | , | |||
double * | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
Stype_t | , | |||
Dtype_t | , | |||
Mtype_t | ||||
) |
void dCreate_Dense_Matrix | ( | SuperMatrix * | , | |
int | , | |||
int | , | |||
double * | , | |||
int | , | |||
Stype_t | , | |||
Dtype_t | , | |||
Mtype_t | ||||
) |
void dCreate_SuperNode_Matrix | ( | SuperMatrix * | , | |
int | , | |||
int | , | |||
int | , | |||
double * | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
Stype_t | , | |||
Dtype_t | , | |||
Mtype_t | ||||
) |
void dCreate_SuperNode_Permuted | ( | SuperMatrix * | , | |
int | , | |||
int | , | |||
int | , | |||
double * | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
Stype_t | , | |||
Dtype_t | , | |||
Mtype_t | ||||
) |
void Destroy_CompCol_Matrix | ( | SuperMatrix * | ) |
void Destroy_CompCol_NCP | ( | SuperMatrix * | ) |
void Destroy_CompCol_Permuted | ( | SuperMatrix * | ) |
A is of type Stype==NCP.
void Destroy_SuperMatrix_Store | ( | SuperMatrix * | ) |
Deallocate the structure pointing to the actual storage of the matrix.
void Destroy_SuperNode_Matrix | ( | SuperMatrix * | ) |
void Destroy_SuperNode_SCP | ( | SuperMatrix * | ) |
void dfill | ( | double * | , | |
int | , | |||
double | ||||
) |
Fills a double precision array with a given value.
void dFillRHS | ( | trans_t | , | |
int | , | |||
double * | , | |||
int | , | |||
SuperMatrix * | , | |||
SuperMatrix * | ||||
) |
Let rhs[i] = sum of i-th row of A, so the solution vector is all 1's.
int dgemm_ | ( | char * | , | |
char * | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
double * | , | |||
double * | , | |||
int * | , | |||
double * | , | |||
int * | , | |||
double * | , | |||
double * | , | |||
int * | ||||
) |
int dgemv_ | ( | char * | , | |
int * | , | |||
int * | , | |||
double * | , | |||
double * | , | |||
int * | , | |||
double * | , | |||
int * | , | |||
double * | , | |||
double * | , | |||
int * | ||||
) |
void dGenXtrue | ( | int | , | |
int | , | |||
double * | , | |||
int | ||||
) |
void dgscon | ( | char * | , | |
SuperMatrix * | , | |||
SuperMatrix * | , | |||
double | , | |||
double * | , | |||
int * | ||||
) |
void dgsequ | ( | SuperMatrix * | A, | |
double * | r, | |||
double * | c, | |||
double * | rowcnd, | |||
double * | colcnd, | |||
double * | amax, | |||
int * | info | |||
) |
Driver related.
Purpose =======
DGSEQU computes row and column scalings intended to equilibrate an M-by-N sparse matrix A and reduce its condition number. R returns the row scale factors and C the column scale factors, chosen to try to make the largest element in each row and column of the matrix B with elements B(i,j)=R(i)*A(i,j)*C(j) have absolute value 1.
R(i) and C(j) are restricted to be between SMLNUM = smallest safe number and BIGNUM = largest safe number. Use of these scaling factors is not guaranteed to reduce the condition number of A but works well in practice.
See supermatrix.h for the definition of 'SuperMatrix' structure.
Arguments =========
A (input) SuperMatrix* The matrix of dimension (A->nrow, A->ncol) whose equilibration factors are to be computed. The type of A can be: Stype = SLU_NC; Dtype = SLU_D; Mtype = SLU_GE.
R (output) double*, size A->nrow If INFO = 0 or INFO > M, R contains the row scale factors for A.
C (output) double*, size A->ncol If INFO = 0, C contains the column scale factors for A.
ROWCND (output) double* If INFO = 0 or INFO > M, ROWCND contains the ratio of the smallest R(i) to the largest R(i). If ROWCND >= 0.1 and AMAX is neither too large nor too small, it is not worth scaling by R.
COLCND (output) double* If INFO = 0, COLCND contains the ratio of the smallest C(i) to the largest C(i). If COLCND >= 0.1, it is not worth scaling by C.
AMAX (output) double* Absolute value of largest matrix element. If AMAX is very close to overflow or very close to underflow, the matrix should be scaled.
INFO (output) int* = 0: successful exit < 0: if INFO = -i, the i-th argument had an illegal value > 0: if INFO = i, and i is <= A->nrow: the i-th row of A is exactly zero > A->ncol: the (i-M)-th column of A is exactly zero
=====================================================================
void dgsrfs | ( | trans_t | , | |
SuperMatrix * | , | |||
SuperMatrix * | , | |||
SuperMatrix * | , | |||
int * | , | |||
int * | , | |||
equed_t | , | |||
double * | , | |||
double * | , | |||
SuperMatrix * | , | |||
SuperMatrix * | , | |||
double * | , | |||
double * | , | |||
Gstat_t * | , | |||
int * | ||||
) |
void dgstrs | ( | trans_t | , | |
SuperMatrix * | , | |||
SuperMatrix * | , | |||
int * | , | |||
int * | , | |||
SuperMatrix * | , | |||
Gstat_t * | , | |||
int * | ||||
) |
void dinf_norm_error | ( | int | , | |
SuperMatrix * | , | |||
double * | ||||
) |
Check the inf-norm of the error vector.
double dlamch_ | ( | ) |
void dlaqgs | ( | SuperMatrix * | , | |
double * | , | |||
double * | , | |||
double | , | |||
double | , | |||
double | , | |||
equed_t * | ||||
) |
void dlsolve | ( | int | ldm, | |
int | ncol, | |||
double * | M, | |||
double * | rhs | |||
) |
Solves a dense UNIT lower triangular system.
The unit lower triangular matrix is stored in a 2D array M(1:nrow,1:ncol). The solution will be returned in the rhs vector.
void dmatvec | ( | int | ldm, | |
int | nrow, | |||
int | ncol, | |||
double * | M, | |||
double * | vec, | |||
double * | Mxvec | |||
) |
Performs a dense matrix-vector multiply: Mxvec = Mxvec + M * vec.
The input matrix is M(1:nrow,1:ncol); The product is returned in Mxvec[].
double* doubleCalloc | ( | int | ) |
double* doubleMalloc | ( | int | ) |
int dParallelInit | ( | int | , | |
pxgstrf_relax_t * | , | |||
superlumt_options_t * | , | |||
pxgstrf_shared_t * | ||||
) |
double dPivotGrowth | ( | int | ncols, | |
SuperMatrix * | A, | |||
int * | perm_c, | |||
SuperMatrix * | L, | |||
SuperMatrix * | U | |||
) |
Purpose =======
Compute the reciprocal pivot growth factor of the leading ncols columns of the matrix, using the formula: min_j ( max_i(abs(A_ij)) / max_i(abs(U_ij)) )
Arguments =========
ncols (input) int The number of columns of matrices A, L and U.
A (input) SuperMatrix* Original matrix A, permuted by columns, of dimension (A->nrow, A->ncol). The type of A can be: Stype = NC; Dtype = SLU_D; Mtype = GE.
L (output) SuperMatrix* The factor L from the factorization Pr*A=L*U; use compressed row subscripts storage for supernodes, i.e., L has type: Stype = SC; Dtype = SLU_D; Mtype = TRLU.
U (output) SuperMatrix* The factor U from the factorization Pr*A*Pc=L*U. Use column-wise storage scheme, i.e., U has types: Stype = NC; Dtype = SLU_D; Mtype = TRU.
int dPresetMap | ( | const | int, | |
SuperMatrix * | , | |||
pxgstrf_relax_t * | , | |||
superlumt_options_t * | , | |||
GlobalLU_t * | ||||
) |
void dPrintPerf | ( | SuperMatrix * | , | |
SuperMatrix * | , | |||
superlu_memusage_t * | , | |||
double | , | |||
double | , | |||
double * | , | |||
double * | , | |||
char * | , | |||
Gstat_t * | ||||
) |
void dreadhb | ( | int * | , | |
int * | , | |||
int * | , | |||
double ** | , | |||
int ** | , | |||
int ** | ||||
) |
void dreadmt | ( | int * | , | |
int * | , | |||
int * | , | |||
double ** | , | |||
int ** | , | |||
int ** | ||||
) |
void dsp_colorder | ( | SuperMatrix * | , | |
int * | , | |||
superlumt_options_t * | , | |||
SuperMatrix * | ||||
) |
void dstat_allocate | ( | int | ) |
int dtrsm_ | ( | char * | , | |
char * | , | |||
char * | , | |||
char * | , | |||
int * | , | |||
int * | , | |||
double * | , | |||
double * | , | |||
int * | , | |||
double * | , | |||
int * | ||||
) |
int dtrsv_ | ( | char * | , | |
char * | , | |||
char * | , | |||
int * | , | |||
double * | , | |||
int * | , | |||
double * | , | |||
int * | ||||
) |
void dusolve | ( | int | ldm, | |
int | ncol, | |||
double * | M, | |||
double * | rhs | |||
) |
Solves a dense upper triangular system.
The upper triangular matrix is stored in a 2-dim array M(1:ldm,1:ncol). The solution will be returned in the rhs vector.
int DynamicSetMap | ( | const | int, | |
const | int, | |||
const | int, | |||
pxgstrf_shared_t * | ||||
) |
int EnqueueDomains | ( | queue_t * | , | |
struct Branch * | , | |||
pxgstrf_shared_t * | ||||
) |
int EnqueueRelaxSnode | ( | queue_t * | , | |
int | , | |||
pxgstrf_relax_t * | , | |||
pxgstrf_shared_t * | ||||
) |
void fixupL | ( | const | int, | |
const int * | , | |||
GlobalLU_t * | ||||
) |
Fix up the data storage lsub for L-subscripts. It removes the subscript sets for structural pruning, and applies permuation to the remaining subscripts.
void get_perm_c | ( | int | ispec, | |
SuperMatrix * | A, | |||
int * | perm_c | |||
) |
Purpose =======
GET_PERM_C obtains a permutation matrix Pc, by applying the multiple minimum degree ordering code by Joseph Liu to matrix A'*A or A+A'. or using approximate minimum degree column ordering by Davis et. al. The LU factorization of A*Pc tends to have less fill than the LU factorization of A.
Arguments =========
ispec (input) int Specifies the type of column ordering to reduce fill: = 1: minimum degree on the structure of A^T * A = 2: minimum degree on the structure of A^T + A = 3: approximate minimum degree for unsymmetric matrices If ispec == 0, the natural ordering (i.e., Pc = I) is returned.
A (input) SuperMatrix* Matrix A in A*X=B, of dimension (A->nrow, A->ncol). The number of the linear equations is A->nrow. Currently, the type of A can be: Stype = NC; Dtype = _D; Mtype = GE. In the future, more general A can be handled.
perm_c (output) int* Column permutation vector of size A->ncol, which defines the permutation matrix Pc; perm_c[i] = j means column i of A is in position j in A*Pc.
int Glu_alloc | ( | const | int, | |
const | int, | |||
const | int, | |||
const | MemType, | |||
int * | , | |||
pxgstrf_shared_t * | ||||
) |
void ifill | ( | int * | , | |
int | , | |||
int | ||||
) |
Fills an integer array with a given value.
int* intCalloc | ( | int | ) |
int* intMalloc | ( | int | ) |
int lockon | ( | int * | ) |
int lsame_ | ( | char * | ca, | |
char * | cb | |||
) |
Purpose =======
LSAME returns .TRUE. if CA is the same letter as CB regardless of case.
Arguments =========
CA (input) CHARACTER*1 CB (input) CHARACTER*1 CA and CB specify the single characters to be compared.
=====================================================================
int memory_usage | ( | ) |
int NewNsuper | ( | const | int, | |
pxgstrf_shared_t * | , | |||
int * | ||||
) |
int ParallelFinalize | ( | ) |
void PartDomains | ( | const | int, | |
const | float, | |||
SuperMatrix * | , | |||
int * | , | |||
int * | ||||
) |
void pdgssv | ( | int | , | |
SuperMatrix * | , | |||
int * | , | |||
int * | , | |||
SuperMatrix * | , | |||
SuperMatrix * | , | |||
SuperMatrix * | , | |||
int * | ||||
) |
void pdgssvx | ( | int | , | |
superlumt_options_t * | , | |||
SuperMatrix * | , | |||
int * | , | |||
int * | , | |||
equed_t * | , | |||
double * | , | |||
double * | , | |||
SuperMatrix * | , | |||
SuperMatrix * | , | |||
SuperMatrix * | , | |||
SuperMatrix * | , | |||
double * | , | |||
double * | , | |||
double * | , | |||
double * | , | |||
superlu_memusage_t * | , | |||
int * | ||||
) |
void pdgstrf | ( | superlumt_options_t * | , | |
SuperMatrix * | , | |||
int * | , | |||
SuperMatrix * | , | |||
SuperMatrix * | , | |||
Gstat_t * | , | |||
int * | ||||
) |
void pdgstrf_bmod1D | ( | const | int, | |
const | int, | |||
const | int, | |||
const | int, | |||
const | int, | |||
const | int, | |||
const | int, | |||
int | , | |||
int | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
double * | , | |||
double * | , | |||
GlobalLU_t * | , | |||
Gstat_t * | ||||
) |
void pdgstrf_bmod1D_mv2 | ( | const | int, | |
const | int, | |||
const | int, | |||
const | int, | |||
const | int, | |||
const | int, | |||
const | int, | |||
int | , | |||
int | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
double * | , | |||
double * | , | |||
GlobalLU_t * | , | |||
Gstat_t * | ||||
) |
void pdgstrf_bmod2D | ( | const | int, | |
const | int, | |||
const | int, | |||
const | int, | |||
const | int, | |||
const | int, | |||
const | int, | |||
int | , | |||
int | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
double * | , | |||
double * | , | |||
GlobalLU_t * | , | |||
Gstat_t * | ||||
) |
void pdgstrf_bmod2D_mv2 | ( | const | int, | |
const | int, | |||
const | int, | |||
const | int, | |||
const | int, | |||
const | int, | |||
const | int, | |||
int | , | |||
int | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
double * | , | |||
double * | , | |||
GlobalLU_t * | , | |||
Gstat_t * | ||||
) |
int pdgstrf_column_bmod | ( | const | int, | |
const | int, | |||
const | int, | |||
const | int, | |||
int * | , | |||
int * | , | |||
double * | , | |||
double * | , | |||
pxgstrf_shared_t * | , | |||
Gstat_t * | ||||
) |
int pdgstrf_column_dfs | ( | const | int, | |
const | int, | |||
const | int, | |||
const | int, | |||
int * | , | |||
int * | , | |||
int * | , | |||
int | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
pxgstrf_shared_t * | ||||
) |
int pdgstrf_copy_to_ucol | ( | const | int, | |
const | int, | |||
const | int, | |||
const int * | , | |||
const int * | , | |||
const int * | , | |||
double * | , | |||
pxgstrf_shared_t * | ||||
) |
int pdgstrf_factor_snode | ( | const | int, | |
const | int, | |||
SuperMatrix * | , | |||
const | double, | |||
yes_no_t * | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
double * | , | |||
double * | , | |||
pxgstrf_shared_t * | , | |||
int * | ||||
) |
void pdgstrf_finalize | ( | superlumt_options_t * | , | |
SuperMatrix * | ||||
) |
void pdgstrf_init | ( | int | , | |
fact_t | , | |||
trans_t | , | |||
yes_no_t | , | |||
int | , | |||
int | , | |||
double | , | |||
yes_no_t | , | |||
double | , | |||
int * | , | |||
int * | , | |||
void * | , | |||
int | , | |||
SuperMatrix * | , | |||
SuperMatrix * | , | |||
superlumt_options_t * | , | |||
Gstat_t * | ||||
) |
float pdgstrf_MemInit | ( | int | , | |
int | , | |||
superlumt_options_t * | , | |||
SuperMatrix * | , | |||
SuperMatrix * | , | |||
GlobalLU_t * | ||||
) |
float pdgstrf_memory_use | ( | const | int, | |
const | int, | |||
const | int | |||
) |
int pdgstrf_MemXpand | ( | int | , | |
int | , | |||
MemType | , | |||
int * | , | |||
GlobalLU_t * | ||||
) |
void pdgstrf_panel_bmod | ( | const | int, | |
const | int, | |||
const | int, | |||
const | int, | |||
const | int, | |||
int * | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
double * | , | |||
double * | , | |||
pxgstrf_shared_t * | ||||
) |
void pdgstrf_panel_dfs | ( | const | int, | |
const | int, | |||
const | int, | |||
const | int, | |||
SuperMatrix * | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
double * | , | |||
GlobalLU_t * | ||||
) |
int pdgstrf_pivotL | ( | const | int, | |
const | int, | |||
const | double, | |||
yes_no_t * | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
GlobalLU_t * | , | |||
Gstat_t * | ||||
) |
void pdgstrf_relax_snode | ( | const | int, | |
superlumt_options_t * | , | |||
pxgstrf_relax_t * | ||||
) |
void pdgstrf_SetRWork | ( | int | , | |
int | , | |||
double * | , | |||
double ** | , | |||
double ** | ||||
) |
int pdgstrf_snode_bmod | ( | const | int, | |
const | int, | |||
const | int, | |||
const | int, | |||
double * | , | |||
double * | , | |||
GlobalLU_t * | , | |||
Gstat_t * | ||||
) |
int pdgstrf_snode_dfs | ( | const | int, | |
const | int, | |||
const | int, | |||
const int * | , | |||
const int * | , | |||
const int * | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
pxgstrf_shared_t * | ||||
) |
void pdgstrf_StackFree | ( | ) |
void pdgstrf_thread_finalize | ( | pdgstrf_threadarg_t * | , | |
pxgstrf_shared_t * | , | |||
SuperMatrix * | , | |||
int * | , | |||
SuperMatrix * | , | |||
SuperMatrix * | ||||
) |
pdgstrf_threadarg_t* pdgstrf_thread_init | ( | SuperMatrix * | , | |
SuperMatrix * | , | |||
SuperMatrix * | , | |||
superlumt_options_t * | , | |||
pxgstrf_shared_t * | , | |||
Gstat_t * | , | |||
int * | ||||
) |
void pdgstrf_WorkFree | ( | int * | , | |
double * | , | |||
GlobalLU_t * | ||||
) |
int pdgstrf_WorkInit | ( | int | , | |
int | , | |||
int ** | , | |||
double ** | ||||
) |
void print_lu_col | ( | int | , | |
char * | , | |||
int | , | |||
int | , | |||
int | , | |||
int * | , | |||
GlobalLU_t * | ||||
) |
void print_panel_seg | ( | int | , | |
int | , | |||
int | , | |||
int | , | |||
int * | , | |||
int * | ||||
) |
Diagnostic print of segment info after panel_dfs().
void PrintSumm | ( | char * | , | |
int | , | |||
int | , | |||
int | ||||
) |
Print a summary of the testing results.
void pxgstrf_finalize | ( | superlumt_options_t * | , | |
SuperMatrix * | ||||
) |
void pxgstrf_mark_busy_descends | ( | int | , | |
int | , | |||
int * | , | |||
pxgstrf_shared_t * | , | |||
int * | , | |||
int * | ||||
) |
void pxgstrf_pruneL | ( | const | int, | |
const int * | , | |||
const | int, | |||
const | int, | |||
const int * | , | |||
const int * | , | |||
int * | , | |||
int * | , | |||
GlobalLU_t * | ||||
) |
void pxgstrf_resetrep_col | ( | const | int, | |
const int * | , | |||
int * | ||||
) |
void pxgstrf_scheduler | ( | const | int, | |
const | int, | |||
const int * | , | |||
int * | , | |||
int * | , | |||
pxgstrf_shared_t * | ||||
) |
void pxgstrf_SetIWork | ( | int | , | |
int | , | |||
int * | , | |||
int ** | , | |||
int ** | , | |||
int ** | , | |||
int ** | , | |||
int ** | , | |||
int ** | , | |||
int ** | ||||
) |
void pxgstrf_super_bnd_dfs | ( | const | int, | |
const | int, | |||
const | int, | |||
const | int, | |||
const | int, | |||
SuperMatrix * | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
pxgstrf_shared_t * | ||||
) |
int qrnzcnt | ( | int | , | |
int | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
int * | , | |||
int * | ||||
) |
int queue_destroy | ( | queue_t * | ) |
int queue_init | ( | queue_t * | , | |
int | ||||
) |
void snode_profile | ( | int | , | |
int * | ||||
) |
int sp_coletree | ( | int * | , | |
int * | , | |||
int * | , | |||
int | , | |||
int | , | |||
int * | ||||
) |
int sp_dgemm | ( | char * | , | |
int | , | |||
int | , | |||
int | , | |||
double | , | |||
SuperMatrix * | , | |||
double * | , | |||
int | , | |||
double | , | |||
double * | , | |||
int | ||||
) |
int sp_dgemv | ( | char * | trans, | |
double | alpha, | |||
SuperMatrix * | A, | |||
double * | x, | |||
int | incx, | |||
double | beta, | |||
double * | y, | |||
int | incy | |||
) |
Performs one of the matrix-vector operations y := alpha*A*x + beta*y, or y := alpha*A'*x + beta*y,.
Purpose =======
sp_dgemv() performs one of the matrix-vector operations y := alpha*A*x + beta*y, or y := alpha*A'*x + beta*y, where alpha and beta are scalars, x and y are vectors and A is a sparse A->nrow by A->ncol matrix.
Parameters ==========
TRANS - (input) char* On entry, TRANS specifies the operation to be performed as follows: TRANS = 'N' or 'n' y := alpha*A*x + beta*y. TRANS = 'T' or 't' y := alpha*A'*x + beta*y. TRANS = 'C' or 'c' y := alpha*A'*x + beta*y.
ALPHA - (input) double On entry, ALPHA specifies the scalar alpha.
A - (input) SuperMatrix* Matrix A with a sparse format, of dimension (A->nrow, A->ncol). Currently, the type of A can be: Stype = NC or NCP; Dtype = SLU_D; Mtype = GE. In the future, more general A can be handled.
X - (input) double*, array of DIMENSION at least ( 1 + ( n - 1 )*abs( INCX ) ) when TRANS = 'N' or 'n' and at least ( 1 + ( m - 1 )*abs( INCX ) ) otherwise. Before entry, the incremented array X must contain the vector x.
INCX - (input) int On entry, INCX specifies the increment for the elements of X. INCX must not be zero.
BETA - (input) double On entry, BETA specifies the scalar beta. When BETA is supplied as zero then Y need not be set on input.
Y - (output) double*, array of DIMENSION at least ( 1 + ( m - 1 )*abs( INCY ) ) when TRANS = 'N' or 'n' and at least ( 1 + ( n - 1 )*abs( INCY ) ) otherwise. Before entry with BETA non-zero, the incremented array Y must contain the vector y. On exit, Y is overwritten by the updated vector y.
INCY - (input) int On entry, INCY specifies the increment for the elements of Y. INCY must not be zero.
==== Sparse Level 2 Blas routine.
int sp_dtrsv | ( | char * | , | |
char * | , | |||
char * | , | |||
SuperMatrix * | , | |||
SuperMatrix * | , | |||
double * | , | |||
int * | ||||
) |
int sp_ienv | ( | int | ispec | ) |
Purpose =======
sp_ienv() is inquired to choose machine-dependent parameters for the local environment. See ISPEC for a description of the parameters.
This version provides a set of parameters which should give good, but not optimal, performance on many of the currently available computers. Users are encouraged to modify this subroutine to set the tuning parameters for their particular machine using the option and problem size information in the arguments.
Arguments =========
ISPEC (input) int Specifies the parameter to be returned as the value of SP_IENV. = 1: the panel size w; a panel consists of w consecutive columns of matrix A in the process of Gaussian elimination. The best value depends on machine's cache characters. = 2: the relaxation parameter relax; if the number of nodes (columns) in a subtree of the elimination tree is less than relax, this subtree is considered as one supernode, regardless of their row structures. = 3: the maximum size for a supernode in complete LU; = 4: the minimum row dimension for 2-D blocking to be used; = 5: the minimum column dimension for 2-D blocking to be used; = 6: the estimated fills factor for L and U, compared with A; = 7: the maximum size for a supernode in ILU.
(SP_IENV) (output) int >= 0: the value of the parameter specified by ISPEC < 0: if SP_IENV = -k, the k-th argument had an illegal value.
=====================================================================
int spcoletree | ( | int * | , | |
int * | , | |||
int * | , | |||
int | , | |||
int | , | |||
int * | ||||
) |
void super_stats | ( | int | , | |
int * | , | |||
int * | ||||
) |
void superlu_abort_and_exit | ( | char * | ) |
Global statistics variale.
int superlu_dQuerySpace | ( | int | , | |
SuperMatrix * | , | |||
SuperMatrix * | , | |||
int | , | |||
superlu_memusage_t * | ||||
) |
double SuperLU_timer_ | ( | ) |
Timer function.
int* TreePostorder | ( | int | , | |
int * | ||||
) |
int xerbla_ | ( | char * | , | |
int * | ||||
) |