CFD/hypre
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

Adding a nodal version of Example 1.

Browse Source
This commit is contained in:
falgout 2009-05-18 20:28:05 +00:00
parent 332c954e36
commit 5328b2857f
11 changed files with 467 additions and 5 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

8
examples/Makefile
View File

@ -67,7 +67,7 @@ LFLAGS90 =
########################################################################
# List of all programs to be compiled
########################################################################
ALLPROGS = ex1 ex2 ex3 ex4 ex5 ex5f ex6 ex7 ex8 ex9 ex10 ex11
ALLPROGS = ex1 ex2 ex3 ex4 ex5 ex5f ex6 ex7 ex8 ex9 ex10 ex11 ex12
BABELPROGS = ex5b ex5b77 ex5bxx ex6b ex6b77
# ... ex5bp and ex5b90 are not in BABELPROGS because they require a
# software environment which many people haven't set up.
@ -187,6 +187,12 @@ ex10: ex10.o
ex11: ex11.o
$(CXX) -o $@ $^ $(LFLAGS)
########################################################################
# Example 12
########################################################################
ex12: ex12.o
$(CC) -o $@ $^ $(LFLAGS)
########################################################################
# Clean up
########################################################################

1
examples/docs/amg.htm
View File

@ -5,6 +5,7 @@
<!--#include file = "ex5.html"-->
<!--#include file = "ex9.html"-->
<!--#include file = "ex12.html"-->
</blockquote>

1
examples/docs/c.htm
View File

@ -14,6 +14,7 @@
<!--#include file = "ex7.html"-->
<!--#include file = "ex8.html"-->
<!--#include file = "ex9.html"-->
<!--#include file = "ex12.html"-->
</blockquote>

10
examples/docs/ex12.htm Normal file
View File

@ -0,0 +1,10 @@
<a href="ex12.c.html"><h2>Example 12</h2></a>
<p>
The grid layout is the same as ex1, but with nodal unknowns. The solver is PCG
preconditioned with either PFMG or BoomerAMG, selected on the command line.
<p>
We recommend viewing the Struct examples before viewing this and the other
SStruct examples. This is one of the simplest SStruct examples, used primarily
to demonstrate how to set up non-cell-centered problems, and to demonstrate how
easy it is to switch between structured solvers (PFMG) and solvers designed for
more general settings (AMG).

1
examples/docs/index.htm
View File

@ -90,6 +90,7 @@ hypre team.
<li> <a href="README_files/ex9.html">Example 9</a>: a SStruct example for the biharmonic problem treated as a system of equations </li>
<li> <a href="README_files/ex10.html">Example 10</a>: FEI example with bilinear finite elements for the 2D Laplace equation </li>
<li> <a href="README_files/ex11.html">Example 11</a>: eigensolver for the 5-pt discretization of the 2D Laplace equation</li>
<li> <a href="README_files/ex12.html">Example 12</a>: nodal version of Example 1</li>
</ul>
<center>

1
examples/docs/pcg-amg.htm
View File

@ -5,6 +5,7 @@
<!--#include file = "ex5.html"-->
<!--#include file = "ex10.html"-->
<!--#include file = "ex12.html"-->
</blockquote>

1
examples/docs/pcg-pfmg.htm
View File

@ -5,6 +5,7 @@
<!--#include file = "ex4.html"-->
<!--#include file = "ex7.html"-->
<!--#include file = "ex12.html"-->
</blockquote>

1
examples/docs/pfmg.htm
View File

@ -5,6 +5,7 @@
<!--#include file = "ex4.html"-->
<!--#include file = "ex7.html"-->
<!--#include file = "ex12.html"-->
</blockquote>

1
examples/docs/sstruct.htm
View File

@ -7,6 +7,7 @@
<!--#include file = "ex7.html"-->
<!--#include file = "ex8.html"-->
<!--#include file = "ex9.html"-->
<!--#include file = "ex12.html"-->
</blockquote>

7
examples/ex1.c
View File

@ -48,7 +48,6 @@ int main (int argc, char *argv[])
return(0);
}
/* 1. Set up a grid. Each processor describes the piece
of the grid that it owns. */
{
@ -184,8 +183,7 @@ int main (int argc, char *argv[])
stencil_indices, values);
}
{
/* values to the right of our box
(that do not border the other box on proc. 0) */
/* values to the right of our box */
int ilower[2]={2,1}, iupper[2]={2,4};
int stencil_indices[1] = {2};
HYPRE_StructMatrixSetBoxValues(A, ilower, iupper, 1,
@ -199,7 +197,8 @@ int main (int argc, char *argv[])
stencil_indices, values);
}
{
/* values to the left of our box */
/* values to the left of our box
(that do not border the other box on proc. 0) */
int ilower[2]={0,3}, iupper[2]={0,4};
int stencil_indices[1] = {1};
HYPRE_StructMatrixSetBoxValues(A, ilower, iupper, 1,

440
examples/ex12.c Normal file
View File

@ -0,0 +1,440 @@
/*
Example 1
Interface: Semi-Structured interface (SStruct)
Compile with: make ex12 (may need to edit HYPRE_DIR in Makefile)
Sample runs: mpirun -np 2 ex12 -pfmg
mpirun -np 2 ex12 -boomeramg
Description: The grid layout is the same as ex1, but with nodal unknowns. The
solver is PCG preconditioned with either PFMG or BoomerAMG,
selected on the command line.
We recommend viewing the Struct examples before viewing this
and the other SStruct examples. This is one of the simplest
SStruct examples, used primarily to demonstrate how to set up
non-cell-centered problems and how easy it is to switch from
structured solvers (PFMG) and to solvers designed for more
general settings (AMG).
*/
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include "HYPRE_sstruct_ls.h"
#include "HYPRE_parcsr_ls.h"
#include "HYPRE_krylov.h"
int main (int argc, char *argv[])
{
int i, j, myid, num_procs;
HYPRE_SStructGrid grid;
HYPRE_SStructGraph graph;
HYPRE_SStructStencil stencil;
HYPRE_SStructMatrix A;
HYPRE_SStructVector b;
HYPRE_SStructVector x;
/* We only have one part and one variable */
int nparts = 1;
int nvars = 1;
int part = 0;
int var = 0;
int precond_id, object_type;
/* Initialize MPI */
MPI_Init(&argc, &argv);
MPI_Comm_rank(MPI_COMM_WORLD, &myid);
MPI_Comm_size(MPI_COMM_WORLD, &num_procs);
if (num_procs != 2)
{
if (myid == 0) printf("Must run with 2 processors!\n");
exit(1);
}
/* Parse the command line to determine the solver */
if (argc != 2)
{
if (myid == 0) printf("Must specify a solver!\n");
exit(1);
}
if ( strcmp(argv[1], "-pfmg") == 0 )
{
precond_id = 1;
object_type = HYPRE_STRUCT;
}
else if ( strcmp(argv[1], "-boomeramg") == 0 )
{
precond_id = 2;
object_type = HYPRE_PARCSR;
}
else
{
if (myid == 0) printf("Invalid solver!\n");
exit(1);
}
/* 1. Set up the grid. Here we use only one part. Each processor describes
the piece of the grid that it owns. */
{
/* Create an empty 2D grid object */
HYPRE_SStructGridCreate(MPI_COMM_WORLD, 2, nparts, &grid);
/* Add boxes to the grid */
if (myid == 0)
{
int ilower[2]={-3,1}, iupper[2]={-1,2};
HYPRE_SStructGridSetExtents(grid, part, ilower, iupper);
}
else if (myid == 1)
{
int ilower[2]={0,1}, iupper[2]={2,4};
HYPRE_SStructGridSetExtents(grid, part, ilower, iupper);
}
/* Set the variable type and number of variables on each part. */
{
HYPRE_SStructVariable vartypes[1] = {HYPRE_SSTRUCT_VARIABLE_NODE};
HYPRE_SStructGridSetVariables(grid, part, nvars, vartypes);
}
/* This is a collective call finalizing the grid assembly.
The grid is now ``ready to be used'' */
HYPRE_SStructGridAssemble(grid);
}
/* 2. Define the discretization stencil */
{
/* Create an empty 2D, 5-pt stencil object */
HYPRE_SStructStencilCreate(2, 5, &stencil);
/* Define the geometry of the stencil. Each represents a relative offset
(in the index space). */
{
int entry;
int offsets[5][2] = {{0,0}, {-1,0}, {1,0}, {0,-1}, {0,1}};
/* Assign numerical values to the offsets so that we can easily refer
to them - the last argument indicates the variable for which we are
assigning this stencil */
for (entry = 0; entry < 5; entry++)
HYPRE_SStructStencilSetEntry(stencil, entry, offsets[entry], var);
}
}
/* 3. Set up the Graph - this determines the non-zero structure of the matrix
and allows non-stencil relationships between the parts */
{
/* Create the graph object */
HYPRE_SStructGraphCreate(MPI_COMM_WORLD, grid, &graph);
/* Now we need to tell the graph which stencil to use for each variable on
each part (we only have one variable and one part) */
HYPRE_SStructGraphSetStencil(graph, part, var, stencil);
/* Here we could establish connections between parts if we had more than
one part using the graph. For example, we could use
HYPRE_GraphAddEntries() routine or HYPRE_GridSetNeighborPart() */
/* Assemble the graph */
HYPRE_SStructGraphAssemble(graph);
}
/* 4. Set up a SStruct Matrix */
{
/* Create an empty matrix object */
HYPRE_SStructMatrixCreate(MPI_COMM_WORLD, graph, &A);
/* Set the object type (by default HYPRE_SSTRUCT). This determines the
data structure used to store the matrix. For PFMG we need to use
HYPRE_STRUCT, and for BoomerAMG we need HYPRE_PARCSR (set above). */
HYPRE_SStructMatrixSetObjectType(A, object_type);
/* Get ready to set values */
HYPRE_SStructMatrixInitialize(A);
/* Set the matrix coefficients. Each processor assigns coefficients for
the boxes in the grid that it owns. Note that the coefficients
associated with each stencil entry may vary from grid point to grid
point if desired. Here, we first set the same stencil entries for each
grid point. Then we make modifications to grid points near the
boundary. Note that the ilower values are different from those used in
ex1 because of the way nodal variables are referenced. Also note that
some of the stencil values are set on both processor 0 and processor 1.
See the User and Reference manuals for more details. */
if (myid == 0)
{
int ilower[2]={-4,0}, iupper[2]={-1,2};
int stencil_indices[5] = {0,1,2,3,4}; /* labels for the stencil entries -
these correspond to the offsets
defined above */
int nentries = 5;
int nvalues = 60; /* 12 grid points, each with 5 stencil entries */
double values[60];
for (i = 0; i < nvalues; i += nentries)
{
values[i] = 4.0;
for (j = 1; j < nentries; j++)
values[i+j] = -1.0;
}
HYPRE_SStructMatrixSetBoxValues(A, part, ilower, iupper, var, nentries,
stencil_indices, values);
}
else if (myid == 1)
{
int ilower[2]={-1,0}, iupper[2]={2,4};
int stencil_indices[5] = {0,1,2,3,4};
int nentries = 5;
int nvalues = 100; /* 20 grid points, each with 5 stencil entries */
double values[100];
for (i = 0; i < nvalues; i += nentries)
{
values[i] = 4.0;
for (j = 1; j < nentries; j++)
values[i+j] = -1.0;
}
HYPRE_SStructMatrixSetBoxValues(A, part, ilower, iupper, var, nentries,
stencil_indices, values);
}
/* Set the coefficients reaching outside of the boundary to 0. Note that
* both ilower *and* iupper may be different from those in ex1. */
if (myid == 0)
{
double values[4];
for (i = 0; i < 4; i++)
values[i] = 0.0;
{
/* values below our box */
int ilower[2]={-4,0}, iupper[2]={-1,0};
int stencil_indices[1] = {3};
HYPRE_SStructMatrixSetBoxValues(A, part, ilower, iupper, var, 1,
stencil_indices, values);
}
{
/* values to the left of our box */
int ilower[2]={-4,0}, iupper[2]={-4,2};
int stencil_indices[1] = {1};
HYPRE_SStructMatrixSetBoxValues(A, part, ilower, iupper, var, 1,
stencil_indices, values);
}
{
/* values above our box */
int ilower[2]={-4,2}, iupper[2]={-2,2};
int stencil_indices[1] = {4};
HYPRE_SStructMatrixSetBoxValues(A, part, ilower, iupper, var, 1,
stencil_indices, values);
}
}
else if (myid == 1)
{
double values[5];
for (i = 0; i < 5; i++)
values[i] = 0.0;
{
/* values below our box */
int ilower[2]={-1,0}, iupper[2]={2,0};
int stencil_indices[1] = {3};
HYPRE_SStructMatrixSetBoxValues(A, part, ilower, iupper, var, 1,
stencil_indices, values);
}
{
/* values to the right of our box */
int ilower[2]={2,0}, iupper[2]={2,4};
int stencil_indices[1] = {2};
HYPRE_SStructMatrixSetBoxValues(A, part, ilower, iupper, var, 1,
stencil_indices, values);
}
{
/* values above our box */
int ilower[2]={-1,4}, iupper[2]={2,4};
int stencil_indices[1] = {4};
HYPRE_SStructMatrixSetBoxValues(A, part, ilower, iupper, var, 1,
stencil_indices, values);
}
{
/* values to the left of our box
(that do not border the other box on proc. 0) */
int ilower[2]={-1,3}, iupper[2]={-1,4};
int stencil_indices[1] = {1};
HYPRE_SStructMatrixSetBoxValues(A, part, ilower, iupper, var, 1,
stencil_indices, values);
}
}
/* This is a collective call finalizing the matrix assembly.
The matrix is now ``ready to be used'' */
HYPRE_SStructMatrixAssemble(A);
}
/* 5. Set up SStruct Vectors for b and x. */
{
/* Create an empty vector object */
HYPRE_SStructVectorCreate(MPI_COMM_WORLD, grid, &b);
HYPRE_SStructVectorCreate(MPI_COMM_WORLD, grid, &x);
/* As with the matrix, set the appropriate object type for the vectors */
HYPRE_SStructVectorSetObjectType(b, object_type);
HYPRE_SStructVectorSetObjectType(x, object_type);
/* Indicate that the vector coefficients are ready to be set */
HYPRE_SStructVectorInitialize(b);
HYPRE_SStructVectorInitialize(x);
/* Set the vector coefficients. Again, note that the ilower values are
different from those used in ex1, and some of the values are set on
both processors. */
if (myid == 0)
{
int ilower[2]={-4,0}, iupper[2]={-1,2};
double values[12]; /* 12 grid points */
for (i = 0; i < 12; i ++)
values[i] = 1.0;
HYPRE_SStructVectorSetBoxValues(b, part, ilower, iupper, var, values);
for (i = 0; i < 12; i ++)
values[i] = 0.0;
HYPRE_SStructVectorSetBoxValues(x, part, ilower, iupper, var, values);
}
else if (myid == 1)
{
int ilower[2]={0,1}, iupper[2]={2,4};
double values[20]; /* 20 grid points */
for (i = 0; i < 20; i ++)
values[i] = 1.0;
HYPRE_SStructVectorSetBoxValues(b, part, ilower, iupper, var, values);
for (i = 0; i < 20; i ++)
values[i] = 0.0;
HYPRE_SStructVectorSetBoxValues(x, part, ilower, iupper, var, values);
}
/* This is a collective call finalizing the vector assembly.
The vectors are now ``ready to be used'' */
HYPRE_SStructVectorAssemble(b);
HYPRE_SStructVectorAssemble(x);
}
/* 6. Set up and use a solver (See the Reference Manual for descriptions
of all of the options.) */
if (precond_id == 1) /* PFMG */
{
HYPRE_StructMatrix sA;
HYPRE_StructVector sb;
HYPRE_StructVector sx;
HYPRE_StructSolver solver;
HYPRE_StructSolver precond;
/* Because we are using a struct solver, we need to get the
object of the matrix and vectors to pass in to the struct solvers */
HYPRE_SStructMatrixGetObject(A, (void **) &sA);
HYPRE_SStructVectorGetObject(b, (void **) &sb);
HYPRE_SStructVectorGetObject(x, (void **) &sx);
/* Create an empty PCG Struct solver */
HYPRE_StructPCGCreate(MPI_COMM_WORLD, &solver);
/* Set PCG parameters */
HYPRE_StructPCGSetTol(solver, 1.0e-06);
HYPRE_StructPCGSetPrintLevel(solver, 2);
HYPRE_StructPCGSetMaxIter(solver, 50);
/* Create the Struct PFMG solver for use as a preconditioner */
HYPRE_StructPFMGCreate(MPI_COMM_WORLD, &precond);
/* Set PFMG parameters */
HYPRE_StructPFMGSetMaxIter(precond, 1);
HYPRE_StructPFMGSetTol(precond, 0.0);
HYPRE_StructPFMGSetZeroGuess(precond);
HYPRE_StructPFMGSetNumPreRelax(precond, 2);
HYPRE_StructPFMGSetNumPostRelax(precond, 2);
/* non-Galerkin coarse grid (more efficient for this problem) */
HYPRE_StructPFMGSetRAPType(precond, 1);
/* R/B Gauss-Seidel */
HYPRE_StructPFMGSetRelaxType(precond, 2);
/* skip relaxation on some levels (more efficient for this problem) */
HYPRE_StructPFMGSetSkipRelax(precond, 1);
/* Set preconditioner and solve */
HYPRE_StructPCGSetPrecond(solver, HYPRE_StructPFMGSolve,
HYPRE_StructPFMGSetup, precond);
HYPRE_StructPCGSetup(solver, sA, sb, sx);
HYPRE_StructPCGSolve(solver, sA, sb, sx);
/* Free memory */
HYPRE_StructPCGDestroy(solver);
HYPRE_StructPFMGDestroy(precond);
}
else if (precond_id == 2) /* BoomerAMG */
{
HYPRE_ParCSRMatrix parA;
HYPRE_ParVector parb;
HYPRE_ParVector parx;
HYPRE_Solver solver;
HYPRE_Solver precond;
/* Because we are using a struct solver, we need to get the
object of the matrix and vectors to pass in to the struct solvers */
HYPRE_SStructMatrixGetObject(A, (void **) &parA);
HYPRE_SStructVectorGetObject(b, (void **) &parb);
HYPRE_SStructVectorGetObject(x, (void **) &parx);
/* Create an empty PCG Struct solver */
HYPRE_ParCSRPCGCreate(MPI_COMM_WORLD, &solver);
/* Set PCG parameters */
HYPRE_ParCSRPCGSetTol(solver, 1.0e-06);
HYPRE_ParCSRPCGSetPrintLevel(solver, 2);
HYPRE_ParCSRPCGSetMaxIter(solver, 50);
/* Create the BoomerAMG solver for use as a preconditioner */
HYPRE_BoomerAMGCreate(&precond);
/* Set BoomerAMG parameters */
HYPRE_BoomerAMGSetMaxIter(precond, 1);
HYPRE_BoomerAMGSetTol(precond, 0.0);
HYPRE_BoomerAMGSetPrintLevel(precond, 1); /* print amg solution info */
HYPRE_BoomerAMGSetCoarsenType(precond, 6);
HYPRE_BoomerAMGSetRelaxType(precond, 6); /* Sym G.S./Jacobi hybrid */
HYPRE_BoomerAMGSetNumSweeps(precond, 1);
/* Set preconditioner and solve */
HYPRE_ParCSRPCGSetPrecond(solver, HYPRE_BoomerAMGSolve,
HYPRE_BoomerAMGSetup, precond);
HYPRE_ParCSRPCGSetup(solver, parA, parb, parx);
HYPRE_ParCSRPCGSolve(solver, parA, parb, parx);
/* Free memory */
HYPRE_ParCSRPCGDestroy(solver);
HYPRE_BoomerAMGDestroy(precond);
}
/* Free memory */
HYPRE_SStructGridDestroy(grid);
HYPRE_SStructStencilDestroy(stencil);
HYPRE_SStructGraphDestroy(graph);
HYPRE_SStructMatrixDestroy(A);
HYPRE_SStructVectorDestroy(b);
HYPRE_SStructVectorDestroy(x);
/* Finalize MPI */
MPI_Finalize();
return (0);
}