Fast Auxiliary Space Preconditioning  1.8.4 Feb/15/2016
interface_superlu.c
Go to the documentation of this file.
1 
8 #include <stdio.h>
9 #include <stdlib.h>
10 #include <time.h>
11 
12 #include "fasp.h"
13 #include "fasp_functs.h"
14 
15 #if WITH_SuperLU
16 #include "slu_ddefs.h"
17 #endif
18 
19 /*---------------------------------*/
20 /*-- Public Functions --*/
21 /*---------------------------------*/
22 
41  dvector *b,
42  dvector *u,
43  const SHORT prtlvl)
44 {
45 
46 #if WITH_SuperLU
47 
48  SuperMatrix A, L, U, B;
49 
50  int *perm_r; /* row permutations from partial pivoting */
51  int *perm_c; /* column permutation vector */
52  int nrhs=1, info, m=ptrA->row, n=ptrA->col, nnz=ptrA->nnz;
53 
54  if ( prtlvl > PRINT_NONE )
55  printf("superlu: nr=%d, nc=%d, nnz=%d\n", m, n, nnz);
56 
57  clock_t LU_start=clock();
58 
59  dCSRmat tempA=fasp_dcsr_create(m,n,nnz); fasp_dcsr_cp(ptrA,&tempA);
60 
61  dvector tempb=fasp_dvec_create(n); fasp_dvec_cp(b, &tempb);
62 
63  /* Create matrix A in the format expected by SuperLU. */
64  dCreate_CompCol_Matrix(&A, m, n, nnz, tempA.val, tempA.JA, tempA.IA, SLU_NR, SLU_D, SLU_GE);
65 
66  /* Create right-hand side B. */
67  dCreate_Dense_Matrix(&B, m, nrhs, tempb.val, m, SLU_DN, SLU_D, SLU_GE);
68 
69  if ( !(perm_r = intMalloc(m)) ) ABORT("Malloc fails for perm_r[].");
70  if ( !(perm_c = intMalloc(n)) ) ABORT("Malloc fails for perm_c[].");
71 
72  /* Set the default input options. */
73  superlu_options_t options;
74  set_default_options(&options);
75  options.ColPerm = COLAMD; //MMD_AT_PLUS_A; MMD_ATA; NATURAL;
76 
77  /* Initialize the statistics variables. */
78  SuperLUStat_t stat;
79  StatInit(&stat);
80 
81  /* SuperLU */
82  dgssv(&options, &A, perm_c, perm_r, &L, &U, &B, &stat, &info);
83 
84  DNformat *BB = (DNformat *) B.Store;
85  u->val = (double *) BB->nzval;
86  u->row = n;
87 
88  if ( prtlvl > PRINT_MIN ) {
89  clock_t LU_end=clock();
90  double LU_time = (double)(LU_end - LU_start)/(double)(CLOCKS_PER_SEC);
91  printf("SuperLU totally costs %f seconds.\n", LU_time);
92  }
93 
94  /* De-allocate storage */
95  SUPERLU_FREE(perm_r);
96  SUPERLU_FREE(perm_c);
97  Destroy_CompCol_Matrix(&A);
98  Destroy_SuperMatrix_Store(&B);
99  Destroy_SuperNode_Matrix(&L);
100  Destroy_CompCol_Matrix(&U);
101  StatFree(&stat);
102 
103  return FASP_SUCCESS;
104 
105 #else
106 
107  printf("### ERROR: SuperLU is not available!\n");
108  return ERROR_SOLVER_EXIT;
109 
110 #endif
111 
112 }
113 
114 /*---------------------------------*/
115 /*-- End of File --*/
116 /*---------------------------------*/
int fasp_solver_superlu(dCSRmat *ptrA, dvector *b, dvector *u, const SHORT prtlvl)
Solve Au=b by SuperLU.
dvector fasp_dvec_create(const INT m)
Create dvector data space of REAL type.
Definition: vec.c:56
dCSRmat fasp_dcsr_create(const INT m, const INT n, const INT nnz)
Create CSR sparse matrix data memory space.
Definition: sparse_csr.c:34
INT * JA
integer array of column indexes, the size is nnz
Definition: fasp.h:163
void fasp_dvec_cp(dvector *x, dvector *y)
Copy dvector x to dvector y.
Definition: vec.c:345
REAL * val
nonzero entries of A
Definition: fasp.h:166
REAL * val
actual vector entries
Definition: fasp.h:348
Vector with n entries of REAL type.
Definition: fasp.h:342
#define FASP_SUCCESS
Definition of return status and error messages.
Definition: fasp_const.h:27
void fasp_dcsr_cp(dCSRmat *A, dCSRmat *B)
copy a dCSRmat to a new one B=A
Definition: sparse_csr.c:723
INT nnz
number of nonzero entries
Definition: fasp.h:157
INT col
column of matrix A, n
Definition: fasp.h:154
Main header file for FASP.
#define PRINT_NONE
Print level for all subroutines – not including DEBUG output.
Definition: fasp_const.h:79
INT row
row number of matrix A, m
Definition: fasp.h:151
Sparse matrix of REAL type in CSR format.
Definition: fasp.h:148
#define ERROR_SOLVER_EXIT
Definition: fasp_const.h:55
INT row
number of rows
Definition: fasp.h:345
#define SHORT
FASP integer and floating point numbers.
Definition: fasp.h:63
INT * IA
integer array of row pointers, the size is m+1
Definition: fasp.h:160
#define PRINT_MIN
Definition: fasp_const.h:80