MFEM  v4.3.0 Finite element discretization library
ex0p.cpp
Go to the documentation of this file.
1 // MFEM Example 0 - Parallel Version
2 //
3 // Compile with: make ex0p
4 //
5 // Sample runs: mpirun -np 4 ex0p
6 // mpirun -np 4 ex0p -m ../data/fichera.mesh
7 // mpirun -np 4 ex0p -m ../data/square-disc.mesh -o 2
8 //
9 // Description: This example code demonstrates the most basic parallel usage of
10 // MFEM to define a simple finite element discretization of the
11 // Laplace problem -Delta u = 1 with zero Dirichlet boundary
12 // conditions. General 2D/3D serial mesh files and finite element
13 // polynomial degrees can be specified by command line options.
14
15 #include "mfem.hpp"
16 #include <fstream>
17 #include <iostream>
18
19 using namespace std;
20 using namespace mfem;
21
22 int main(int argc, char *argv[])
23 {
24  // 1. Initialize MPI
25  MPI_Session mpi(argc, argv);
26
27  // 2. Parse command line options
28  const char *mesh_file = "../data/star.mesh";
29  int order = 1;
30
31  OptionsParser args(argc, argv);
32  args.AddOption(&mesh_file, "-m", "--mesh", "Mesh file to use.");
33  args.AddOption(&order, "-o", "--order", "Finite element polynomial degree");
34  args.ParseCheck();
35
36  // 3. Read the serial mesh from the given mesh file.
37  Mesh serial_mesh(mesh_file);
38
39  // 4. Define a parallel mesh by a partitioning of the serial mesh. Refine
40  // this mesh once in parallel to increase the resolution.
41  ParMesh mesh(MPI_COMM_WORLD, serial_mesh);
42  serial_mesh.Clear(); // the serial mesh is no longer needed
43  mesh.UniformRefinement();
44
45  // 5. Define a finite element space on the mesh. Here we use H1 continuous
46  // high-order Lagrange finite elements of the given order.
47  H1_FECollection fec(order, mesh.Dimension());
48  ParFiniteElementSpace fespace(&mesh, &fec);
49  HYPRE_BigInt total_num_dofs = fespace.GlobalTrueVSize();
50  if (mpi.Root()) { cout << "Number of unknowns: " << total_num_dofs << endl; }
51
52  // 6. Extract the list of all the boundary DOFs. These will be marked as
53  // Dirichlet in order to enforce zero boundary conditions.
54  Array<int> boundary_dofs;
55  fespace.GetBoundaryTrueDofs(boundary_dofs);
56
57  // 7. Define the solution x as a finite element grid function in fespace. Set
58  // the initial guess to zero, which also sets the boundary conditions.
59  ParGridFunction x(&fespace);
60  x = 0.0;
61
62  // 8. Set up the linear form b(.) corresponding to the right-hand side.
63  ConstantCoefficient one(1.0);
64  ParLinearForm b(&fespace);
66  b.Assemble();
67
68  // 9. Set up the bilinear form a(.,.) corresponding to the -Delta operator.
69  ParBilinearForm a(&fespace);
71  a.Assemble();
72
73  // 10. Form the linear system A X = B. This includes eliminating boundary
74  // conditions, applying AMR constraints, parallel assembly, etc.
76  Vector B, X;
77  a.FormLinearSystem(boundary_dofs, x, b, A, X, B);
78
79  // 11. Solve the system using PCG with hypre's BoomerAMG preconditioner.
80  HypreBoomerAMG M(A);
81  CGSolver cg(MPI_COMM_WORLD);
82  cg.SetRelTol(1e-12);
83  cg.SetMaxIter(2000);
84  cg.SetPrintLevel(1);
85  cg.SetPreconditioner(M);
86  cg.SetOperator(A);
87  cg.Mult(B, X);
88
89  // 12. Recover the solution x as a grid function and save to file. The output
90  // can be viewed using GLVis as follows: "glvis -np <np> -m mesh -g sol"
91  a.RecoverFEMSolution(X, b, x);
92  x.Save("sol");
93  mesh.Save("mesh");
94
95  return 0;
96 }
Class for domain integration L(v) := (f, v)
Definition: lininteg.hpp:97
Definition: solvers.hpp:316
A coefficient that is constant across space and time.
Definition: coefficient.hpp:78
virtual void FormLinearSystem(const Array< int > &ess_tdof_list, Vector &x, Vector &b, OperatorHandle &A, Vector &X, Vector &B, int copy_interior=0)
Form the linear system A X = B, corresponding to this bilinear form and the linear form b(...
virtual void Mult(const Vector &b, Vector &x) const
Operator application: y=A(x).
Definition: solvers.cpp:616
virtual void Save(std::ostream &out) const
Definition: pgridfunc.cpp:841
Abstract parallel finite element space.
Definition: pfespace.hpp:28
virtual void Save(const char *fname, int precision=16) const
Definition: pmesh.cpp:4499
The BoomerAMG solver in hypre.
Definition: hypre.hpp:1387
Class for parallel linear form.
Definition: plinearform.hpp:26
A simple convenience class that calls MPI_Init() at construction and MPI_Finalize() at destruction...
void SetPrintLevel(int print_lvl)
Definition: solvers.cpp:71
double b
Definition: lissajous.cpp:42
void UniformRefinement(int i, const DSTable &, int *, int *, int *)
Definition: mesh.cpp:9143
void SetMaxIter(int max_it)
Definition: solvers.hpp:100
void Assemble(int skip_zeros=1)
Assemble the local matrix.
bool Root() const
Return true if WorldRank() == 0.
int Dimension() const
Definition: mesh.hpp:911
Adds new Domain Integrator. Assumes ownership of lfi.
Definition: linearform.cpp:39
void SetRelTol(double rtol)
Definition: solvers.hpp:98
void AddOption(bool *var, const char *enable_short_name, const char *enable_long_name, const char *disable_short_name, const char *disable_long_name, const char *description, bool required=false)
Add a boolean option and set &#39;var&#39; to receive the value. Enable/disable tags are used to set the bool...
Definition: optparser.hpp:82
HYPRE_Int HYPRE_BigInt
virtual void RecoverFEMSolution(const Vector &X, const Vector &b, Vector &x)
double a
Definition: lissajous.cpp:41
Adds new Domain Integrator. Assumes ownership of bfi.
Class for parallel bilinear form.
void Clear()
Clear the contents of the Mesh.
Definition: mesh.hpp:828
virtual void SetOperator(const Operator &op)
Also calls SetOperator for the preconditioner if there is one.
Definition: solvers.hpp:330
Vector data type.
Definition: vector.hpp:60
virtual void SetPreconditioner(Solver &pr)
This should be called before SetOperator.
Definition: solvers.cpp:92
Arbitrary order H1-conforming (continuous) finite elements.
Definition: fe_coll.hpp:216
Class for parallel grid function.
Definition: pgridfunc.hpp:32
Wrapper for hypre&#39;s ParCSR matrix class.
Definition: hypre.hpp:277
Class for parallel meshes.
Definition: pmesh.hpp:32
int main()
void ParseCheck(std::ostream &out=mfem::out)
Definition: optparser.cpp:251