32int main(
int argc,
char *argv[])
34 MPI_Init(&argc, &argv);
38 MPI_Comm_size(MPI_COMM_WORLD, &num_procs);
39 MPI_Comm_rank(MPI_COMM_WORLD, &rank);
44 const char *source_mesh_file =
"../../data/inline-tri.mesh";
45 const char *destination_mesh_file =
"../../data/inline-quad.mesh";
47 int src_n_refinements = 0;
48 int dest_n_refinements = 0;
51 int source_fe_order = 1;
52 int dest_fe_order = 0;
53 bool visualization =
true;
55 int max_iterations = 30000;
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",
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");
80 if (source_fe_order == 0 && dest_fe_order != 0)
83 "Source fe order should not be 0 unless destination fe order is also 0!\n";
86 return MPI_Finalize();
89 ifstream imesh(source_mesh_file);
90 shared_ptr<Mesh> src_mesh, dest_mesh;
93 src_mesh = make_shared<Mesh>(imesh, 1, 1);
99 mfem::err <<
"WARNING: Source mesh file not found: " << source_mesh_file
101 <<
"Using default 2D triangle mesh.";
105 imesh.open(destination_mesh_file);
108 dest_mesh = make_shared<Mesh>(imesh, 1, 1);
114 mfem::err <<
"WARNING: Destination mesh file not found: "
115 << destination_mesh_file <<
"\n"
116 <<
"Using default 2D quad mesh.";
122 for (
int i = 0; i < src_n_refinements; ++i)
124 src_mesh->UniformRefinement();
127 for (
int i = 0; i < dest_n_refinements; ++i)
129 dest_mesh->UniformRefinement();
132 src_mesh->EnsureNCMesh();
133 dest_mesh->EnsureNCMesh();
136 for (
int l = 0; l < 4; l++)
138 src_mesh->RandomRefinement(0.1);
143 for (
int l = 0; l < 4; l++)
145 dest_mesh->RandomRefinement(0.1);
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);
153 make_shared<DG_FECollection>(source_fe_order, p_src_mesh->Dimension());
155 make_shared<ParFiniteElementSpace>(p_src_mesh.get(), src_fe_coll.get());
158 make_shared<DG_FECollection>(dest_fe_order, p_dest_mesh->Dimension());
160 make_shared<ParFiniteElementSpace>(p_dest_mesh.get(), dest_fe_coll.get());
175 if (assembler.
Transfer(src_fun, dest_fun))
185 mfem::out <<
"l2 error: src: " << src_err <<
", dest: " << dest_err
189 plot(*p_src_mesh, src_fun,
"source");
190 plot(*p_dest_mesh, dest_fun,
"destination");
195 mfem::out <<
"Transfer failed! Use --verbose option for diagnostic!" <<
201 return MPI_Finalize();
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...
Class for parallel grid function.
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).
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...
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)
void destination_transform(const Vector &x, Vector &x_new)
OutStream out(std::cout)
Global stream used by the library for standard output. Initially it uses the same std::streambuf as s...
OutStream err(std::cerr)
Global stream used by the library for standard error output. Initially it uses the same std::streambu...
void InitTransfer(int argc, char *argv[])
Initializes the par_moonolith library. It also calls MPI_Init.
int FinalizeTransfer()
Finalize the par_moonolith library.