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 1
10 
11    Interface:    Structured interface (Struct)
12 
13    Compile with: make ex1 (may need to edit HYPRE_DIR in Makefile)
14 
15    Sample run:   mpirun -np 2 ex1
16 
17    Description:  This is a two processor example.  Each processor owns one
18                  box in the grid.  For reference, the two grid boxes are those
19                  in the example diagram in the struct interface chapter
20                  of the User's Manual. Note that in this example code, we have
21                  used the two boxes shown in the diagram as belonging
22                  to processor 0 (and given one box to each processor). The
23                  solver is PCG with no preconditioner.
24 
25                  We recommend viewing examples 1-4 sequentially for
26                  a nice overview/tutorial of the struct interface.
27 */
28 
29 /* Struct linear solvers header */
30 #include "HYPRE_struct_ls.h"
31 
32 #ifdef HYPRE_FORTRAN
33 #include "fortran.h"
34 #include "hypre_struct_fortran_test.h"
35 #endif
36 
main(HYPRE_Int argc,char * argv[])37 HYPRE_Int main (HYPRE_Int argc, char *argv[])
38 {
39    HYPRE_Int i, j, myid;
40 
41 #ifdef HYPRE_FORTRAN
42    hypre_F90_Obj grid;
43    hypre_F90_Obj stencil;
44    hypre_F90_Obj A;
45    hypre_F90_Obj b;
46    hypre_F90_Obj x;
47    hypre_F90_Obj solver;
48         HYPRE_Int temp_COMM;
49         HYPRE_Int one = 1;
50         HYPRE_Int two = 2;
51         HYPRE_Int five = 5;
52      HYPRE_Real tol = 1.e-6;
53 #else
54    HYPRE_StructGrid     grid;
55    HYPRE_StructStencil  stencil;
56    HYPRE_StructMatrix   A;
57    HYPRE_StructVector   b;
58    HYPRE_StructVector   x;
59    HYPRE_StructSolver   solver;
60 #endif
61 
62    /* Initialize MPI */
63    hypre_MPI_Init(&argc, &argv);
64    hypre_MPI_Comm_rank(hypre_MPI_COMM_WORLD, &myid);
65 
66    /* 1. Set up a grid. Each processor describes the piece
67       of the grid that it owns. */
68    {
69       /* Create an empty 2D grid object */
70 #ifdef HYPRE_FORTRAN
71       temp_COMM = (HYPRE_Int) hypre_MPI_COMM_WORLD;
72       HYPRE_StructGridCreate(&temp_COMM, &two, &grid);
73 #else
74       HYPRE_StructGridCreate(hypre_MPI_COMM_WORLD, 2, &grid);
75 #endif
76 
77       /* Add boxes to the grid */
78       if (myid == 0)
79       {
80          HYPRE_Int ilower[2]={-3,1}, iupper[2]={-1,2};
81 #ifdef HYPRE_FORTRAN
82          HYPRE_StructGridSetExtents(&grid, &ilower[0], &iupper[0]);
83 #else
84          HYPRE_StructGridSetExtents(grid, ilower, iupper);
85 #endif
86       }
87       else if (myid == 1)
88       {
89          HYPRE_Int ilower[2]={0,1}, iupper[2]={2,4};
90 #ifdef HYPRE_FORTRAN
91          HYPRE_StructGridSetExtents(&grid, &ilower[0], &iupper[0]);
92 #else
93          HYPRE_StructGridSetExtents(grid, ilower, iupper);
94 #endif
95       }
96 
97       /* This is a collective call finalizing the grid assembly.
98          The grid is now ``ready to be used'' */
99 #ifdef HYPRE_FORTRAN
100       HYPRE_StructGridAssemble(&grid);
101 #else
102       HYPRE_StructGridAssemble(grid);
103 #endif
104    }
105 
106    /* 2. Define the discretization stencil */
107    {
108       /* Create an empty 2D, 5-pt stencil object */
109 #ifdef HYPRE_FORTRAN
110       HYPRE_StructStencilCreate(&two, &five, &stencil);
111 #else
112       HYPRE_StructStencilCreate(2, 5, &stencil);
113 #endif
114 
115       /* Define the geometry of the stencil. Each represents a
116          relative offset (in the index space). */
117       {
118          HYPRE_Int entry;
119          HYPRE_Int offsets[5][2] = {{0,0}, {-1,0}, {1,0}, {0,-1}, {0,1}};
120 
121          /* Assign each of the 5 stencil entries */
122 #ifdef HYPRE_FORTRAN
123          for (entry = 0; entry < 5; entry++)
124             HYPRE_StructStencilSetElement(&stencil, &entry, offsets[entry]);
125 #else
126          for (entry = 0; entry < 5; entry++)
127             HYPRE_StructStencilSetElement(stencil, entry, offsets[entry]);
128 #endif
129       }
130    }
131 
132    /* 3. Set up a Struct Matrix */
133    {
134       /* Create an empty matrix object */
135 #ifdef HYPRE_FORTRAN
136       HYPRE_StructMatrixCreate(&temp_COMM, &grid, &stencil, &A);
137 #else
138       HYPRE_StructMatrixCreate(hypre_MPI_COMM_WORLD, grid, stencil, &A);
139 #endif
140 
141       /* Indicate that the matrix coefficients are ready to be set */
142 #ifdef HYPRE_FORTRAN
143       HYPRE_StructMatrixInitialize(&A);
144 #else
145       HYPRE_StructMatrixInitialize(A);
146 #endif
147 
148       /* Set the matrix coefficients.  Each processor assigns coefficients
149          for the boxes in the grid that it owns. Note that the coefficients
150          associated with each stencil entry may vary from grid point to grid
151          point if desired.  Here, we first set the same stencil entries for
152          each grid point.  Then we make modifications to grid points near
153          the boundary. */
154       if (myid == 0)
155       {
156          HYPRE_Int ilower[2]={-3,1}, iupper[2]={-1,2};
157          HYPRE_Int stencil_indices[5] = {0,1,2,3,4}; /* labels for the stencil entries -
158                                                   these correspond to the offsets
159                                                   defined above */
160          HYPRE_Int nentries = 5;
161          HYPRE_Int nvalues  = 30; /* 6 grid points, each with 5 stencil entries */
162          HYPRE_Real values[30];
163 
164          /* We have 6 grid points, each with 5 stencil entries */
165          for (i = 0; i < nvalues; i += nentries)
166          {
167             values[i] = 4.0;
168             for (j = 1; j < nentries; j++)
169                values[i+j] = -1.0;
170          }
171 
172 #ifdef HYPRE_FORTRAN
173          HYPRE_StructMatrixSetBoxValues(&A, &ilower[0], &iupper[0], &nentries,
174                                         &stencil_indices[0], &values[0]);
175 #else
176          HYPRE_StructMatrixSetBoxValues(A, ilower, iupper, nentries,
177                                         stencil_indices, values);
178 #endif
179       }
180       else if (myid == 1)
181       {
182          HYPRE_Int ilower[2]={0,1}, iupper[2]={2,4};
183          HYPRE_Int stencil_indices[5] = {0,1,2,3,4};
184          HYPRE_Int nentries = 5;
185          HYPRE_Int nvalues  = 60; /* 12 grid points, each with 5 stencil entries */
186          HYPRE_Real values[60];
187 
188          for (i = 0; i < nvalues; i += nentries)
189          {
190             values[i] = 4.0;
191             for (j = 1; j < nentries; j++)
192                values[i+j] = -1.0;
193          }
194 
195 #ifdef HYPRE_FORTRAN
196          HYPRE_StructMatrixSetBoxValues(&A, &ilower[0], &iupper[0], &nentries,
197                                         &stencil_indices[0], &values[0]);
198 #else
199          HYPRE_StructMatrixSetBoxValues(A, ilower, iupper, nentries,
200                                         stencil_indices, values);
201 #endif
202       }
203 
204       /* Set the coefficients reaching outside of the boundary to 0 */
205       if (myid == 0)
206       {
207          HYPRE_Real values[3];
208          for (i = 0; i < 3; i++)
209             values[i] = 0.0;
210          {
211             /* values below our box */
212             HYPRE_Int ilower[2]={-3,1}, iupper[2]={-1,1};
213             HYPRE_Int stencil_indices[1] = {3};
214 #ifdef HYPRE_FORTRAN
215             HYPRE_StructMatrixSetBoxValues(&A, &ilower[0], &iupper[0], &one,
216                                         &stencil_indices[0], &values[0]);
217 #else
218             HYPRE_StructMatrixSetBoxValues(A, ilower, iupper, 1,
219                                            stencil_indices, values);
220 #endif
221          }
222          {
223             /* values to the left of our box */
224             HYPRE_Int ilower[2]={-3,1}, iupper[2]={-3,2};
225             HYPRE_Int stencil_indices[1] = {1};
226 #ifdef HYPRE_FORTRAN
227             HYPRE_StructMatrixSetBoxValues(&A, &ilower[0], &iupper[0], &one,
228                                         &stencil_indices[0], &values[0]);
229 #else
230             HYPRE_StructMatrixSetBoxValues(A, ilower, iupper, 1,
231                                            stencil_indices, values);
232 #endif
233          }
234          {
235             /* values above our box */
236             HYPRE_Int ilower[2]={-3,2}, iupper[2]={-1,2};
237             HYPRE_Int stencil_indices[1] = {4};
238 #ifdef HYPRE_FORTRAN
239             HYPRE_StructMatrixSetBoxValues(&A, &ilower[0], &iupper[0], &one,
240                                         &stencil_indices[0], &values[0]);
241 #else
242             HYPRE_StructMatrixSetBoxValues(A, ilower, iupper, 1,
243                                            stencil_indices, values);
244 #endif
245          }
246       }
247       else if (myid == 1)
248       {
249          HYPRE_Real values[4];
250          for (i = 0; i < 4; i++)
251             values[i] = 0.0;
252          {
253             /* values below our box */
254             HYPRE_Int ilower[2]={0,1}, iupper[2]={2,1};
255             HYPRE_Int stencil_indices[1] = {3};
256 #ifdef HYPRE_FORTRAN
257             HYPRE_StructMatrixSetBoxValues(&A, &ilower[0], &iupper[0], &one,
258                                         &stencil_indices[0], &values[0]);
259 #else
260             HYPRE_StructMatrixSetBoxValues(A, ilower, iupper, 1,
261                                            stencil_indices, values);
262 #endif
263          }
264          {
265             /* values to the right of our box
266                (that do not border the other box on proc. 0) */
267             HYPRE_Int ilower[2]={2,1}, iupper[2]={2,4};
268             HYPRE_Int stencil_indices[1] = {2};
269 #ifdef HYPRE_FORTRAN
270             HYPRE_StructMatrixSetBoxValues(&A, &ilower[0], &iupper[0], &one,
271                                         &stencil_indices[0], &values[0]);
272 #else
273             HYPRE_StructMatrixSetBoxValues(A, ilower, iupper, 1,
274                                            stencil_indices, values);
275 #endif
276          }
277          {
278             /* values above our box */
279             HYPRE_Int ilower[2]={0,4}, iupper[2]={2,4};
280             HYPRE_Int stencil_indices[1] = {4};
281 #ifdef HYPRE_FORTRAN
282             HYPRE_StructMatrixSetBoxValues(&A, &ilower[0], &iupper[0], &one,
283                                         &stencil_indices[0], &values[0]);
284 #else
285             HYPRE_StructMatrixSetBoxValues(A, ilower, iupper, 1,
286                                            stencil_indices, values);
287 #endif
288          }
289          {
290             /* values to the left of our box */
291             HYPRE_Int ilower[2]={0,3}, iupper[2]={0,4};
292             HYPRE_Int stencil_indices[1] = {1};
293 #ifdef HYPRE_FORTRAN
294             HYPRE_StructMatrixSetBoxValues(&A, &ilower[0], &iupper[0], &one,
295                                         &stencil_indices[0], &values[0]);
296 #else
297             HYPRE_StructMatrixSetBoxValues(A, ilower, iupper, 1,
298                                            stencil_indices, values);
299 #endif
300          }
301       }
302 
303       /* This is a collective call finalizing the matrix assembly.
304          The matrix is now ``ready to be used'' */
305 #ifdef HYPRE_FORTRAN
306       HYPRE_StructMatrixAssemble(&A);
307 #else
308       HYPRE_StructMatrixAssemble(A);
309 #endif
310    }
311 
312    /* 4. Set up Struct Vectors for b and x.  Each processor sets the vectors
313       corresponding to its boxes. */
314    {
315       /* Create an empty vector object */
316 #ifdef HYPRE_FORTRAN
317       HYPRE_StructVectorCreate(&temp_COMM, &grid, &b);
318       HYPRE_StructVectorCreate(&temp_COMM, &grid, &x);
319 #else
320       HYPRE_StructVectorCreate(hypre_MPI_COMM_WORLD, grid, &b);
321       HYPRE_StructVectorCreate(hypre_MPI_COMM_WORLD, grid, &x);
322 #endif
323 
324       /* Indicate that the vector coefficients are ready to be set */
325 #ifdef HYPRE_FORTRAN
326       HYPRE_StructVectorInitialize(&b);
327       HYPRE_StructVectorInitialize(&x);
328 #else
329       HYPRE_StructVectorInitialize(b);
330       HYPRE_StructVectorInitialize(x);
331 #endif
332 
333       /* Set the vector coefficients */
334       if (myid == 0)
335       {
336          HYPRE_Int ilower[2]={-3,1}, iupper[2]={-1,2};
337          HYPRE_Real values[6]; /* 6 grid points */
338 
339          for (i = 0; i < 6; i ++)
340             values[i] = 1.0;
341 #ifdef HYPRE_FORTRAN
342          HYPRE_StructVectorSetBoxValues(&b, &ilower[0], &iupper[0], &values[0]);
343 #else
344          HYPRE_StructVectorSetBoxValues(b, ilower, iupper, values);
345 #endif
346 
347          for (i = 0; i < 6; i ++)
348             values[i] = 0.0;
349 #ifdef HYPRE_FORTRAN
350          HYPRE_StructVectorSetBoxValues(&x, &ilower[0], &iupper[0], &values[0]);
351 #else
352          HYPRE_StructVectorSetBoxValues(x, ilower, iupper, values);
353 #endif
354       }
355       else if (myid == 1)
356       {
357          HYPRE_Int ilower[2]={0,1}, iupper[2]={2,4};
358          HYPRE_Real values[12]; /* 12 grid points */
359 
360          for (i = 0; i < 12; i ++)
361             values[i] = 1.0;
362 #ifdef HYPRE_FORTRAN
363          HYPRE_StructVectorSetBoxValues(&b, &ilower[0], &iupper[0], &values[0]);
364 #else
365          HYPRE_StructVectorSetBoxValues(b, ilower, iupper, values);
366 #endif
367 
368          for (i = 0; i < 12; i ++)
369             values[i] = 0.0;
370 #ifdef HYPRE_FORTRAN
371          HYPRE_StructVectorSetBoxValues(&x, &ilower[0], &iupper[0], &values[0]);
372 #else
373          HYPRE_StructVectorSetBoxValues(x, ilower, iupper, values);
374 #endif
375       }
376 
377       /* This is a collective call finalizing the vector assembly.
378          The vectors are now ``ready to be used'' */
379 #ifdef HYPRE_FORTRAN
380       HYPRE_StructVectorAssemble(&b);
381       HYPRE_StructVectorAssemble(&x);
382 #else
383       HYPRE_StructVectorAssemble(b);
384       HYPRE_StructVectorAssemble(x);
385 #endif
386    }
387 
388    /* 5. Set up and use a solver (See the Reference Manual for descriptions
389       of all of the options.) */
390    {
391       /* Create an empty PCG Struct solver */
392 #ifdef HYPRE_FORTRAN
393       HYPRE_StructPCGCreate(&temp_COMM, &solver);
394 #else
395       HYPRE_StructPCGCreate(hypre_MPI_COMM_WORLD, &solver);
396 #endif
397 
398       /* Set some parameters */
399 #ifdef HYPRE_FORTRAN
400       HYPRE_StructPCGSetTol(&solver, &tol); /* convergence tolerance */
401       HYPRE_StructPCGSetPrintLevel(&solver, &two); /* amount of info. printed */
402 #else
403       HYPRE_StructPCGSetTol(solver, 1.0e-06); /* convergence tolerance */
404       HYPRE_StructPCGSetPrintLevel(solver, 2); /* amount of info. printed */
405 #endif
406 
407       /* Setup and solve */
408 #ifdef HYPRE_FORTRAN
409       HYPRE_StructPCGSetup(&solver, &A, &b, &x);
410       HYPRE_StructPCGSolve(&solver, &A, &b, &x);
411 #else
412       HYPRE_StructPCGSetup(solver, A, b, x);
413       HYPRE_StructPCGSolve(solver, A, b, x);
414 #endif
415    }
416 
417    /* Free memory */
418 #ifdef HYPRE_FORTRAN
419    HYPRE_StructGridDestroy(&grid);
420    HYPRE_StructStencilDestroy(&stencil);
421    HYPRE_StructMatrixDestroy(&A);
422    HYPRE_StructVectorDestroy(&b);
423    HYPRE_StructVectorDestroy(&x);
424    HYPRE_StructPCGDestroy(&solver);
425 #else
426    HYPRE_StructGridDestroy(grid);
427    HYPRE_StructStencilDestroy(stencil);
428    HYPRE_StructMatrixDestroy(A);
429    HYPRE_StructVectorDestroy(b);
430    HYPRE_StructVectorDestroy(x);
431    HYPRE_StructPCGDestroy(solver);
432 #endif
433 
434    /* Finalize MPI */
435    hypre_MPI_Finalize();
436 
437    return (0);
438 }
439