SuperLU 6.0.1
slu_ddefs.h
Go to the documentation of this file.
1
73#ifndef __SUPERLU_dSP_DEFS /* allow multiple inclusions */
74#define __SUPERLU_dSP_DEFS
75
76/*
77 * File name: dsp_defs.h
78 * Purpose: Sparse matrix types and function prototypes
79 * History:
80 */
81
82#ifdef _CRAY
83#include <fortran.h>
84#endif
85
86#include <math.h>
87#include <limits.h>
88#include <stdio.h>
89#include <stdlib.h>
90#include <stdint.h>
91#include <string.h>
92#include "slu_Cnames.h"
93#include "superlu_config.h"
94#include "supermatrix.h"
95#include "slu_util.h"
96
97
98/* -------- Prototypes -------- */
99
100#ifdef __cplusplus
101extern "C" {
102#endif
103
105extern void
106dgssv(superlu_options_t *, SuperMatrix *, int *, int *, SuperMatrix *,
108extern void
109dgssvx(superlu_options_t *, SuperMatrix *, int *, int *, int *,
110 char *, double *, double *, SuperMatrix *, SuperMatrix *,
111 void *, int_t lwork, SuperMatrix *, SuperMatrix *,
112 double *, double *, double *, double *,
114 /* ILU */
115extern void
117 SuperMatrix *, SuperMatrix *, SuperLUStat_t *, int *);
118extern void
119dgsisx(superlu_options_t *options, SuperMatrix *A, int *perm_c, int *perm_r,
120 int *etree, char *equed, double *R, double *C,
121 SuperMatrix *L, SuperMatrix *U, void *work, int_t lwork,
122 SuperMatrix *B, SuperMatrix *X, double *recip_pivot_growth, double *rcond,
123 GlobalLU_t *Glu, mem_usage_t *mem_usage, SuperLUStat_t *stat, int_t *info);
124
125
127extern void
128dCreate_CompCol_Matrix(SuperMatrix *, int, int, int_t, double *,
130extern void
131dCreate_CompRow_Matrix(SuperMatrix *, int, int, int_t, double *,
133extern void dCompRow_to_CompCol(int, int, int_t, double*, int_t*, int_t*,
134 double **, int_t **, int_t **);
135extern void
137extern void
138dCreate_Dense_Matrix(SuperMatrix *, int, int, double *, int,
140extern void
141dCreate_SuperNode_Matrix(SuperMatrix *, int, int, int_t, double *,
142 int_t *, int_t *, int_t *, int *, int *,
144extern void
145dCopy_Dense_Matrix(int, int, double *, int, double *, int);
146
147extern void dallocateA (int, int_t, double **, int_t **, int_t **);
148extern void dgstrf (superlu_options_t*, SuperMatrix*,
149 int, int, int*, void *, int_t, int *, int *,
151 SuperLUStat_t*, int_t *info);
152extern int_t dsnode_dfs (const int, const int, const int_t *, const int_t *,
153 const int_t *, int_t *, int *, GlobalLU_t *);
154extern int dsnode_bmod (const int, const int, const int, double *,
155 double *, GlobalLU_t *, SuperLUStat_t*);
156extern void dpanel_dfs (const int, const int, const int, SuperMatrix *,
157 int *, int *, double *, int *, int *, int *,
158 int_t *, int *, int *, int_t *, GlobalLU_t *);
159extern void dpanel_bmod (const int, const int, const int, const int,
160 double *, double *, int *, int *,
162extern int dcolumn_dfs (const int, const int, int *, int *, int *, int *,
163 int *, int_t *, int *, int *, int_t *, GlobalLU_t *);
164extern int dcolumn_bmod (const int, const int, double *,
165 double *, int *, int *, int,
167extern int dcopy_to_ucol (int, int, int *, int *, int *,
168 double *, GlobalLU_t *);
169extern int dpivotL (const int, const double, int *, int *,
170 int *, int *, int *, GlobalLU_t *, SuperLUStat_t*);
171extern void dpruneL (const int, const int *, const int, const int,
172 const int *, const int *, int_t *, GlobalLU_t *);
173extern void dreadmt (int *, int *, int_t *, double **, int_t **, int_t **);
174extern void dGenXtrue (int, int, double *, int);
175extern void dFillRHS (trans_t, int, double *, int, SuperMatrix *,
176 SuperMatrix *);
177extern void dgstrs (trans_t, SuperMatrix *, SuperMatrix *, int *, int *,
178 SuperMatrix *, SuperLUStat_t*, int *);
179/* ILU */
180extern void dgsitrf (superlu_options_t*, SuperMatrix*, int, int, int*,
181 void *, int_t, int *, int *, SuperMatrix *, SuperMatrix *,
182 GlobalLU_t *, SuperLUStat_t*, int_t *info);
183extern int dldperm(int, int, int_t, int_t [], int_t [], double [],
184 int [], double [], double []);
185extern int ilu_dsnode_dfs (const int, const int, const int_t *, const int_t *,
186 const int_t *, int *, GlobalLU_t *);
187extern void ilu_dpanel_dfs (const int, const int, const int, SuperMatrix *,
188 int *, int *, double *, double *, int *, int *,
189 int *, int *, int *, int_t *, GlobalLU_t *);
190extern int ilu_dcolumn_dfs (const int, const int, int *, int *, int *,
191 int *, int *, int *, int *, int_t *, GlobalLU_t *);
192extern int ilu_dcopy_to_ucol (int, int, int *, int *, int *,
193 double *, int, milu_t, double, int,
194 double *, int *, GlobalLU_t *, double *);
195extern int ilu_dpivotL (const int, const double, int *, int *, int, int *,
196 int *, int *, int *, double, milu_t,
197 double, GlobalLU_t *, SuperLUStat_t*);
198extern int ilu_ddrop_row (superlu_options_t *, int, int, double,
199 int, int *, double *, GlobalLU_t *,
200 double *, double *, int);
201
202
205extern void dgsequ (SuperMatrix *, double *, double *, double *,
206 double *, double *, int *);
207extern void dlaqgs (SuperMatrix *, double *, double *, double,
208 double, double, char *);
209extern void dgscon (char *, SuperMatrix *, SuperMatrix *,
210 double, double *, SuperLUStat_t*, int *);
211extern double dPivotGrowth(int, SuperMatrix *, int *,
213extern void dgsrfs (trans_t, SuperMatrix *, SuperMatrix *,
214 SuperMatrix *, int *, int *, char *, double *,
215 double *, SuperMatrix *, SuperMatrix *,
216 double *, double *, SuperLUStat_t*, int *);
217
218extern int sp_dtrsv (char *, char *, char *, SuperMatrix *,
219 SuperMatrix *, double *, SuperLUStat_t*, int *);
220extern int sp_dgemv (char *, double, SuperMatrix *, double *,
221 int, double, double *, int);
222
223extern int sp_dgemm (char *, char *, int, int, int, double,
224 SuperMatrix *, double *, int, double,
225 double *, int);
226extern double dmach(char *); /* from C99 standard, in float.h */
227
229extern int_t dLUMemInit (fact_t, void *, int_t, int, int, int_t, int,
230 double, SuperMatrix *, SuperMatrix *,
231 GlobalLU_t *, int **, double **);
232extern void dSetRWork (int, int, double *, double **, double **);
233extern void dLUWorkFree (int *, double *, GlobalLU_t *);
234extern int_t dLUMemXpand (int, int_t, MemType, int_t *, GlobalLU_t *);
235
236extern double *doubleMalloc(size_t);
237extern double *doubleCalloc(size_t);
238extern int_t dmemory_usage(const int_t, const int_t, const int_t, const int);
241
243extern void dreadhb(FILE *, int *, int *, int_t *, double **, int_t **, int_t **);
244extern void dreadrb(int *, int *, int_t *, double **, int_t **, int_t **);
245extern void dreadtriple(int *, int *, int_t *, double **, int_t **, int_t **);
246extern void dreadtriple_noheader(int *, int *, int_t *, double **, int_t **, int_t **);
247extern void dreadMM(FILE *, int *, int *, int_t *, double **, int_t **, int_t **);
248extern void dfill (double *, int, double);
249extern void dinf_norm_error (int, SuperMatrix *, double *);
250extern double dqselect(int, double *, int);
251
252
254extern void dPrint_CompCol_Matrix(char *, SuperMatrix *);
255extern void dPrint_SuperNode_Matrix(char *, SuperMatrix *);
256extern void dPrint_Dense_Matrix(char *, SuperMatrix *);
257extern void dprint_lu_col(char *, int, int, int_t *, GlobalLU_t *);
258extern int print_double_vec(char *, int, double *);
259extern void dcheck_tempv(int, double *);
260
263extern int dgemm_(const char*, const char*, const int*, const int*, const int*,
264 const double*, const double*, const int*, const double*,
265 const int*, const double*, double*, const int*);
266extern int dtrsv_(char*, char*, char*, int*, double*, int*,
267 double*, int*);
268extern int dtrsm_(char*, char*, char*, char*, int*, int*,
269 double*, double*, int*, double*, int*);
270extern int dgemv_(char *, int *, int *, double *, double *a, int *,
271 double *, int *, double *, double *, int *);
272
273extern void dusolve(int, int, double*, double*);
274extern void dlsolve(int, int, double*, double*);
275extern void dmatvec(int, int, int, double*, double*, double*);
276
277#ifdef __cplusplus
278 }
279#endif
280
281#endif /* __SUPERLU_dSP_DEFS */
282
#define X(I)
#define A(I, J)
#define U(I)
Macros defining how C routines will be called.
void dgsisv(superlu_options_t *, SuperMatrix *, int *, int *, SuperMatrix *, SuperMatrix *, SuperMatrix *, SuperLUStat_t *, int *)
void ilu_dpanel_dfs(const int, const int, const int, SuperMatrix *, int *, int *, double *, double *, int *, int *, int *, int *, int *, int_t *, GlobalLU_t *)
Definition: ilu_dpanel_dfs.c:56
double * doubleMalloc(size_t)
Definition: dmemory.c:679
void dCreate_Dense_Matrix(SuperMatrix *, int, int, double *, int, Stype_t, Dtype_t, Mtype_t)
Definition: dutil.c:103
int dpivotL(const int, const double, int *, int *, int *, int *, int *, GlobalLU_t *, SuperLUStat_t *)
Definition: dpivotL.c:66
void dreadtriple_noheader(int *, int *, int_t *, double **, int_t **, int_t **)
Read matrix in triplet format from stdin.
Definition: dreadtriple_noheader.c:38
int ilu_dpivotL(const int, const double, int *, int *, int, int *, int *, int *, int *, double, milu_t, double, GlobalLU_t *, SuperLUStat_t *)
Definition: ilu_dpivotL.c:56
void dgscon(char *, SuperMatrix *, SuperMatrix *, double, double *, SuperLUStat_t *, int *)
Definition: dgscon.c:84
void dCopy_CompCol_Matrix(SuperMatrix *, SuperMatrix *)
Copy matrix A into matrix B.
Definition: dutil.c:82
void dpanel_bmod(const int, const int, const int, const int, double *, double *, int *, int *, GlobalLU_t *, SuperLUStat_t *)
Definition: dpanel_bmod.c:61
void dreadtriple(int *, int *, int_t *, double **, int_t **, int_t **)
Definition: dreadtriple.c:26
void dCreate_CompCol_Matrix(SuperMatrix *, int, int, int_t, double *, int_t *, int_t *, Stype_t, Dtype_t, Mtype_t)
Supernodal LU factor related.
Definition: dutil.c:39
int dsnode_bmod(const int, const int, const int, double *, double *, GlobalLU_t *, SuperLUStat_t *)
Performs numeric block updates within the relaxed snode.
Definition: dsnode_bmod.c:41
double dqselect(int, double *, int)
double dPivotGrowth(int, SuperMatrix *, int *, SuperMatrix *, SuperMatrix *)
Definition: dpivotgrowth.c:59
void dpanel_dfs(const int, const int, const int, SuperMatrix *, int *, int *, double *, int *, int *, int *, int_t *, int *, int *, int_t *, GlobalLU_t *)
Definition: dpanel_dfs.c:69
void dreadrb(int *, int *, int_t *, double **, int_t **, int_t **)
Definition: dreadrb.c:281
int dcolumn_bmod(const int, const int, double *, double *, int *, int *, int, GlobalLU_t *, SuperLUStat_t *)
Definition: dcolumn_bmod.c:52
void dreadhb(FILE *, int *, int *, int_t *, double **, int_t **, int_t **)
Auxiliary routines.
Definition: dreadhb.c:283
int_t dsnode_dfs(const int, const int, const int_t *, const int_t *, const int_t *, int_t *, int *, GlobalLU_t *)
Definition: dsnode_dfs.c:55
void dgsrfs(trans_t, SuperMatrix *, SuperMatrix *, SuperMatrix *, int *, int *, char *, double *, double *, SuperMatrix *, SuperMatrix *, double *, double *, SuperLUStat_t *, int *)
Definition: dgsrfs.c:141
void dgssvx(superlu_options_t *, SuperMatrix *, int *, int *, int *, char *, double *, double *, SuperMatrix *, SuperMatrix *, void *, int_t lwork, SuperMatrix *, SuperMatrix *, double *, double *, double *, double *, GlobalLU_t *, mem_usage_t *, SuperLUStat_t *, int_t *info)
Definition: dgssvx.c:357
int dldperm(int, int, int_t, int_t[], int_t[], double[], int[], double[], double[])
void dPrint_CompCol_Matrix(char *, SuperMatrix *)
Routines for debugging.
Definition: dutil.c:202
void dallocateA(int, int_t, double **, int_t **, int_t **)
Allocate storage for original matrix A.
Definition: dmemory.c:671
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,
Definition: dsp_blas2.c:374
double dmach(char *)
Definition: dmach.c:18
void dlsolve(int, int, double *, double *)
Solves a dense UNIT lower triangular system.
Definition: dmyblas2.c:38
void dgssv(superlu_options_t *, SuperMatrix *, int *, int *, SuperMatrix *, SuperMatrix *, SuperMatrix *, SuperLUStat_t *, int_t *info)
Driver routines.
Definition: dgssv.c:143
void dgstrf(superlu_options_t *, SuperMatrix *, int, int, int *, void *, int_t, int *, int *, SuperMatrix *, SuperMatrix *, GlobalLU_t *, SuperLUStat_t *, int_t *info)
Definition: dgstrf.c:200
void dSetRWork(int, int, double *, double **, double **)
Set up pointers for real working arrays.
Definition: dmemory.c:391
void dmatvec(int, int, int, double *, double *, double *)
Performs a dense matrix-vector multiply: Mxvec = Mxvec + M * vec.
Definition: dmyblas2.c:161
int ilu_dcopy_to_ucol(int, int, int *, int *, int *, double *, int, milu_t, double, int, double *, int *, GlobalLU_t *, double *)
Definition: ilu_dcopy_to_ucol.c:44
int dcolumn_dfs(const int, const int, int *, int *, int *, int *, int *, int_t *, int *, int *, int_t *, GlobalLU_t *)
Definition: dcolumn_dfs.c:76
int ilu_dQuerySpace(SuperMatrix *, SuperMatrix *, mem_usage_t *)
Definition: dmemory.c:147
void dlaqgs(SuperMatrix *, double *, double *, double, double, double, char *)
Definition: dlaqgs.c:92
void dreadMM(FILE *, int *, int *, int_t *, double **, int_t **, int_t **)
Definition: dreadMM.c:35
void dgsitrf(superlu_options_t *, SuperMatrix *, int, int, int *, void *, int_t, int *, int *, SuperMatrix *, SuperMatrix *, GlobalLU_t *, SuperLUStat_t *, int_t *info)
Definition: dgsitrf.c:187
int ilu_dsnode_dfs(const int, const int, const int_t *, const int_t *, const int_t *, int *, GlobalLU_t *)
Definition: ilu_dsnode_dfs.c:42
int dtrsv_(char *, char *, char *, int *, double *, int *, double *, int *)
void dgsisx(superlu_options_t *options, SuperMatrix *A, int *perm_c, int *perm_r, int *etree, char *equed, double *R, double *C, SuperMatrix *L, SuperMatrix *U, void *work, int_t lwork, SuperMatrix *B, SuperMatrix *X, double *recip_pivot_growth, double *rcond, GlobalLU_t *Glu, mem_usage_t *mem_usage, SuperLUStat_t *stat, int_t *info)
Definition: dgsisx.c:404
void dgstrs(trans_t, SuperMatrix *, SuperMatrix *, int *, int *, SuperMatrix *, SuperLUStat_t *, int *)
Definition: dgstrs.c:93
int dtrsm_(char *, char *, char *, char *, int *, int *, double *, double *, int *, double *, int *)
void dCopy_Dense_Matrix(int, int, double *, int, double *, int)
Definition: dutil.c:121
int dcopy_to_ucol(int, int, int *, int *, int *, double *, GlobalLU_t *)
Definition: dcopy_to_ucol.c:36
int dQuerySpace(SuperMatrix *, SuperMatrix *, mem_usage_t *)
Definition: dmemory.c:110
void dpruneL(const int, const int *, const int, const int, const int *, const int *, int_t *, GlobalLU_t *)
Definition: dpruneL.c:48
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.
Definition: dutil.c:366
int print_double_vec(char *, int, double *)
Definition: dutil.c:473
void dfill(double *, int, double)
Fills a double precision array with a given value.
Definition: dutil.c:393
int sp_dgemm(char *, char *, int, int, int, double, SuperMatrix *, double *, int, double, double *, int)
Definition: dsp_blas3.c:126
int dgemv_(char *, int *, int *, double *, double *a, int *, double *, int *, double *, double *, int *)
int ilu_ddrop_row(superlu_options_t *, int, int, double, int, int *, double *, GlobalLU_t *, double *, double *, int)
void dCompRow_to_CompCol(int, int, int_t, double *, int_t *, int_t *, double **, int_t **, int_t **)
Convert a row compressed storage into a column compressed storage.
Definition: dutil.c:164
void dusolve(int, int, double *, double *)
Solves a dense upper triangular system.
Definition: dmyblas2.c:136
int dgemm_(const char *, const char *, const int *, const int *, const int *, const double *, const double *, const int *, const double *, const int *, const double *, double *, const int *)
BLAS.
double * doubleCalloc(size_t)
Definition: dmemory.c:689
void dcheck_tempv(int, double *)
Check whether tempv[] == 0. This should be true before and after calling any numeric routines,...
Definition: dutil.c:339
void dreadmt(int *, int *, int_t *, double **, int_t **, int_t **)
void dprint_lu_col(char *, int, int, int_t *, GlobalLU_t *)
Diagnostic print of column "jcol" in the U/L factor.
Definition: dutil.c:298
void dPrint_Dense_Matrix(char *, SuperMatrix *)
Definition: dutil.c:276
void dPrint_SuperNode_Matrix(char *, SuperMatrix *)
Definition: dutil.c:226
void dinf_norm_error(int, SuperMatrix *, double *)
Check the inf-norm of the error vector.
Definition: dutil.c:403
int_t dmemory_usage(const int_t, const int_t, const int_t, const int)
Definition: dmemory.c:703
void dLUWorkFree(int *, double *, GlobalLU_t *)
Free the working storage used by factor routines.
Definition: dmemory.c:406
void dCreate_CompRow_Matrix(SuperMatrix *, int, int, int_t, double *, int_t *, int_t *, Stype_t, Dtype_t, Mtype_t)
Definition: dutil.c:60
void dGenXtrue(int, int, double *, int)
Definition: dutil.c:354
int_t dLUMemXpand(int, int_t, MemType, int_t *, GlobalLU_t *)
Expand the data structures for L and U during the factorization.
Definition: dmemory.c:429
int_t dLUMemInit(fact_t, void *, int_t, int, int, int_t, int, double, SuperMatrix *, SuperMatrix *, GlobalLU_t *, int **, double **)
Memory-related.
Definition: dmemory.c:188
int ilu_dcolumn_dfs(const int, const int, int *, int *, int *, int *, int *, int *, int *, int_t *, GlobalLU_t *)
Definition: ilu_dcolumn_dfs.c:61
int sp_dtrsv(char *, char *, char *, SuperMatrix *, SuperMatrix *, double *, SuperLUStat_t *, int *)
Solves one of the systems of equations A*x = b, or A'*x = b.
Definition: dsp_blas2.c:86
void dgsequ(SuperMatrix *, double *, double *, double *, double *, double *, int *)
Driver related.
Definition: dgsequ.c:94
void dCreate_SuperNode_Matrix(SuperMatrix *, int, int, int_t, double *, int_t *, int_t *, int_t *, int *, int *, Stype_t, Dtype_t, Mtype_t)
Definition: dutil.c:134
Utility header file.
Definition: slu_util.h:336
Definition: slu_util.h:321
Definition: supermatrix.h:54
Definition: slu_util.h:330
Definition: slu_util.h:277
int64_t int_t
Definition: superlu_config.h:18
trans_t
Definition: superlu_enum_consts.h:34
milu_t
Definition: superlu_enum_consts.h:46
MemType
Definition: superlu_enum_consts.h:38
fact_t
Definition: superlu_enum_consts.h:30
Matrix type definitions.
Mtype_t
Definition: supermatrix.h:42
Dtype_t
Definition: supermatrix.h:35
Stype_t
Definition: supermatrix.h:22