1 /******************************************************************************
2 * Copyright 1998-2019 Lawrence Livermore National Security, LLC and other
3 * HYPRE Project Developers. See the top-level COPYRIGHT file for details.
4 *
5 * SPDX-License-Identifier: (Apache-2.0 OR MIT)
6 ******************************************************************************/
7
8 /*
9 Example 5
10
11 Interface: Linear-Algebraic (IJ)
12
13 Compile with: make ex5
14
15 Sample run: mpirun -np 4 ex5
16
17 Description: This example solves the 2-D Laplacian problem with zero boundary
18 conditions on an n x n grid. The number of unknowns is N=n^2.
19 The standard 5-point stencil is used, and we solve for the
20 interior nodes only.
21
22 This example solves the same problem as Example 3. Available
23 solvers are AMG, PCG, and PCG with AMG or Parasails
24 preconditioners. */
25
26 #include <stdio.h>
27 #include <stdlib.h>
28 #include <string.h>
29 #include <math.h>
30 #include "HYPRE_krylov.h"
31 #include "HYPRE.h"
32 #include "HYPRE_parcsr_ls.h"
33 #include "ex.h"
34
35 #ifdef HYPRE_EXVIS
36 #include "vis.c"
37 #endif
38
39 int hypre_FlexGMRESModifyPCAMGExample(void *precond_data, int iterations,
40 double rel_residual_norm);
41
42 #define my_min(a,b) (((a)<(b)) ? (a) : (b))
43
main(int argc,char * argv[])44 int main (int argc, char *argv[])
45 {
46 int i;
47 int myid, num_procs;
48 int N, n;
49
50 int ilower, iupper;
51 int local_size, extra;
52
53 int solver_id;
54 int vis, print_system;
55
56 double h, h2;
57
58 HYPRE_IJMatrix A;
59 HYPRE_ParCSRMatrix parcsr_A;
60 HYPRE_IJVector b;
61 HYPRE_ParVector par_b;
62 HYPRE_IJVector x;
63 HYPRE_ParVector par_x;
64
65 HYPRE_Solver solver, precond;
66
67 /* Initialize MPI */
68 MPI_Init(&argc, &argv);
69 MPI_Comm_rank(MPI_COMM_WORLD, &myid);
70 MPI_Comm_size(MPI_COMM_WORLD, &num_procs);
71
72 /* Initialize HYPRE */
73 HYPRE_Init();
74
75 /* Print GPU info */
76 /* HYPRE_PrintDeviceInfo(); */
77
78 /* Default problem parameters */
79 n = 33;
80 solver_id = 0;
81 vis = 0;
82 print_system = 0;
83
84
85 /* Parse command line */
86 {
87 int arg_index = 0;
88 int print_usage = 0;
89
90 while (arg_index < argc)
91 {
92 if ( strcmp(argv[arg_index], "-n") == 0 )
93 {
94 arg_index++;
95 n = atoi(argv[arg_index++]);
96 }
97 else if ( strcmp(argv[arg_index], "-solver") == 0 )
98 {
99 arg_index++;
100 solver_id = atoi(argv[arg_index++]);
101 }
102 else if ( strcmp(argv[arg_index], "-vis") == 0 )
103 {
104 arg_index++;
105 vis = 1;
106 }
107 else if ( strcmp(argv[arg_index], "-print_system") == 0 )
108 {
109 arg_index++;
110 print_system = 1;
111 }
112 else if ( strcmp(argv[arg_index], "-help") == 0 )
113 {
114 print_usage = 1;
115 break;
116 }
117 else
118 {
119 arg_index++;
120 }
121 }
122
123 if ((print_usage) && (myid == 0))
124 {
125 printf("\n");
126 printf("Usage: %s [<options>]\n", argv[0]);
127 printf("\n");
128 printf(" -n <n> : problem size in each direction (default: 33)\n");
129 printf(" -solver <ID> : solver ID\n");
130 printf(" 0 - AMG (default) \n");
131 printf(" 1 - AMG-PCG\n");
132 printf(" 8 - ParaSails-PCG\n");
133 printf(" 50 - PCG\n");
134 printf(" 61 - AMG-FlexGMRES\n");
135 printf(" -vis : save the solution for GLVis visualization\n");
136 printf(" -print_system : print the matrix and rhs\n");
137 printf("\n");
138 }
139
140 if (print_usage)
141 {
142 MPI_Finalize();
143 return (0);
144 }
145 }
146
147 /* Preliminaries: want at least one processor per row */
148 if (n*n < num_procs) n = sqrt(num_procs) + 1;
149 N = n*n; /* global number of rows */
150 h = 1.0/(n+1); /* mesh size*/
151 h2 = h*h;
152
153 /* Each processor knows only of its own rows - the range is denoted by ilower
154 and upper. Here we partition the rows. We account for the fact that
155 N may not divide evenly by the number of processors. */
156 local_size = N/num_procs;
157 extra = N - local_size*num_procs;
158
159 ilower = local_size*myid;
160 ilower += my_min(myid, extra);
161
162 iupper = local_size*(myid+1);
163 iupper += my_min(myid+1, extra);
164 iupper = iupper - 1;
165
166 /* How many rows do I have? */
167 local_size = iupper - ilower + 1;
168
169 /* Create the matrix.
170 Note that this is a square matrix, so we indicate the row partition
171 size twice (since number of rows = number of cols) */
172 HYPRE_IJMatrixCreate(MPI_COMM_WORLD, ilower, iupper, ilower, iupper, &A);
173
174 /* Choose a parallel csr format storage (see the User's Manual) */
175 HYPRE_IJMatrixSetObjectType(A, HYPRE_PARCSR);
176
177 /* Initialize before setting coefficients */
178 HYPRE_IJMatrixInitialize(A);
179
180 /* Now go through my local rows and set the matrix entries.
181 Each row has at most 5 entries. For example, if n=3:
182
183 A = [M -I 0; -I M -I; 0 -I M]
184 M = [4 -1 0; -1 4 -1; 0 -1 4]
185
186 Note that here we are setting one row at a time, though
187 one could set all the rows together (see the User's Manual).
188 */
189 {
190 int nnz;
191 /* OK to use constant-length arrays for CPUs
192 double values[5];
193 int cols[5];
194 */
195 double *values = (double *) malloc(5*sizeof(double));
196 int *cols = (int *) malloc(5*sizeof(int));
197 int *tmp = (int *) malloc(2*sizeof(int));
198
199 for (i = ilower; i <= iupper; i++)
200 {
201 nnz = 0;
202
203 /* The left identity block:position i-n */
204 if ((i-n)>=0)
205 {
206 cols[nnz] = i-n;
207 values[nnz] = -1.0;
208 nnz++;
209 }
210
211 /* The left -1: position i-1 */
212 if (i%n)
213 {
214 cols[nnz] = i-1;
215 values[nnz] = -1.0;
216 nnz++;
217 }
218
219 /* Set the diagonal: position i */
220 cols[nnz] = i;
221 values[nnz] = 4.0;
222 nnz++;
223
224 /* The right -1: position i+1 */
225 if ((i+1)%n)
226 {
227 cols[nnz] = i+1;
228 values[nnz] = -1.0;
229 nnz++;
230 }
231
232 /* The right identity block:position i+n */
233 if ((i+n)< N)
234 {
235 cols[nnz] = i+n;
236 values[nnz] = -1.0;
237 nnz++;
238 }
239
240 /* Set the values for row i */
241 tmp[0] = nnz;
242 tmp[1] = i;
243 HYPRE_IJMatrixSetValues(A, 1, &tmp[0], &tmp[1], cols, values);
244 }
245
246 free(values);
247 free(cols);
248 free(tmp);
249 }
250
251 /* Assemble after setting the coefficients */
252 HYPRE_IJMatrixAssemble(A);
253
254 /* Note: for the testing of small problems, one may wish to read
255 in a matrix in IJ format (for the format, see the output files
256 from the -print_system option).
257 In this case, one would use the following routine:
258 HYPRE_IJMatrixRead( <filename>, MPI_COMM_WORLD,
259 HYPRE_PARCSR, &A );
260 <filename> = IJ.A.out to read in what has been printed out
261 by -print_system (processor numbers are omitted).
262 A call to HYPRE_IJMatrixRead is an *alternative* to the
263 following sequence of HYPRE_IJMatrix calls:
264 Create, SetObjectType, Initialize, SetValues, and Assemble
265 */
266
267
268 /* Get the parcsr matrix object to use */
269 HYPRE_IJMatrixGetObject(A, (void**) &parcsr_A);
270
271
272 /* Create the rhs and solution */
273 HYPRE_IJVectorCreate(MPI_COMM_WORLD, ilower, iupper,&b);
274 HYPRE_IJVectorSetObjectType(b, HYPRE_PARCSR);
275 HYPRE_IJVectorInitialize(b);
276
277 HYPRE_IJVectorCreate(MPI_COMM_WORLD, ilower, iupper,&x);
278 HYPRE_IJVectorSetObjectType(x, HYPRE_PARCSR);
279 HYPRE_IJVectorInitialize(x);
280
281 /* Set the rhs values to h^2 and the solution to zero */
282 {
283 double *rhs_values, *x_values;
284 int *rows;
285
286 rhs_values = (double*) calloc(local_size, sizeof(double));
287 x_values = (double*) calloc(local_size, sizeof(double));
288 rows = (int*) calloc(local_size, sizeof(int));
289
290 for (i=0; i<local_size; i++)
291 {
292 rhs_values[i] = h2;
293 x_values[i] = 0.0;
294 rows[i] = ilower + i;
295 }
296
297 HYPRE_IJVectorSetValues(b, local_size, rows, rhs_values);
298 HYPRE_IJVectorSetValues(x, local_size, rows, x_values);
299
300 free(x_values);
301 free(rhs_values);
302 free(rows);
303 }
304
305
306 HYPRE_IJVectorAssemble(b);
307 /* As with the matrix, for testing purposes, one may wish to read in a rhs:
308 HYPRE_IJVectorRead( <filename>, MPI_COMM_WORLD,
309 HYPRE_PARCSR, &b );
310 as an alternative to the
311 following sequence of HYPRE_IJVectors calls:
312 Create, SetObjectType, Initialize, SetValues, and Assemble
313 */
314 HYPRE_IJVectorGetObject(b, (void **) &par_b);
315
316 HYPRE_IJVectorAssemble(x);
317 HYPRE_IJVectorGetObject(x, (void **) &par_x);
318
319
320 /* Print out the system - files names will be IJ.out.A.XXXXX
321 and IJ.out.b.XXXXX, where XXXXX = processor id */
322 if (print_system)
323 {
324 HYPRE_IJMatrixPrint(A, "IJ.out.A");
325 HYPRE_IJVectorPrint(b, "IJ.out.b");
326 }
327
328
329 /* Choose a solver and solve the system */
330
331 /* AMG */
332 if (solver_id == 0)
333 {
334 int num_iterations;
335 double final_res_norm;
336
337 /* Create solver */
338 HYPRE_BoomerAMGCreate(&solver);
339
340 /* Set some parameters (See Reference Manual for more parameters) */
341 HYPRE_BoomerAMGSetPrintLevel(solver, 3); /* print solve info + parameters */
342 HYPRE_BoomerAMGSetOldDefault(solver); /* Falgout coarsening with modified classical interpolaiton */
343 HYPRE_BoomerAMGSetRelaxType(solver, 3); /* G-S/Jacobi hybrid relaxation */
344 HYPRE_BoomerAMGSetRelaxOrder(solver, 1); /* uses C/F relaxation */
345 HYPRE_BoomerAMGSetNumSweeps(solver, 1); /* Sweeeps on each level */
346 HYPRE_BoomerAMGSetMaxLevels(solver, 20); /* maximum number of levels */
347 HYPRE_BoomerAMGSetTol(solver, 1e-7); /* conv. tolerance */
348
349 /* Now setup and solve! */
350 HYPRE_BoomerAMGSetup(solver, parcsr_A, par_b, par_x);
351 HYPRE_BoomerAMGSolve(solver, parcsr_A, par_b, par_x);
352
353 /* Run info - needed logging turned on */
354 HYPRE_BoomerAMGGetNumIterations(solver, &num_iterations);
355 HYPRE_BoomerAMGGetFinalRelativeResidualNorm(solver, &final_res_norm);
356 if (myid == 0)
357 {
358 printf("\n");
359 printf("Iterations = %d\n", num_iterations);
360 printf("Final Relative Residual Norm = %e\n", final_res_norm);
361 printf("\n");
362 }
363
364 /* Destroy solver */
365 HYPRE_BoomerAMGDestroy(solver);
366 }
367 /* PCG */
368 else if (solver_id == 50)
369 {
370 int num_iterations;
371 double final_res_norm;
372
373 /* Create solver */
374 HYPRE_ParCSRPCGCreate(MPI_COMM_WORLD, &solver);
375
376 /* Set some parameters (See Reference Manual for more parameters) */
377 HYPRE_PCGSetMaxIter(solver, 1000); /* max iterations */
378 HYPRE_PCGSetTol(solver, 1e-7); /* conv. tolerance */
379 HYPRE_PCGSetTwoNorm(solver, 1); /* use the two norm as the stopping criteria */
380 HYPRE_PCGSetPrintLevel(solver, 2); /* prints out the iteration info */
381 HYPRE_PCGSetLogging(solver, 1); /* needed to get run info later */
382
383 /* Now setup and solve! */
384 HYPRE_ParCSRPCGSetup(solver, parcsr_A, par_b, par_x);
385 HYPRE_ParCSRPCGSolve(solver, parcsr_A, par_b, par_x);
386
387 /* Run info - needed logging turned on */
388 HYPRE_PCGGetNumIterations(solver, &num_iterations);
389 HYPRE_PCGGetFinalRelativeResidualNorm(solver, &final_res_norm);
390 if (myid == 0)
391 {
392 printf("\n");
393 printf("Iterations = %d\n", num_iterations);
394 printf("Final Relative Residual Norm = %e\n", final_res_norm);
395 printf("\n");
396 }
397
398 /* Destroy solver */
399 HYPRE_ParCSRPCGDestroy(solver);
400 }
401 /* PCG with AMG preconditioner */
402 else if (solver_id == 1)
403 {
404 int num_iterations;
405 double final_res_norm;
406
407 /* Create solver */
408 HYPRE_ParCSRPCGCreate(MPI_COMM_WORLD, &solver);
409
410 /* Set some parameters (See Reference Manual for more parameters) */
411 HYPRE_PCGSetMaxIter(solver, 1000); /* max iterations */
412 HYPRE_PCGSetTol(solver, 1e-7); /* conv. tolerance */
413 HYPRE_PCGSetTwoNorm(solver, 1); /* use the two norm as the stopping criteria */
414 HYPRE_PCGSetPrintLevel(solver, 2); /* print solve info */
415 HYPRE_PCGSetLogging(solver, 1); /* needed to get run info later */
416
417 /* Now set up the AMG preconditioner and specify any parameters */
418 HYPRE_BoomerAMGCreate(&precond);
419 HYPRE_BoomerAMGSetPrintLevel(precond, 1); /* print amg solution info */
420 HYPRE_BoomerAMGSetCoarsenType(precond, 6);
421 HYPRE_BoomerAMGSetOldDefault(precond);
422 HYPRE_BoomerAMGSetRelaxType(precond, 6); /* Sym G.S./Jacobi hybrid */
423 HYPRE_BoomerAMGSetNumSweeps(precond, 1);
424 HYPRE_BoomerAMGSetTol(precond, 0.0); /* conv. tolerance zero */
425 HYPRE_BoomerAMGSetMaxIter(precond, 1); /* do only one iteration! */
426
427 /* Set the PCG preconditioner */
428 HYPRE_PCGSetPrecond(solver, (HYPRE_PtrToSolverFcn) HYPRE_BoomerAMGSolve,
429 (HYPRE_PtrToSolverFcn) HYPRE_BoomerAMGSetup, precond);
430
431 /* Now setup and solve! */
432 HYPRE_ParCSRPCGSetup(solver, parcsr_A, par_b, par_x);
433 HYPRE_ParCSRPCGSolve(solver, parcsr_A, par_b, par_x);
434
435 /* Run info - needed logging turned on */
436 HYPRE_PCGGetNumIterations(solver, &num_iterations);
437 HYPRE_PCGGetFinalRelativeResidualNorm(solver, &final_res_norm);
438 if (myid == 0)
439 {
440 printf("\n");
441 printf("Iterations = %d\n", num_iterations);
442 printf("Final Relative Residual Norm = %e\n", final_res_norm);
443 printf("\n");
444 }
445
446 /* Destroy solver and preconditioner */
447 HYPRE_ParCSRPCGDestroy(solver);
448 HYPRE_BoomerAMGDestroy(precond);
449 }
450 /* PCG with Parasails Preconditioner */
451 else if (solver_id == 8)
452 {
453 int num_iterations;
454 double final_res_norm;
455
456 int sai_max_levels = 1;
457 double sai_threshold = 0.1;
458 double sai_filter = 0.05;
459 int sai_sym = 1;
460
461 /* Create solver */
462 HYPRE_ParCSRPCGCreate(MPI_COMM_WORLD, &solver);
463
464 /* Set some parameters (See Reference Manual for more parameters) */
465 HYPRE_PCGSetMaxIter(solver, 1000); /* max iterations */
466 HYPRE_PCGSetTol(solver, 1e-7); /* conv. tolerance */
467 HYPRE_PCGSetTwoNorm(solver, 1); /* use the two norm as the stopping criteria */
468 HYPRE_PCGSetPrintLevel(solver, 2); /* print solve info */
469 HYPRE_PCGSetLogging(solver, 1); /* needed to get run info later */
470
471 /* Now set up the ParaSails preconditioner and specify any parameters */
472 HYPRE_ParaSailsCreate(MPI_COMM_WORLD, &precond);
473
474 /* Set some parameters (See Reference Manual for more parameters) */
475 HYPRE_ParaSailsSetParams(precond, sai_threshold, sai_max_levels);
476 HYPRE_ParaSailsSetFilter(precond, sai_filter);
477 HYPRE_ParaSailsSetSym(precond, sai_sym);
478 HYPRE_ParaSailsSetLogging(precond, 3);
479
480 /* Set the PCG preconditioner */
481 HYPRE_PCGSetPrecond(solver, (HYPRE_PtrToSolverFcn) HYPRE_ParaSailsSolve,
482 (HYPRE_PtrToSolverFcn) HYPRE_ParaSailsSetup, precond);
483
484 /* Now setup and solve! */
485 HYPRE_ParCSRPCGSetup(solver, parcsr_A, par_b, par_x);
486 HYPRE_ParCSRPCGSolve(solver, parcsr_A, par_b, par_x);
487
488
489 /* Run info - needed logging turned on */
490 HYPRE_PCGGetNumIterations(solver, &num_iterations);
491 HYPRE_PCGGetFinalRelativeResidualNorm(solver, &final_res_norm);
492 if (myid == 0)
493 {
494 printf("\n");
495 printf("Iterations = %d\n", num_iterations);
496 printf("Final Relative Residual Norm = %e\n", final_res_norm);
497 printf("\n");
498 }
499
500 /* Destory solver and preconditioner */
501 HYPRE_ParCSRPCGDestroy(solver);
502 HYPRE_ParaSailsDestroy(precond);
503 }
504 /* Flexible GMRES with AMG Preconditioner */
505 else if (solver_id == 61)
506 {
507 int num_iterations;
508 double final_res_norm;
509 int restart = 30;
510 int modify = 1;
511
512
513 /* Create solver */
514 HYPRE_ParCSRFlexGMRESCreate(MPI_COMM_WORLD, &solver);
515
516 /* Set some parameters (See Reference Manual for more parameters) */
517 HYPRE_FlexGMRESSetKDim(solver, restart);
518 HYPRE_FlexGMRESSetMaxIter(solver, 1000); /* max iterations */
519 HYPRE_FlexGMRESSetTol(solver, 1e-7); /* conv. tolerance */
520 HYPRE_FlexGMRESSetPrintLevel(solver, 2); /* print solve info */
521 HYPRE_FlexGMRESSetLogging(solver, 1); /* needed to get run info later */
522
523
524 /* Now set up the AMG preconditioner and specify any parameters */
525 HYPRE_BoomerAMGCreate(&precond);
526 HYPRE_BoomerAMGSetPrintLevel(precond, 1); /* print amg solution info */
527 HYPRE_BoomerAMGSetCoarsenType(precond, 6);
528 HYPRE_BoomerAMGSetOldDefault(precond);
529 HYPRE_BoomerAMGSetRelaxType(precond, 6); /* Sym G.S./Jacobi hybrid */
530 HYPRE_BoomerAMGSetNumSweeps(precond, 1);
531 HYPRE_BoomerAMGSetTol(precond, 0.0); /* conv. tolerance zero */
532 HYPRE_BoomerAMGSetMaxIter(precond, 1); /* do only one iteration! */
533
534 /* Set the FlexGMRES preconditioner */
535 HYPRE_FlexGMRESSetPrecond(solver, (HYPRE_PtrToSolverFcn) HYPRE_BoomerAMGSolve,
536 (HYPRE_PtrToSolverFcn) HYPRE_BoomerAMGSetup, precond);
537
538
539 if (modify)
540 /* this is an optional call - if you don't call it, hypre_FlexGMRESModifyPCDefault
541 is used - which does nothing. Otherwise, you can define your own, similar to
542 the one used here */
543 HYPRE_FlexGMRESSetModifyPC( solver,
544 (HYPRE_PtrToModifyPCFcn) hypre_FlexGMRESModifyPCAMGExample);
545
546
547 /* Now setup and solve! */
548 HYPRE_ParCSRFlexGMRESSetup(solver, parcsr_A, par_b, par_x);
549 HYPRE_ParCSRFlexGMRESSolve(solver, parcsr_A, par_b, par_x);
550
551 /* Run info - needed logging turned on */
552 HYPRE_FlexGMRESGetNumIterations(solver, &num_iterations);
553 HYPRE_FlexGMRESGetFinalRelativeResidualNorm(solver, &final_res_norm);
554 if (myid == 0)
555 {
556 printf("\n");
557 printf("Iterations = %d\n", num_iterations);
558 printf("Final Relative Residual Norm = %e\n", final_res_norm);
559 printf("\n");
560 }
561
562 /* Destory solver and preconditioner */
563 HYPRE_ParCSRFlexGMRESDestroy(solver);
564 HYPRE_BoomerAMGDestroy(precond);
565
566 }
567 else
568 {
569 if (myid ==0) printf("Invalid solver id specified.\n");
570 }
571
572 /* Save the solution for GLVis visualization, see vis/glvis-ex5.sh */
573 if (vis)
574 {
575 #ifdef HYPRE_EXVIS
576 FILE *file;
577 char filename[255];
578
579 int nvalues = local_size;
580 int *rows = (int*) calloc(nvalues, sizeof(int));
581 double *values = (double*) calloc(nvalues, sizeof(double));
582
583 for (i = 0; i < nvalues; i++)
584 rows[i] = ilower + i;
585
586 /* get the local solution */
587 HYPRE_IJVectorGetValues(x, nvalues, rows, values);
588
589 sprintf(filename, "%s.%06d", "vis/ex5.sol", myid);
590 if ((file = fopen(filename, "w")) == NULL)
591 {
592 printf("Error: can't open output file %s\n", filename);
593 MPI_Finalize();
594 exit(1);
595 }
596
597 /* save solution */
598 for (i = 0; i < nvalues; i++)
599 fprintf(file, "%.14e\n", values[i]);
600
601 fflush(file);
602 fclose(file);
603
604 free(rows);
605 free(values);
606
607 /* save global finite element mesh */
608 if (myid == 0)
609 GLVis_PrintGlobalSquareMesh("vis/ex5.mesh", n-1);
610 #endif
611 }
612
613 /* Clean up */
614 HYPRE_IJMatrixDestroy(A);
615 HYPRE_IJVectorDestroy(b);
616 HYPRE_IJVectorDestroy(x);
617
618 /* Finalize HYPRE */
619 HYPRE_Finalize();
620
621 /* Finalize MPI*/
622 MPI_Finalize();
623
624 return(0);
625 }
626
627 /*--------------------------------------------------------------------------
628 hypre_FlexGMRESModifyPCAMGExample -
629
630 This is an example (not recommended)
631 of how we can modify things about AMG that
632 affect the solve phase based on how FlexGMRES is doing...For
633 another preconditioner it may make sense to modify the tolerance..
634
635 *--------------------------------------------------------------------------*/
636
hypre_FlexGMRESModifyPCAMGExample(void * precond_data,int iterations,double rel_residual_norm)637 int hypre_FlexGMRESModifyPCAMGExample(void *precond_data, int iterations,
638 double rel_residual_norm)
639 {
640
641
642 if (rel_residual_norm > .1)
643 {
644 HYPRE_BoomerAMGSetNumSweeps((HYPRE_Solver)precond_data, 10);
645 }
646 else
647 {
648 HYPRE_BoomerAMGSetNumSweeps((HYPRE_Solver)precond_data, 1);
649 }
650
651
652 return 0;
653 }
654