download the original source code.
1 /*
2 Example 9
3
4 Interface: Semi-Structured interface (SStruct)
5
6 Compile with: make ex9
7
8 Sample run: mpirun -np 16 ex9 -n 33 -solver 0 -v 1 1
9
10 To see options: ex9 -help
11
12 Description: This code solves a system corresponding to a discretization
13 of the biharmonic problem treated as a system of equations
14 on the unit square. Specifically, instead of solving
15 Delta^2(u) = f with zero boundary conditions for u and
16 Delta(u), we solve the system A x = b, where
17
18 A = [ Delta -I ; 0 Delta], x = [ u ; v] and b = [ 0 ; f]
19
20 The corresponding boundary conditions are u = 0 and v = 0.
21
22 The domain is split into an N x N processor grid. Thus, the
23 given number of processors should be a perfect square.
24 Each processor's piece of the grid has n x n cells with n x n
25 nodes. We use cell-centered variables, and, therefore, the
26 nodes are not shared. Note that we have two variables, u and
27 v, and need only one part to describe the domain. We use the
28 standard 5-point stencil to discretize the Laplace operators.
29 The boundary conditions are incorporated as in Example 3.
30
31 We recommend viewing Examples 3, 6 and 7 before this example.
32 */
33
34 #include <math.h>
35 #include "_hypre_utilities.h"
36 #include "HYPRE_sstruct_ls.h"
37 #include "HYPRE_krylov.h"
38
39 int main (int argc, char *argv[])
40 {
41 int i, j;
42
43 int myid, num_procs;
44
45 int n, N, pi, pj;
46 double h, h2;
47 int ilower[2], iupper[2];
48
49 int solver_id;
50 int n_pre, n_post;
51
52 int print_solution;
53 int object_type;
54
55 HYPRE_SStructGrid grid;
56 HYPRE_SStructGraph graph;
57 HYPRE_SStructStencil stencil_v;
58 HYPRE_SStructStencil stencil_u;
59 HYPRE_SStructMatrix A;
60 HYPRE_SStructVector b;
61 HYPRE_SStructVector x;
62
63 /* sstruct solvers */
64 HYPRE_SStructSolver solver;
65 HYPRE_SStructSolver precond;
66
67 /* parcsr solvers */
68 HYPRE_Solver par_solver;
69 HYPRE_Solver par_precond;
70
71 /* Initialize MPI */
72 MPI_Init(&argc, &argv);
73 MPI_Comm_rank(MPI_COMM_WORLD, &myid);
74 MPI_Comm_size(MPI_COMM_WORLD, &num_procs);
75
76 /* Set defaults */
77 n = 33;
78 solver_id = 0;
79 n_pre = 1;
80 n_post = 1;
81 print_solution = 0;
82
83 /* Parse command line */
84 {
85 int arg_index = 0;
86 int print_usage = 0;
87
88 while (arg_index < argc)
89 {
90 if ( strcmp(argv[arg_index], "-n") == 0 )
91 {
92 arg_index++;
93 n = atoi(argv[arg_index++]);
94 }
95 else if ( strcmp(argv[arg_index], "-solver") == 0 )
96 {
97 arg_index++;
98 solver_id = atoi(argv[arg_index++]);
99 }
100 else if ( strcmp(argv[arg_index], "-v") == 0 )
101 {
102 arg_index++;
103 n_pre = atoi(argv[arg_index++]);
104 n_post = atoi(argv[arg_index++]);
105 }
106 else if ( strcmp(argv[arg_index], "-print_solution") == 0 )
107 {
108 arg_index++;
109 print_solution = 1;
110 }
111 else if ( strcmp(argv[arg_index], "-help") == 0 )
112 {
113 print_usage = 1;
114 break;
115 }
116 else
117 {
118 arg_index++;
119 }
120 }
121
122 if ((print_usage) && (myid == 0))
123 {
124 printf("\n");
125 printf("Usage: %s [<options>]\n", argv[0]);
126 printf("\n");
127 printf(" -n <n> : problem size per processor (default: 33)\n");
128 printf(" -solver <ID> : solver ID\n");
129 printf(" 0 - GMRES with sysPFMG precond (default)\n");
130 printf(" 1 - sysPFMG\n");
131 printf(" 2 - GMRES with AMG precond\n");
132 printf(" 3 - AMG\n");
133 printf(" -v <n_pre> <n_post> : number of pre and post relaxations for SysPFMG (default: 1 1)\n");
134 printf(" -print_solution : print the solution vector\n");
135 printf("\n");
136 }
137
138 if (print_usage)
139 {
140 MPI_Finalize();
141 return (0);
142 }
143 }
144
145 /* Figure out the processor grid (N x N). The local problem
146 size for the interior nodes is indicated by n (n x n).
147 pi and pj indicate position in the processor grid. */
148 N = sqrt(num_procs);
149 h = 1.0 / (N*n+1); /* note that when calculating h we must
150 remember to count the boundary nodes */
151 h2 = h*h;
152 pj = myid / N;
153 pi = myid - pj*N;
154
155 /* Figure out the extents of each processor's piece of the grid. */
156 ilower[0] = pi*n;
157 ilower[1] = pj*n;
158
159 iupper[0] = ilower[0] + n-1;
160 iupper[1] = ilower[1] + n-1;
161
162 /* 1. Set up a grid - we have one part and two variables */
163 {
164 int nparts = 1;
165 int part = 0;
166 int ndim = 2;
167
168 /* Create an empty 2D grid object */
169 HYPRE_SStructGridCreate(MPI_COMM_WORLD, ndim, nparts, &grid);
170
171 /* Add a new box to the grid */
172 HYPRE_SStructGridSetExtents(grid, part, ilower, iupper);
173
174 /* Set the variable type and number of variables on each part.*/
175 {
176 int i;
177 int nvars = 2;
178 HYPRE_SStructVariable vartypes[2] = {HYPRE_SSTRUCT_VARIABLE_CELL,
179 HYPRE_SSTRUCT_VARIABLE_CELL };
180
181 for (i = 0; i< nparts; i++)
182 HYPRE_SStructGridSetVariables(grid, i, nvars, vartypes);
183 }
184
185 /* This is a collective call finalizing the grid assembly.
186 The grid is now ``ready to be used'' */
187 HYPRE_SStructGridAssemble(grid);
188 }
189
190 /* 2. Define the discretization stencils */
191 {
192 int entry;
193 int stencil_size;
194 int var;
195 int ndim = 2;
196
197 /* Stencil object for variable u (labeled as variable 0) */
198 {
199 int offsets[6][2] = {{0,0}, {-1,0}, {1,0}, {0,-1}, {0,1}, {0,0}};
200 stencil_size = 6;
201
202 HYPRE_SStructStencilCreate(ndim, stencil_size, &stencil_u);
203
204 /* The first 5 entries are for the u-u connections */
205 var = 0; /* connect to variable 0 */
206 for (entry = 0; entry < stencil_size-1 ; entry++)
207 HYPRE_SStructStencilSetEntry(stencil_u, entry, offsets[entry], var);
208
209 /* The last entry is for the u-v connection */
210 var = 1; /* connect to variable 1 */
211 entry = 5;
212 HYPRE_SStructStencilSetEntry(stencil_u, entry, offsets[entry], var);
213 }
214
215 /* Stencil object for variable v (variable 1) */
216 {
217 int offsets[5][2] = {{0,0}, {-1,0}, {1,0}, {0,-1}, {0,1}};
218 stencil_size = 5;
219
220 HYPRE_SStructStencilCreate(ndim, stencil_size, &stencil_v);
221
222 /* These are all v-v connections */
223 var = 1; /* Connect to variable 1 */
224 for (entry = 0; entry < stencil_size; entry++)
225 HYPRE_SStructStencilSetEntry(stencil_v, entry, offsets[entry], var);
226 }
227 }
228
229 /* 3. Set up the Graph - this determines the non-zero structure
230 of the matrix and allows non-stencil relationships between the parts. */
231 {
232 int var;
233 int part = 0;
234
235 /* Create the graph object */
236 HYPRE_SStructGraphCreate(MPI_COMM_WORLD, grid, &graph);
237
238 /* Assign the u-stencil we created to variable u (variable 0) */
239 var = 0;
240 HYPRE_SStructGraphSetStencil(graph, part, var, stencil_u);
241
242 /* Assign the v-stencil we created to variable v (variable 1) */
243 var = 1;
244 HYPRE_SStructGraphSetStencil(graph, part, var, stencil_v);
245
246 /* Assemble the graph */
247 HYPRE_SStructGraphAssemble(graph);
248 }
249
250 /* 4. Set up the SStruct Matrix */
251 {
252 int nentries;
253 int nvalues;
254 int var;
255 int part = 0;
256
257 /* Create an empty matrix object */
258 HYPRE_SStructMatrixCreate(MPI_COMM_WORLD, graph, &A);
259
260 /* Set the object type (by default HYPRE_SSTRUCT). This determines the
261 data structure used to store the matrix. If you want to use
262 unstructured solvers, e.g. BoomerAMG, the object type should be
263 HYPRE_PARCSR. If the problem is purely structured (with one part), you
264 may want to use HYPRE_STRUCT to access the structured solvers. */
265 if (solver_id > 1 && solver_id < 4)
266 {
267 object_type = HYPRE_PARCSR;
268 }
269 else
270 {
271 object_type = HYPRE_SSTRUCT;
272 }
273 HYPRE_SStructMatrixSetObjectType(A, object_type);
274
275 /* Indicate that the matrix coefficients are ready to be set */
276 HYPRE_SStructMatrixInitialize(A);
277
278 /* Each processor must set the stencil values for their boxes on each part.
279 In this example, we only set stencil entries and therefore use
280 HYPRE_SStructMatrixSetBoxValues. If we need to set non-stencil entries,
281 we have to use HYPRE_SStructMatrixSetValues. */
282
283 /* First set the u-stencil entries. Note that
284 HYPRE_SStructMatrixSetBoxValues can only set values corresponding
285 to stencil entries for the same variable. Therefore, we must set the
286 entries for each variable within a stencil with separate function calls.
287 For example, below the u-u connections and u-v connections are handled
288 in separate calls. */
289 {
290 int i, j;
291 double *u_values;
292 int u_v_indices[1] = {5};
293 int u_u_indices[5] = {0, 1, 2, 3, 4};
294
295 var = 0; /* Set values for the u connections */
296
297 /* First the u-u connections */
298 nentries = 5;
299 nvalues = nentries*n*n;
300 u_values = calloc(nvalues, sizeof(double));
301
302 for (i = 0; i < nvalues; i += nentries)
303 {
304 u_values[i] = 4.0;
305 for (j = 1; j < nentries; j++)
306 u_values[i+j] = -1.0;
307 }
308
309 HYPRE_SStructMatrixSetBoxValues(A, part, ilower, iupper,
310 var, nentries,
311 u_u_indices, u_values);
312 free(u_values);
313
314 /* Next the u-v connections */
315 nentries = 1;
316 nvalues = nentries*n*n;
317 u_values = calloc(nvalues, sizeof(double));
318
319 for (i = 0; i < nvalues; i++)
320 {
321 u_values[i] = -h2;
322 }
323
324 HYPRE_SStructMatrixSetBoxValues(A, part, ilower, iupper,
325 var, nentries,
326 u_v_indices, u_values);
327
328 free(u_values);
329 }
330
331 /* Now set the v-stencil entries */
332 {
333 int i, j;
334 double *v_values;
335 int v_v_indices[5] = {0, 1, 2, 3, 4};
336
337 var = 1; /* the v connections */
338
339 /* the v-v connections */
340 nentries = 5;
341 nvalues = nentries*n*n;
342 v_values = calloc(nvalues, sizeof(double));
343
344 for (i = 0; i < nvalues; i += nentries)
345 {
346 v_values[i] = 4.0;
347 for (j = 1; j < nentries; j++)
348 v_values[i+j] = -1.0;
349 }
350
351 HYPRE_SStructMatrixSetBoxValues(A, part, ilower, iupper,
352 var, nentries,
353 v_v_indices, v_values);
354
355 free(v_values);
356
357 /* There are no v-u connections to set */
358 }
359 }
360
361 /* 5. Incorporate the zero boundary conditions: go along each edge of
362 the domain and set the stencil entry that reaches to the boundary
363 to zero.*/
364 {
365 int bc_ilower[2];
366 int bc_iupper[2];
367 int nentries = 1;
368 int nvalues = nentries*n; /* number of stencil entries times the length
369 of one side of my grid box */
370 int var;
371 double *values;
372 int stencil_indices[1];
373
374 int part = 0;
375
376 values = calloc(nvalues, sizeof(double));
377 for (j = 0; j < nvalues; j++)
378 values[j] = 0.0;
379
380 /* Recall: pi and pj describe position in the processor grid */
381 if (pj == 0)
382 {
383 /* Bottom row of grid points */
384 bc_ilower[0] = pi*n;
385 bc_ilower[1] = pj*n;
386
387 bc_iupper[0] = bc_ilower[0] + n-1;
388 bc_iupper[1] = bc_ilower[1];
389
390 stencil_indices[0] = 3;
391
392 /* Need to do this for u and for v */
393 var = 0;
394 HYPRE_SStructMatrixSetBoxValues(A, part, bc_ilower, bc_iupper,
395 var, nentries,
396 stencil_indices, values);
397
398 var = 1;
399 HYPRE_SStructMatrixSetBoxValues(A, part, bc_ilower, bc_iupper,
400 var, nentries,
401 stencil_indices, values);
402 }
403
404 if (pj == N-1)
405 {
406 /* upper row of grid points */
407 bc_ilower[0] = pi*n;
408 bc_ilower[1] = pj*n + n-1;
409
410 bc_iupper[0] = bc_ilower[0] + n-1;
411 bc_iupper[1] = bc_ilower[1];
412
413 stencil_indices[0] = 4;
414
415 /* Need to do this for u and for v */
416 var = 0;
417 HYPRE_SStructMatrixSetBoxValues(A, part, bc_ilower, bc_iupper,
418 var, nentries,
419 stencil_indices, values);
420
421 var = 1;
422 HYPRE_SStructMatrixSetBoxValues(A, part, bc_ilower, bc_iupper,
423 var, nentries,
424 stencil_indices, values);
425
426 }
427
428 if (pi == 0)
429 {
430 /* Left row of grid points */
431 bc_ilower[0] = pi*n;
432 bc_ilower[1] = pj*n;
433
434 bc_iupper[0] = bc_ilower[0];
435 bc_iupper[1] = bc_ilower[1] + n-1;
436
437 stencil_indices[0] = 1;
438
439 /* Need to do this for u and for v */
440 var = 0;
441 HYPRE_SStructMatrixSetBoxValues(A, part, bc_ilower, bc_iupper,
442 var, nentries,
443 stencil_indices, values);
444
445 var = 1;
446 HYPRE_SStructMatrixSetBoxValues(A, part, bc_ilower, bc_iupper,
447 var, nentries,
448 stencil_indices, values);
449 }
450
451 if (pi == N-1)
452 {
453 /* Right row of grid points */
454 bc_ilower[0] = pi*n + n-1;
455 bc_ilower[1] = pj*n;
456
457 bc_iupper[0] = bc_ilower[0];
458 bc_iupper[1] = bc_ilower[1] + n-1;
459
460 stencil_indices[0] = 2;
461
462 /* Need to do this for u and for v */
463 var = 0;
464 HYPRE_SStructMatrixSetBoxValues(A, part, bc_ilower, bc_iupper,
465 var, nentries,
466 stencil_indices, values);
467
468 var = 1;
469 HYPRE_SStructMatrixSetBoxValues(A, part, bc_ilower, bc_iupper,
470 var, nentries,
471 stencil_indices, values);
472 }
473
474 free(values);
475 }
476
477 /* This is a collective call finalizing the matrix assembly.
478 The matrix is now ``ready to be used'' */
479 HYPRE_SStructMatrixAssemble(A);
480
481 /* 5. Set up SStruct Vectors for b and x */
482 {
483 int nvalues = n*n;
484 double *values;
485 int part = 0;
486 int var;
487
488 values = calloc(nvalues, sizeof(double));
489
490 /* Create an empty vector object */
491 HYPRE_SStructVectorCreate(MPI_COMM_WORLD, grid, &b);
492 HYPRE_SStructVectorCreate(MPI_COMM_WORLD, grid, &x);
493
494 /* Set the object type for the vectors
495 to be the same as was already set for the matrix */
496 HYPRE_SStructVectorSetObjectType(b, object_type);
497 HYPRE_SStructVectorSetObjectType(x, object_type);
498
499 /* Indicate that the vector coefficients are ready to be set */
500 HYPRE_SStructVectorInitialize(b);
501 HYPRE_SStructVectorInitialize(x);
502
503 /* Set the values for b */
504 for (i = 0; i < nvalues; i ++)
505 values[i] = h2;
506 var = 1;
507 HYPRE_SStructVectorSetBoxValues(b, part, ilower, iupper, var, values);
508
509 for (i = 0; i < nvalues; i ++)
510 values[i] = 0.0;
511 var = 0;
512 HYPRE_SStructVectorSetBoxValues(b, part, ilower, iupper, var, values);
513
514 /* Set the values for the initial guess */
515 var = 0;
516 HYPRE_SStructVectorSetBoxValues(x, part, ilower, iupper, var, values);
517
518 var = 1;
519 HYPRE_SStructVectorSetBoxValues(x, part, ilower, iupper, var, values);
520
521 free(values);
522
523 /* This is a collective call finalizing the vector assembly.
524 The vector is now ``ready to be used'' */
525 HYPRE_SStructVectorAssemble(b);
526 HYPRE_SStructVectorAssemble(x);
527 }
528
529 /* 6. Set up and use a solver
530 (Solver options can be found in the Reference Manual.) */
531 {
532 double final_res_norm;
533 int its;
534
535 HYPRE_ParCSRMatrix par_A;
536 HYPRE_ParVector par_b;
537 HYPRE_ParVector par_x;
538
539 /* If we are using a parcsr solver, we need to get the object for the
540 matrix and vectors. */
541 if (object_type == HYPRE_PARCSR)
542 {
543 HYPRE_SStructMatrixGetObject(A, (void **) &par_A);
544 HYPRE_SStructVectorGetObject(b, (void **) &par_b);
545 HYPRE_SStructVectorGetObject(x, (void **) &par_x);
546 }
547
548 if (solver_id ==0 ) /* GMRES with SysPFMG - the default*/
549 {
550 HYPRE_SStructGMRESCreate(MPI_COMM_WORLD, &solver);
551
552 /* GMRES parameters */
553 HYPRE_SStructGMRESSetMaxIter(solver, 50 );
554 HYPRE_SStructGMRESSetTol(solver, 1.0e-06 );
555 HYPRE_SStructGMRESSetPrintLevel(solver, 2 ); /* print each GMRES
556 iteration */
557 HYPRE_SStructGMRESSetLogging(solver, 1);
558
559 /* use SysPFMG as precondititioner */
560 HYPRE_SStructSysPFMGCreate(MPI_COMM_WORLD, &precond);
561
562 /* Set sysPFMG parameters */
563 HYPRE_SStructSysPFMGSetTol(precond, 0.0);
564 HYPRE_SStructSysPFMGSetMaxIter(precond, 1);
565 HYPRE_SStructSysPFMGSetNumPreRelax(precond, n_pre);
566 HYPRE_SStructSysPFMGSetNumPostRelax(precond, n_post);
567 HYPRE_SStructSysPFMGSetPrintLevel(precond, 0);
568 HYPRE_SStructSysPFMGSetZeroGuess(precond);
569
570 /* Set the preconditioner*/
571 HYPRE_SStructGMRESSetPrecond(solver, HYPRE_SStructSysPFMGSolve,
572 HYPRE_SStructSysPFMGSetup, precond);
573 /* do the setup */
574 HYPRE_SStructGMRESSetup(solver, A, b, x);
575
576 /* do the solve */
577 HYPRE_SStructGMRESSolve(solver, A, b, x);
578
579 /* get some info */
580 HYPRE_SStructGMRESGetFinalRelativeResidualNorm(solver,
581 &final_res_norm);
582 HYPRE_SStructGMRESGetNumIterations(solver, &its);
583
584 /* clean up */
585 HYPRE_SStructGMRESDestroy(solver);
586 }
587 else if (solver_id == 1) /* SysPFMG */
588 {
589 HYPRE_SStructSysPFMGCreate(MPI_COMM_WORLD, &solver);
590
591 /* Set sysPFMG parameters */
592 HYPRE_SStructSysPFMGSetTol(solver, 1.0e-6);
593 HYPRE_SStructSysPFMGSetMaxIter(solver, 50);
594 HYPRE_SStructSysPFMGSetNumPreRelax(solver, n_pre);
595 HYPRE_SStructSysPFMGSetNumPostRelax(solver, n_post);
596 HYPRE_SStructSysPFMGSetPrintLevel(solver, 0);
597 HYPRE_SStructSysPFMGSetLogging(solver, 1);
598
599 /* do the setup */
600 HYPRE_SStructSysPFMGSetup(solver, A, b, x);
601
602 /* do the solve */
603 HYPRE_SStructSysPFMGSolve(solver, A, b, x);
604
605 /* get some info */
606 HYPRE_SStructSysPFMGGetFinalRelativeResidualNorm(solver,
607 &final_res_norm);
608 HYPRE_SStructSysPFMGGetNumIterations(solver, &its);
609
610 /* clean up */
611 HYPRE_SStructSysPFMGDestroy(solver);
612 }
613 else if (solver_id == 2) /* GMRES with AMG */
614 {
615 HYPRE_ParCSRGMRESCreate(MPI_COMM_WORLD, &par_solver);
616
617 /* set the GMRES paramaters */
618 HYPRE_GMRESSetKDim(par_solver, 5);
619 HYPRE_GMRESSetMaxIter(par_solver, 100);
620 HYPRE_GMRESSetTol(par_solver, 1.0e-06);
621 HYPRE_GMRESSetPrintLevel(par_solver, 2);
622 HYPRE_GMRESSetLogging(par_solver, 1);
623
624 /* use BoomerAMG as preconditioner */
625 HYPRE_BoomerAMGCreate(&par_precond);
626 HYPRE_BoomerAMGSetCoarsenType(par_precond, 6);
627 HYPRE_BoomerAMGSetStrongThreshold(par_precond, 0.25);
628 HYPRE_BoomerAMGSetTol(par_precond, 0.0);
629 HYPRE_BoomerAMGSetPrintLevel(par_precond, 1);
630 HYPRE_BoomerAMGSetPrintFileName(par_precond, "ex9.out.log");
631 HYPRE_BoomerAMGSetMaxIter(par_precond, 1);
632
633 /* set the preconditioner */
634 HYPRE_ParCSRGMRESSetPrecond(par_solver,
635 HYPRE_BoomerAMGSolve,
636 HYPRE_BoomerAMGSetup,
637 par_precond);
638
639 /* do the setup */
640 HYPRE_ParCSRGMRESSetup(par_solver, par_A, par_b, par_x);
641
642 /* do the solve */
643 HYPRE_ParCSRGMRESSolve(par_solver, par_A, par_b, par_x);
644
645 /* get some info */
646 HYPRE_GMRESGetNumIterations(par_solver, &its);
647 HYPRE_GMRESGetFinalRelativeResidualNorm(par_solver,
648 &final_res_norm);
649 /* clean up */
650 HYPRE_ParCSRGMRESDestroy(par_solver);
651 }
652 else if (solver_id == 3) /* AMG */
653 {
654 HYPRE_BoomerAMGCreate(&par_solver);
655 HYPRE_BoomerAMGSetCoarsenType(par_solver, 6);
656 HYPRE_BoomerAMGSetStrongThreshold(par_solver, 0.25);
657 HYPRE_BoomerAMGSetTol(par_solver, 1.9e-6);
658 HYPRE_BoomerAMGSetPrintLevel(par_solver, 1);
659 HYPRE_BoomerAMGSetPrintFileName(par_solver, "ex9.out.log");
660 HYPRE_BoomerAMGSetMaxIter(par_solver, 50);
661
662 /* do the setup */
663 HYPRE_BoomerAMGSetup(par_solver, par_A, par_b, par_x);
664
665 /* do the solve */
666 HYPRE_BoomerAMGSolve(par_solver, par_A, par_b, par_x);
667
668 /* get some info */
669 HYPRE_BoomerAMGGetNumIterations(par_solver, &its);
670 HYPRE_BoomerAMGGetFinalRelativeResidualNorm(par_solver,
671 &final_res_norm);
672 /* clean up */
673 HYPRE_BoomerAMGDestroy(par_solver);
674 }
675 else
676 {
677 if (myid ==0) printf("\n ERROR: Invalid solver id specified.\n");
678 }
679
680 /* Gather the solution vector. This needs to be done if:
681 (1) the object type is parcsr OR
682 (2) any one of the variables is NOT cell-centered */
683 if (object_type == HYPRE_PARCSR)
684 {
685 HYPRE_SStructVectorGather(x);
686 }
687
688 /* Print the solution and other info */
689 if (print_solution)
690 {
691 HYPRE_SStructVectorPrint("sstruct.out.x", x, 0);
692 }
693
694 if (myid == 0)
695 {
696 printf("\n");
697 printf("Iterations = %d\n", its);
698 printf("Final Relative Residual Norm = %g\n", final_res_norm);
699 printf("\n");
700 }
701 }
702
703 /* Free memory */
704 HYPRE_SStructGridDestroy(grid);
705 HYPRE_SStructStencilDestroy(stencil_v);
706 HYPRE_SStructStencilDestroy(stencil_u);
707 HYPRE_SStructGraphDestroy(graph);
708 HYPRE_SStructMatrixDestroy(A);
709 HYPRE_SStructVectorDestroy(b);
710 HYPRE_SStructVectorDestroy(x);
711
712 /* Finalize MPI */
713 MPI_Finalize();
714
715 return (0);
716 }
syntax highlighted by Code2HTML, v. 0.9.1