42int main(
int argc,
char *argv[])
51 const char *mesh_file =
"../data/star.mesh";
53 bool visualization = 1;
56 args.
AddOption(&mesh_file,
"-m",
"--mesh",
59 "Finite element order (polynomial degree).");
60 args.
AddOption(&visualization,
"-vis",
"--visualization",
"-no-vis",
62 "Enable or disable GLVis visualization.");
80 Mesh *mesh =
new Mesh(mesh_file, 1, 1);
89 (int)floor(log(10000./mesh->
GetNE())/log(2.)/
dim);
90 for (
int l = 0; l < ref_levels; l++)
102 int par_ref_levels = 1;
103 for (
int l = 0; l < par_ref_levels; l++)
117 unsigned int trial_order = order;
118 unsigned int trace_order = order - 1;
119 unsigned int test_order = order;
121 if (
dim == 2 && (order%2 == 0 || (pmesh->
MeshGenerator() & 2 && order > 1)))
125 if (test_order < trial_order)
129 cerr <<
"Warning, test space not enriched enough to handle primal"
151 cout <<
"\nNumber of Unknowns:\n"
152 <<
" Trial space, X0 : " << glob_true_s0
153 <<
" (order " << trial_order <<
")\n"
154 <<
" Interface space, Xhat : " << glob_true_s1
155 <<
" (order " << trace_order <<
")\n"
156 <<
" Test space, Y : " << glob_true_s_test
157 <<
" (order " << test_order <<
")\n\n";
213 enum {x0_var, xhat_var, NVAR};
217 int true_s_test = test_space->
TrueVSize();
221 true_offsets[1] = true_s0;
222 true_offsets[2] = true_s0+true_s1;
225 true_offsets_test[0] = 0;
226 true_offsets_test[1] = true_s_test;
244 matSinv->
Mult(*trueF, SinvF);
263 if (
dim == 2) { Shatinv =
new HypreAMS(*Shat, xhat_space); }
264 else { Shatinv =
new HypreADS(*Shat, xhat_space); }
285 matSinv->
Mult(LSres, tmp);
289 cout <<
"\n|| B0*x0 + Bhat*xhat - F ||_{S^-1} = " << res << endl;
298 ostringstream mesh_name, sol_name;
299 mesh_name <<
"mesh." << setfill(
'0') << setw(6) << myid;
300 sol_name <<
"sol." << setfill(
'0') << setw(6) << myid;
302 ofstream mesh_ofs(mesh_name.str().c_str());
303 mesh_ofs.precision(8);
304 pmesh->
Print(mesh_ofs);
306 ofstream sol_ofs(sol_name.str().c_str());
307 sol_ofs.precision(8);
317 sol_sock <<
"parallel " << num_procs <<
" " << myid <<
"\n";
318 sol_sock.precision(8);
319 sol_sock <<
"solution\n" << *pmesh << *x0 << flush;
T Max() const
Find the maximal element in the array, using the comparison operator < for class T.
A class to handle Block diagonal preconditioners in a matrix-free implementation.
void SetDiagonalBlock(int iblock, Operator *op)
Add a square block op in the block-entry (iblock, iblock).
A class to handle Block systems in a matrix-free implementation.
void SetBlock(int iRow, int iCol, Operator *op, real_t c=1.0)
Add a block op in the block-entry (iblock, jblock).
virtual void Mult(const Vector &x, Vector &y) const
Operator application.
virtual void MultTranspose(const Vector &x, Vector &y) const
Action of the transpose operator.
A class to handle Vectors in a block fashion.
Conjugate gradient method.
virtual void SetOperator(const Operator &op)
Also calls SetOperator for the preconditioner if there is one.
virtual void Mult(const Vector &b, Vector &x) const
Iterative solution of the linear system using the Conjugate Gradient method.
A coefficient that is constant across space and time.
Class for domain integration .
Collection of finite elements from the same family in multiple dimensions. This class is used to matc...
Arbitrary order H1-conforming (continuous) finite elements.
The Auxiliary-space Divergence Solver in hypre.
The Auxiliary-space Maxwell Solver in hypre.
The BoomerAMG solver in hypre.
void SetPrintLevel(int print_level)
Wrapper for hypre's ParCSR matrix class.
HYPRE_Int Mult(HypreParVector &x, HypreParVector &y, real_t alpha=1.0, real_t beta=0.0) const
Computes y = alpha * A * x + beta * y.
Wrapper for hypre's parallel vector class.
Abstract class for hypre's solvers and preconditioners.
static void Init()
Initialize hypre by calling HYPRE_Init() and set default options. After calling Hypre::Init(),...
Integrator that inverts the matrix assembled by another integrator.
void SetRelTol(real_t rtol)
virtual void SetPreconditioner(Solver &pr)
This should be called before SetOperator.
virtual void SetPrintLevel(int print_lvl)
Legacy method to set the level of verbosity of the solver output.
void SetMaxIter(int max_it)
Arbitrary order "L2-conforming" discontinuous finite elements.
Array< int > bdr_attributes
A list of all unique boundary attributes used by the Mesh.
int GetNE() const
Returns number of elements.
int Dimension() const
Dimension of the reference space used within the elements.
int MeshGenerator() const
Get the mesh generator/type.
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).
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
void GetEssentialVDofs(const Array< int > &bdr_attr_is_ess, Array< int > &ess_dofs, int component=-1) const override
Determine the boundary degrees of freedom.
int TrueVSize() const
Obsolete, kept for backward compatibility.
Class for parallel grid function.
void Save(std::ostream &out) const override
void Distribute(const Vector *tv)
Class for parallel meshes.
void Print(std::ostream &out=mfem::out, const std::string &comments="") const override
The operator x -> R*A*P*x constructed through the actions of R^T, A and P.
Arbitrary order "H^{-1/2}-conforming" face finite elements defined on the interface between mesh elem...
Integrator defining a sum of multiple Integrators.
real_t InnerProduct(HypreParVector *x, HypreParVector *y)
void RAP(const DenseMatrix &A, const DenseMatrix &P, DenseMatrix &RAP)