download the original source code.
1 /*
2 Example 8
3
4 Interface: Semi-Structured interface (SStruct)
5
6 Compile with: make ex8
7
8 Sample run: mpirun -np 2 ex8
9
10 Description: This is a two processor example which solves a similar
11 problem to the one in Example 2, and Example 6 (The grid
12 boxes are exactly those in the example diagram in the
13 struct interface chapter of the User's Manual.)
14
15 The difference with the previous examples is that we use
16 three parts, two with a 5-point and one with a 9-point
17 discretization stencil. The solver is PCG with split-SMG
18 preconditioner.
19 */
20
21 #include <stdio.h>
22
23 /* SStruct linear solvers headers */
24 #include "HYPRE_sstruct_ls.h"
25
26 #include "vis.c"
27
28 int main (int argc, char *argv[])
29 {
30 int myid, num_procs;
31
32 int vis = 0;
33
34 HYPRE_SStructGrid grid;
35 HYPRE_SStructGraph graph;
36 HYPRE_SStructStencil stencil_5pt;
37 HYPRE_SStructStencil stencil_9pt;
38 HYPRE_SStructMatrix A;
39 HYPRE_SStructVector b;
40 HYPRE_SStructVector x;
41 HYPRE_SStructSolver solver;
42 HYPRE_SStructSolver precond;
43
44 int object_type;
45
46 /* Initialize MPI */
47 MPI_Init(&argc, &argv);
48 MPI_Comm_rank(MPI_COMM_WORLD, &myid);
49 MPI_Comm_size(MPI_COMM_WORLD, &num_procs);
50
51 if (num_procs != 2)
52 {
53 if (myid ==0) printf("Must run with 2 processors!\n");
54 MPI_Finalize();
55
56 return(0);
57 }
58
59 /* Parse command line */
60 {
61 int arg_index = 0;
62 int print_usage = 0;
63
64 while (arg_index < argc)
65 {
66 if ( strcmp(argv[arg_index], "-vis") == 0 )
67 {
68 arg_index++;
69 vis = 1;
70 }
71 else if ( strcmp(argv[arg_index], "-help") == 0 )
72 {
73 print_usage = 1;
74 break;
75 }
76 else
77 {
78 arg_index++;
79 }
80 }
81
82 if ((print_usage) && (myid == 0))
83 {
84 printf("\n");
85 printf("Usage: %s [<options>]\n", argv[0]);
86 printf("\n");
87 printf(" -vis : save the solution for GLVis visualization\n");
88 printf("\n");
89 }
90
91 if (print_usage)
92 {
93 MPI_Finalize();
94 return (0);
95 }
96 }
97
98 /* 1. Set up the 2D grid. This gives the index space in each part.
99 We have one variable in each part. */
100 {
101 int ndim = 2;
102 int nparts = 3;
103 int part;
104
105 /* Create an empty 2D grid object */
106 HYPRE_SStructGridCreate(MPI_COMM_WORLD, ndim, nparts, &grid);
107
108 /* Set the extents of the grid - each processor sets its grid
109 boxes. Each part has its own relative index space numbering. */
110
111 /* Processor 0 owns two boxes - one in part 0 and one in part 1. */
112 if (myid == 0)
113 {
114 /* Add the first box to the grid in part 0 */
115 {
116 int ilower[2] = {-3, 1};
117 int iupper[2] = {-1, 2};
118
119 part = 0;
120 HYPRE_SStructGridSetExtents(grid, part, ilower, iupper);
121 }
122
123 /* Add the second box to the grid in part 1 */
124 {
125 /* For convenience we use the same index space across all
126 parts, but this is not a requirement. For example, on this
127 part we could have used ilower=[23,24] and iupper=[25,27]. */
128 int ilower[2] = {0, 1};
129 int iupper[2] = {2, 4};
130
131 part = 1;
132 HYPRE_SStructGridSetExtents(grid, part, ilower, iupper);
133 }
134 }
135
136 /* Processor 1 owns one box in part 2. */
137 else if (myid == 1)
138 {
139 /* Add a new box to the grid in part 2 */
140 {
141 int ilower[2] = {3, 1};
142 int iupper[2] = {6, 4};
143
144 part = 2;
145 HYPRE_SStructGridSetExtents(grid, part, ilower, iupper);
146 }
147 }
148
149 /* Set the variable type and number of variables on each part. */
150 {
151 int i;
152 int nvars = 1;
153 HYPRE_SStructVariable vartypes[1] = {HYPRE_SSTRUCT_VARIABLE_CELL};
154
155 for (i = 0; i< nparts; i++)
156 HYPRE_SStructGridSetVariables(grid, i, nvars, vartypes);
157 }
158
159 /* Now we need to set the spatial relation between each of the parts.
160 Since we have the same types of variables on both parts, we can
161 use HYPRE_GridSetNeighborPart(). Each processor calls this function
162 for each part on which it owns boxes that border a different part. */
163
164 if (myid == 0)
165 {
166 /* Relation between part 0 and part 1 on processor 0 */
167 {
168 int part = 0;
169 int nbor_part = 1;
170 /* Cells just outside of the boundary of part 0 in
171 its coordinates */
172 int b_ilower[2] = {0,1}, b_iupper[2] = {0,2};
173 /* The same cells in part 1's coordinates. Since we use the same
174 index space across all parts, the coordinates coincide. */
175 int nbor_ilower[2] = {0,1}, nbor_iupper[2] = {0,2};
176 /* These parts have the same orientation, so no
177 rotation is necessary */
178 int index_map[2] = {0,1};
179 /* These parts map increasing values to increasing values
180 for both variables (note: if decreasing maps to increasing, use -1)*/
181 int index_dir[2] = {1,1};
182
183 HYPRE_SStructGridSetNeighborPart(grid, part, b_ilower, b_iupper,
184 nbor_part, nbor_ilower, nbor_iupper,
185 index_map, index_dir);
186 }
187
188 /* Relation between part 1 and part 0 on processor 0 */
189 {
190 int part = 1;
191 int nbor_part = 0;
192 /* Cells just outside of the boundary of part 1 in
193 its coordinates */
194 int b_ilower[2] = {-1,1}, b_iupper[2] = {-1,2};
195 /* The same cells in part 0's coordinates. Since we use the same
196 index space across all parts, the coordinates coincide. */
197 int nbor_ilower[2] = {-1,1}, nbor_iupper[2] = {-1,2};
198 /* These parts have the same orientation, so no
199 rotation is necessary */
200 int index_map[2] = {0,1};
201 /* These parts map increasing values to increasing values
202 for both variables (note: if decreasing maps to increasing, use -1)*/
203 int index_dir[2] = {1,1};
204
205 HYPRE_SStructGridSetNeighborPart(grid, part, b_ilower, b_iupper,
206 nbor_part, nbor_ilower, nbor_iupper,
207 index_map, index_dir);
208 }
209
210 /* Relation between part 1 and part 2 on processor 0 */
211 {
212 int part = 1;
213 int nbor_part = 2;
214 /* Cells just outside of the boundary of part 1 in
215 its coordinates */
216 int b_ilower[2] = {3,1}, b_iupper[2] = {3,4};
217 /* The same cells in part 2's coordinates. Since we use the same
218 index space across all parts, the coordinates coincide. */
219 int nbor_ilower[2] = {3,1}, nbor_iupper[2] = {3,4};
220 /* These parts have the same orientation, so no
221 rotation is necessary */
222 int index_map[2] = {0,1};
223 /* These parts map increasing values to increasing values
224 for both variables (note: if decreasing maps to increasing, use -1)*/
225 int index_dir[2] = {1,1};
226
227 HYPRE_SStructGridSetNeighborPart(grid, part, b_ilower, b_iupper,
228 nbor_part, nbor_ilower, nbor_iupper,
229 index_map, index_dir);
230 }
231 }
232 else if (myid == 1)
233 {
234 /* Relation between part 2 and part 1 on processor 1 */
235 {
236 int part = 2;
237 int nbor_part = 1;
238 /* Cells just outside of the boundary of part 2 in
239 its coordinates */
240 int b_ilower[2] = {2,1}, b_iupper[2] = {2,4};
241 /* The same cells in part 1's coordinates. Since we use the same
242 index space across all parts, the coordinates coincide. */
243 int nbor_ilower[2] = {2,1}, nbor_iupper[2] = {2,4};
244 /* These parts have the same orientation, so no
245 rotation is necessary */
246 int index_map[2] = {0,1};
247 /* These parts map increasing values to increasing values
248 for both variables (note: if decreasing maps to increasing, use -1)*/
249 int index_dir[2] = {1,1};
250
251 HYPRE_SStructGridSetNeighborPart(grid, part, b_ilower, b_iupper,
252 nbor_part, nbor_ilower, nbor_iupper,
253 index_map, index_dir);
254 }
255 }
256
257 /* Now the grid is ready to use */
258 HYPRE_SStructGridAssemble(grid);
259 }
260
261 /* 2. Define the discretization stencils */
262 {
263 int ndim = 2;
264 int var = 0;
265 int entry;
266
267 /* the 5-pt stencil in 2D */
268 {
269 int offsets[5][2] = {{0,0}, {-1,0}, {1,0}, {0,-1}, {0,1}};
270 int stencil_size = 5;
271
272 HYPRE_SStructStencilCreate(ndim, stencil_size, &stencil_5pt);
273
274 for (entry = 0; entry < 5; entry++)
275 HYPRE_SStructStencilSetEntry(stencil_5pt, entry, offsets[entry], var);
276 }
277
278 /* the 9-pt stencil in 2D */
279 {
280 int offsets[9][2] = {{0,0}, {-1,0}, {1,0}, {0,-1}, {0,1},
281 {-1,-1}, {1,-1}, {1,1}, {-1,1}};
282 int stencil_size = 9;
283 HYPRE_SStructStencilCreate(ndim, stencil_size, &stencil_9pt);
284
285 for (entry = 0; entry < stencil_size; entry++)
286 HYPRE_SStructStencilSetEntry(stencil_9pt, entry, offsets[entry], var);
287 }
288 }
289
290 /* 3. Set up the Graph - this determines the non-zero structure
291 of the matrix and allows non-stencil relationships between the parts */
292 {
293 int var = 0;
294 int part;
295
296 /* Create the graph object */
297 HYPRE_SStructGraphCreate(MPI_COMM_WORLD, grid, &graph);
298
299 /* See MatrixSetObjectType below */
300 object_type = HYPRE_SSTRUCT;
301 HYPRE_SStructGraphSetObjectType(graph, object_type);
302
303 /* Use the 5-pt stencil on part 0 */
304 part = 0;
305 HYPRE_SStructGraphSetStencil(graph, part, var, stencil_5pt);
306
307 /* Use the 9-pt stencil on part 1 */
308 part = 1;
309 HYPRE_SStructGraphSetStencil(graph, part, var, stencil_9pt);
310
311 /* Use the 5-pt stencil on part 2 */
312 part = 2;
313 HYPRE_SStructGraphSetStencil(graph, part, var, stencil_5pt);
314
315 /* Since we have only stencil connections between parts, we don't need to
316 call HYPRE_SStructGraphAddEntries. */
317
318 /* Assemble the graph */
319 HYPRE_SStructGraphAssemble(graph);
320 }
321
322 /* 4. Set up a SStruct Matrix */
323 {
324 int i,j;
325 int part;
326 int var = 0;
327
328 /* Create the empty matrix object */
329 HYPRE_SStructMatrixCreate(MPI_COMM_WORLD, graph, &A);
330
331 /* Set the object type (by default HYPRE_SSTRUCT). This determines the
332 data structure used to store the matrix. If you want to use unstructured
333 solvers, e.g. BoomerAMG, the object type should be HYPRE_PARCSR.
334 If the problem is purely structured (with one part), you may want to use
335 HYPRE_STRUCT to access the structured solvers. Since we have two parts
336 with different stencils, we set the object type to HYPRE_SSTRUCT. */
337 object_type = HYPRE_SSTRUCT;
338 HYPRE_SStructMatrixSetObjectType(A, object_type);
339
340 /* Get ready to set values */
341 HYPRE_SStructMatrixInitialize(A);
342
343 /* Each processor must set the stencil values for their boxes on each part.
344 In this example, we only set stencil entries and therefore use
345 HYPRE_SStructMatrixSetBoxValues. If we need to set non-stencil entries,
346 we have to use HYPRE_SStructMatrixSetValues. */
347
348 if (myid == 0)
349 {
350 /* Set the matrix coefficients for some set of stencil entries
351 over all the gridpoints in my first box (account for boundary
352 grid points later) */
353 {
354 int ilower[2] = {-3, 1};
355 int iupper[2] = {-1, 2};
356
357 int nentries = 5;
358 int nvalues = 30; /* 6 grid points, each with 5 stencil entries */
359 double values[30];
360
361 int stencil_indices[5];
362 for (j = 0; j < nentries; j++) /* label the stencil indices -
363 these correspond to the offsets
364 defined above */
365 stencil_indices[j] = j;
366
367 for (i = 0; i < nvalues; i += nentries)
368 {
369 values[i] = 4.0;
370 for (j = 1; j < nentries; j++)
371 values[i+j] = -1.0;
372 }
373
374 part = 0;
375 HYPRE_SStructMatrixSetBoxValues(A, part, ilower, iupper,
376 var, nentries,
377 stencil_indices, values);
378 }
379
380 /* Set the matrix coefficients for some set of stencil entries
381 over the gridpoints in my second box */
382 {
383 int ilower[2] = {0, 1};
384 int iupper[2] = {2, 4};
385
386 int nentries = 9;
387 int nvalues = 108; /* 12 grid points, each with 5 stencil entries */
388 double values[108];
389
390 int stencil_indices[9];
391 for (j = 0; j < nentries; j++)
392 stencil_indices[j] = j;
393
394 for (i = 0; i < nvalues; i += nentries)
395 {
396 values[i] = 8./3.;
397 for (j = 1; j < nentries; j++)
398 values[i+j] = -1./3.;
399 }
400
401 part = 1;
402 HYPRE_SStructMatrixSetBoxValues(A, part, ilower, iupper,
403 var, nentries,
404 stencil_indices, values);
405 }
406 }
407 else if (myid == 1)
408 {
409 /* Set the matrix coefficients for some set of stencil entries
410 over the gridpoints in my box */
411 {
412 int ilower[2] = {3, 1};
413 int iupper[2] = {6, 4};
414
415 int nentries = 5;
416 int nvalues = 80; /* 16 grid points, each with 5 stencil entries */
417 double values[80];
418
419 int stencil_indices[5];
420 for (j = 0; j < nentries; j++)
421 stencil_indices[j] = j;
422
423 for (i = 0; i < nvalues; i += nentries)
424 {
425 values[i] = 4.0;
426 for (j = 1; j < nentries; j++)
427 values[i+j] = -1.0;
428 }
429
430 part = 2;
431 HYPRE_SStructMatrixSetBoxValues(A, part, ilower, iupper,
432 var, nentries,
433 stencil_indices, values);
434 }
435 }
436
437 /* Modify the 9-pt stencil on the boundary between parts to ensure
438 symmetry and good global approximation. */
439 if (myid == 0)
440 {
441 int nentries = 6;
442 int nvalues = 24; /* 4 grid points, each with 6 stencil entries */
443 double values[24];
444
445 part = 1;
446
447 for (i = 0; i < nvalues; i += nentries)
448 {
449 values[i] = 10./3.;
450 values[i+1] = -1.;
451 values[i+2] = -2./3.;
452 values[i+3] = -2./3.;
453 values[i+4] = 0.0;
454 values[i+5] = 0.0;
455 }
456
457 {
458 /* Values to the right of the second box */
459 int ilower[2] = { 2, 1};
460 int iupper[2] = { 2, 4};
461
462 int stencil_indices[6] = {0,2,3,4,6,7};
463
464 HYPRE_SStructMatrixSetBoxValues(A, part, ilower, iupper,
465 var, nentries,
466 stencil_indices, values);
467 }
468
469 {
470 /* Values to the left of the second box */
471 int ilower[2] = { 0, 1};
472 int iupper[2] = { 0, 4};
473
474 int stencil_indices[6] = {0,1,3,4,5,8};
475
476 HYPRE_SStructMatrixSetBoxValues(A, part, ilower, iupper,
477 var, nentries,
478 stencil_indices, values);
479 }
480 }
481
482 /* For each box, set any coefficients that reach ouside of the
483 boundary to 0 */
484 if (myid == 0)
485 {
486 int maxnvalues = 9;
487 double values[9];
488
489 for (i = 0; i < maxnvalues; i++)
490 values[i] = 0.0;
491
492 part = 0;
493
494 {
495 /* Values below our first box */
496 int ilower[2] = {-3, 1};
497 int iupper[2] = {-1, 1};
498
499 int stencil_indices[1] = {3};
500
501 HYPRE_SStructMatrixSetBoxValues(A, part, ilower, iupper,
502 var, 1,
503 stencil_indices, values);
504 }
505
506 {
507 /* Values to the left of our first box */
508 int ilower[2] = {-3, 1};
509 int iupper[2] = {-3, 2};
510
511 int stencil_indices[1] = {1};
512
513 HYPRE_SStructMatrixSetBoxValues(A, part, ilower, iupper,
514 var, 1,
515 stencil_indices, values);
516 }
517
518 {
519 /* Values above our first box */
520 int ilower[2] = {-3, 2};
521 int iupper[2] = {-1, 2};
522
523 int stencil_indices[1] = {4};
524
525 HYPRE_SStructMatrixSetBoxValues(A, part, ilower, iupper,
526 var, 1,
527 stencil_indices, values);
528 }
529
530 part = 1;
531
532 {
533 /* Values below our second box */
534 int ilower[2] = { 0, 1};
535 int iupper[2] = { 2, 1};
536
537 int stencil_indices[3] = {3,5,6};
538
539 HYPRE_SStructMatrixSetBoxValues(A, part, ilower, iupper,
540 var, 3,
541 stencil_indices, values);
542 }
543
544 {
545 /* Values to the left of our second box (that do not border the
546 first box). */
547 int ilower[2] = { 0, 3};
548 int iupper[2] = { 0, 4};
549
550 int stencil_indices[3] = {1,5,8};
551
552 HYPRE_SStructMatrixSetBoxValues(A, part, ilower, iupper,
553 var, 3,
554 stencil_indices, values);
555 }
556
557 {
558 /* Values above our second box */
559 int ilower[2] = { 0, 4};
560 int iupper[2] = { 2, 4};
561
562 int stencil_indices[3] = {4,7,8};
563
564 HYPRE_SStructMatrixSetBoxValues(A, part, ilower, iupper,
565 var, 3,
566 stencil_indices, values);
567 }
568 }
569 else if (myid == 1)
570 {
571 int maxnvalues = 4;
572 double values[4];
573
574 for (i = 0; i < maxnvalues; i++)
575 values[i] = 0.0;
576
577 part = 2;
578
579 {
580 /* Values below our box */
581 int ilower[2] = { 3, 1};
582 int iupper[2] = { 6, 1};
583
584 int stencil_indices[1] = {3};
585
586 HYPRE_SStructMatrixSetBoxValues(A, part, ilower, iupper,
587 var, 1,
588 stencil_indices, values);
589 }
590
591 {
592 /* Values to the right of our box */
593 int ilower[2] = { 6, 1};
594 int iupper[2] = { 6, 4};
595
596 int stencil_indices[1] = {2};
597
598 HYPRE_SStructMatrixSetBoxValues(A, part, ilower, iupper,
599 var, 1,
600 stencil_indices, values);
601 }
602
603 {
604 /* Values above our box */
605 int ilower[2] = { 3, 4};
606 int iupper[2] = { 6, 4};
607
608 int stencil_indices[1] = {4};
609
610 HYPRE_SStructMatrixSetBoxValues(A, part, ilower, iupper,
611 var, 1,
612 stencil_indices, values);
613 }
614 }
615
616 /* This is a collective call finalizing the matrix assembly.
617 The matrix is now ``ready to be used'' */
618 HYPRE_SStructMatrixAssemble(A);
619 }
620
621 /* 5. Set up SStruct Vectors for b and x */
622 {
623 int i;
624 int part;
625 int var = 0;
626
627 /* Create an empty vector object */
628 HYPRE_SStructVectorCreate(MPI_COMM_WORLD, grid, &b);
629 HYPRE_SStructVectorCreate(MPI_COMM_WORLD, grid, &x);
630
631 /* As with the matrix, set the object type for the vectors
632 to be the sstruct type */
633 object_type = HYPRE_SSTRUCT;
634 HYPRE_SStructVectorSetObjectType(b, object_type);
635 HYPRE_SStructVectorSetObjectType(x, object_type);
636
637 /* Indicate that the vector coefficients are ready to be set */
638 HYPRE_SStructVectorInitialize(b);
639 HYPRE_SStructVectorInitialize(x);
640
641 if (myid == 0)
642 {
643 /* Set the vector coefficients over the gridpoints in my first box */
644 {
645 int ilower[2] = {-3, 1};
646 int iupper[2] = {-1, 2};
647
648 int nvalues = 6; /* 6 grid points */
649 double values[6];
650
651 part = 0;
652
653 for (i = 0; i < nvalues; i ++)
654 values[i] = 1.0;
655 HYPRE_SStructVectorSetBoxValues(b, part, ilower, iupper, var, values);
656
657 for (i = 0; i < nvalues; i ++)
658 values[i] = 0.0;
659 HYPRE_SStructVectorSetBoxValues(x, part, ilower, iupper, var, values);
660 }
661
662 /* Set the vector coefficients over the gridpoints in my second box */
663 {
664 int ilower[2] = { 0, 1};
665 int iupper[2] = { 2, 4};
666
667 int nvalues = 12; /* 12 grid points */
668 double values[12];
669
670 part = 1;
671
672 for (i = 0; i < nvalues; i ++)
673 values[i] = 1.0;
674 HYPRE_SStructVectorSetBoxValues(b, part, ilower, iupper, var, values);
675
676 for (i = 0; i < nvalues; i ++)
677 values[i] = 0.0;
678 HYPRE_SStructVectorSetBoxValues(x, part, ilower, iupper, var, values);
679 }
680 }
681 else if (myid == 1)
682 {
683 /* Set the vector coefficients over the gridpoints in my box */
684 {
685 int ilower[2] = { 3, 1};
686 int iupper[2] = { 6, 4};
687
688 int nvalues = 16; /* 16 grid points */
689 double values[16];
690
691 part = 2;
692
693 for (i = 0; i < nvalues; i ++)
694 values[i] = 1.0;
695 HYPRE_SStructVectorSetBoxValues(b, part, ilower, iupper, var, values);
696
697 for (i = 0; i < nvalues; i ++)
698 values[i] = 0.0;
699 HYPRE_SStructVectorSetBoxValues(x, part, ilower, iupper, var, values);
700 }
701 }
702
703 /* This is a collective call finalizing the vector assembly.
704 The vectors are now ``ready to be used'' */
705 HYPRE_SStructVectorAssemble(b);
706 HYPRE_SStructVectorAssemble(x);
707 }
708
709 /* 6. Set up and use a solver (See the Reference Manual for descriptions
710 of all of the options.) */
711 {
712 /* Create an empty PCG Struct solver */
713 HYPRE_SStructPCGCreate(MPI_COMM_WORLD, &solver);
714
715 /* Set PCG parameters */
716 HYPRE_SStructPCGSetTol(solver, 1.0e-6 );
717 HYPRE_SStructPCGSetPrintLevel(solver, 2);
718 HYPRE_SStructPCGSetMaxIter(solver, 50);
719
720 /* Create a split SStruct solver for use as a preconditioner */
721 HYPRE_SStructSplitCreate(MPI_COMM_WORLD, &precond);
722 HYPRE_SStructSplitSetMaxIter(precond, 1);
723 HYPRE_SStructSplitSetTol(precond, 0.0);
724 HYPRE_SStructSplitSetZeroGuess(precond);
725
726 /* Set the preconditioner type to split-SMG */
727 HYPRE_SStructSplitSetStructSolver(precond, HYPRE_SMG);
728
729 /* Set preconditioner and solve */
730 HYPRE_SStructPCGSetPrecond(solver, HYPRE_SStructSplitSolve,
731 HYPRE_SStructSplitSetup, precond);
732 HYPRE_SStructPCGSetup(solver, A, b, x);
733 HYPRE_SStructPCGSolve(solver, A, b, x);
734 }
735
736 /* Save the solution for GLVis visualization, see vis/glvis-ex8.sh */
737 if (vis)
738 {
739 GLVis_PrintSStructGrid(grid, "vis/ex8.mesh", myid, NULL, NULL);
740 GLVis_PrintSStructVector(x, 0, "vis/ex8.sol", myid);
741 GLVis_PrintData("vis/ex8.data", myid, num_procs);
742 }
743
744 /* Free memory */
745 HYPRE_SStructGridDestroy(grid);
746 HYPRE_SStructStencilDestroy(stencil_5pt);
747 HYPRE_SStructStencilDestroy(stencil_9pt);
748 HYPRE_SStructGraphDestroy(graph);
749 HYPRE_SStructMatrixDestroy(A);
750 HYPRE_SStructVectorDestroy(b);
751 HYPRE_SStructVectorDestroy(x);
752
753 HYPRE_SStructPCGDestroy(solver);
754 HYPRE_SStructSplitDestroy(precond);
755
756 /* Finalize MPI */
757 MPI_Finalize();
758
759 return (0);
760 }
syntax highlighted by Code2HTML, v. 0.9.1