Files
mfem/examples/ex0.cpp
T

98 lines
3.1 KiB
C++

// MFEM Example 0
//
// Compile with: make ex0
//
// Sample runs: ex0
// ex0 -m ../data/fichera.mesh
// ex0 -m ../data/square-disc.mesh -o 2
//
// Description: This example code demonstrates the most basic usage of MFEM
// to define a simple finite element discretization of the Laplace
// problem -Delta u = 1 with homogeneous Dirichlet boundary
// conditions. The mesh file and finite element polynomial degree
// are given by command line options.
#include "mfem.hpp"
#include <fstream>
#include <iostream>
using namespace std;
using namespace mfem;
int main(int argc, char *argv[])
{
// 1. Parse command line options
const char *mesh_file = "../data/star.mesh";
int order = 1;
OptionsParser args(argc, argv);
args.AddOption(&mesh_file, "-m", "--mesh", "Mesh file to use.");
args.AddOption(&order, "-o", "--order",
"Finite element order (polynomial degree)");
args.Parse();
if (!args.Good())
{
args.PrintUsage(cout);
return 1;
}
args.PrintOptions(cout);
// 2. Read the mesh from the given mesh file, and refine once uniformly.
Mesh mesh(mesh_file);
mesh.UniformRefinement();
// 3. Define a finite element space on the mesh. Here we use continuous
// Lagrange finite elements of the specified order.
H1_FECollection fec(order, mesh.Dimension());
FiniteElementSpace fespace(&mesh, &fec);
cout << "Number of unknowns: " << fespace.GetTrueVSize() << endl;
// 4. Get a list of all the boundary DOFs. These will be marked as essential
// in order to enforce Dirichlet boundary conditions.
Array<int> boundary_dofs;
fespace.GetBoundaryTrueDofs(boundary_dofs);
// 5. Define the solution vector x as a finite element grid function
// corresponding to fespace. Initialize x with initial guess of zero,
// which also determines the boundary conditions.
GridFunction x(&fespace);
x = 0.0;
// 6. Set up the linear form b(.) corresponding to the right-hand side.
ConstantCoefficient one(1.0);
LinearForm b(&fespace);
b.AddDomainIntegrator(new DomainLFIntegrator(one));
b.Assemble();
// 7. Set up the bilinear form a(.,.) corresponding to the Laplacian -Delta.
BilinearForm a(&fespace);
a.AddDomainIntegrator(new DiffusionIntegrator);
a.Assemble();
// 8. Form the linear system A X = B. This includes eliminating boundary
// conditions, applying conforming constraints for non-conforming AMR,
// etc.
SparseMatrix A;
Vector B, X;
a.FormLinearSystem(boundary_dofs, x, b, A, X, B);
// 9. Solve using preconditioned CG with a symmetric Gauss-Seidel
// preconditioner.
GSSmoother M(A);
PCG(A, M, B, X, 1, 200, 1e-12, 0.0);
// 10. Recover the solution as a grid function and save to file. The output
// can be viewed using GLVis with the command:
// glvis -m mesh.mesh -g sol.gf
a.RecoverFEMSolution(X, b, x);
ofstream mesh_ofs("mesh.mesh");
mesh_ofs.precision(8);
mesh.Print(mesh_ofs);
ofstream sol_ofs("sol.gf");
sol_ofs.precision(8);
x.Save(sol_ofs);
return 0;
}