39 int main(
int argc,
char *argv[])
42 const char *mesh_file =
"../data/star.mesh";
44 bool visualization = 1;
47 args.
AddOption(&mesh_file,
"-m",
"--mesh",
50 "Finite element order (polynomial degree).");
51 args.
AddOption(&visualization,
"-vis",
"--visualization",
"-no-vis",
53 "Enable or disable GLVis visualization.");
65 Mesh *mesh =
new Mesh(mesh_file, 1, 1);
74 (int)floor(log(10000./mesh->
GetNE())/log(2.)/
dim);
75 for (
int l = 0; l < ref_levels; l++)
89 int trial_order = order;
90 int trace_order = order - 1;
91 int test_order = order;
92 if (dim == 2 && (order%2 == 0 || (mesh->
MeshGenerator() & 2 && order > 1)))
96 if (test_order < trial_order)
97 cerr <<
"Warning, test space not enriched enough to handle primal"
113 enum {x0_var, xhat_var, NVAR};
117 int s_test = test_space->
GetVSize();
126 offsets_test[1] = s_test;
128 std::cout <<
"\nNumber of Unknowns:\n"
129 <<
" Trial space, X0 : " << s0
130 <<
" (order " << trial_order <<
")\n"
131 <<
" Interface space, Xhat : " << s1
132 <<
" (order " << trace_order <<
")\n"
133 <<
" Test space, Y : " << s_test
134 <<
" (order " << test_order <<
")\n\n";
193 matSinv.
Mult(F,SinvF);
205 #ifndef MFEM_USE_SUITESPARSE
206 const double prec_rtol = 1e-3;
207 const int prec_maxit = 200;
234 PCG(A, P, b, x, 1, 200, 1e-12, 0.0);
241 cout <<
"\n|| B0*x0 + Bhat*xhat - F ||_{S^-1} = " << res << endl;
245 x0.
MakeRef(x0_space, x.GetBlock(x0_var), 0);
250 ofstream mesh_ofs(
"refined.mesh");
251 mesh_ofs.precision(8);
252 mesh->
Print(mesh_ofs);
253 ofstream sol_ofs(
"sol.gf");
254 sol_ofs.precision(8);
261 char vishost[] =
"localhost";
264 sol_sock.precision(8);
265 sol_sock <<
"solution\n" << *mesh << x0 << flush;
Class for domain integration L(v) := (f, v)
Conjugate gradient method.
Class for grid function - Vector with associated FE space.
Integrator defining a sum of multiple Integrators.
Subclass constant coefficient.
virtual void MultTranspose(const Vector &x, Vector &y) const
Action of the transpose operator.
int GetNE() const
Returns number of elements.
void AddIntegrator(BilinearFormIntegrator *integ)
void SetBlock(int iRow, int iCol, Operator *op)
Add a block op in the block-entry (iblock, jblock).
bool iterative_mode
If true, use the second argument of Mult() as an initial guess.
int main(int argc, char *argv[])
Direct sparse solver using UMFPACK.
Integrator that inverts the matrix assembled by another integrator.
virtual void Save(std::ostream &out) const
Save the GridFunction to an output stream.
void SetPrintLevel(int print_lvl)
A class to handle Block diagonal preconditioners in a matrix-free implementation. ...
void UniformRefinement(int i, const DSTable &, int *, int *, int *)
void SetMaxIter(int max_it)
HypreParMatrix * RAP(HypreParMatrix *A, HypreParMatrix *P)
Returns the matrix P^t * A * P.
int MeshGenerator()
Get the mesh generator/type.
T Max() const
Find the maximal element in the array, using the comparison operator < for class T.
The operator x -> R*A*P*x.
void MakeRef(FiniteElementSpace *f, double *v)
Make the GridFunction reference external data on a new FiniteElementSpace.
void PCG(const Operator &A, Solver &B, const Vector &b, Vector &x, int print_iter, int max_num_iter, double RTOLERANCE, double ATOLERANCE)
Preconditioned conjugate gradient method. (tolerances are squared)
void PrintUsage(std::ostream &out) const
virtual void Mult(const Vector &x, Vector &y) const
Operator application.
Array< int > bdr_attributes
A list of all unique boundary attributes used by the Mesh.
void SetRelTol(double rtol)
virtual void Mult(const Vector &x, Vector &y) const
Matrix vector multiplication.
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)
void PrintOptions(std::ostream &out) const
virtual void SetOperator(const Operator &op)
Also calls SetOperator for the preconditioner if there is one.
Arbitrary order H1-conforming (continuous) finite elements.
A class to handle Block systems in a matrix-free implementation.
double InnerProduct(const Vector &x, const Vector &y) const
Compute y^t A x.
void SetDiagonalBlock(int iblock, Operator *op)
Add a square block op in the block-entry (iblock, iblock).
Arbitrary order "L2-conforming" discontinuous finite elements.
virtual void Print(std::ostream &out=std::cout) const