Fast Auxiliary Space Preconditioning  1.8.4 Feb/15/2016
amg_setup_cr.c
Go to the documentation of this file.
1 
13 #include <math.h>
14 #include <time.h>
15 
16 #include "fasp.h"
17 #include "fasp_functs.h"
18 
19 /*---------------------------------*/
20 /*-- Public Functions --*/
21 /*---------------------------------*/
22 
39  AMG_param *param)
40 {
41  const SHORT prtlvl = param->print_level;
42  const SHORT min_cdof = MAX(param->coarse_dof,50);
43  const INT m = mgl[0].A.row;
44 
45  // local variables
46  INT i_0 = 0, i_n;
47  SHORT level = 0, status = FASP_SUCCESS;
48  SHORT max_levels = param->max_levels;
49  REAL setup_start, setup_end;
50 
51  // The variable vertices stores level info (fine: 0; coarse: 1)
52  ivector vertices = fasp_ivec_create(m);
53 
54  fasp_gettime(&setup_start);
55 
56 #if DEBUG_MODE > 0
57  printf("### DEBUG: %s ...... [Start]\n", __FUNCTION__);
58  printf("### DEBUG: nr=%d, nc=%d, nnz=%d\n",
59  mgl[0].A.row, mgl[0].A.col, mgl[0].A.nnz);
60 #endif
61 
62 #if DIAGONAL_PREF
63  fasp_dcsr_diagpref(&mgl[0].A); // reorder each row to make diagonal appear first
64 #endif
65 
66  // Main AMG setup loop
67  while ( (mgl[level].A.row > min_cdof) && (level < max_levels-1) ) {
68 
69  /*-- Coarsen and form the structure of interpolation --*/
70  i_n = mgl[level].A.row-1;
71 
72  fasp_amg_coarsening_cr(i_0,i_n,&mgl[level].A, &vertices, param);
73 
74  /*-- Form interpolation --*/
75  /* 1. SPARSITY -- Form ip and jp */
76  /* First a symbolic one
77  then gather the list */
78  /* 2. COEFFICIENTS -- Form P */
79  // energymin(mgl[level].A, &vertices[level], mgl[level].P, param);
80  // fasp_mem_free(vertices[level].val);
81 
82  /*-- Form coarse level stiffness matrix --*/
83  // fasp_dcsr_trans(mgl[level].P, mgl[level].R);
84 
85  /*-- Form coarse level stiffness matrix --*/
86  //fasp_blas_dcsr_rap(mgl[level].R, mgl[level].A, mgl[level].P, mgl[level+1].A);
87 
88  ++level;
89 
90 #if DIAGONAL_PREF
91  fasp_dcsr_diagpref(&mgl[level].A); // reorder each row to make diagonal appear first
92 #endif
93  }
94 
95  // setup total level number and current level
96  mgl[0].num_levels = max_levels = level+1;
97  mgl[0].w = fasp_dvec_create(m);
98 
99  for ( level = 1; level < max_levels; ++level ) {
100  INT mm = mgl[level].A.row;
101  mgl[level].num_levels = max_levels;
102  mgl[level].b = fasp_dvec_create(mm);
103  mgl[level].x = fasp_dvec_create(mm);
104  mgl[level].w = fasp_dvec_create(mm);
105  }
106 
107  if ( prtlvl > PRINT_NONE ) {
108  fasp_gettime(&setup_end);
109  print_amgcomplexity(mgl,prtlvl);
110  print_cputime("Compatible relaxation setup", setup_end - setup_start);
111  }
112 
113  fasp_ivec_free(&vertices);
114 
115 #if DEBUG_MODE > 0
116  printf("### DEBUG: %s ...... [Finish]\n", __FUNCTION__);
117 #endif
118 
119  return status;
120 }
121 
122 /*---------------------------------*/
123 /*-- End of File --*/
124 /*---------------------------------*/
dvector fasp_dvec_create(const INT m)
Create dvector data space of REAL type.
Definition: vec.c:56
Parameters for AMG solver.
Definition: fasp.h:583
void fasp_dcsr_diagpref(dCSRmat *A)
Re-order the column and data arrays of a CSR matrix, so that the first entry in each row is the diago...
Definition: sparse_csr.c:553
INT coarse_dof
max number of coarsest level DOF
Definition: fasp.h:601
#define REAL
Definition: fasp.h:67
Vector with n entries of INT type.
Definition: fasp.h:356
SHORT num_levels
number of levels in use <= max_levels
Definition: fasp.h:730
dvector b
pointer to the right-hand side at level level_num
Definition: fasp.h:744
#define INT
Definition: fasp.h:64
#define FASP_SUCCESS
Definition of return status and error messages.
Definition: fasp_const.h:27
void fasp_ivec_free(ivector *u)
Free vector data space of INT type.
Definition: vec.c:159
dvector x
pointer to the iterative solution at level level_num
Definition: fasp.h:747
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.
SHORT max_levels
max number of levels of AMG
Definition: fasp.h:598
void fasp_gettime(REAL *time)
Get system time.
Definition: timing.c:28
void print_amgcomplexity(AMG_data *mgl, const SHORT prtlvl)
Print complexities of AMG method.
Definition: message.c:79
#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
SHORT fasp_amg_setup_cr(AMG_data *mgl, AMG_param *param)
Set up phase of Brannick Falgout CR coarsening for classic AMG.
Definition: amg_setup_cr.c:38
#define SHORT
FASP integer and floating point numbers.
Definition: fasp.h:63
dCSRmat A
pointer to the matrix at level level_num
Definition: fasp.h:735
void print_cputime(const char *message, const REAL cputime)
Print CPU walltime.
Definition: message.c:165
SHORT print_level
print level for AMG
Definition: fasp.h:589
Data for AMG solvers.
Definition: fasp.h:722
ivector fasp_ivec_create(const INT m)
Create vector data space of INT type.
Definition: vec.c:78
dvector w
Temporary work space.
Definition: fasp.h:781
INT fasp_amg_coarsening_cr(const INT i_0, const INT i_n, dCSRmat *A, ivector *vertices, AMG_param *param)
CR coarsening.
Definition: coarsening_cr.c:42
#define MAX(a, b)
Definition of max, min, abs.
Definition: fasp.h:72