SuperLU_DIST  4.0
superlu_dist on CPU and GPU clusters
 All Classes Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
superlu_ddefs.h
Go to the documentation of this file.
1 
2 
13 #ifndef __SUPERLU_dDEFS /* allow multiple inclusions */
14 #define __SUPERLU_dDEFS
15 
16 /*
17  * File name: superlu_ddefs.h
18  * Purpose: Distributed SuperLU data types and function prototypes
19  * History:
20  */
21 
22 #include "superlu_defs.h"
23 
24 /*-- Auxiliary data type used in PxGSTRS/PxGSTRS1. */
25 typedef struct {
26  int_t lbnum; /* Row block number (local). */
27  int_t indpos; /* Starting position in Uindex[]. */
28 } Ucb_indptr_t;
29 
30 /*
31  * On each processor, the blocks in L are stored in compressed block
32  * column format, the blocks in U are stored in compressed block row format.
33  */
34 #define MAX_LOOKAHEADS 50
35 typedef struct {
36  int_t **Lrowind_bc_ptr; /* size ceil(NSUPERS/Pc) */
37  double **Lnzval_bc_ptr; /* size ceil(NSUPERS/Pc) */
38  int_t **Ufstnz_br_ptr; /* size ceil(NSUPERS/Pr) */
39  double **Unzval_br_ptr; /* size ceil(NSUPERS/Pr) */
40 #if 0
41  int_t *Lsub_buf; /* Buffer for the remote subscripts of L */
42  double *Lval_buf; /* Buffer for the remote nonzeros of L */
43  int_t *Usub_buf; /* Buffer for the remote subscripts of U */
44  double *Uval_buf; /* Buffer for the remote nonzeros of U */
45 #endif
46  int_t *Lsub_buf_2[MAX_LOOKAHEADS]; /* Buffers for the remote subscripts of L*/
47  double *Lval_buf_2[MAX_LOOKAHEADS]; /* Buffers for the remote nonzeros of L */
48  int_t *Usub_buf_2[MAX_LOOKAHEADS]; /* Buffer for the remote subscripts of U */
49  double *Uval_buf_2[MAX_LOOKAHEADS]; /* Buffer for the remote nonzeros of U */
50  double *ujrow; /* used in panel factorization. */
51  int_t bufmax[NBUFFERS]; /* Buffer size; 5 entries
52  * 0 : size of Lsub_buf[]
53  * 1 : size of Lval_buf[]
54  * 2 : size of Usub_buf[]
55  * 3 : size of Uval_buf[]
56  * 4 : size of tempv[LDA]
57  */
58 
59  /*-- Record communication schedule for factorization. --*/
60  int *ToRecv; /* Recv from no one (0), left (1), and up (2).*/
61  int *ToSendD; /* Whether need to send down block row. */
62  int **ToSendR; /* List of processes to send right block col. */
63 
64  /*-- Record communication schedule for forward/back solves. --*/
65  int_t *fmod; /* Modification count for L-solve */
66  int_t **fsendx_plist; /* Column process list to send down Xk */
67  int_t *frecv; /* Modifications to be recv'd in proc row */
68  int_t nfrecvx; /* Number of Xk I will receive in L-solve */
69  int_t nfsendx; /* Number of Xk I will send in L-solve */
70  int_t *bmod; /* Modification count for U-solve */
71  int_t **bsendx_plist; /* Column process list to send down Xk */
72  int_t *brecv; /* Modifications to be recv'd in proc row */
73  int_t nbrecvx; /* Number of Xk I will receive in U-solve */
74  int_t nbsendx; /* Number of Xk I will send in U-solve */
75  int_t *mod_bit; /* Flag contribution from each row blocks */
76 
77  /*-- Auxiliary arrays used for forward/back solves. --*/
78  int_t *ilsum; /* Starting position of each supernode in lsum
79  (local) */
80  int_t ldalsum; /* LDA of lsum (local) */
81  int_t SolveMsgSent; /* Number of actual messages sent in LU-solve */
82  int_t SolveMsgVol; /* Volume of messages sent in the solve phase */
83 
84 
85  /*********************/
86  /* The following variables are used in the hybrid solver */
87 
88  /*-- Counts to be used in U^{-T} triangular solve. -- */
92  int_t ut_ldalsum; /* LDA of lsum (local) */
93  int_t *ut_ilsum; /* ilsum in column-wise */
94  int_t *utmod; /* Modification count for Ut-solve. */
95  int_t **ut_sendx_plist; /* Row process list to send down Xk */
96  int_t *utrecv; /* Modifications to be recev'd in proc column. */
97  int_t n_utsendx; /* Number of Xk I will receive */
98  int_t n_utrecvx; /* Number of Xk I will send */
103  Ucb_indptr_t **Ucb_indptr;/* Vertical linked list pointing to Uindex[] */
104  int_t **Ucb_valptr; /* Vertical linked list pointing to Unzval[] */
105 
106  /* some additional counters for L solve */
110 } LocalLU_t;
111 
112 
113 typedef struct {
117 } LUstruct_t;
118 
119 
120 /*-- Data structure for communication during matrix-vector multiplication. */
121 typedef struct {
123  int_t *ind_tosend; /* X indeices to be sent to other processes */
124  int_t *ind_torecv; /* X indeices to be received from other processes */
125  int_t *ptr_ind_tosend;/* Printers to ind_tosend[] (Size procs)
126  (also point to val_torecv) */
127  int_t *ptr_ind_torecv;/* Printers to ind_torecv[] (Size procs)
128  (also point to val_tosend) */
129  int *SendCounts; /* Numbers of X indices to be sent
130  (also numbers of X values to be received) */
131  int *RecvCounts; /* Numbers of X indices to be received
132  (also numbers of X values to be sent) */
133  double *val_tosend; /* X values to be sent to other processes */
134  double *val_torecv; /* X values to be received from other processes */
135  int_t TotalIndSend; /* Total number of indices to be sent
136  (also total number of values to be received) */
137  int_t TotalValSend; /* Total number of values to be sent.
138  (also total number of indices to be received) */
139 } pdgsmv_comm_t;
140 
141 /*-- Data structure for redistribution of B and X --*/
142 typedef struct {
145  int *ptr_to_ibuf, *ptr_to_dbuf;
146 
147  /* the following are needed in the hybrid solver */
150  int *disp_ibuf;
152  void *send_dbuf;
153 
154  int_t x2b, b2x;
157  void *send_dbuf2;
158  void *recv_dbuf2;
160 
161 /*-- Data structure holding the information for the solution phase --*/
162 typedef struct {
165  int_t num_diag_procs, *diag_procs, *diag_len;
168  int_t *A_colind_gsmv; /* After pdgsmv_init(), the global column
169  indices of A are translated into the relative
170  positions in the gathered x-vector.
171  This is re-used in repeated calls to pdgsmv() */
173 } SOLVEstruct_t;
174 
175 
176 /***********************************************************************
177  * Function prototypes
178  ***********************************************************************/
179 
180 #ifdef __cplusplus
181 extern "C" {
182 #endif
183 
184 
185 /* Supernodal LU factor related */
186 extern void
188  int_t *, int_t *, Stype_t, Dtype_t, Mtype_t);
189 extern void
191  int_t, double *, int_t *, int_t *,
193 extern void
195  double **, int_t **, int_t **);
196 extern int
198  SuperMatrix *);
199 extern void
201 extern void
204 extern void
206  int_t *, int_t *, int_t *, int_t *, int_t *,
208 extern void
210  double *, int_t);
211 
212 extern void dallocateA_dist (int_t, int_t, double **, int_t **, int_t **);
213 extern void dGenXtrue_dist (int_t, int_t, double *, int_t);
214 extern void dFillRHS_dist (char *, int_t, double *, int_t,
215  SuperMatrix *, double *, int_t);
216 extern int dcreate_matrix(SuperMatrix *, int, double **, int *,
217  double **, int *, FILE *, gridinfo_t *);
218 extern int dcreate_matrix_rb(SuperMatrix *, int, double **, int *,
219  double **, int *, FILE *, gridinfo_t *);
220 extern int dcreate_matrix_dat(SuperMatrix *, int, double **, int *,
221  double **, int *, FILE *, gridinfo_t *);
222 
223 /* Driver related */
224 extern void dgsequ_dist (SuperMatrix *, double *, double *, double *,
225  double *, double *, int_t *);
226 extern double dlangs_dist (char *, SuperMatrix *);
227 extern void dlaqgs_dist (SuperMatrix *, double *, double *, double,
228  double, double, char *);
229 extern void pdgsequ (SuperMatrix *, double *, double *, double *,
230  double *, double *, int_t *, gridinfo_t *);
231 extern double pdlangs (char *, SuperMatrix *, gridinfo_t *);
232 extern void pdlaqgs (SuperMatrix *, double *, double *, double,
233  double, double, char *);
234 extern int pdPermute_Dense_Matrix(int_t, int_t, int_t [], int_t[],
235  double [], int, double [], int, int,
236  gridinfo_t *);
237 
238 extern int sp_dtrsv_dist (char *, char *, char *, SuperMatrix *,
239  SuperMatrix *, double *, int *);
240 extern int sp_dgemv_dist (char *, double, SuperMatrix *, double *,
241  int, double, double *, int);
242 extern int sp_dgemm_dist (char *, int, double, SuperMatrix *,
243  double *, int, double, double *, int);
244 
245 extern float ddistribute(fact_t, int_t, SuperMatrix *, Glu_freeable_t *,
246  LUstruct_t *, gridinfo_t *);
248  ScalePermstruct_t *, double *,
249  int, int, gridinfo_t *, LUstruct_t *, double *,
250  SuperLUStat_t *, int *);
251 extern float pddistribute(fact_t, int_t, SuperMatrix *,
253  LUstruct_t *, gridinfo_t *);
254 extern void pdgssvx(superlu_options_t *, SuperMatrix *,
255  ScalePermstruct_t *, double *,
256  int, int, gridinfo_t *, LUstruct_t *,
257  SOLVEstruct_t *, double *, SuperLUStat_t *, int *);
258 extern int dSolveInit(superlu_options_t *, SuperMatrix *, int_t [], int_t [],
261  int_t [], int_t [], gridinfo_t *grid,
263 extern void pxgstrs_finalize(pxgstrs_comm_t *);
265 extern void dldperm_dist(int_t, int_t, int_t, int_t [], int_t [],
266  double [], int_t *, double [], double []);
267 extern int static_schedule(superlu_options_t *, int, int,
269  int_t *, int_t *, int *);
270 
271 /* #define GPU_PROF
272 #define IPM_PROF */
273 
274 extern int_t pdgstrf(superlu_options_t *, int, int, double,
275  LUstruct_t*, gridinfo_t*, SuperLUStat_t*, int*);
276 extern void pdgstrs_Bglobal(int_t, LUstruct_t *, gridinfo_t *,
277  double *, int_t, int, SuperLUStat_t *, int *);
279  double *, int_t, int_t, int_t, int, SOLVEstruct_t *,
280  SuperLUStat_t *, int *);
281 extern void dlsum_fmod(double *, double *, double *, double *,
282  int, int, int_t , int_t *, int_t, int_t, int_t,
283  int_t *, gridinfo_t *, LocalLU_t *,
284  MPI_Request [], SuperLUStat_t *);
285 extern void dlsum_bmod(double *, double *, double *,
286  int, int_t, int_t *, int_t *, Ucb_indptr_t **,
287  int_t **, int_t *, gridinfo_t *, LocalLU_t *,
288  MPI_Request [], SuperLUStat_t *);
289 extern void pdgsrfs(int_t, SuperMatrix *, double, LUstruct_t *,
291  double [], int_t, double [], int_t, int,
292  SOLVEstruct_t *, double *, SuperLUStat_t *, int *);
293 extern void pdgsrfs_ABXglobal(int_t, SuperMatrix *, double, LUstruct_t *,
294  gridinfo_t *, double *, int_t, double *, int_t,
295  int, double *, SuperLUStat_t *, int *);
297  gridinfo_t *, int_t *, int_t *[],
298  double *[], int_t *[], int_t []);
299 extern int pdgsmv_AXglobal(int_t, int_t [], double [], int_t [],
300  double [], double []);
301 extern int pdgsmv_AXglobal_abs(int_t, int_t [], double [], int_t [],
302  double [], double []);
303 extern void pdgsmv_init(SuperMatrix *, int_t *, gridinfo_t *,
304  pdgsmv_comm_t *);
305 extern void pdgsmv(int_t, SuperMatrix *, gridinfo_t *, pdgsmv_comm_t *,
306  double x[], double ax[]);
307 extern void pdgsmv_finalize(pdgsmv_comm_t *);
308 
309 /* Memory-related */
310 extern double *doubleMalloc_dist(int_t);
311 extern double *doubleCalloc_dist(int_t);
312 extern void *duser_malloc_dist (int_t, int_t);
313 extern void duser_free_dist (int_t, int_t);
316 extern void Destroy_LU(int_t, gridinfo_t *, LUstruct_t *);
317 extern void LUstructInit(const int_t, LUstruct_t *);
318 extern void LUstructFree(LUstruct_t *);
319 
320 /* Auxiliary routines */
321 extern void dfill_dist (double *, int_t, double);
322 extern void dinf_norm_error_dist (int_t, int_t, double*, int_t,
323  double*, int_t, gridinfo_t*);
324 extern void pdinf_norm_error(int, int_t, int_t, double [], int_t,
325  double [], int_t , gridinfo_t *);
326 extern void dreadhb_dist (int, FILE *, int_t *, int_t *, int_t *,
327  double **, int_t **, int_t **);
328 extern void dreadtriple(FILE *, int_t *, int_t *, int_t *,
329  double **, int_t **, int_t **);
330 extern void dreadrb_dist(FILE *, int *, int *, int *,
331  double **, int **, int **);
332 
333 /* Distribute the data for numerical factorization */
334 extern float ddist_psymbtonum(fact_t, int_t, SuperMatrix *,
336  LUstruct_t *, gridinfo_t *);
337 
338 
339 /* Routines for debugging */
340 extern void dPrintLblocks(int, int_t, gridinfo_t *, Glu_persist_t *,
341  LocalLU_t *);
342 extern void dPrintUblocks(int, int_t, gridinfo_t *, Glu_persist_t *,
343  LocalLU_t *);
347 extern int file_PrintDouble5(FILE *, char *, int_t, double *);
348 
349 
350 /* BLAS */
351 
352 #ifdef USE_VENDOR_BLAS
353 extern int dgemm_(const char*, const char*, const int*, const int*, const int*,
354  const double*, const double*, const int*, const double*,
355  const int*, const double*, double*, const int*, int, int);
356 extern int dtrsv_(char*, char*, char*, int*, double*, int*,
357  double*, int*, int, int, int);
358 extern int dtrsm_(char*, char*, char*, char*, int*, int*,
359  double*, double*, int*, double*,
360  int*, int, int, int, int);
361 extern int dgemv_(char *, int *, int *, double *, double *a, int *,
362  double *, int *, double *, double *, int *, int);
363 #else
364 extern int dgemm_(const char*, const char*, const int*, const int*, const int*,
365  const double*, const double*, const int*,
366  const double*, const int*, const double*,
367  double*, const int*);
368 extern int dtrsv_(char*, char*, char*, int*, double*, int*,
369  double*, int*);
370 extern int dtrsm_(char*, char*, char*, char*, int*, int*,
371  double*, double*, int*, double*, int*);
372 extern int dgemv_(char *, int *, int *, double *, double *a, int *,
373  double *, int *, double *, double *, int *);
374 #endif
375 
376 extern int dger_(int*, int*, double*, double*, int*,
377  double*, int*, double*, int*);
378 
379 
380 #ifdef __cplusplus
381  }
382 #endif
383 
384 #endif /* __SUPERLU_dDEFS */
385 
void * duser_malloc_dist(int_t bytes, int_t which_end)
Definition: dmemory.c:20
int dSolveInit(superlu_options_t *options, SuperMatrix *A, int_t perm_r[], int_t perm_c[], int_t nrhs, LUstruct_t *LUstruct, gridinfo_t *grid, SOLVEstruct_t *SOLVEstruct)
Initialize the data structure for the solution phase.
Definition: pdutil.c:393
void pdgsmv(int_t abs, SuperMatrix *A_internal, gridinfo_t *grid, pdgsmv_comm_t *gsmv_comm, double x[], double ax[])
Definition: pdgsmv.c:225
int_t * extern_start
Definition: superlu_ddefs.h:122
int_t n_utsendx
Definition: superlu_ddefs.h:97
int * ptr_to_ibuf
Definition: superlu_ddefs.h:145
void pdgsequ(SuperMatrix *A, double *r, double *c, double *rowcnd, double *colcnd, double *amax, int_t *info, gridinfo_t *grid)
Definition: pdgsequ.c:76
void pdinf_norm_error(int iam, int_t n, int_t nrhs, double x[], int_t ldx, double xtrue[], int_t ldxtrue, gridinfo_t *grid)
Check the inf-norm of the error vector.
Definition: pdutil.c:499
void dCopy_CompCol_Matrix_dist(SuperMatrix *A, SuperMatrix *B)
Copy matrix A into matrix B.
Definition: dutil.c:113
int_t nfrecvmod
Definition: superlu_ddefs.h:109
void dFillRHS_dist(char *trans, int_t nrhs, double *x, int_t ldx, SuperMatrix *A, double *rhs, int_t ldb)
Let rhs[i] = sum of i-th row of A, so the solution vector is all 1's.
Definition: dutil.c:298
void dCompRow_to_CompCol_dist(int_t m, int_t n, int_t nnz, double *a, int_t *colind, int_t *rowptr, double **at, int_t **rowind, int_t **colptr)
Convert a row compressed storage into a column compressed storage.
Definition: dutil.c:75
int_t ** Ucb_valptr
Definition: superlu_ddefs.h:104
Definition: superlu_defs.h:250
Definition: psymbfact.h:47
int_t lbnum
Definition: superlu_ddefs.h:26
int_t n_utrecvx
Definition: superlu_ddefs.h:98
void * send_dbuf
Definition: superlu_ddefs.h:152
void pdlaqgs(SuperMatrix *A, double *r, double *c, double rowcnd, double colcnd, double amax, char *equed)
Definition: pdlaqgs.c:75
int * disp_ibuf
Definition: superlu_ddefs.h:150
pdgsmv_comm_t * gsmv_comm
Definition: superlu_ddefs.h:166
Definition: superlu_defs.h:501
int_t * fmod
Definition: superlu_ddefs.h:65
int ** ToSendR
Definition: superlu_ddefs.h:62
int_t nbrecvx
Definition: superlu_ddefs.h:73
int sp_dgemv_dist(char *trans, double alpha, SuperMatrix *A, double *x, int incx, double beta, double *y, int incy)
Definition: dsp_blas2.c:382
int * B_to_X_SendCnt
Definition: superlu_ddefs.h:143
double * ujrow
Definition: superlu_ddefs.h:50
int_t nfrecvx
Definition: superlu_ddefs.h:68
int_t SolveMsgVol
Definition: superlu_ddefs.h:82
void dreadrb_dist(FILE *fp, int_t *nrow, int_t *ncol, int_t *nonz, double **nzval, int **rowind, int **colptr)
Definition: dreadrb.c:273
Definition: superlu_ddefs.h:121
int dgemv_(char *, int *, int *, double *, double *a, int *, double *, int *, double *, double *, int *)
Definition: dgemv.c:9
void pdgsrfs_ABXglobal(int_t n, SuperMatrix *A, double anorm, LUstruct_t *LUstruct, gridinfo_t *grid, double *B, int_t ldb, double *X, int_t ldx, int nrhs, double *berr, SuperLUStat_t *stat, int *info)
Definition: pdgsrfs_ABXglobal.c:112
double dlangs_dist(char *norm, SuperMatrix *A)
Definition: dlangs.c:52
Definition: util_dist.h:72
double * doubleMalloc_dist(int_t n)
Definition: dmemory.c:142
int_t n_utrecvmod
Definition: superlu_ddefs.h:99
int_t * row_to_proc
Definition: superlu_ddefs.h:163
Definition: supermatrix.h:44
void pdgsmv_init(SuperMatrix *A, int_t *row_to_proc, gridinfo_t *grid, pdgsmv_comm_t *gsmv_comm)
Definition: pdgsmv.c:17
int_t * bmod
Definition: superlu_ddefs.h:70
fact_t
Definition: superlu_enum_consts.h:17
int_t * etree
Definition: superlu_ddefs.h:114
int_t TotalValSend
Definition: superlu_ddefs.h:137
int_t indpos
Definition: superlu_ddefs.h:27
Stype_t
Definition: supermatrix.h:12
int_t * ilsum
Definition: superlu_ddefs.h:78
Definition: superlu_ddefs.h:25
void pdgsrfs(int_t n, SuperMatrix *A, double anorm, LUstruct_t *LUstruct, ScalePermstruct_t *ScalePermstruct, gridinfo_t *grid, double *B, int_t ldb, double *X, int_t ldx, int nrhs, SOLVEstruct_t *SOLVEstruct, double *berr, SuperLUStat_t *stat, int *info)
Definition: pdgsrfs.c:107
void LUstructFree(LUstruct_t *)
Deallocate LUstruct.
Definition: util.c:205
int pdPermute_Dense_Matrix(int_t fst_row, int_t m_loc, int_t row_to_proc[], int_t perm[], double X[], int ldx, double B[], int ldb, int nrhs, gridinfo_t *grid)
Permute the distributed dense matrix: B <= perm(X). perm[i] = j means the i-th row of X is in the j-t...
Definition: pdutil.c:281
int_t * Urbs
Definition: superlu_ddefs.h:102
void Destroy_LU(int_t, gridinfo_t *, LUstruct_t *)
Destroy distributed L & U matrices.
Definition: util.c:101
void dldperm_dist(int_t job, int_t n, int_t nnz, int_t colptr[], int_t adjncy[], double nzval[], int_t *perm, double u[], double v[])
Definition: dldperm_dist.c:85
Definition: superlu_ddefs.h:162
double * val_torecv
Definition: superlu_ddefs.h:134
void dfill_dist(double *a, int_t alen, double dval)
Fills a double precision array with a given value.
Definition: dutil.c:311
int_t ** fsendx_plist
Definition: superlu_ddefs.h:66
int * ToSendD
Definition: superlu_ddefs.h:61
Glu_persist_t * Glu_persist
Definition: superlu_ddefs.h:115
int_t nfsendx
Definition: superlu_ddefs.h:69
void dGenXtrue_dist(int_t n, int_t nrhs, double *x, int_t ldx)
Definition: dutil.c:285
int_t UT_SOLVE
Definition: superlu_ddefs.h:89
double pdlangs(char *norm, SuperMatrix *A, gridinfo_t *grid)
Definition: pdlangs.c:55
int dtrsm_(char *, char *, char *, char *, int *, int *, double *, double *, int *, double *, int *)
Definition: dtrsm.c:9
float ddistribute(fact_t fact, int_t n, SuperMatrix *A, Glu_freeable_t *Glu_freeable, LUstruct_t *LUstruct, gridinfo_t *grid)
Definition: ddistribute.c:51
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 *)
#define NBUFFERS
Definition: superlu_defs.h:100
Definition: superlu_defs.h:355
int dcreate_matrix_dat(SuperMatrix *, int, double **, int *, double **, int *, FILE *, gridinfo_t *)
void dPrintLblocks(int iam, int_t nsupers, gridinfo_t *grid, Glu_persist_t *Glu_persist, LocalLU_t *Llu)
Print the blocks in the factored matrix L.
Definition: dutil.c:369
int_t * send_ibuf
Definition: superlu_ddefs.h:151
void dPrint_CompCol_Matrix_dist(SuperMatrix *A)
Definition: dutil.c:133
void dPrintUblocks(int iam, int_t nsupers, gridinfo_t *grid, Glu_persist_t *Glu_persist, LocalLU_t *Llu)
Print the blocks in the factored matrix U.
Definition: dutil.c:418
int pdgsmv_AXglobal_setup(SuperMatrix *A, Glu_persist_t *Glu_persist, gridinfo_t *grid, int_t *m, int_t *update[], double *val[], int_t *bindx[], int_t *mv_sup_to_proc)
Definition: pdgsmv_AXglobal.c:23
void dreadhb_dist(int iam, FILE *fp, int_t *nrow, int_t *ncol, int_t *nonz, double **nzval, int_t **rowind, int_t **colptr)
Definition: dreadhb.c:97
void dPrint_Dense_Matrix_dist(SuperMatrix *A)
Definition: dutil.c:157
void pdgstrs_Bglobal(int_t n, LUstruct_t *LUstruct, gridinfo_t *grid, double *B, int_t ldb, int nrhs, SuperLUStat_t *stat, int *info)
Definition: pdgstrs_Bglobal.c:94
Definition: superlu_defs.h:227
void dinf_norm_error_dist(int_t n, int_t nrhs, double *x, int_t ldx, double *xtrue, int_t ldxtrue, gridinfo_t *grid)
Check the inf-norm of the error vector.
Definition: dutil.c:321
void dlsum_fmod(double *lsum, double *x, double *xk, double *rtemp, int nrhs, int knsupc, int_t k, int_t *fmod, int_t nlb, int_t lptr, int_t luptr, int_t *xsup, gridinfo_t *grid, LocalLU_t *Llu, MPI_Request send_req[], SuperLUStat_t *stat)
Definition: pdgstrs_lsum.c:45
int_t * A_colind_gsmv
Definition: superlu_ddefs.h:168
int_t ldalsum
Definition: superlu_ddefs.h:80
double ** Unzval_br_ptr
Definition: superlu_ddefs.h:39
void dallocateA_dist(int_t n, int_t nnz, double **a, int_t **asub, int_t **xa)
Definition: dmemory.c:134
int_t FRECV
Definition: superlu_ddefs.h:91
int_t ** Lrowind_bc_ptr
Definition: superlu_ddefs.h:36
void dlaqgs_dist(SuperMatrix *A, double *r, double *c, double rowcnd, double colcnd, double amax, char *equed)
Definition: dlaqgs.c:71
int_t * utrecv
Definition: superlu_ddefs.h:96
int * RecvCounts
Definition: superlu_ddefs.h:131
int_t * mod_bit
Definition: superlu_ddefs.h:75
int_t * frecv
Definition: superlu_ddefs.h:67
Definition: superlu_defs.h:305
int_t nleaf
Definition: superlu_ddefs.h:108
int_t num_diag_procs
Definition: superlu_ddefs.h:165
void pdgssvx_ABglobal(superlu_options_t *options, SuperMatrix *A, ScalePermstruct_t *ScalePermstruct, double B[], int ldb, int nrhs, gridinfo_t *grid, LUstruct_t *LUstruct, double *berr, SuperLUStat_t *stat, int *info)
Definition: pdgssvx_ABglobal.c:455
void duser_free_dist(int_t bytes, int_t which_end)
Definition: dmemory.c:39
int_t * recv_ibuf2
Definition: superlu_ddefs.h:156
void pdgstrs(int_t n, LUstruct_t *LUstruct, ScalePermstruct_t *ScalePermstruct, gridinfo_t *grid, double *B, int_t m_loc, int_t fst_row, int_t ldb, int nrhs, SOLVEstruct_t *SOLVEstruct, SuperLUStat_t *stat, int *info)
Definition: pdgstrs.c:456
int_t x2b
Definition: superlu_ddefs.h:154
double ** Lnzval_bc_ptr
Definition: superlu_ddefs.h:37
int dPrint_CompRowLoc_Matrix_dist(SuperMatrix *A)
Definition: dutil.c:174
Definition: superlu_ddefs.h:142
Mtype_t
Definition: supermatrix.h:32
void dSolveFinalize(superlu_options_t *options, SOLVEstruct_t *SOLVEstruct)
Release the resources used for the solution phase.
Definition: pdutil.c:480
int static_schedule(superlu_options_t *options, int m, int n, LUstruct_t *LUstruct, gridinfo_t *grid, SuperLUStat_t *stat, int_t *perm_c_supno, int_t *iperm_c_supno, int *info)
Definition: static_schedule.c:36
int_t * send_ibuf2
Definition: superlu_ddefs.h:155
void dCopy_Dense_Matrix_dist(int_t M, int_t N, double *X, int_t ldx, double *Y, int_t ldy)
Definition: dutil.c:237
int_t * utmod
Definition: superlu_ddefs.h:94
int * ToRecv
Definition: superlu_ddefs.h:60
LocalLU_t * Llu
Definition: superlu_ddefs.h:116
int_t pdgstrf(superlu_options_t *options, int m, int n, double anorm, LUstruct_t *LUstruct, gridinfo_t *grid, SuperLUStat_t *stat, int *info)
Definition: pdgstrf_irecv.c:183
double * val_tosend
Definition: superlu_ddefs.h:133
int pdgsmv_AXglobal(int_t m, int_t update[], double val[], int_t bindx[], double X[], double ax[])
Definition: pdgsmv_AXglobal.c:249
void dCreate_SuperNode_Matrix_dist(SuperMatrix *L, int_t m, int_t n, int_t nnz, double *nzval, int_t *nzval_colptr, int_t *rowind, int_t *rowind_colptr, int_t *col_to_sup, int_t *sup_to_col, Stype_t stype, Dtype_t dtype, Mtype_t mtype)
Definition: dutil.c:257
int_t dQuerySpace_dist(int_t n, LUstruct_t *LUstruct, gridinfo_t *grid, SuperLUStat_t *stat, mem_usage_t *mem_usage)
Definition: dmemory.c:63
int int_t
Definition: superlu_defs.h:37
int_t * xrow_to_proc
Definition: superlu_ddefs.h:172
int dcreate_matrix(SuperMatrix *, int, double **, int *, double **, int *, FILE *, gridinfo_t *)
Definition: dcreate_matrix.c:56
int_t pxgstrs_init(int_t, int_t, int_t, int_t, int_t[], int_t[], gridinfo_t *grid, Glu_persist_t *, SOLVEstruct_t *)
Definition: util.c:415
int_t * ptr_ind_tosend
Definition: superlu_ddefs.h:125
int_t ** Ufstnz_br_ptr
Definition: superlu_ddefs.h:38
int sp_dgemm_dist(char *transa, int n, double alpha, SuperMatrix *A, double *b, int ldb, double beta, double *c, int ldc)
Definition: dsp_blas3.c:114
int_t * ind_tosend
Definition: superlu_ddefs.h:123
void * send_dbuf2
Definition: superlu_ddefs.h:157
void pdgsmv_finalize(pdgsmv_comm_t *gsmv_comm)
Definition: pdgsmv.c:361
int_t * ind_torecv
Definition: superlu_ddefs.h:124
float pddistribute(fact_t fact, int_t n, SuperMatrix *A, ScalePermstruct_t *ScalePermstruct, Glu_freeable_t *Glu_freeable, LUstruct_t *LUstruct, gridinfo_t *grid)
Definition: pddistribute.c:310
int_t L_SOLVE
Definition: superlu_ddefs.h:90
int_t nroot
Definition: superlu_ddefs.h:100
int * X_to_B_vSendCnt
Definition: superlu_ddefs.h:149
void LUstructInit(const int_t, LUstruct_t *)
Allocate storage in LUstruct.
Definition: util.c:192
int_t TotalIndSend
Definition: superlu_ddefs.h:135
Definition: superlu_ddefs.h:35
int * X_to_B_iSendCnt
Definition: superlu_ddefs.h:148
Dtype_t
Definition: supermatrix.h:25
int * X_to_B_SendCnt
Definition: superlu_ddefs.h:144
int sp_dtrsv_dist(char *uplo, char *trans, char *diag, SuperMatrix *L, SuperMatrix *U, double *x, int *info)
Definition: dsp_blas2.c:83
void dgsequ_dist(SuperMatrix *A, double *r, double *c, double *rowcnd, double *colcnd, double *amax, int_t *info)
Definition: dgsequ.c:74
Definitions which are precision-neutral.
void dCreate_Dense_Matrix_dist(SuperMatrix *X, int_t m, int_t n, double *x, int_t ldx, Stype_t stype, Dtype_t dtype, Mtype_t mtype)
Definition: dutil.c:218
int pdCompRow_loc_to_CompCol_global(int_t need_value, SuperMatrix *A, gridinfo_t *grid, SuperMatrix *GA)
Gather A from the distributed compressed row format to global A in compressed column format...
Definition: pdutil.c:19
int dger_(int *, int *, double *, double *, int *, double *, int *, double *, int *)
Definition: dger.c:9
int_t ** ut_sendx_plist
Definition: superlu_ddefs.h:95
Definition: superlu_defs.h:531
int_t * inv_perm_c
Definition: superlu_ddefs.h:164
void pxgstrs_finalize(pxgstrs_comm_t *)
Definition: util.c:542
int file_PrintDouble5(FILE *fp, char *name, int_t len, double *x)
Definition: dutil.c:354
int_t ut_ldalsum
Definition: superlu_ddefs.h:92
Definition: superlu_ddefs.h:113
void dlsum_bmod(double *lsum, double *x, double *xk, int nrhs, int_t k, int_t *bmod, int_t *Urbs, Ucb_indptr_t **Ucb_indptr, int_t **Ucb_valptr, int_t *xsup, gridinfo_t *grid, LocalLU_t *Llu, MPI_Request send_req[], SuperLUStat_t *stat)
Definition: pdgstrs_lsum.c:208
int pdgsmv_AXglobal_abs(int_t m, int_t update[], double val[], int_t bindx[], double X[], double ax[])
Definition: pdgsmv_AXglobal.c:275
pxgstrs_comm_t * gstrs_comm
Definition: superlu_ddefs.h:167
void dreadtriple(FILE *fp, int_t *m, int_t *n, int_t *nonz, double **nzval, int_t **rowind, int_t **colptr)
Definition: dreadtriple.c:25
double * doubleCalloc_dist(int_t n)
Definition: dmemory.c:149
int_t * ut_modbit
Definition: superlu_ddefs.h:101
int dcreate_matrix_rb(SuperMatrix *, int, double **, int *, double **, int *, FILE *, gridinfo_t *)
void * recv_dbuf2
Definition: superlu_ddefs.h:158
void dCreate_CompCol_Matrix_dist(SuperMatrix *A, int_t m, int_t n, int_t nnz, double *nzval, int_t *rowind, int_t *colptr, Stype_t stype, Dtype_t dtype, Mtype_t mtype)
Definition: dutil.c:28
void pdgssvx(superlu_options_t *options, SuperMatrix *A, ScalePermstruct_t *ScalePermstruct, double B[], int ldb, int nrhs, gridinfo_t *grid, LUstruct_t *LUstruct, SOLVEstruct_t *SOLVEstruct, double *berr, SuperLUStat_t *stat, int *info)
Definition: pdgssvx.c:486
Ucb_indptr_t ** Ucb_indptr
Definition: superlu_ddefs.h:103
void dCreate_CompRowLoc_Matrix_dist(SuperMatrix *A, int_t m, int_t n, int_t nnz_loc, int_t m_loc, int_t fst_row, double *nzval, int_t *colind, int_t *rowptr, Stype_t stype, Dtype_t dtype, Mtype_t mtype)
Definition: dutil.c:49
int dtrsv_(char *, char *, char *, int *, double *, int *, double *, int *)
Definition: dtrsv.c:9
int_t nbsendx
Definition: superlu_ddefs.h:74
int_t ** bsendx_plist
Definition: superlu_ddefs.h:71
int_t * brecv
Definition: superlu_ddefs.h:72
int_t * ut_ilsum
Definition: superlu_ddefs.h:93
#define MAX_LOOKAHEADS
Definition: superlu_ddefs.h:34
int_t SolveMsgSent
Definition: superlu_ddefs.h:81
int_t n
Definition: superlu_ddefs.h:107
int * SendCounts
Definition: superlu_ddefs.h:129
float ddist_psymbtonum(fact_t fact, int_t n, SuperMatrix *A, ScalePermstruct_t *ScalePermstruct, Pslu_freeable_t *Pslu_freeable, LUstruct_t *LUstruct, gridinfo_t *grid)
Definition: pdsymbfact_distdata.c:1177
int_t * ptr_ind_torecv
Definition: superlu_ddefs.h:127