MFEM v4.7.0
Finite element discretization library
Loading...
Searching...
No Matches
ex2p.cpp
Go to the documentation of this file.
1// MFEM + Moonolith Example 2 (parallel version)
2//
3// Compile with: make ex2p
4//
5// Moonolith sample runs:
6// mpirun -np 4 ex2p
7// mpirun -np 4 ex2p --source_refinements 1 --dest_refinements 2
8// mpirun -np 4 ex2p -s ../../data/inline-hex.mesh -d ../../data/inline-tet.mesh
9//
10// Description: This example code demonstrates the use of MFEM for transferring
11// discrete fields from one finite element mesh to another. The
12// meshes can be of arbitrary shape and completely unrelated with
13// each other. This feature can be used for implementing immersed
14// domain methods for fluid-structure interaction or general
15// multi-physics applications.
16//
17// This particular example concerns discontinuous Galerkin FEM with
18// adaptive mesh refinement for parallel runtimes.
19
20#include "example_utils.hpp"
21#include "mfem.hpp"
22
23using namespace mfem;
24using namespace std;
25
26void destination_transform(const Vector &x, Vector &x_new)
27{
28 x_new = x;
29 // x_new *= .5;
30}
31
32int main(int argc, char *argv[])
33{
34 MPI_Init(&argc, &argv);
35
36 int num_procs, rank;
37
38 MPI_Comm_size(MPI_COMM_WORLD, &num_procs);
39 MPI_Comm_rank(MPI_COMM_WORLD, &rank);
40
41 // Init transfer library context, with MPI handled outside the library
42 InitTransfer(argc, argv, MPI_COMM_WORLD);
43
44 const char *source_mesh_file = "../../data/inline-tri.mesh";
45 const char *destination_mesh_file = "../../data/inline-quad.mesh";
46
47 int src_n_refinements = 0;
48 int dest_n_refinements = 0;
49
50 // Source fe order has to be greater or equal than destination order
51 int source_fe_order = 1;
52 int dest_fe_order = 0;
53 bool visualization = true;
54 bool verbose = false;
55 int max_iterations = 30000;
56
57 OptionsParser args(argc, argv);
58 args.AddOption(&source_mesh_file, "-s", "--source_mesh",
59 "Mesh file to use for src.");
60 args.AddOption(&destination_mesh_file, "-d", "--destination_mesh",
61 "Mesh file to use for dest.");
62 args.AddOption(&src_n_refinements, "-sr", "--source_refinements",
63 "Number of src refinements");
64 args.AddOption(&dest_n_refinements, "-dr", "--dest_refinements",
65 "Number of dest refinements");
66 args.AddOption(&visualization, "-vis", "--visualization", "-no-vis",
67 "--no-visualization",
68 "Enable or disable GLVis visualization.");
69 args.AddOption(&source_fe_order, "-so", "--source_fe_order",
70 "Order of the src finite elements");
71 args.AddOption(&dest_fe_order, "-do", "--dest_fe_order",
72 "Order of the dest finite elements");
73 args.AddOption(&verbose, "-verb", "--verbose", "--no-verb", "--no-verbose",
74 "Enable/Disable verbose output");
75 args.AddOption(&max_iterations, "-m", "--max_iterations",
76 "Max number of solver iterations");
77 args.Parse();
78 check_options(args);
79
80 if (source_fe_order == 0 && dest_fe_order != 0)
81 {
82 mfem::out <<
83 "Source fe order should not be 0 unless destination fe order is also 0!\n";
84
86 return MPI_Finalize();
87 }
88
89 ifstream imesh(source_mesh_file);
90 shared_ptr<Mesh> src_mesh, dest_mesh;
91 if (imesh)
92 {
93 src_mesh = make_shared<Mesh>(imesh, 1, 1);
94 imesh.close();
95 }
96 else
97 {
98 if (rank == 0)
99 mfem::err << "WARNING: Source mesh file not found: " << source_mesh_file
100 << "\n"
101 << "Using default 2D triangle mesh.";
102 src_mesh = make_shared<Mesh>(4, 4, Element::TRIANGLE);
103 }
104
105 imesh.open(destination_mesh_file);
106 if (imesh)
107 {
108 dest_mesh = make_shared<Mesh>(imesh, 1, 1);
109 imesh.close();
110 }
111 else
112 {
113 if (rank == 0)
114 mfem::err << "WARNING: Destination mesh file not found: "
115 << destination_mesh_file << "\n"
116 << "Using default 2D quad mesh.";
117 dest_mesh = make_shared<Mesh>(4, 4, Element::QUADRILATERAL);
118 }
119
120 dest_mesh->Transform(&destination_transform);
121
122 for (int i = 0; i < src_n_refinements; ++i)
123 {
124 src_mesh->UniformRefinement();
125 }
126
127 for (int i = 0; i < dest_n_refinements; ++i)
128 {
129 dest_mesh->UniformRefinement();
130 }
131
132 src_mesh->EnsureNCMesh();
133 dest_mesh->EnsureNCMesh();
134
135 {
136 for (int l = 0; l < 4; l++)
137 {
138 src_mesh->RandomRefinement(0.1); // 10% probability
139 }
140 }
141
142 {
143 for (int l = 0; l < 4; l++)
144 {
145 dest_mesh->RandomRefinement(0.1); // 10% probability
146 }
147 }
148
149 auto p_src_mesh = make_shared<ParMesh>(MPI_COMM_WORLD, *src_mesh);
150 auto p_dest_mesh = make_shared<ParMesh>(MPI_COMM_WORLD, *dest_mesh);
151
152 auto src_fe_coll =
153 make_shared<DG_FECollection>(source_fe_order, p_src_mesh->Dimension());
154 auto src_fe =
155 make_shared<ParFiniteElementSpace>(p_src_mesh.get(), src_fe_coll.get());
156
157 auto dest_fe_coll =
158 make_shared<DG_FECollection>(dest_fe_order, p_dest_mesh->Dimension());
159 auto dest_fe =
160 make_shared<ParFiniteElementSpace>(p_dest_mesh.get(), dest_fe_coll.get());
161
162 ParGridFunction src_fun(src_fe.get());
164 make_fun(*src_fe, coeff, src_fun);
165
166 ParGridFunction dest_fun(dest_fe.get());
167 dest_fun = 0.0;
168 dest_fun.Update();
169
170 ParMortarAssembler assembler(src_fe, dest_fe);
171 assembler.SetVerbose(verbose);
172 assembler.SetMaxSolverIterations(max_iterations);
173
174 assembler.AddMortarIntegrator(make_shared<L2MortarIntegrator>());
175 if (assembler.Transfer(src_fun, dest_fun))
176 {
177 if (visualization)
178 {
179
180 const double src_err = src_fun.ComputeL2Error(coeff);
181 const double dest_err = dest_fun.ComputeL2Error(coeff);
182
183 if (rank == 0)
184 {
185 mfem::out << "l2 error: src: " << src_err << ", dest: " << dest_err
186 << std::endl;
187 }
188
189 plot(*p_src_mesh, src_fun, "source");
190 plot(*p_dest_mesh, dest_fun, "destination");
191 }
192 }
193 else
194 {
195 mfem::out << "Transfer failed! Use --verbose option for diagnostic!" <<
196 std::endl;
197 }
198
199 // Finalize transfer library context
201 return MPI_Finalize();
202}
A general function coefficient.
void Parse()
Parse the command-line options. Note that this function expects all the options provided through the ...
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 'var' to receive the value. Enable/disable tags are used to set the bool...
Definition optparser.hpp:82
Class for parallel grid function.
Definition pgridfunc.hpp:33
real_t ComputeL2Error(Coefficient *exsol[], const IntegrationRule *irs[]=NULL, const Array< int > *elems=NULL) const override
void Update() override
Transform by the Space UpdateMatrix (e.g., on Mesh change).
Definition pgridfunc.cpp:90
This class implements the parallel variational transfer between finite element spaces....
void SetMaxSolverIterations(const int max_solver_iterations)
Control the maximum numbers of conjugate gradients steps for mass matrix inversion.
bool Transfer(const ParGridFunction &src_fun, ParGridFunction &dest_fun)
transfer a function from source to destination. if the transfer is to be performed multiple times use...
void SetVerbose(const bool verbose)
Expose process details with verbose output.
void AddMortarIntegrator(const std::shared_ptr< MortarIntegrator > &integrator)
This method must be called before Assemble or Transfer. It will assemble the operator in all intersec...
Vector data type.
Definition vector.hpp:80
void plot(mfem::Mesh &mesh, mfem::GridFunction &x, std::string title)
void make_fun(mfem::FiniteElementSpace &fe, mfem::Coefficient &c, mfem::GridFunction &f)
double example_fun(const mfem::Vector &x)
void check_options(mfem::OptionsParser &args)
int main()
void destination_transform(const Vector &x, Vector &x_new)
Definition ex2p.cpp:26
OutStream out(std::cout)
Global stream used by the library for standard output. Initially it uses the same std::streambuf as s...
Definition globals.hpp:66
OutStream err(std::cerr)
Global stream used by the library for standard error output. Initially it uses the same std::streambu...
Definition globals.hpp:71
void InitTransfer(int argc, char *argv[])
Initializes the par_moonolith library. It also calls MPI_Init.
Definition transfer.cpp:19
int FinalizeTransfer()
Finalize the par_moonolith library.
Definition transfer.cpp:23