source: CIVL/examples/mpi-omp/AMG2013/parcsr_ls/par_amg_setup.c

main
Last change on this file was ea777aa, checked in by Alex Wilton <awilton@…>, 3 years ago

Moved examples, include, build_default.properties, common.xml, and README out from dev.civl.com into the root of the repo.

git-svn-id: ​svn://vsl.cis.udel.edu/civl/trunk@5704 fb995dde-84ed-4084-dfe6-e5aef3e2452c
  • Property mode set to 100644
File size: 73.2 KB
Line 
1/*BHEADER**********************************************************************
2 * Copyright (c) 2008, Lawrence Livermore National Security, LLC.
3 * Produced at the Lawrence Livermore National Laboratory.
4 * This file is part of HYPRE. See file COPYRIGHT for details.
5 *
6 * HYPRE is free software; you can redistribute it and/or modify it under the
7 * terms of the GNU Lesser General Public License (as published by the Free
8 * Software Foundation) version 2.1 dated February 1999.
9 *
10 * $Revision: 2.4 $
11 ***********************************************************************EHEADER*/
12
13
14#include "headers.h"
15#include "par_amg.h"
16/*#include "par_csr_block_matrix.h"*/
17
18#define DEBUG 0
19
20#define THETA 0
21
22#define MPIP_SETUP_ON 0
23
24/*****************************************************************************
25 *
26 * Routine for driving the setup phase of AMG
27 *
28 *****************************************************************************/
29
30/*****************************************************************************
31 * hypre_BoomerAMGSetup
32 *****************************************************************************/
33
34int
35hypre_BoomerAMGSetup( void *amg_vdata,
36 hypre_ParCSRMatrix *A,
37 hypre_ParVector *f,
38 hypre_ParVector *u )
39{
40 MPI_Comm comm = hypre_ParCSRMatrixComm(A);
41
42 hypre_ParAMGData *amg_data = amg_vdata;
43
44 /* Data Structure variables */
45
46 hypre_ParCSRMatrix **A_array;
47 hypre_ParVector **F_array;
48 hypre_ParVector **U_array;
49 hypre_ParVector *Vtemp;
50 hypre_ParVector *Rtemp;
51 hypre_ParVector *Ptemp;
52 hypre_ParVector *Ztemp;
53 hypre_ParCSRMatrix **P_array;
54 hypre_ParVector *Residual_array;
55 int **CF_marker_array;
56 int **dof_func_array;
57 int *dof_func;
58 int *dof_func1;
59 int *col_offd_S_to_A;
60 int *col_offd_SN_to_AN;
61 double *relax_weight;
62 double *omega;
63 double strong_threshold;
64 double max_row_sum;
65 double trunc_factor, jacobi_trunc_threshold;
66 double agg_trunc_factor;
67 double S_commpkg_switch;
68 /*double CR_rate;*/
69
70 int max_levels;
71 int amg_logging;
72 int amg_print_level;
73 int debug_flag;
74 int local_num_vars;
75 int P_max_elmts;
76 int P_max1;
77 int P_max2;
78 int agg_P_max_elmts;
79 /*int IS_type;
80 int num_CR_relax_steps;
81 int CR_use_CG; */
82
83
84 /*hypre_ParCSRBlockMatrix **A_block_array, **P_block_array;*/
85
86 /* Local variables */
87 int *CF_marker;
88 int *CFN_marker;
89 int *CF2_marker;
90 hypre_ParCSRMatrix *S;
91 hypre_ParCSRMatrix *S2;
92 hypre_ParCSRMatrix *SN;
93 hypre_ParCSRMatrix *P;
94 hypre_ParCSRMatrix *P1;
95 hypre_ParCSRMatrix *P2;
96 hypre_ParCSRMatrix *PN;
97 hypre_ParCSRMatrix *A_H;
98 hypre_ParCSRMatrix *AN;
99 double **l1_norms = NULL;
100 /*double *SmoothVecs = NULL;*/
101
102 int old_num_levels, num_levels;
103 int level;
104 int local_size, i;
105 HYPRE_BigInt first_local_row, num_fun;
106 HYPRE_BigInt coarse_size;
107 int coarsen_type;
108 int measure_type;
109 int setup_type;
110 HYPRE_BigInt fine_size;
111 int rest, tms, indx;
112 double size;
113 int not_finished_coarsening = 1;
114 int Setup_err_flag = 0;
115 int seq_threshold = hypre_ParAMGDataSeqThreshold(amg_data);
116 int coarse_threshold = hypre_ParAMGDataMaxCoarseSize(amg_data);
117 int j, k;
118 int num_procs,my_id;
119 int num_threads;
120 int *grid_relax_type = hypre_ParAMGDataGridRelaxType(amg_data);
121 int num_functions = hypre_ParAMGDataNumFunctions(amg_data);
122 int nodal = hypre_ParAMGDataNodal(amg_data);
123 int num_paths = hypre_ParAMGDataNumPaths(amg_data);
124 int agg_num_levels = hypre_ParAMGDataAggNumLevels(amg_data);
125 int *coarse_dof_func = NULL;
126 HYPRE_BigInt *coarse_pnts_global = NULL;
127 HYPRE_BigInt *coarse_pnts_global1 = NULL;
128 int num_cg_sweeps;
129
130 double *max_eig_est = NULL;
131 double *min_eig_est = NULL;
132
133 HYPRE_Solver *smoother = NULL;
134 HYPRE_Solver *smoother_prec = NULL;
135 /*
136 int smooth_type = hypre_ParAMGDataSmoothType(amg_data);
137 int smooth_num_levels = hypre_ParAMGDataSmoothNumLevels(amg_data);
138 int sym;
139 int nlevel;
140 double thresh;
141 double filter;
142 double drop_tol;
143 int max_nz_per_row;
144 char *euclidfile;*/
145
146 int interp_type;
147 int agg_interp_type;
148 int post_interp_type; /* what to do after computing the interpolation matrix
149 0 for nothing, 1 for a Jacobi step */
150
151 /*hypre_ParCSRBlockMatrix *A_H_block;*/
152
153 int block_mode = 0;
154
155 double wall_time; /* for debugging instrumentation */
156
157 MPI_Comm_size(comm, &num_procs);
158 MPI_Comm_rank(comm,&my_id);
159 num_threads = hypre_NumThreads();
160
161 old_num_levels = hypre_ParAMGDataNumLevels(amg_data);
162 max_levels = hypre_ParAMGDataMaxLevels(amg_data);
163 amg_logging = hypre_ParAMGDataLogging(amg_data);
164 amg_print_level = hypre_ParAMGDataPrintLevel(amg_data);
165 coarsen_type = hypre_ParAMGDataCoarsenType(amg_data);
166 measure_type = hypre_ParAMGDataMeasureType(amg_data);
167 setup_type = hypre_ParAMGDataSetupType(amg_data);
168 debug_flag = hypre_ParAMGDataDebugFlag(amg_data);
169 relax_weight = hypre_ParAMGDataRelaxWeight(amg_data);
170 omega = hypre_ParAMGDataOmega(amg_data);
171 dof_func = hypre_ParAMGDataDofFunc(amg_data);
172 /*sym = hypre_ParAMGDataSym(amg_data);
173 nlevel = hypre_ParAMGDataLevel(amg_data);
174 filter = hypre_ParAMGDataFilter(amg_data);
175 thresh = hypre_ParAMGDataThreshold(amg_data);
176 drop_tol = hypre_ParAMGDataDropTol(amg_data);
177 max_nz_per_row = hypre_ParAMGDataMaxNzPerRow(amg_data);
178 euclidfile = hypre_ParAMGDataEuclidFile(amg_data);*/
179 agg_interp_type = hypre_ParAMGDataAggInterpType(amg_data);
180 interp_type = hypre_ParAMGDataInterpType(amg_data);
181 post_interp_type = hypre_ParAMGDataPostInterpType(amg_data);
182 /*IS_type = hypre_ParAMGDataISType(amg_data);
183 num_CR_relax_steps = hypre_ParAMGDataNumCRRelaxSteps(amg_data);
184 CR_rate = hypre_ParAMGDataCRRate(amg_data);
185 CR_use_CG = hypre_ParAMGDataCRUseCG(amg_data);*/
186
187 hypre_ParCSRMatrixSetNumNonzeros(A);
188 hypre_ParCSRMatrixSetDNumNonzeros(A);
189 hypre_ParAMGDataNumVariables(amg_data) = hypre_ParCSRMatrixNumRows(A);
190
191 if (setup_type == 0) return Setup_err_flag;
192
193 S = NULL;
194
195 A_array = hypre_ParAMGDataAArray(amg_data);
196 P_array = hypre_ParAMGDataPArray(amg_data);
197 CF_marker_array = hypre_ParAMGDataCFMarkerArray(amg_data);
198 dof_func_array = hypre_ParAMGDataDofFuncArray(amg_data);
199 local_size = hypre_CSRMatrixNumRows(hypre_ParCSRMatrixDiag(A));
200
201
202 /*A_block_array = hypre_ParAMGDataABlockArray(amg_data);
203 P_block_array = hypre_ParAMGDataPBlockArray(amg_data);*/
204
205 grid_relax_type[3] = hypre_ParAMGDataUserCoarseRelaxType(amg_data);
206
207
208
209 hypre_ParAMGDataBlockMode(amg_data) = block_mode;
210 /* end of systems checks */
211
212
213
214 /*if (A_array || A_block_array || P_array || P_block_array || CF_marker_array || dof_func_array)*/
215 if (A_array || P_array || CF_marker_array || dof_func_array)
216 {
217 for (j = 1; j < old_num_levels; j++)
218 {
219 if (A_array[j])
220 {
221 hypre_ParCSRMatrixDestroy(A_array[j]);
222 A_array[j] = NULL;
223 }
224
225 /*if (A_block_array[j])
226 {
227 hypre_ParCSRBlockMatrixDestroy(A_block_array[j]);
228 A_block_array[j] = NULL;
229 }*/
230
231
232
233 if (dof_func_array[j])
234 {
235 hypre_TFree(dof_func_array[j]);
236 dof_func_array[j] = NULL;
237 }
238 }
239
240 for (j = 0; j < old_num_levels-1; j++)
241 {
242 if (P_array[j])
243 {
244 hypre_ParCSRMatrixDestroy(P_array[j]);
245 P_array[j] = NULL;
246 }
247
248 /*if (P_block_array[j])
249 {
250 hypre_ParCSRBlockMatrixDestroy(P_block_array[j]);
251 P_array[j] = NULL;
252 }*/
253
254 }
255
256/* Special case use of CF_marker_array when old_num_levels == 1
257 requires us to attempt this deallocation every time */
258 if (CF_marker_array[0])
259 {
260 hypre_TFree(CF_marker_array[0]);
261 CF_marker_array[0] = NULL;
262 }
263
264 for (j = 1; j < old_num_levels-1; j++)
265 {
266 if (CF_marker_array[j])
267 {
268 hypre_TFree(CF_marker_array[j]);
269 CF_marker_array[j] = NULL;
270 }
271 }
272 }
273
274 if (A_array == NULL)
275 A_array = hypre_CTAlloc(hypre_ParCSRMatrix*, max_levels);
276 /*if (A_block_array == NULL)
277 A_block_array = hypre_CTAlloc(hypre_ParCSRBlockMatrix*, max_levels);*/
278
279
280 if (P_array == NULL && max_levels > 1)
281 P_array = hypre_CTAlloc(hypre_ParCSRMatrix*, max_levels-1);
282 /*if (P_block_array == NULL && max_levels > 1)
283 P_block_array = hypre_CTAlloc(hypre_ParCSRBlockMatrix*, max_levels-1);*/
284
285
286 if (CF_marker_array == NULL)
287 CF_marker_array = hypre_CTAlloc(int*, max_levels);
288 if (dof_func_array == NULL)
289 dof_func_array = hypre_CTAlloc(int*, max_levels);
290 if (num_functions > 1 && dof_func == NULL)
291 {
292 first_local_row = hypre_ParCSRMatrixFirstRowIndex(A);
293 dof_func = hypre_CTAlloc(int,local_size);
294 num_fun = (HYPRE_BigInt) num_functions;
295 rest = (int)(first_local_row-((first_local_row/num_fun)*num_fun));
296 indx = num_functions-rest;
297 if (rest == 0) indx = 0;
298 k = num_functions - 1;
299 for (j = indx-1; j > -1; j--)
300 dof_func[j] = k--;
301 tms = local_size/num_functions;
302 if (tms*num_functions+indx > local_size) tms--;
303 for (j=0; j < tms; j++)
304 {
305 for (k=0; k < num_functions; k++)
306 dof_func[indx++] = k;
307 }
308 k = 0;
309 while (indx < local_size)
310 dof_func[indx++] = k++;
311 hypre_ParAMGDataDofFunc(amg_data) = dof_func;
312 }
313
314 A_array[0] = A;
315
316
317
318 dof_func_array[0] = dof_func;
319 hypre_ParAMGDataCFMarkerArray(amg_data) = CF_marker_array;
320 hypre_ParAMGDataDofFuncArray(amg_data) = dof_func_array;
321 hypre_ParAMGDataAArray(amg_data) = A_array;
322 hypre_ParAMGDataPArray(amg_data) = P_array;
323 hypre_ParAMGDataRArray(amg_data) = P_array;
324
325
326 Vtemp = hypre_ParAMGDataVtemp(amg_data);
327
328 if (Vtemp != NULL)
329 {
330 hypre_ParVectorDestroy(Vtemp);
331 Vtemp = NULL;
332 }
333
334 Vtemp = hypre_ParVectorCreate(hypre_ParCSRMatrixComm(A_array[0]),
335 hypre_ParCSRMatrixGlobalNumRows(A_array[0]),
336 hypre_ParCSRMatrixRowStarts(A_array[0]));
337 hypre_ParVectorInitialize(Vtemp);
338 hypre_ParVectorSetPartitioningOwner(Vtemp,0);
339 hypre_ParAMGDataVtemp(amg_data) = Vtemp;
340
341 /*if ((smooth_num_levels > 0 && smooth_type > 9)
342 || relax_weight[0] < 0 || omega[0] < 0 ||
343 hypre_ParAMGDataSchwarzRlxWeight(amg_data) < 0)*/
344 if ( relax_weight[0] < 0 || omega[0] < 0 )
345 {
346 Ptemp = hypre_ParVectorCreate(hypre_ParCSRMatrixComm(A_array[0]),
347 hypre_ParCSRMatrixGlobalNumRows(A_array[0]),
348 hypre_ParCSRMatrixRowStarts(A_array[0]));
349 hypre_ParVectorInitialize(Ptemp);
350 hypre_ParVectorSetPartitioningOwner(Ptemp,0);
351 hypre_ParAMGDataPtemp(amg_data) = Ptemp;
352 Rtemp = hypre_ParVectorCreate(hypre_ParCSRMatrixComm(A_array[0]),
353 hypre_ParCSRMatrixGlobalNumRows(A_array[0]),
354 hypre_ParCSRMatrixRowStarts(A_array[0]));
355 hypre_ParVectorInitialize(Rtemp);
356 hypre_ParVectorSetPartitioningOwner(Rtemp,0);
357 hypre_ParAMGDataRtemp(amg_data) = Rtemp;
358 }
359 /*if ((smooth_num_levels > 0 && smooth_type > 6)
360 || relax_weight[0] < 0 || omega[0] < 0 ||
361 hypre_ParAMGDataSchwarzRlxWeight(amg_data) < 0)*/
362 if ( relax_weight[0] < 0 || omega[0] < 0 )
363 {
364 Ztemp = hypre_ParVectorCreate(hypre_ParCSRMatrixComm(A_array[0]),
365 hypre_ParCSRMatrixGlobalNumRows(A_array[0]),
366 hypre_ParCSRMatrixRowStarts(A_array[0]));
367 hypre_ParVectorInitialize(Ztemp);
368 hypre_ParVectorSetPartitioningOwner(Ztemp,0);
369 hypre_ParAMGDataZtemp(amg_data) = Ztemp;
370 }
371 F_array = hypre_ParAMGDataFArray(amg_data);
372 U_array = hypre_ParAMGDataUArray(amg_data);
373
374
375
376 if (F_array != NULL || U_array != NULL)
377 {
378 for (j = 1; j < old_num_levels; j++)
379 {
380 if (F_array[j] != NULL)
381 {
382 hypre_ParVectorDestroy(F_array[j]);
383 F_array[j] = NULL;
384 }
385 if (U_array[j] != NULL)
386 {
387 hypre_ParVectorDestroy(U_array[j]);
388 U_array[j] = NULL;
389 }
390 }
391 }
392
393 if (F_array == NULL)
394 F_array = hypre_CTAlloc(hypre_ParVector*, max_levels);
395 if (U_array == NULL)
396 U_array = hypre_CTAlloc(hypre_ParVector*, max_levels);
397
398 F_array[0] = f;
399 U_array[0] = u;
400
401 hypre_ParAMGDataFArray(amg_data) = F_array;
402 hypre_ParAMGDataUArray(amg_data) = U_array;
403
404 /*----------------------------------------------------------
405 * Initialize hypre_ParAMGData
406 *----------------------------------------------------------*/
407
408 not_finished_coarsening = 1;
409 level = 0;
410
411 strong_threshold = hypre_ParAMGDataStrongThreshold(amg_data);
412 max_row_sum = hypre_ParAMGDataMaxRowSum(amg_data);
413 trunc_factor = hypre_ParAMGDataTruncFactor(amg_data);
414 agg_trunc_factor = hypre_ParAMGDataAggTruncFactor(amg_data);
415 P_max_elmts = hypre_ParAMGDataPMaxElmts(amg_data);
416 P_max1 = hypre_ParAMGDataPMax1(amg_data);
417 P_max2 = hypre_ParAMGDataPMax2(amg_data);
418 agg_P_max_elmts = hypre_ParAMGDataAggPMaxElmts(amg_data);
419 jacobi_trunc_threshold = hypre_ParAMGDataJacobiTruncThreshold(amg_data);
420 S_commpkg_switch = hypre_ParAMGDataSCommPkgSwitch(amg_data);
421 /*if (smooth_num_levels > level)
422 {
423 smoother = hypre_CTAlloc(HYPRE_Solver, smooth_num_levels);
424 hypre_ParAMGDataSmoother(amg_data) = smoother;
425 }*/
426
427 /*-----------------------------------------------------
428 * Enter Coarsening Loop
429 *-----------------------------------------------------*/
430
431 while (not_finished_coarsening)
432 {
433
434#if MPIP_SETUP_ON
435 MPI_Pcontrol(2);
436 MPI_Pcontrol(1);
437#endif
438
439 {
440 fine_size = hypre_ParCSRMatrixGlobalNumRows(A_array[level]);
441 }
442
443 if (level > 0)
444 {
445 {
446 F_array[level] =
447 hypre_ParVectorCreate(hypre_ParCSRMatrixComm(A_array[level]),
448 hypre_ParCSRMatrixGlobalNumRows(A_array[level]),
449 hypre_ParCSRMatrixRowStarts(A_array[level]));
450 hypre_ParVectorInitialize(F_array[level]);
451 hypre_ParVectorSetPartitioningOwner(F_array[level],0);
452
453 U_array[level] =
454 hypre_ParVectorCreate(hypre_ParCSRMatrixComm(A_array[level]),
455 hypre_ParCSRMatrixGlobalNumRows(A_array[level]),
456 hypre_ParCSRMatrixRowStarts(A_array[level]));
457 hypre_ParVectorInitialize(U_array[level]);
458 hypre_ParVectorSetPartitioningOwner(U_array[level],0);
459 }
460
461 }
462
463 /*-------------------------------------------------------------
464 * Select coarse-grid points on 'level' : returns CF_marker
465 * for the level. Returns strength matrix, S
466 *--------------------------------------------------------------*/
467
468 if (debug_flag==1) wall_time = time_getWallclockSeconds();
469 if (debug_flag==3)
470 {
471 printf("\n ===== Proc = %d Level = %d =====\n",
472 my_id, level);
473 fflush(NULL);
474 }
475
476 /*if ((smooth_type == 6 || smooth_type == 16) && smooth_num_levels > level)
477 {
478
479 schwarz_relax_wt = hypre_ParAMGDataSchwarzRlxWeight(amg_data);
480 HYPRE_SchwarzCreate(&smoother[level]);
481 HYPRE_SchwarzSetNumFunctions(smoother[level],num_functions);
482 HYPRE_SchwarzSetVariant(smoother[level],
483 hypre_ParAMGDataVariant(amg_data));
484 HYPRE_SchwarzSetOverlap(smoother[level],
485 hypre_ParAMGDataOverlap(amg_data));
486 HYPRE_SchwarzSetDomainType(smoother[level],
487 hypre_ParAMGDataDomainType(amg_data));
488 if (schwarz_relax_wt > 0)
489 HYPRE_SchwarzSetRelaxWeight(smoother[level],schwarz_relax_wt);
490 HYPRE_SchwarzSetup(smoother[level],
491 (HYPRE_ParCSRMatrix) A_array[level],
492 (HYPRE_ParVector) f,
493 (HYPRE_ParVector) u);
494 }*/
495
496 if (max_levels > 1)
497 {
498 {
499 local_num_vars =
500 hypre_CSRMatrixNumRows(hypre_ParCSRMatrixDiag(A_array[level]));
501 }
502
503
504 /**** Get the Strength Matrix ****/
505
506 if (hypre_ParAMGDataGSMG(amg_data) == 0)
507 {
508 if (nodal) /* if we are solving systems and
509 not using the unknown approach then we need to
510 convert A to a nodal matrix - values that represent the
511 blocks - before getting the strength matrix*/
512 {
513
514 {
515 hypre_BoomerAMGCreateNodalA(A_array[level],num_functions,
516 dof_func_array[level], abs(nodal), &AN);
517 }
518
519 /* dof array not needed for creating S because we pass in that
520 the number of functions is 1 */
521 if (nodal == 3 || nodal == -3 || nodal == 7) /* option 3 may have negative entries in AN -
522 all other options are pos numbers only */
523 hypre_BoomerAMGCreateS(AN, strong_threshold, max_row_sum,
524 1, NULL,&SN);
525 else
526 hypre_BoomerAMGCreateSabs(AN, strong_threshold, max_row_sum,
527 1, NULL,&SN);
528#if 1
529/* TEMP - to compare with serial */
530 hypre_BoomerAMGCreateS(AN, strong_threshold, max_row_sum,
531 1, NULL,&SN);
532#endif
533
534
535 col_offd_S_to_A = NULL;
536 col_offd_SN_to_AN = NULL;
537 if (strong_threshold > S_commpkg_switch)
538 hypre_BoomerAMGCreateSCommPkg(AN,SN,&col_offd_SN_to_AN);
539 }
540 else
541 {
542 hypre_BoomerAMGCreateS(A_array[level],
543 strong_threshold, max_row_sum,
544 num_functions, dof_func_array[level],&S);
545 col_offd_S_to_A = NULL;
546 if (strong_threshold > S_commpkg_switch)
547 hypre_BoomerAMGCreateSCommPkg(A_array[level],S,
548 &col_offd_S_to_A);
549 }
550 }
551 /*else
552 {
553 hypre_BoomerAMGCreateSmoothDirs(amg_data, A_array[level],
554 SmoothVecs, strong_threshold,
555 num_functions, dof_func_array[level], &S);
556 }*/
557
558
559 /**** Do the appropriate coarsening ****/
560
561
562 if (coarsen_type == 6) /* falgout */
563 {
564 if (nodal == 0) /* nonsystems or unknown approach for systems*/
565 {
566 hypre_BoomerAMGCoarsenFalgout(S, A_array[level], measure_type,
567 debug_flag, &CF_marker);
568
569 if (level < agg_num_levels && (agg_interp_type < 30 || agg_interp_type > 32))
570 {
571 /* set num_functions 1 in CoarseParms, since coarse_dof_func
572 is not needed here */
573 hypre_BoomerAMGCoarseParms(comm, local_num_vars,
574 1, dof_func_array[level], CF_marker,
575 &coarse_dof_func,&coarse_pnts_global1);
576 hypre_BoomerAMGCreate2ndS (S, CF_marker, num_paths,
577 coarse_pnts_global1, &S2);
578 hypre_BoomerAMGCoarsenFalgout(S2, S2, measure_type,
579 debug_flag, &CFN_marker);
580 hypre_ParCSRMatrixDestroy(S2);
581 hypre_BoomerAMGCorrectCFMarker (CF_marker, local_num_vars, CFN_marker);
582 hypre_TFree(coarse_pnts_global1);
583 hypre_TFree(CFN_marker);
584 }
585 }
586 else if (nodal < 0 || block_mode ) /* nodal interpolation
587 or if nodal < 0 then we
588 build interpolation normally
589 using the nodal matrix */
590 {
591 hypre_BoomerAMGCoarsenFalgout(SN, SN, measure_type,
592 debug_flag, &CF_marker);
593 }
594
595 else if (nodal > 0) /* nodal = 1,2,3: here we convert our nodal coarsening
596 so that we can perform regular interpolation */
597 {
598 hypre_BoomerAMGCoarsenFalgout(SN, SN, measure_type,
599 debug_flag, &CFN_marker);
600
601 if (level < agg_num_levels)
602 {
603 /* set num_functions 1 in CoarseParms, since coarse_dof_func
604 is not needed here */
605 hypre_BoomerAMGCoarseParms(comm, local_num_vars,
606 1, dof_func_array[level], CFN_marker,
607 &coarse_dof_func,&coarse_pnts_global1);
608 hypre_BoomerAMGCreate2ndS (SN, CFN_marker, num_paths,
609 coarse_pnts_global1, &S2);
610 hypre_BoomerAMGCoarsenFalgout(S2, S2, measure_type,
611 debug_flag, &CF2_marker);
612 hypre_ParCSRMatrixDestroy(S2);
613 hypre_BoomerAMGCorrectCFMarker (CFN_marker, local_num_vars, CF2_marker);
614 hypre_TFree(coarse_pnts_global1);
615 hypre_TFree(CF2_marker);
616 }
617
618 col_offd_S_to_A = NULL;
619 hypre_BoomerAMGCreateScalarCFS(SN, CFN_marker, col_offd_SN_to_AN,
620 num_functions, nodal, 0, NULL, &CF_marker,
621 &col_offd_S_to_A, &S);
622 if (col_offd_SN_to_AN == NULL)
623 col_offd_S_to_A = NULL;
624 hypre_TFree(CFN_marker);
625 hypre_TFree(col_offd_SN_to_AN);
626 hypre_ParCSRMatrixDestroy(AN);
627 hypre_ParCSRMatrixDestroy(SN);
628
629 }
630 }
631 else if (coarsen_type == 7) /* cljp1 */
632 {
633 if (nodal == 0) /* nonsystems or unknown approach for systems*/
634 {
635 hypre_BoomerAMGCoarsen(S, A_array[level], 2,
636 debug_flag, &CF_marker);
637 if (level < agg_num_levels && (agg_interp_type < 30 || agg_interp_type > 32))
638 {
639 /* set num_functions 1 in CoarseParms, since coarse_dof_func
640 is not needed here */
641 hypre_BoomerAMGCoarseParms(comm, local_num_vars,
642 1, dof_func_array[level], CF_marker,
643 &coarse_dof_func,&coarse_pnts_global1);
644 hypre_BoomerAMGCreate2ndS (S, CF_marker, num_paths,
645 coarse_pnts_global1, &S2);
646 hypre_BoomerAMGCoarsen(S2, S2, 2, debug_flag, &CFN_marker);
647 hypre_ParCSRMatrixDestroy(S2);
648 hypre_BoomerAMGCorrectCFMarker (CF_marker, local_num_vars, CFN_marker);
649 hypre_TFree(CFN_marker);
650 hypre_TFree(coarse_pnts_global1);
651 }
652 }
653 else if (nodal < 0 || block_mode ) /* nodal interpolation
654 or if nodal < 0 then we
655 build interpolation normally
656 using the nodal matrix */
657 {
658 hypre_BoomerAMGCoarsen(SN, SN, 2, debug_flag, &CF_marker);
659 }
660
661 else if (nodal > 0) /* nodal = 1,2,3: here we convert our nodal coarsening
662 so that we can perform regualr interpolation */
663 {
664 hypre_BoomerAMGCoarsen(SN, SN, 2, debug_flag, &CFN_marker);
665
666 col_offd_S_to_A = NULL;
667 hypre_BoomerAMGCreateScalarCFS(SN, CFN_marker, col_offd_SN_to_AN,
668 num_functions, nodal, 0, NULL, &CF_marker,
669 &col_offd_S_to_A, &S);
670 if (col_offd_SN_to_AN == NULL)
671 col_offd_S_to_A = NULL;
672 hypre_TFree(CFN_marker);
673 hypre_TFree(col_offd_SN_to_AN);
674
675 hypre_ParCSRMatrixDestroy(AN);
676 hypre_ParCSRMatrixDestroy(SN);
677
678 }
679 }
680
681 else if (coarsen_type == 8) /* pmis */
682 {
683 if (nodal == 0) /* nonsystems or unknown approach for systems*/
684 {
685 hypre_BoomerAMGCoarsenPMIS(S, A_array[level], 0,
686 debug_flag, &CF_marker);
687
688 if (level < agg_num_levels && (agg_interp_type < 30 || agg_interp_type > 32))
689 {
690 /* set num_functions 1 in CoarseParms, since coarse_dof_func
691 is not needed here */
692 hypre_BoomerAMGCoarseParms(comm, local_num_vars,
693 1, dof_func_array[level], CF_marker,
694 &coarse_dof_func,&coarse_pnts_global1);
695 hypre_BoomerAMGCreate2ndS (S, CF_marker, num_paths,
696 coarse_pnts_global1, &S2);
697 hypre_BoomerAMGCoarsenPMIS(S2, S2, 0,
698 debug_flag, &CFN_marker);
699 hypre_ParCSRMatrixDestroy(S2);
700 hypre_BoomerAMGCorrectCFMarker (CF_marker, local_num_vars, CFN_marker);
701 hypre_TFree(CFN_marker);
702 hypre_TFree(coarse_pnts_global1);
703 }
704 }
705 else if (nodal < 0 || block_mode ) /* nodal interpolation
706 or if nodal < 0 then we
707 build interpolation normally
708 using the nodal matrix */
709 {
710 hypre_BoomerAMGCoarsenPMIS(SN, SN, 0,
711 debug_flag, &CF_marker);
712 }
713 else if (nodal > 0) /* nodal = 1,2,3: here we convert our nodal coarsening
714 so that we can perform regualr interpolation */
715 {
716 hypre_BoomerAMGCoarsenPMIS(SN, SN, 0,
717 debug_flag, &CFN_marker);
718 /*hypre_BoomerAMGCreateScalarCFS(SN, CFN_marker,
719 num_functions, nodal, 0, NULL, &CF_marker, &S);*/
720 col_offd_S_to_A = NULL;
721 hypre_BoomerAMGCreateScalarCFS(SN, CFN_marker, col_offd_SN_to_AN,
722 num_functions, nodal, 0, NULL, &CF_marker,
723 &col_offd_S_to_A, &S);
724 if (col_offd_SN_to_AN == NULL)
725 col_offd_S_to_A = NULL;
726 hypre_TFree(CFN_marker);
727 hypre_TFree(col_offd_SN_to_AN);
728 hypre_ParCSRMatrixDestroy(AN);
729 hypre_ParCSRMatrixDestroy(SN);
730 }
731 }
732
733 else if (coarsen_type == 9) /* pmis1 */
734 {
735 if (nodal == 0) /* nonsystems or unknown approach for systems*/
736 {
737 hypre_BoomerAMGCoarsenPMIS(S, A_array[level], 2,
738 debug_flag, &CF_marker);
739 if (level < agg_num_levels && (agg_interp_type < 30 || agg_interp_type > 32))
740 {
741 /* set num_functions 1 in CoarseParms, since coarse_dof_func
742 is not needed here */
743 hypre_BoomerAMGCoarseParms(comm, local_num_vars,
744 1, dof_func_array[level], CF_marker,
745 &coarse_dof_func,&coarse_pnts_global1);
746 hypre_BoomerAMGCreate2ndS (S, CF_marker, num_paths,
747 coarse_pnts_global1, &S2);
748 hypre_BoomerAMGCoarsenPMIS(S2, S2, 2,
749 debug_flag, &CFN_marker);
750 hypre_ParCSRMatrixDestroy(S2);
751 hypre_BoomerAMGCorrectCFMarker (CF_marker, local_num_vars, CFN_marker);
752 hypre_TFree(CFN_marker);
753 hypre_TFree(coarse_pnts_global1);
754 }
755 }
756
757 else if (nodal < 0 || block_mode ) /* nodal interpolation
758 or if nodal < 0 then we
759 build interpolation normally
760 using the nodal matrix */
761 {
762 hypre_BoomerAMGCoarsenPMIS(SN, SN, 2,
763 debug_flag, &CF_marker);
764
765 }
766 else if (nodal > 0) /* nodal = 1,2,3: here we convert our nodal coarsening
767 so that we can perform regualr interpolation */
768 {
769 hypre_BoomerAMGCoarsenPMIS(SN, SN, 2,
770 debug_flag, &CFN_marker);
771 /*hypre_BoomerAMGCreateScalarCFS(SN, CFN_marker,
772 num_functions, nodal, 0, NULL, &CF_marker, &S);*/
773 col_offd_S_to_A = NULL;
774 hypre_BoomerAMGCreateScalarCFS(SN, CFN_marker, col_offd_SN_to_AN,
775 num_functions, nodal, 0, NULL, &CF_marker,
776 &col_offd_S_to_A, &S);
777 if (col_offd_SN_to_AN == NULL)
778 col_offd_S_to_A = NULL;
779 hypre_TFree(CFN_marker);
780 hypre_TFree(col_offd_SN_to_AN);
781 hypre_ParCSRMatrixDestroy(AN);
782 hypre_ParCSRMatrixDestroy(SN);
783 }
784
785 }
786 else if (coarsen_type == 10) /*hmis */
787 {
788 if (nodal == 0) /* nonsystems or unknown approach for systems*/
789 {
790 hypre_BoomerAMGCoarsenHMIS(S, A_array[level], measure_type,
791 debug_flag, &CF_marker);
792 if (level < agg_num_levels && (agg_interp_type < 30 || agg_interp_type > 32))
793 {
794 /* set num_functions 1 in CoarseParms, since coarse_dof_func
795 is not needed here */
796 hypre_BoomerAMGCoarseParms(comm, local_num_vars,
797 1, dof_func_array[level], CF_marker,
798 &coarse_dof_func,&coarse_pnts_global1);
799 hypre_BoomerAMGCreate2ndS (S, CF_marker, num_paths,
800 coarse_pnts_global1, &S2);
801 hypre_BoomerAMGCoarsenHMIS(S2, S2, measure_type,
802 debug_flag, &CFN_marker);
803 hypre_ParCSRMatrixDestroy(S2);
804 hypre_BoomerAMGCorrectCFMarker (CF_marker, local_num_vars, CFN_marker);
805 hypre_TFree(CFN_marker);
806 hypre_TFree(coarse_pnts_global1);
807 }
808 }
809 else if (nodal < 0 || block_mode ) /* nodal interpolation
810 or if nodal < 0 then we
811 build interpolation normally
812 using the nodal matrix */
813 {
814 hypre_BoomerAMGCoarsenHMIS(SN, SN, measure_type,
815 debug_flag, &CF_marker);
816 }
817 else if (nodal > 0) /* nodal = 1,2,3: here we convert our nodal coarsening
818 so that we can perform regular interpolation */
819 {
820 hypre_BoomerAMGCoarsenHMIS(SN, SN, measure_type,
821 debug_flag, &CFN_marker);
822 /*hypre_BoomerAMGCreateScalarCFS(SN, CFN_marker,
823 num_functions, nodal, 0, NULL, &CF_marker, &S);*/
824 if (level < agg_num_levels)
825 {
826 /* set num_functions 1 in CoarseParms, since coarse_dof_func
827 is not needed here */
828 hypre_BoomerAMGCoarseParms(comm, local_num_vars,
829 1, dof_func_array[level], CFN_marker,
830 &coarse_dof_func,&coarse_pnts_global1);
831 hypre_BoomerAMGCreate2ndS (SN, CFN_marker, num_paths,
832 coarse_pnts_global1, &S2);
833 hypre_BoomerAMGCoarsenHMIS(S2, S2, measure_type,
834 debug_flag, &CF2_marker);
835 hypre_ParCSRMatrixDestroy(S2);
836 hypre_BoomerAMGCorrectCFMarker (CFN_marker, local_num_vars, CF2_marker);
837 hypre_TFree(CF2_marker);
838 hypre_TFree(coarse_pnts_global1);
839 }
840 col_offd_S_to_A = NULL;
841 hypre_BoomerAMGCreateScalarCFS(SN, CFN_marker, col_offd_SN_to_AN,
842 num_functions, nodal, 0, NULL, &CF_marker,
843 &col_offd_S_to_A, &S);
844 if (col_offd_SN_to_AN == NULL)
845 col_offd_S_to_A = NULL;
846 hypre_TFree(CFN_marker);
847 hypre_TFree(col_offd_SN_to_AN);
848 hypre_ParCSRMatrixDestroy(AN);
849 hypre_ParCSRMatrixDestroy(SN);
850 }
851
852 }
853
854 else if (coarsen_type) /* ruge, ruge1p, ruge2b, ruge3, ruge3c, ruge1x */
855 {
856 if (nodal == 0) /* nonsystems or unknown approach for systems*/
857 {
858 hypre_BoomerAMGCoarsenRuge(S, A_array[level],
859 measure_type, coarsen_type, debug_flag,
860 &CF_marker);
861 if (level < agg_num_levels && (agg_interp_type < 30 || agg_interp_type > 32))
862 {
863 /* set num_functions 1 in CoarseParms, since coarse_dof_func
864 is not needed here */
865 hypre_BoomerAMGCoarseParms(comm, local_num_vars,
866 1, dof_func_array[level], CF_marker,
867 &coarse_dof_func,&coarse_pnts_global1);
868 hypre_BoomerAMGCreate2ndS (S, CF_marker, num_paths,
869 coarse_pnts_global1, &S2);
870 hypre_BoomerAMGCoarsenRuge(S2, S2, measure_type, coarsen_type,
871 debug_flag, &CFN_marker);
872 hypre_ParCSRMatrixDestroy(S2);
873 hypre_BoomerAMGCorrectCFMarker (CF_marker, local_num_vars, CFN_marker);
874 hypre_TFree(CFN_marker);
875 hypre_TFree(coarse_pnts_global1);
876 }
877 }
878 else if (nodal < 0 || block_mode ) /* nodal interpolation
879 or if nodal < 0 then we
880 build interpolation normally
881 using the nodal matrix */
882 {
883 hypre_BoomerAMGCoarsenRuge(SN, SN,
884 measure_type, coarsen_type, debug_flag,
885 &CF_marker);
886 }
887 else if (nodal > 0) /* nodal = 1,2,3: here we convert our nodal coarsening
888 so that we can perform regualr interpolation */
889 {
890 hypre_BoomerAMGCoarsenRuge(SN, SN,
891 measure_type, coarsen_type, debug_flag,
892 &CFN_marker);
893 /*hypre_BoomerAMGCreateScalarCFS(SN, CFN_marker,
894 num_functions, 0, &CF_marker, &S);*/
895 /*hypre_BoomerAMGCreateScalarCFS(SN, CFN_marker,
896 num_functions, nodal, 0, NULL, &CF_marker, &S);*/
897 col_offd_S_to_A = NULL;
898 hypre_BoomerAMGCreateScalarCFS(SN, CFN_marker, col_offd_SN_to_AN,
899 num_functions, nodal, 0, NULL, &CF_marker,
900 &col_offd_S_to_A, &S);
901 if (col_offd_SN_to_AN == NULL)
902 col_offd_S_to_A = NULL;
903 hypre_TFree(CFN_marker);
904 hypre_TFree(col_offd_SN_to_AN);
905 hypre_ParCSRMatrixDestroy(AN);
906 hypre_ParCSRMatrixDestroy(SN);
907 }
908 }
909 else /* coarsen_type = 0 (or anything negative) - cljp */
910 {
911
912 if (nodal == 0) /* nonsystems or unknown approach for systems*/
913 {
914
915 hypre_BoomerAMGCoarsen(S, A_array[level], 0,
916 debug_flag, &CF_marker);
917 if (level < agg_num_levels && (agg_interp_type < 30 || agg_interp_type > 32))
918 {
919 /* set num_functions 1 in CoarseParms, since coarse_dof_func
920 is not needed here */
921 hypre_BoomerAMGCoarseParms(comm, local_num_vars,
922 1, dof_func_array[level], CF_marker,
923 &coarse_dof_func,&coarse_pnts_global1);
924 hypre_BoomerAMGCreate2ndS (S, CF_marker, num_paths,
925 coarse_pnts_global1, &S2);
926 hypre_BoomerAMGCoarsen(S2, S2, 0, debug_flag, &CFN_marker);
927 hypre_ParCSRMatrixDestroy(S2);
928 hypre_BoomerAMGCorrectCFMarker (CF_marker, local_num_vars, CFN_marker);
929 hypre_TFree(CFN_marker);
930 hypre_TFree(coarse_pnts_global1);
931 }
932 }
933 else if (nodal < 0 || block_mode ) /* nodal interpolation
934 or if nodal < 0 then we
935 build interpolation normally
936 using the nodal matrix */
937 {
938 hypre_BoomerAMGCoarsen(SN, SN, 0, debug_flag, &CF_marker);
939 }
940 else if (nodal > 0) /* nodal = 1,2,3: here we convert our nodal coarsening
941 so that we can perform regualr interpolation */
942 {
943 hypre_BoomerAMGCoarsen(SN, SN, 0, debug_flag, &CFN_marker);
944
945 col_offd_S_to_A = NULL;
946 hypre_BoomerAMGCreateScalarCFS(SN, CFN_marker, col_offd_SN_to_AN,
947 num_functions, nodal, 0, NULL, &CF_marker,
948 &col_offd_S_to_A, &S);
949 if (col_offd_SN_to_AN == NULL)
950 col_offd_S_to_A = NULL;
951 hypre_TFree(CFN_marker);
952 hypre_TFree(col_offd_SN_to_AN);
953
954 hypre_ParCSRMatrixDestroy(AN);
955 hypre_ParCSRMatrixDestroy(SN);
956 }
957 }
958
959 /** end of coarsen type options **/
960
961 /* store the CF array */
962 CF_marker_array[level] = CF_marker;
963
964
965 if (relax_weight[level] == 0.0)
966 {
967 hypre_ParCSRMatrixScaledNorm(A_array[level], &relax_weight[level]);
968 if (relax_weight[level] != 0.0)
969 relax_weight[level] = 4.0/3.0/relax_weight[level];
970 else
971 printf (" Warning ! Matrix norm is zero !!!");
972 }
973 if (relax_weight[level] < 0 )
974 {
975 num_cg_sweeps = (int) (-relax_weight[level]);
976 hypre_BoomerAMGCGRelaxWt(amg_data, level, num_cg_sweeps,
977 &relax_weight[level]);
978 }
979 if (omega[level] < 0 )
980 {
981 num_cg_sweeps = (int) (-omega[level]);
982 hypre_BoomerAMGCGRelaxWt(amg_data, level, num_cg_sweeps,
983 &omega[level]);
984 }
985 /*if (schwarz_relax_wt < 0 )
986 {
987 num_cg_sweeps = (int) (-schwarz_relax_wt);
988 hypre_BoomerAMGCGRelaxWt(amg_data, level, num_cg_sweeps,
989 &schwarz_relax_wt);*/
990 /*printf (" schwarz weight %f \n", schwarz_relax_wt);*/
991 /*HYPRE_SchwarzSetRelaxWeight(smoother[level], schwarz_relax_wt);
992 if (hypre_ParAMGDataVariant(amg_data) > 0)
993 {
994 local_size = hypre_CSRMatrixNumRows
995 (hypre_ParCSRMatrixDiag(A_array[level]));
996 hypre_SchwarzReScale(smoother[level], local_size,
997 schwarz_relax_wt);
998 }
999 schwarz_relax_wt = 1;
1000 }*/
1001 if (debug_flag==1)
1002 {
1003 wall_time = time_getWallclockSeconds() - wall_time;
1004 printf("Proc = %d Level = %d Coarsen Time = %f\n",
1005 my_id,level, wall_time);
1006 fflush(NULL);
1007 }
1008
1009 /**** Get the coarse parameters ****/
1010 if (nodal < 0 || block_mode )
1011 {
1012 /* here we will determine interpolation using a nodal matrix */
1013 hypre_BoomerAMGCoarseParms(comm,
1014 hypre_CSRMatrixNumRows(hypre_ParCSRMatrixDiag(AN)),
1015 1, NULL, CF_marker, NULL, &coarse_pnts_global);
1016 }
1017 else
1018 {
1019 hypre_BoomerAMGCoarseParms(comm, local_num_vars,
1020 num_functions, dof_func_array[level], CF_marker,
1021 &coarse_dof_func,&coarse_pnts_global);
1022 }
1023 dof_func_array[level+1] = NULL;
1024 if (num_functions > 1 && nodal > -1 && (!block_mode) )
1025 dof_func_array[level+1] = coarse_dof_func;
1026
1027#ifdef HYPRE_NO_GLOBAL_PARTITION
1028 if (my_id == (num_procs -1)) coarse_size = coarse_pnts_global[1];
1029 MPI_Bcast(&coarse_size, 1, MPI_HYPRE_BIG_INT, num_procs-1, comm);
1030#else
1031 coarse_size = coarse_pnts_global[num_procs];
1032#endif
1033
1034 }
1035 else
1036 {
1037 S = NULL;
1038 coarse_pnts_global = NULL;
1039 CF_marker = hypre_CTAlloc(int, local_size );
1040 for (i=0; i < local_size ; i++)
1041 CF_marker[i] = 1;
1042 /*CF_marker_array = hypre_CTAlloc(int*, 1);*/
1043 CF_marker_array[level] = CF_marker;
1044 coarse_size = fine_size;
1045 }
1046
1047 /* if no coarse-grid, stop coarsening, and set the
1048 * coarsest solve to be a single sweep of Jacobi */
1049 if ((coarse_size == 0) ||
1050 (coarse_size == fine_size))
1051 {
1052 int *num_grid_sweeps =
1053 hypre_ParAMGDataNumGridSweeps(amg_data);
1054 int **grid_relax_points =
1055 hypre_ParAMGDataGridRelaxPoints(amg_data);
1056 if (grid_relax_type[3] == 9)
1057 {
1058 grid_relax_type[3] = grid_relax_type[0];
1059 num_grid_sweeps[3] = 1;
1060 if (grid_relax_points) grid_relax_points[3][0] = 0;
1061 }
1062 if (S)
1063 hypre_ParCSRMatrixDestroy(S);
1064 hypre_TFree(coarse_pnts_global);
1065 if (level > 0)
1066 {
1067 /* note special case treatment of CF_marker is necessary
1068 * to do CF relaxation correctly when num_levels = 1 */
1069 hypre_TFree(CF_marker_array[level]);
1070 hypre_ParVectorDestroy(F_array[level]);
1071 hypre_ParVectorDestroy(U_array[level]);
1072 }
1073
1074#if MPIP_SETUP_ON
1075 MPI_Pcontrol(3);
1076 MPI_Pcontrol(0);
1077#endif
1078
1079 break;
1080 }
1081
1082 /*-------------------------------------------------------------
1083 * Build prolongation matrix, P, and place in P_array[level]
1084 *--------------------------------------------------------------*/
1085
1086 if (debug_flag==1) wall_time = time_getWallclockSeconds();
1087
1088 if (hypre_ParAMGDataAggInterpType(amg_data) == 30 && (level < agg_num_levels
1089 && (!block_mode)))
1090 {
1091 hypre_BoomerAMGBuildExtPIInterp(A_array[level], CF_marker_array[level],
1092 S, coarse_pnts_global, num_functions, dof_func_array[level],
1093 debug_flag, 0, P_max1, col_offd_S_to_A, &P1);
1094 hypre_TFree(col_offd_S_to_A);
1095 /*if (my_id == 0 && debug_flag==1)
1096 {
1097 wall_time = time_getWallclockSeconds() - wall_time;
1098 printf("Proc = %d Level = %d Interp P1 Time = %f\n",
1099 my_id,level, wall_time);
1100 fflush(NULL);
1101 }
1102 if (my_id == 0 && debug_flag==1) wall_time = time_getWallclockSeconds();
1103 hypre_BoomerAMGCreate2ndS (S, CF_marker, num_paths,
1104 coarse_pnts_global, &S2);
1105 if (my_id == 0 && debug_flag==1)
1106 {
1107 wall_time = time_getWallclockSeconds() - wall_time;
1108 printf("Proc = %d Level = %d Create S2 Time = %f\n",
1109 my_id,level, wall_time);
1110 fflush(NULL);
1111 }
1112 if (my_id == 0 && debug_flag==1) wall_time = time_getWallclockSeconds();*/
1113 hypre_BoomerAMGCreate2ndS (S, CF_marker, num_paths,
1114 coarse_pnts_global, &S2);
1115 if (coarsen_type == 10)
1116 hypre_BoomerAMGCoarsenHMIS(S2, S2, measure_type,
1117 debug_flag, &CFN_marker);
1118 else if (coarsen_type == 8)
1119 hypre_BoomerAMGCoarsenPMIS(S2, S2, 3,
1120 debug_flag, &CFN_marker);
1121 else if (coarsen_type == 9)
1122 hypre_BoomerAMGCoarsenHMIS(S2, S2, measure_type+3,
1123 debug_flag, &CFN_marker);
1124 else if (coarsen_type == 11)
1125 hypre_BoomerAMGCoarsenPMIS(S2, S2, 4,
1126 debug_flag, &CFN_marker);
1127 else if (coarsen_type == 6)
1128 hypre_BoomerAMGCoarsenFalgout(S2, S2, measure_type,
1129 debug_flag, &CFN_marker);
1130 else if (coarsen_type == 0)
1131 hypre_BoomerAMGCoarsen(S2, S2, 0, debug_flag, &CFN_marker);
1132 hypre_ParCSRMatrixDestroy(S2);
1133 hypre_BoomerAMGCorrectCFMarker2 (CF_marker, local_num_vars, CFN_marker);
1134 hypre_TFree(CFN_marker);
1135 hypre_TFree(coarse_dof_func);
1136 coarse_dof_func = NULL;
1137 hypre_BoomerAMGCoarseParms(comm, local_num_vars,
1138 num_functions, dof_func_array[level], CF_marker,
1139 &coarse_dof_func,&coarse_pnts_global1);
1140 if (num_functions > 1 && nodal > -1 && (!block_mode) )
1141 dof_func_array[level+1] = coarse_dof_func;
1142 /*if (my_id == 0 && debug_flag==1)
1143 {
1144 wall_time = time_getWallclockSeconds() - wall_time;
1145 printf("Proc = %d Level = %d Coarsen 2 Time = %f\n",
1146 my_id,level, wall_time);
1147 fflush(NULL);
1148 }
1149 if (my_id == 0 && debug_flag==1) wall_time = time_getWallclockSeconds();*/
1150 hypre_BoomerAMGBuildPartialExtPIInterp(A_array[level], CF_marker_array[level],
1151 S, coarse_pnts_global1, coarse_pnts_global, num_functions, dof_func_array[level],
1152 debug_flag, 0, P_max2, col_offd_S_to_A, &P2);
1153 /*if (my_id == 0 && debug_flag==1)
1154 {
1155 wall_time = time_getWallclockSeconds() - wall_time;
1156 printf("Proc = %d Level = %d Interp P2 Time = %f\n",
1157 my_id,level, wall_time);
1158 fflush(NULL);
1159 }
1160 if (my_id == 0 && debug_flag==1) wall_time = time_getWallclockSeconds(); */
1161 P = hypre_ParMatmul(P1,P2);
1162 /*if (my_id == 0 && debug_flag==1)
1163 {
1164 wall_time = time_getWallclockSeconds() - wall_time;
1165 printf("Proc = %d Level = %d P1*P2 Time = %f\n",
1166 my_id,level, wall_time);
1167 fflush(NULL);
1168 }
1169 if (my_id == 0 && debug_flag==1) wall_time = time_getWallclockSeconds();*/
1170 hypre_BoomerAMGInterpTruncation(P,agg_trunc_factor,agg_P_max_elmts);
1171 hypre_ParCSRMatrixDestroy(P1);
1172 hypre_ParCSRMatrixOwnsColStarts(P2) = 0;
1173 hypre_ParCSRMatrixDestroy(P2);
1174 hypre_ParCSRMatrixOwnsColStarts(P) = 1;
1175 coarse_pnts_global = coarse_pnts_global1;
1176#ifdef HYPRE_NO_GLOBAL_PARTITION
1177 hypre_NewCommPkgCreate(P);
1178 if (my_id == (num_procs -1)) coarse_size = coarse_pnts_global[1];
1179 MPI_Bcast(&coarse_size, 1, MPI_HYPRE_BIG_INT, num_procs-1, comm);
1180#else
1181 hypre_MatvecCommPkgCreate(P);
1182 coarse_size = coarse_pnts_global[num_procs];
1183#endif
1184
1185 }
1186 else if (hypre_ParAMGDataAggInterpType(amg_data) == 31 && (level < agg_num_levels
1187 && (!block_mode)))
1188 {
1189 hypre_BoomerAMGBuildExtInterp(A_array[level], CF_marker_array[level],
1190 S, coarse_pnts_global, num_functions, dof_func_array[level],
1191 debug_flag, 0, P_max1, col_offd_S_to_A, &P1);
1192 hypre_TFree(col_offd_S_to_A);
1193 hypre_BoomerAMGCreate2ndS (S, CF_marker, num_paths,
1194 coarse_pnts_global, &S2);
1195 if (coarsen_type == 10)
1196 hypre_BoomerAMGCoarsenHMIS(S2, S2, measure_type,
1197 debug_flag, &CFN_marker);
1198 else if (coarsen_type == 8)
1199 hypre_BoomerAMGCoarsenPMIS(S2, S2, 3,
1200 debug_flag, &CFN_marker);
1201 else if (coarsen_type == 9)
1202 hypre_BoomerAMGCoarsenHMIS(S2, S2, measure_type+3,
1203 debug_flag, &CFN_marker);
1204 else if (coarsen_type == 11)
1205 hypre_BoomerAMGCoarsenPMIS(S2, S2, 4,
1206 debug_flag, &CFN_marker);
1207 else if (coarsen_type == 6)
1208 hypre_BoomerAMGCoarsenFalgout(S2, S2, measure_type,
1209 debug_flag, &CFN_marker);
1210 else if (coarsen_type == 0)
1211 hypre_BoomerAMGCoarsen(S2, S2, 0, debug_flag, &CFN_marker);
1212 hypre_ParCSRMatrixDestroy(S2);
1213 hypre_BoomerAMGCorrectCFMarker2 (CF_marker, local_num_vars, CFN_marker);
1214 hypre_TFree(CFN_marker);
1215 hypre_TFree(coarse_dof_func);
1216 coarse_dof_func = NULL;
1217 hypre_BoomerAMGCoarseParms(comm, local_num_vars,
1218 num_functions, dof_func_array[level], CF_marker,
1219 &coarse_dof_func,&coarse_pnts_global1);
1220 if (num_functions > 1 && nodal > -1 && (!block_mode) )
1221 dof_func_array[level+1] = coarse_dof_func;
1222 /*if (my_id == 0 && debug_flag==1)
1223 {
1224 wall_time = time_getWallclockSeconds() - wall_time;
1225 printf("Proc = %d Level = %d Coarsen 2 Time = %f\n",
1226 my_id,level, wall_time);
1227 fflush(NULL);
1228 }
1229 if (my_id == 0 && debug_flag==1) wall_time = time_getWallclockSeconds();*/
1230 hypre_BoomerAMGBuildPartialExtInterp(A_array[level], CF_marker_array[level],
1231 S, coarse_pnts_global1, coarse_pnts_global, num_functions, dof_func_array[level],
1232 debug_flag, 0, P_max2, col_offd_S_to_A, &P2);
1233 /*if (my_id == 0 && debug_flag==1)
1234 {
1235 wall_time = time_getWallclockSeconds() - wall_time;
1236 printf("Proc = %d Level = %d Interp P2 Time = %f\n",
1237 my_id,level, wall_time);
1238 fflush(NULL);
1239 }
1240 if (my_id == 0 && debug_flag==1) wall_time = time_getWallclockSeconds(); */
1241 P = hypre_ParMatmul(P1,P2);
1242 /*if (my_id == 0 && debug_flag==1)
1243 {
1244 wall_time = time_getWallclockSeconds() - wall_time;
1245 printf("Proc = %d Level = %d P1*P2 Time = %f\n",
1246 my_id,level, wall_time);
1247 fflush(NULL);
1248 }
1249 if (my_id == 0 && debug_flag==1) wall_time = time_getWallclockSeconds();*/
1250 hypre_BoomerAMGInterpTruncation(P,agg_trunc_factor,agg_P_max_elmts);
1251 hypre_ParCSRMatrixDestroy(P1);
1252 hypre_ParCSRMatrixOwnsColStarts(P2) = 0;
1253 hypre_ParCSRMatrixDestroy(P2);
1254 hypre_ParCSRMatrixOwnsColStarts(P) = 1;
1255 coarse_pnts_global = coarse_pnts_global1;
1256#ifdef HYPRE_NO_GLOBAL_PARTITION
1257 hypre_NewCommPkgCreate(P);
1258 if (my_id == (num_procs -1)) coarse_size = coarse_pnts_global[1];
1259 MPI_Bcast(&coarse_size, 1, MPI_HYPRE_BIG_INT, num_procs-1, comm);
1260#else
1261 hypre_MatvecCommPkgCreate(P);
1262 coarse_size = coarse_pnts_global[num_procs];
1263#endif
1264
1265 }
1266 else if (hypre_ParAMGDataAggInterpType(amg_data) == 32 && (level < agg_num_levels
1267 && (!block_mode)))
1268 {
1269 hypre_BoomerAMGBuildStdInterp(A_array[level], CF_marker_array[level],
1270 S, coarse_pnts_global, num_functions, dof_func_array[level],
1271 debug_flag, 0, P_max1, 0, col_offd_S_to_A, &P1);
1272 hypre_TFree(col_offd_S_to_A);
1273 hypre_BoomerAMGCreate2ndS (S, CF_marker, num_paths,
1274 coarse_pnts_global, &S2);
1275 if (coarsen_type == 10)
1276 hypre_BoomerAMGCoarsenHMIS(S2, S2, measure_type,
1277 debug_flag, &CFN_marker);
1278 else if (coarsen_type == 8)
1279 hypre_BoomerAMGCoarsenPMIS(S2, S2, 3,
1280 debug_flag, &CFN_marker);
1281 else if (coarsen_type == 9)
1282 hypre_BoomerAMGCoarsenHMIS(S2, S2, measure_type+3,
1283 debug_flag, &CFN_marker);
1284 else if (coarsen_type == 11)
1285 hypre_BoomerAMGCoarsenPMIS(S2, S2, 4,
1286 debug_flag, &CFN_marker);
1287 else if (coarsen_type == 6)
1288 hypre_BoomerAMGCoarsenFalgout(S2, S2, measure_type,
1289 debug_flag, &CFN_marker);
1290 else if (coarsen_type == 0)
1291 hypre_BoomerAMGCoarsen(S2, S2, 0, debug_flag, &CFN_marker);
1292 hypre_ParCSRMatrixDestroy(S2);
1293 hypre_BoomerAMGCorrectCFMarker2 (CF_marker, local_num_vars, CFN_marker);
1294 hypre_TFree(CFN_marker);
1295 hypre_TFree(coarse_dof_func);
1296 coarse_dof_func = NULL;
1297 hypre_BoomerAMGCoarseParms(comm, local_num_vars,
1298 num_functions, dof_func_array[level], CF_marker,
1299 &coarse_dof_func,&coarse_pnts_global1);
1300 if (num_functions > 1 && nodal > -1 && (!block_mode) )
1301 dof_func_array[level+1] = coarse_dof_func;
1302 hypre_BoomerAMGBuildPartialStdInterp(A_array[level], CF_marker_array[level],
1303 S, coarse_pnts_global1, coarse_pnts_global, num_functions, dof_func_array[level],
1304 debug_flag, 0, P_max2, 0, col_offd_S_to_A, &P2);
1305 P = hypre_ParMatmul(P1,P2);
1306 hypre_BoomerAMGInterpTruncation(P,agg_trunc_factor,agg_P_max_elmts);
1307 hypre_ParCSRMatrixDestroy(P1);
1308 hypre_ParCSRMatrixOwnsColStarts(P2) = 0;
1309 hypre_ParCSRMatrixDestroy(P2);
1310 hypre_ParCSRMatrixOwnsColStarts(P) = 1;
1311 coarse_pnts_global = coarse_pnts_global1;
1312#ifdef HYPRE_NO_GLOBAL_PARTITION
1313 hypre_NewCommPkgCreate(P);
1314 if (my_id == (num_procs -1)) coarse_size = coarse_pnts_global[1];
1315 MPI_Bcast(&coarse_size, 1, MPI_HYPRE_BIG_INT, num_procs-1, comm);
1316#else
1317 hypre_MatvecCommPkgCreate(P);
1318 coarse_size = coarse_pnts_global[num_procs];
1319#endif
1320
1321 }
1322 else if ((hypre_ParAMGDataAggInterpType(amg_data) == 4 &&
1323 (level < agg_num_levels && (!block_mode)))||
1324 (hypre_ParAMGDataInterpType(amg_data) == 4 &&
1325 level >= agg_num_levels))
1326 {
1327 if (nodal > -1)
1328 {
1329 hypre_BoomerAMGBuildMultipass(A_array[level], CF_marker_array[level],
1330 S, coarse_pnts_global, num_functions, dof_func_array[level],
1331 debug_flag, agg_trunc_factor, agg_P_max_elmts, 0, col_offd_S_to_A, &P);
1332 hypre_TFree(col_offd_S_to_A);
1333 }
1334 else
1335 {
1336 CFN_marker = CF_marker_array[level];
1337 hypre_BoomerAMGBuildMultipass(AN, CFN_marker,
1338 SN, coarse_pnts_global, 1, NULL,
1339 debug_flag, agg_trunc_factor, agg_P_max_elmts, 0, col_offd_SN_to_AN, &PN);
1340 hypre_TFree(col_offd_SN_to_AN);
1341 hypre_BoomerAMGCreateScalarCFS(PN, CFN_marker, NULL,
1342 num_functions, nodal, 1, &dof_func1,
1343 &CF_marker, NULL, &P);
1344 dof_func_array[level+1] = dof_func1;
1345 CF_marker_array[level] = CF_marker;
1346 hypre_TFree (CFN_marker);
1347 hypre_ParCSRMatrixDestroy(AN);
1348 hypre_ParCSRMatrixDestroy(PN);
1349 hypre_ParCSRMatrixDestroy(SN);
1350 }
1351 }
1352 else if (hypre_ParAMGDataInterpType(amg_data) == 5)
1353 {
1354 if (nodal > -1)
1355 {
1356 hypre_BoomerAMGBuildMultipass(A_array[level], CF_marker_array[level],
1357 S, coarse_pnts_global, num_functions, dof_func_array[level],
1358 debug_flag, trunc_factor, P_max_elmts, 1, col_offd_S_to_A, &P);
1359 hypre_TFree(col_offd_S_to_A);
1360 }
1361 else
1362 {
1363 CFN_marker = CF_marker_array[level];
1364 hypre_BoomerAMGBuildMultipass(AN, CFN_marker,
1365 SN, coarse_pnts_global, 1, NULL,
1366 debug_flag, trunc_factor, P_max_elmts, 1, col_offd_SN_to_AN, &PN);
1367 hypre_TFree(col_offd_SN_to_AN);
1368 hypre_BoomerAMGCreateScalarCFS(PN, CFN_marker, NULL,
1369 num_functions, nodal, 1, &dof_func1,
1370 &CF_marker, NULL, &P);
1371 dof_func_array[level+1] = dof_func1;
1372 CF_marker_array[level] = CF_marker;
1373 hypre_TFree (CFN_marker);
1374 hypre_ParCSRMatrixDestroy(AN);
1375 hypre_ParCSRMatrixDestroy(PN);
1376 hypre_ParCSRMatrixDestroy(SN);
1377 }
1378 }
1379 else if (hypre_ParAMGDataInterpType(amg_data) == 6) /*Extended+i interpolation */
1380 {
1381 hypre_BoomerAMGBuildExtPIInterp(A_array[level], CF_marker_array[level],
1382 S, coarse_pnts_global, num_functions, dof_func_array[level],
1383 debug_flag, trunc_factor, P_max_elmts, col_offd_S_to_A, &P);
1384 hypre_TFree(col_offd_S_to_A);
1385 }
1386 else if (hypre_ParAMGDataInterpType(amg_data) == 7) /*Extended+i interpolation (common C only)*/
1387 {
1388 hypre_BoomerAMGBuildExtPICCInterp(A_array[level], CF_marker_array[level],
1389 S, coarse_pnts_global, num_functions, dof_func_array[level],
1390 debug_flag, trunc_factor, P_max_elmts, col_offd_S_to_A, &P);
1391 hypre_TFree(col_offd_S_to_A);
1392 }
1393 /*else if (hypre_ParAMGDataInterpType(amg_data) == 12)*/ /*FF interpolation */
1394 /*{
1395 hypre_BoomerAMGBuildFFInterp(A_array[level], CF_marker_array[level],
1396 S, coarse_pnts_global, num_functions, dof_func_array[level],
1397 debug_flag, trunc_factor, P_max_elmts, col_offd_S_to_A, &P);
1398 hypre_TFree(col_offd_S_to_A);
1399 }*/
1400 else if (hypre_ParAMGDataInterpType(amg_data) == 13) /*FF1 interpolation */
1401 {
1402 hypre_BoomerAMGBuildFF1Interp(A_array[level], CF_marker_array[level],
1403 S, coarse_pnts_global, num_functions, dof_func_array[level],
1404 debug_flag, trunc_factor, P_max_elmts, col_offd_S_to_A, &P);
1405 hypre_TFree(col_offd_S_to_A);
1406 }
1407 /*else if (hypre_ParAMGDataInterpType(amg_data) == 8)*/ /*Standard interpolation */
1408 /*{
1409 hypre_BoomerAMGBuildStdInterp(A_array[level], CF_marker_array[level],
1410 S, coarse_pnts_global, num_functions, dof_func_array[level],
1411 debug_flag, trunc_factor, P_max_elmts, 0, col_offd_S_to_A, &P);
1412 hypre_TFree(col_offd_S_to_A);
1413 }
1414 else if (hypre_ParAMGDataInterpType(amg_data) == 9) */ /*Standard interpolation with separation of
1415 negative and positive weights*/
1416 /*{
1417 hypre_BoomerAMGBuildStdInterp(A_array[level], CF_marker_array[level],
1418 S, coarse_pnts_global, num_functions, dof_func_array[level],
1419 debug_flag, trunc_factor, P_max_elmts, 1, col_offd_S_to_A, &P);
1420 hypre_TFree(col_offd_S_to_A);
1421 }*/
1422 else if (hypre_ParAMGDataGSMG(amg_data) == 0)
1423 {
1424 /* no nodal interp OR nodal =0 */
1425 {
1426 if (nodal > -1) /* 0, 1, 2, 3 - unknown approach for interpolation */
1427 {
1428
1429 hypre_BoomerAMGBuildInterp(A_array[level], CF_marker_array[level],
1430 S, coarse_pnts_global, num_functions,
1431 dof_func_array[level],
1432 debug_flag, trunc_factor, P_max_elmts, col_offd_S_to_A, &P);
1433
1434 hypre_TFree(col_offd_S_to_A);
1435 }
1436 else /* -1, -2, -3: here we build interp using the nodal matrix and then convert
1437 to regular size */
1438 {
1439 CFN_marker = CF_marker_array[level];
1440 hypre_BoomerAMGBuildInterp(AN, CFN_marker,
1441 SN, coarse_pnts_global, 1, NULL,
1442 debug_flag, trunc_factor, P_max_elmts, col_offd_SN_to_AN, &PN);
1443 hypre_TFree(col_offd_SN_to_AN);
1444 hypre_BoomerAMGCreateScalarCFS(PN, CFN_marker, NULL,
1445 num_functions, nodal, 1, &dof_func1,
1446 &CF_marker, NULL, &P);
1447 dof_func_array[level+1] = dof_func1;
1448 CF_marker_array[level] = CF_marker;
1449 hypre_TFree (CFN_marker);
1450 hypre_ParCSRMatrixDestroy(AN);
1451 hypre_ParCSRMatrixDestroy(PN);
1452 hypre_ParCSRMatrixDestroy(SN);
1453 }
1454 }
1455 }
1456
1457 /*if ( post_interp_type>=1 && level < agg_num_levels)*/
1458 for (i=0; i < post_interp_type; i++)
1459 /* Improve on P with Jacobi interpolation */
1460 hypre_BoomerAMGJacobiInterp( A_array[level], &P, S,
1461 num_functions, dof_func,
1462 CF_marker_array[level],
1463 level, jacobi_trunc_threshold, 0.5*jacobi_trunc_threshold );
1464
1465
1466 if (debug_flag==1)
1467 {
1468 wall_time = time_getWallclockSeconds() - wall_time;
1469 printf("Proc = %d Level = %d Build Interp Time = %f\n",
1470 my_id,level, wall_time);
1471 fflush(NULL);
1472 }
1473
1474 if (!block_mode)
1475 {
1476 P_array[level] = P;
1477 }
1478
1479 if (S) hypre_ParCSRMatrixDestroy(S);
1480 S = NULL;
1481
1482
1483 /*-------------------------------------------------------------
1484 * Build coarse-grid operator, A_array[level+1] by R*A*P
1485 *--------------------------------------------------------------*/
1486
1487 if (debug_flag==1) wall_time = time_getWallclockSeconds();
1488
1489
1490 hypre_BoomerAMGBuildCoarseOperator(P_array[level], A_array[level] ,
1491 P_array[level], &A_H);
1492
1493 if (debug_flag==1)
1494 {
1495 wall_time = time_getWallclockSeconds() - wall_time;
1496 printf("Proc = %d Level = %d Build Coarse Operator Time = %f\n",
1497 my_id,level, wall_time);
1498 fflush(NULL);
1499 }
1500
1501 ++level;
1502
1503 if (!block_mode)
1504 {
1505 hypre_ParCSRMatrixSetNumNonzeros(A_H);
1506 hypre_ParCSRMatrixSetDNumNonzeros(A_H);
1507 A_array[level] = A_H;
1508 }
1509
1510 size = ((double) fine_size )*.75;
1511 if (coarsen_type > 0 && coarse_size >= (HYPRE_BigInt) size)
1512 {
1513 coarsen_type = 0;
1514 }
1515
1516 {
1517
1518 int thresh = hypre_max(coarse_threshold, seq_threshold);
1519
1520 //stop coarsening ?
1521 if ( (level == max_levels-1) || (coarse_size <= (HYPRE_BigInt) thresh) )
1522 {
1523 not_finished_coarsening = 0;
1524 }
1525
1526 }
1527
1528#if MPIP_SETUP_ON
1529 MPI_Pcontrol(3);
1530 MPI_Pcontrol(0);
1531#endif
1532
1533 } /* end of coarsening */
1534
1535 //more seq coarsening?
1536 if ( (seq_threshold > coarse_threshold) && (coarse_size > coarse_threshold) && (level != max_levels-1))
1537 {
1538 hypre_seqAMGSetup( amg_data, level, coarse_threshold);
1539 }
1540
1541
1542
1543 if (level > 0)
1544 {
1545 F_array[level] =
1546 hypre_ParVectorCreate(hypre_ParCSRMatrixComm(A_array[level]),
1547 hypre_ParCSRMatrixGlobalNumRows(A_array[level]),
1548 hypre_ParCSRMatrixRowStarts(A_array[level]));
1549 hypre_ParVectorInitialize(F_array[level]);
1550 hypre_ParVectorSetPartitioningOwner(F_array[level],0);
1551
1552 U_array[level] =
1553 hypre_ParVectorCreate(hypre_ParCSRMatrixComm(A_array[level]),
1554 hypre_ParCSRMatrixGlobalNumRows(A_array[level]),
1555 hypre_ParCSRMatrixRowStarts(A_array[level]));
1556 hypre_ParVectorInitialize(U_array[level]);
1557 hypre_ParVectorSetPartitioningOwner(U_array[level],0);
1558 }
1559
1560 /*-----------------------------------------------------------------------
1561 * enter all the stuff created, A[level], P[level], CF_marker[level],
1562 * for levels 1 through coarsest, into amg_data data structure
1563 *-----------------------------------------------------------------------*/
1564
1565 num_levels = level+1;
1566 hypre_ParAMGDataNumLevels(amg_data) = num_levels;
1567
1568
1569 /*-----------------------------------------------------------------------
1570 * Setup F and U arrays
1571 *-----------------------------------------------------------------------*/
1572
1573 if (grid_relax_type[1] == 8 || grid_relax_type[1] == 19 || grid_relax_type[1] == 20 )
1574 {
1575 l1_norms = hypre_CTAlloc(double *, num_levels);
1576 hypre_ParAMGDataL1Norms(amg_data) = l1_norms;
1577
1578 Ztemp = hypre_ParVectorCreate(hypre_ParCSRMatrixComm(A_array[0]),
1579 hypre_ParCSRMatrixGlobalNumRows(A_array[0]),
1580 hypre_ParCSRMatrixRowStarts(A_array[0]));
1581 hypre_ParVectorInitialize(Ztemp);
1582 hypre_ParVectorSetPartitioningOwner(Ztemp,0);
1583 hypre_ParAMGDataZtemp(amg_data) = Ztemp;
1584 }
1585 else if (grid_relax_type[1] == 18)
1586 {
1587 l1_norms = hypre_CTAlloc(double *, num_levels);
1588 hypre_ParAMGDataL1Norms(amg_data) = l1_norms;
1589 }
1590
1591 if (grid_relax_type[1] == 11 || grid_relax_type[1] == 15
1592 || grid_relax_type[1] == 16 || grid_relax_type[1] == 17 ) /* Chebyshev */
1593 {
1594 max_eig_est = hypre_CTAlloc(double, num_levels);
1595 min_eig_est = hypre_CTAlloc(double, num_levels);
1596 hypre_ParAMGDataMaxEigEst(amg_data) = max_eig_est;
1597 hypre_ParAMGDataMinEigEst(amg_data) = min_eig_est;
1598
1599 Ztemp = hypre_ParVectorCreate(hypre_ParCSRMatrixComm(A_array[0]),
1600 hypre_ParCSRMatrixGlobalNumRows(A_array[0]),
1601 hypre_ParCSRMatrixRowStarts(A_array[0]));
1602 hypre_ParVectorInitialize(Ztemp);
1603 hypre_ParVectorSetPartitioningOwner(Ztemp,0);
1604 hypre_ParAMGDataZtemp(amg_data) = Ztemp;
1605
1606 }
1607 if (grid_relax_type[1] == 13 || grid_relax_type[3] == 13 ||
1608 grid_relax_type[1] == 14 || grid_relax_type[3] == 14 ) /* CG */
1609 {
1610 smoother = hypre_CTAlloc(HYPRE_Solver, num_levels);
1611 hypre_ParAMGDataSmoother(amg_data) = smoother;
1612 if (grid_relax_type[1] == 14 || grid_relax_type[3] == 14)
1613 {
1614 smoother_prec = hypre_CTAlloc(HYPRE_Solver, num_levels);
1615 hypre_ParAMGDataSmootherPrec(amg_data) = smoother_prec;
1616
1617 }
1618 }
1619 if ( grid_relax_type[1] == 3 ||
1620 grid_relax_type[1] == 6) /* GS (3 and 6) all need a temp vector*/
1621 {
1622 Ztemp = hypre_ParVectorCreate(hypre_ParCSRMatrixComm(A_array[0]),
1623 hypre_ParCSRMatrixGlobalNumRows(A_array[0]),
1624 hypre_ParCSRMatrixRowStarts(A_array[0]));
1625 hypre_ParVectorInitialize(Ztemp);
1626 hypre_ParVectorSetPartitioningOwner(Ztemp,0);
1627 hypre_ParAMGDataZtemp(amg_data) = Ztemp;
1628 }
1629 for (j = 0; j < num_levels; j++)
1630 {
1631 if (num_threads == 1)
1632 {
1633 if ((grid_relax_type[1] == 8 ) && j < num_levels-1)
1634 {
1635 hypre_ParCSRComputeL1Norms(A_array[j], 4, NULL, &l1_norms[j]);
1636 }
1637 else if ((grid_relax_type[3] == 8 ) && j == num_levels-1)
1638 {
1639 hypre_ParCSRComputeL1Norms(A_array[j], 4, NULL, &l1_norms[j]);
1640 }
1641 if (grid_relax_type[1] == 18 && j < num_levels-1)
1642 {
1643 if (hypre_ParAMGDataRelaxOrder(amg_data))
1644 hypre_ParCSRComputeL1Norms(A_array[j], 4, CF_marker_array[j], &l1_norms[j]);
1645 else
1646 hypre_ParCSRComputeL1Norms(A_array[j], 1, NULL, &l1_norms[j]);
1647 }
1648 else if (grid_relax_type[3] == 18 && j == num_levels-1)
1649 {
1650 hypre_ParCSRComputeL1Norms(A_array[j], 1, NULL, &l1_norms[j]);
1651 }
1652 }
1653 else
1654 {
1655 if ((grid_relax_type[1] == 8 ) && j < num_levels-1)
1656 {
1657 hypre_ParCSRComputeL1NormsThreads(A_array[j], 4, num_threads, NULL, &l1_norms[j]);
1658 }
1659 else if ((grid_relax_type[3] == 8 ) && j == num_levels-1)
1660 {
1661 hypre_ParCSRComputeL1NormsThreads(A_array[j], 4, num_threads, NULL, &l1_norms[j]);
1662 }
1663 if (grid_relax_type[1] == 18 && j < num_levels-1)
1664 {
1665 if (hypre_ParAMGDataRelaxOrder(amg_data))
1666 hypre_ParCSRComputeL1NormsThreads(A_array[j], 4, num_threads, CF_marker_array[j], &l1_norms[j]);
1667 else
1668 hypre_ParCSRComputeL1NormsThreads(A_array[j], 1, num_threads, NULL, &l1_norms[j]);
1669 }
1670 else if (grid_relax_type[3] == 18 && j == num_levels-1)
1671 {
1672 hypre_ParCSRComputeL1NormsThreads(A_array[j], 1, num_threads, NULL, &l1_norms[j]);
1673 }
1674 }
1675 if (grid_relax_type[1] == 11 || grid_relax_type[1] == 15
1676 || grid_relax_type[1] == 16 || grid_relax_type[1] == 17 )
1677 {
1678 int scale = 0;
1679 double temp_d, temp_d2;
1680 /* hypre_ParCSRMaxEigEstimate(A_array[j], 0, &temp_d);*/
1681 if (grid_relax_type[1] == 16 || grid_relax_type[1] == 17)
1682 scale = 1;
1683
1684 hypre_ParCSRMaxEigEstimateCG(A_array[j], scale, 10, &temp_d, &temp_d2);
1685 /* if (my_id ==0) printf("cg est: max = %g, min = %g\n", temp_d, temp_d2); */
1686
1687 max_eig_est[j] = temp_d;
1688 min_eig_est[j] = temp_d2;
1689 }
1690
1691 if (grid_relax_type[1] == 13 || (grid_relax_type[3] == 13 && j == (num_levels-1)) )
1692 {
1693
1694 HYPRE_ParCSRPCGCreate(comm, &smoother[j]);
1695 HYPRE_ParCSRPCGSetup(smoother[j],
1696 (HYPRE_ParCSRMatrix) A_array[j],
1697 (HYPRE_ParVector) F_array[j],
1698 (HYPRE_ParVector) U_array[j]);
1699
1700 HYPRE_PCGSetTol(smoother[j], 1e-9); /* do a fixed number of iterations, so the
1701 conv. tolerance is small*/
1702 HYPRE_PCGSetTwoNorm(smoother[j], 1); /* use 2-norm*/
1703
1704 /* HYPRE_PCGSetLogging(smoother[j], 1);*/ /* needed to get run info later */
1705
1706 /* Do we want diagonal scaling? */
1707 /*
1708 HYPRE_PCGSetPrecond(smoother[j],
1709 (HYPRE_PtrToSolverFcn) HYPRE_ParCSRDiagScale,
1710 (HYPRE_PtrToSolverFcn) HYPRE_ParCSRDiagScaleSetup,
1711 NULL);
1712 */
1713
1714
1715
1716
1717
1718
1719 HYPRE_ParCSRPCGSetup(smoother[j],
1720 (HYPRE_ParCSRMatrix) A_array[j],
1721 (HYPRE_ParVector) F_array[j],
1722 (HYPRE_ParVector) U_array[j]);
1723
1724
1725 }
1726
1727 if (grid_relax_type[1] == 14 || (grid_relax_type[3] == 14 && j == (num_levels-1)) )
1728 {
1729
1730 HYPRE_ParCSRPCGCreate(comm, &smoother[j]);
1731 HYPRE_ParCSRPCGSetup(smoother[j],
1732 (HYPRE_ParCSRMatrix) A_array[j],
1733 (HYPRE_ParVector) F_array[j],
1734 (HYPRE_ParVector) U_array[j]);
1735
1736 HYPRE_PCGSetTol(smoother[j], 1e-9); /* want to do a fixed number of iterations, so the
1737 conv. tolerance is very small*/
1738 HYPRE_PCGSetTwoNorm(smoother[j], 1); /* use 2-norm*/
1739
1740 HYPRE_BoomerAMGCreate(&smoother_prec[j]);
1741 HYPRE_BoomerAMGSetCoarsenType(smoother_prec[j], 10); /* hmis*/
1742 HYPRE_BoomerAMGSetInterpType(smoother_prec[j], 6); /* ext+i (5) */
1743 HYPRE_BoomerAMGSetPMaxElmts(smoother_prec[j], 5);
1744 HYPRE_BoomerAMGSetRelaxType(smoother_prec[j], 8); /* L1-SGS */
1745 HYPRE_BoomerAMGSetMaxIter(smoother_prec[j], 1);
1746
1747 HYPRE_BoomerAMGSetup(smoother_prec[j],
1748 (HYPRE_ParCSRMatrix) A_array[j],
1749 (HYPRE_ParVector) F_array[j],
1750 (HYPRE_ParVector) U_array[j]);
1751
1752
1753 HYPRE_PCGSetPrecond( smoother[j],
1754 (HYPRE_PtrToSolverFcn) HYPRE_BoomerAMGSolve,
1755 (HYPRE_PtrToSolverFcn) HYPRE_BoomerAMGSetup,
1756 smoother_prec[j]);
1757
1758
1759 HYPRE_ParCSRPCGSetup(smoother[j],
1760 (HYPRE_ParCSRMatrix) A_array[j],
1761 (HYPRE_ParVector) F_array[j],
1762 (HYPRE_ParVector) U_array[j]);
1763
1764
1765 }
1766
1767
1768 }
1769
1770 if ( amg_logging > 1 ) {
1771
1772 Residual_array=
1773 hypre_ParVectorCreate(hypre_ParCSRMatrixComm(A_array[0]),
1774 hypre_ParCSRMatrixGlobalNumRows(A_array[0]),
1775 hypre_ParCSRMatrixRowStarts(A_array[0]) );
1776 hypre_ParVectorInitialize(Residual_array);
1777 hypre_ParVectorSetPartitioningOwner(Residual_array,0);
1778 hypre_ParAMGDataResidual(amg_data) = Residual_array;
1779 }
1780 else
1781 hypre_ParAMGDataResidual(amg_data) = NULL;
1782
1783 /*-----------------------------------------------------------------------
1784 * Print some stuff
1785 *-----------------------------------------------------------------------*/
1786
1787 if (amg_print_level == 1 || amg_print_level == 3)
1788 hypre_BoomerAMGSetupStats(amg_data,A);
1789
1790
1791/* print out matrices on all levels */
1792#if 0
1793{
1794 char filename[256];
1795
1796 for (level = 0; level < num_levels; level++)
1797 {
1798 sprintf(filename, "BoomerAMG.out.A.%02d", level);
1799 hypre_ParCSRMatrixPrint(A_array[level], filename);
1800
1801 }
1802}
1803
1804#endif
1805
1806
1807
1808 return(Setup_err_flag);
1809}
Note: See TracBrowser for help on using the repository browser.