44int main(
int argc,
char *argv[])
53 const char *mesh_file =
"../data/beam-tet.mesh";
54 int ser_ref_levels = 2;
55 int par_ref_levels = 1;
59 bool visualization = 1;
60 const char *device_config =
"cpu";
63 args.
AddOption(&mesh_file,
"-m",
"--mesh",
65 args.
AddOption(&ser_ref_levels,
"-rs",
"--refine-serial",
66 "Number of times to refine the mesh uniformly in serial.");
67 args.
AddOption(&par_ref_levels,
"-rp",
"--refine-parallel",
68 "Number of times to refine the mesh uniformly in parallel.");
70 "Finite element order (polynomial degree) or -1 for"
71 " isoparametric space.");
73 "Number of desired eigenmodes.");
74 args.
AddOption(&nc,
"-nc",
"--non-conforming",
"-c",
76 "Mark the mesh as nonconforming before partitioning.");
77 args.
AddOption(&visualization,
"-vis",
"--visualization",
"-no-vis",
79 "Enable or disable GLVis visualization.");
80 args.
AddOption(&device_config,
"-d",
"--device",
81 "Device configuration string, see Device::Configure().");
98 Device device(device_config);
99 if (myid == 0) { device.
Print(); }
104 Mesh *mesh =
new Mesh(mesh_file, 1, 1);
114 for (
int lev = 0; lev < ser_ref_levels; lev++)
125 for (
int lev = 0; lev < par_ref_levels; lev++)
137 cout <<
"Number of unknowns: " << size << endl;
166 a->EliminateEssentialBCDiag(ess_bdr, 1.0);
209 ostringstream mesh_name, mode_name;
210 mesh_name <<
"mesh." << setfill(
'0') << setw(6) << myid;
212 ofstream mesh_ofs(mesh_name.str().c_str());
213 mesh_ofs.precision(8);
214 pmesh->
Print(mesh_ofs);
216 for (
int i=0; i<nev; i++)
221 mode_name <<
"mode_" << setfill(
'0') << setw(2) << i <<
"."
222 << setfill(
'0') << setw(6) << myid;
224 ofstream mode_ofs(mode_name.str().c_str());
225 mode_ofs.precision(8);
237 mode_sock.precision(8);
239 for (
int i=0; i<nev; i++)
243 cout <<
"Eigenmode " << i+1 <<
'/' << nev
244 <<
", Lambda = " << eigenvalues[i] << endl;
250 mode_sock <<
"parallel " << num_procs <<
" " << myid <<
"\n"
251 <<
"solution\n" << *pmesh << x << flush
252 <<
"window_title 'Eigenmode " << i+1 <<
'/' << nev
253 <<
", Lambda = " << eigenvalues[i] <<
"'" << endl;
258 cout <<
"press (q)uit or (c)ontinue --> " << flush;
261 MPI_Bcast(&c, 1, MPI_CHAR, 0, MPI_COMM_WORLD);
T Max() const
Find the maximal element in the array, using the comparison operator < for class T.
void SetSize(int nsize)
Change the logical size of the array, keep existing entries.
int Size() const
Return the logical size of the array.
A coefficient that is constant across space and time.
Integrator for for Nedelec elements.
The MFEM Device class abstracts hardware devices such as GPUs, as well as programming models such as ...
void Print(std::ostream &out=mfem::out)
Print the configuration of the MFEM virtual device object.
Collection of finite elements from the same family in multiple dimensions. This class is used to matc...
void SetNumModes(int num_eigs)
void SetPreconditioner(HypreSolver &precond)
void SetPrintLevel(int logging)
void SetMassMatrix(const HypreParMatrix &M)
const HypreParVector & GetEigenvector(unsigned int i) const
Extract a single eigenvector.
void GetEigenvalues(Array< real_t > &eigenvalues) const
Collect the converged eigenvalues.
void SetOperator(const HypreParMatrix &A)
void Solve()
Solve the eigenproblem.
void SetMaxIter(int max_iter)
The Auxiliary-space Maxwell Solver in hypre.
void SetPrintLevel(int print_lvl)
void SetSingularProblem()
Set this option when solving a curl-curl problem with zero mass term.
Wrapper for hypre's ParCSR matrix class.
static void Init()
Initialize hypre by calling HYPRE_Init() and set default options. After calling Hypre::Init(),...
Array< int > bdr_attributes
A list of all unique boundary attributes used by the Mesh.
int Dimension() const
Dimension of the reference space used within the elements.
void EnsureNCMesh(bool simplices_nonconforming=false)
void UniformRefinement(int i, const DSTable &, int *, int *, int *)
static int WorldRank()
Return the MPI rank in MPI_COMM_WORLD.
static int WorldSize()
Return the size of MPI_COMM_WORLD.
static void Init(int &argc, char **&argv, int required=default_thread_required, int *provided=nullptr)
Singleton creation with Mpi::Init(argc, argv).
Arbitrary order H(curl)-conforming Nedelec finite elements.
void Parse()
Parse the command-line options. Note that this function expects all the options provided through the ...
void PrintUsage(std::ostream &out) const
Print the usage message.
void PrintOptions(std::ostream &out) const
Print the options.
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...
bool Good() const
Return true if the command line options were parsed successfully.
Abstract parallel finite element space.
HYPRE_BigInt GlobalTrueVSize() const
Class for parallel grid function.
void Save(std::ostream &out) const override
Class for parallel meshes.
void Print(std::ostream &out=mfem::out, const std::string &comments="") const override
int close()
Close the socketstream.