MFEM v4.8.0
Finite element discretization library
Loading...
Searching...
No Matches
phpref.cpp
Go to the documentation of this file.
1// Parallel hp-refinement example
2//
3// Compile with: make phpref
4//
5// Sample runs: mpirun -np 4 phpref -dim 2 -n 1000
6// mpirun -np 8 phpref -dim 3 -n 200
7//
8// Description: This example demonstrates h- and p-refinement in a parallel
9// finite element discretization of the Poisson problem (cf. ex1p)
10// -Delta u = 1 with homogeneous Dirichlet boundary conditions.
11// Refinements are performed iteratively, each iteration having h-
12// or p-refinements on all MPI processes. For simplicity, we
13// randomly choose the elements and the type of refinement, for
14// each iteration. In practice, these choices may be made in a
15// problem-dependent way, but this example serves only to
16// illustrate the capabilities of hp-refinement in parallel.
17//
18// We recommend viewing Example 1 before viewing this example.
19
20#include "mfem.hpp"
21#include <fstream>
22#include <iostream>
23
24using namespace std;
25using namespace mfem;
26
28
29// Deterministic function for "random" integers.
30int DetRand(int & seed)
31{
32 seed++;
33 return int(std::abs(1.0e5 * sin(seed * 1.1234 * M_PI)));
34}
35
36void f_exact(const Vector &x, Vector &f);
37
38int main(int argc, char *argv[])
39{
40 // 1. Initialize MPI and HYPRE.
41 Mpi::Init();
42 const int num_procs = Mpi::WorldSize();
43 const int myid = Mpi::WorldRank();
45
46 // 2. Parse command-line options.
47 int order = 1;
48 const char *device_config = "cpu";
49 bool visualization = true;
50 int numIter = 0;
51 int dim = 2;
52 bool deterministic = true;
53 bool projectSolution = false;
54
55 OptionsParser args(argc, argv);
56 args.AddOption(&order, "-o", "--order",
57 "Finite element order (polynomial degree) or -1 for"
58 " isoparametric space.");
59 args.AddOption(&device_config, "-d", "--device",
60 "Device configuration string, see Device::Configure().");
61 args.AddOption(&visualization, "-vis", "--visualization", "-no-vis",
62 "--no-visualization",
63 "Enable or disable GLVis visualization.");
64 args.AddOption(&numIter, "-n", "--num-iter", "Number of hp-ref iterations");
65 args.AddOption(&dim, "-dim", "--dim", "Mesh dimension (2 or 3)");
66 args.AddOption(&deterministic, "-det", "--deterministic", "-not-det",
67 "--not-deterministic",
68 "Whether to use deterministic random refinements");
69 args.AddOption(&projectSolution, "-proj", "--project-solution", "-no-proj",
70 "--no-project",
71 "Whether to project a coefficient to solution");
72
73 args.Parse();
74 if (!args.Good())
75 {
76 if (myid == 0)
77 {
78 args.PrintUsage(cout);
79 }
80 return 1;
81 }
82 if (myid == 0)
83 {
84 args.PrintOptions(cout);
85 }
86
87 // 3. Enable hardware devices such as GPUs, and programming models such as
88 // CUDA, OCCA, RAJA and OpenMP based on command line options.
89 Device device(device_config);
90 if (myid == 0) { device.Print(); }
91
92 // 4. Construct a uniform coarse mesh on all processors.
93 Mesh mesh;
94 if (dim == 3)
95 {
97 }
98 else
99 {
101 }
102
103 mesh.EnsureNCMesh();
104
105 // 5. Define a parallel mesh by a partitioning of the serial mesh. Refine
106 // this mesh further in parallel to increase the resolution. Once the
107 // parallel mesh is defined, the serial mesh can be deleted.
108 ParMesh pmesh(MPI_COMM_WORLD, mesh);
109 mesh.Clear();
110 {
111 int par_ref_levels = 0;
112 for (int l = 0; l < par_ref_levels; l++)
113 {
114 pmesh.UniformRefinement();
115 }
116 }
117
118 // 6. Define a parallel finite element space on the parallel mesh. Here we
119 // use continuous Lagrange finite elements of the specified order. If
120 // order < 1, we instead use an isoparametric/isogeometric space.
122 bool delete_fec;
123 if (order > 0)
124 {
125 fec = new H1_FECollection(order, dim);
126 delete_fec = true;
127 }
128 else if (pmesh.GetNodes())
129 {
130 fec = pmesh.GetNodes()->OwnFEC();
131 delete_fec = false;
132 if (myid == 0)
133 {
134 cout << "Using isoparametric FEs: " << fec->Name() << endl;
135 }
136 }
137 else
138 {
139 fec = new H1_FECollection(order = 1, dim);
140 delete_fec = true;
141 }
142
143 const int fespaceDim = projectSolution ? dim : 1;
144 ParFiniteElementSpace fespace(&pmesh, fec, fespaceDim);
145
146 // 7. Iteratively perform h- and p-refinements.
147
148 int numH = 0;
149 int numP = 0;
150 int seed = myid;
151
152 const std::vector<char> hp_char = {'h', 'p'};
153
154 for (int iter=0; iter<numIter; ++iter)
155 {
156 const int r1 = deterministic ? DetRand(seed) : rand();
157 const int r2 = deterministic ? DetRand(seed) : rand();
158 const int elem = r1 % pmesh.GetNE();
159 int hp = r2 % 2;
160 MPI_Bcast(&hp, 1, MPI_INT, 0, MPI_COMM_WORLD);
161
162 if (myid == 0)
163 cout << "hp-refinement iteration " << iter << ": "
164 << hp_char[hp] << "-refinement" << endl;
165
166 if (hp == 1)
167 {
168 // p-ref
170 refs.Append(pRefinement(elem, 1)); // Increase the element order by 1
171 fespace.PRefineAndUpdate(refs);
172 numP++;
173 }
174 else
175 {
176 // h-ref
178 refs.Append(Refinement(elem));
179 pmesh.GeneralRefinement(refs);
180 fespace.Update(false);
181 numH++;
182 }
183 }
184
185 const HYPRE_BigInt size = fespace.GlobalTrueVSize();
186 const int maxP = fespace.GetMaxElementOrder();
187 if (myid == 0)
188 {
189 cout << "Number of finite element unknowns: " << size << endl;
190 cout << "Total number of h-refinements: " << numH
191 << "\nTotal number of p-refinements: " << numP
192 << "\nMaximum order " << maxP << "\n";
193 }
194
195 ParGridFunction x(&fespace);
196 Vector X;
197
198 if (projectSolution)
199 {
201 x.ProjectCoefficient(vec_coef);
202
203 X.SetSize(fespace.GetTrueVSize());
204
205 fespace.GetRestrictionMatrix()->Mult(x, X);
206 fespace.GetProlongationMatrix()->Mult(X, x);
207
208 // Compute and print the L^2 norm of the error.
209 const real_t error = x.ComputeL2Error(vec_coef);
210 if (myid == 0)
211 {
212 cout << "\n|| E_h - E ||_{L^2} = " << error << '\n' << endl;
213 }
214 }
215 else
216 {
217 // 8. Determine the list of true (i.e. parallel conforming) essential
218 // boundary dofs. In this example, the boundary conditions are defined
219 // by marking all the boundary attributes from the mesh as essential
220 // (Dirichlet) and converting them to a list of true dofs.
221 Array<int> ess_tdof_list;
222 if (pmesh.bdr_attributes.Size())
223 {
224 Array<int> ess_bdr(pmesh.bdr_attributes.Max());
225 ess_bdr = 1;
226 fespace.GetEssentialTrueDofs(ess_bdr, ess_tdof_list);
227 }
228
229 // 9. Set up the parallel linear form b(.) which corresponds to the
230 // right-hand side of the FEM linear system, which in this case is
231 // (1,phi_i) where phi_i are the basis functions in fespace.
232 ParLinearForm b(&fespace);
233 ConstantCoefficient one(1.0);
234 b.AddDomainIntegrator(new DomainLFIntegrator(one));
235 b.Assemble();
236
237 // 10. Define the solution vector x as a parallel finite element grid
238 // function corresponding to fespace. Initialize x with initial guess of
239 // zero, which satisfies the boundary conditions.
240 x = 0.0;
241
242 // 11. Set up the parallel bilinear form a(.,.) on the finite element space
243 // corresponding to the Laplacian operator -Delta, by adding the
244 // diffusion domain integrator.
245 ParBilinearForm a(&fespace);
246 a.AddDomainIntegrator(new DiffusionIntegrator(one));
247
248 // 12. Assemble the parallel bilinear form and the corresponding linear
249 // system, applying any necessary transformations such as: parallel
250 // assembly, eliminating boundary conditions, applying conforming
251 // constraints for non-conforming AMR, static condensation, etc.
252 a.Assemble();
253
254 OperatorPtr A;
255 Vector B;
256 a.FormLinearSystem(ess_tdof_list, x, b, A, X, B);
257
258 // 13. Solve the linear system A X = B.
259 // * With full assembly, use the BoomerAMG preconditioner from hypre.
260 // * With partial assembly, use Jacobi smoothing, for now.
261 Solver *prec = new HypreBoomerAMG;
262 CGSolver cg(MPI_COMM_WORLD);
263 cg.SetRelTol(1e-12);
264 cg.SetMaxIter(2000);
265 cg.SetPrintLevel(1);
266 if (prec) { cg.SetPreconditioner(*prec); }
267 cg.SetOperator(*A);
268 cg.Mult(B, X);
269 delete prec;
270
271 // 14. Recover the parallel grid function corresponding to X. This is the
272 // local finite element solution on each processor.
273 a.RecoverFEMSolution(X, b, x);
274 }
275
276 if (fespaceDim == 1)
277 {
278 const real_t h1error = CheckH1Continuity(x);
279 cout << myid << ": H1 continuity error " << h1error << endl;
280 MFEM_VERIFY(h1error < 1.0e-12, "H1 continuity is not satisfied");
281 }
282
283 L2_FECollection fecL2(0, dim);
284 ParFiniteElementSpace l2fespace(&pmesh, &fecL2);
285 ParGridFunction xo(&l2fespace);
286 xo = 0.0;
287
288 for (int e=0; e<pmesh.GetNE(); ++e)
289 {
290 const int p_elem = fespace.GetElementOrder(e);
291 Array<int> dofs;
292 l2fespace.GetElementDofs(e, dofs);
293 MFEM_VERIFY(dofs.Size() == 1, "");
294 xo[dofs[0]] = p_elem;
295 }
296
297 // 15. Save the refined mesh and the solution in parallel. This output can
298 // be viewed later using GLVis: "glvis -np <np> -m mesh -g sol".
299 std::unique_ptr<GridFunction> vis_x = x.ProlongateToMaxOrder();
300 {
301 ostringstream mesh_name, sol_name, order_name;
302 mesh_name << "mesh." << setfill('0') << setw(6) << myid;
303 sol_name << "sol." << setfill('0') << setw(6) << myid;
304 order_name << "order." << setfill('0') << setw(6) << myid;
305
306 ofstream mesh_ofs(mesh_name.str().c_str());
307 mesh_ofs.precision(8);
308 pmesh.ParPrint(mesh_ofs);
309
310 ofstream sol_ofs(sol_name.str().c_str());
311 sol_ofs.precision(8);
312
313 vis_x->Save(sol_ofs);
314
315 ofstream order_ofs(order_name.str().c_str());
316 order_ofs.precision(8);
317 xo.Save(order_ofs);
318 }
319
320 // 16. Send the solution by socket to a GLVis server.
321 if (visualization)
322 {
323 char vishost[] = "localhost";
324 int visport = 19916;
325 socketstream sol_sock(vishost, visport);
326 sol_sock << "parallel " << num_procs << " " << myid << "\n";
327 sol_sock.precision(8);
328 sol_sock << "solution\n" << pmesh << *vis_x << flush;
329 }
330
331 // 17. Free the used memory.
332 if (delete_fec)
333 {
334 delete fec;
335 }
336
337 return 0;
338}
339
341{
343
344 const ParFiniteElementSpace *fes = x.ParFESpace();
345 ParMesh *mesh = fes->GetParMesh();
346
347 const int dim = mesh->Dimension();
348
349 // Following the example of KellyErrorEstimator::ComputeEstimates(),
350 // we loop over interior faces and then shared faces.
351
352 // Compute error contribution from local interior faces
353 real_t errorMax = 0.0;
354 for (int f = 0; f < mesh->GetNumFaces(); f++)
355 {
356 if (mesh->FaceIsInterior(f))
357 {
358 int Inf1, Inf2, NCFace;
359 mesh->GetFaceInfos(f, &Inf1, &Inf2, &NCFace);
360
361 auto FT = mesh->GetFaceElementTransformations(f);
362
363 const int faceOrder = dim == 3 ? fes->GetFaceOrder(f) :
364 fes->GetEdgeOrder(f);
365 auto &int_rule = IntRules.Get(FT->FaceGeom, 2 * faceOrder);
366 const auto nip = int_rule.GetNPoints();
367
368 // Convention
369 // * Conforming face: Face side with smaller element id handles
370 // the integration
371 // * Non-conforming face: The slave handles the integration.
372 // See FaceInfo documentation for details.
373 bool isNCSlave = FT->Elem2No >= 0 && NCFace >= 0;
374 bool isConforming = FT->Elem2No >= 0 && NCFace == -1;
375 if ((FT->Elem1No < FT->Elem2No && isConforming) || isNCSlave)
376 {
377 for (int i = 0; i < nip; i++)
378 {
379 const auto &fip = int_rule.IntPoint(i);
381
382 FT->Loc1.Transform(fip, ip);
383 const real_t v1 = x.GetValue(FT->Elem1No, ip);
384
385 FT->Loc2.Transform(fip, ip);
386 const real_t v2 = x.GetValue(FT->Elem2No, ip);
387
388 const real_t err_i = std::abs(v1 - v2);
389 errorMax = std::max(errorMax, err_i);
390 }
391 }
392 }
393 }
394
395 // Compute error contribution from shared interior faces
396 for (int sf = 0; sf < mesh->GetNSharedFaces(); sf++)
397 {
398 const int f = mesh->GetSharedFace(sf);
399 const bool trueInterior = mesh->FaceIsTrueInterior(f);
400 if (!trueInterior) { continue; }
401
402 auto FT = mesh->GetSharedFaceTransformations(sf, true);
403 const int faceOrder = dim == 3 ? fes->GetFaceOrder(f) : fes->GetEdgeOrder(f);
404 const auto &int_rule = IntRules.Get(FT->FaceGeom, 2 * faceOrder);
405 const auto nip = int_rule.GetNPoints();
406
407 for (int i = 0; i < nip; i++)
408 {
409 const auto &fip = int_rule.IntPoint(i);
411
412 FT->Loc1.Transform(fip, ip);
413 const real_t v1 = x.GetValue(FT->Elem1No, ip);
414
415 FT->Loc2.Transform(fip, ip);
416 const real_t v2 = x.GetValue(FT->Elem2No, ip);
417
418 const real_t err_i = std::abs(v1 - v2);
419 errorMax = std::max(errorMax, err_i);
420 }
421 }
422
423 return errorMax;
424}
425
426void f_exact(const Vector &x, Vector &f)
427{
428 constexpr real_t freq = 1.0;
429 constexpr real_t kappa = freq * M_PI;
430
431 if (x.Size() == 3)
432 {
433 f(0) = (1. + kappa * kappa) * sin(kappa * x(1));
434 f(1) = (1. + kappa * kappa) * sin(kappa * x(2));
435 f(2) = (1. + kappa * kappa) * sin(kappa * x(0));
436 }
437 else
438 {
439 f(0) = (1. + kappa * kappa) * sin(kappa * x(1));
440 f(1) = (1. + kappa * kappa) * sin(kappa * x(0));
441 if (x.Size() == 3) { f(2) = 0.0; }
442 }
443}
T Max() const
Find the maximal element in the array, using the comparison operator < for class T.
Definition array.cpp:68
int Size() const
Return the logical size of the array.
Definition array.hpp:147
int Append(const T &el)
Append element 'el' to array, resize if necessary.
Definition array.hpp:830
Conjugate gradient method.
Definition solvers.hpp:538
void Mult(const Vector &b, Vector &x) const override
Iterative solution of the linear system using the Conjugate Gradient method.
Definition solvers.cpp:751
void SetOperator(const Operator &op) override
Set/update the solver for the given operator.
Definition solvers.hpp:551
A coefficient that is constant across space and time.
The MFEM Device class abstracts hardware devices such as GPUs, as well as programming models such as ...
Definition device.hpp:123
void Print(std::ostream &out=mfem::out)
Print the configuration of the MFEM virtual device object.
Definition device.cpp:297
Class for domain integration .
Definition lininteg.hpp:106
Collection of finite elements from the same family in multiple dimensions. This class is used to matc...
Definition fe_coll.hpp:27
virtual const char * Name() const
Definition fe_coll.hpp:79
int GetEdgeOrder(int edge, int variant=0) const
Definition fespace.cpp:3356
int GetFaceOrder(int face, int variant=0) const
Returns the polynomial degree of the i'th face finite element.
Definition fespace.cpp:3372
int GetElementOrder(int i) const
Returns the order of the i'th finite element.
Definition fespace.cpp:221
Arbitrary order H1-conforming (continuous) finite elements.
Definition fe_coll.hpp:275
The BoomerAMG solver in hypre.
Definition hypre.hpp:1717
static void Init()
Initialize hypre by calling HYPRE_Init() and set default options. After calling Hypre::Init(),...
Definition hypre.cpp:33
Class for integration point with weight.
Definition intrules.hpp:35
int GetNPoints() const
Returns the number of the points in the integration rule.
Definition intrules.hpp:256
const IntegrationRule & Get(int GeomType, int Order)
Returns an integration rule for given GeomType and Order.
void SetRelTol(real_t rtol)
Definition solvers.hpp:229
virtual void SetPreconditioner(Solver &pr)
This should be called before SetOperator.
Definition solvers.cpp:174
virtual void SetPrintLevel(int print_lvl)
Legacy method to set the level of verbosity of the solver output.
Definition solvers.cpp:72
void SetMaxIter(int max_it)
Definition solvers.hpp:231
Arbitrary order "L2-conforming" discontinuous finite elements.
Definition fe_coll.hpp:346
Mesh data type.
Definition mesh.hpp:64
Array< int > bdr_attributes
A list of all unique boundary attributes used by the Mesh.
Definition mesh.hpp:290
void GetFaceInfos(int Face, int *Inf1, int *Inf2) const
Definition mesh.cpp:1463
int GetNumFaces() const
Return the number of faces (3D), edges (2D) or vertices (1D).
Definition mesh.cpp:6531
void GeneralRefinement(const Array< Refinement > &refinements, int nonconforming=-1, int nc_limit=0)
Definition mesh.cpp:10883
void Clear()
Clear the contents of the Mesh.
Definition mesh.hpp:761
int GetNE() const
Returns number of elements.
Definition mesh.hpp:1282
int Dimension() const
Dimension of the reference space used within the elements.
Definition mesh.hpp:1216
bool FaceIsInterior(int FaceNo) const
Return true if the given face is interior.
Definition mesh.hpp:1462
static Mesh MakeCartesian3D(int nx, int ny, int nz, Element::Type type, real_t sx=1.0, real_t sy=1.0, real_t sz=1.0, bool sfc_ordering=true)
Creates a mesh for the parallelepiped [0,sx]x[0,sy]x[0,sz], divided into nx*ny*nz hexahedra if type =...
Definition mesh.cpp:4481
void GetNodes(Vector &node_coord) const
Definition mesh.cpp:9294
void EnsureNCMesh(bool simplices_nonconforming=false)
Definition mesh.cpp:10951
static Mesh MakeCartesian2D(int nx, int ny, Element::Type type, bool generate_edges=false, real_t sx=1.0, real_t sy=1.0, bool sfc_ordering=true)
Creates mesh for the rectangle [0,sx]x[0,sy], divided into nx*ny quadrilaterals if type = QUADRILATER...
Definition mesh.cpp:4471
void UniformRefinement(int i, const DSTable &, int *, int *, int *)
Definition mesh.cpp:11295
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).
Pointer to an Operator of a specified type.
Definition handle.hpp:34
virtual void Mult(const Vector &x, Vector &y) const =0
Operator application: y=A(x).
virtual void RecoverFEMSolution(const Vector &X, const Vector &b, Vector &x)
Reconstruct a solution vector x (e.g. a GridFunction) from the solution X of a constrained linear sys...
Definition operator.cpp:148
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...
Definition optparser.hpp:82
bool Good() const
Return true if the command line options were parsed successfully.
Class for parallel bilinear form.
Abstract parallel finite element space.
Definition pfespace.hpp:29
int GetMaxElementOrder() const override
Returns the maximum polynomial order over all elements globally.
void GetEssentialTrueDofs(const Array< int > &bdr_attr_is_ess, Array< int > &ess_tdof_list, int component=-1) const override
HYPRE_BigInt GlobalTrueVSize() const
Definition pfespace.hpp:346
int GetTrueVSize() const override
Return the number of local vector true dofs.
Definition pfespace.hpp:350
const Operator * GetProlongationMatrix() const override
void GetElementDofs(int i, Array< int > &dofs, DofTransformation &doftrans) const override
The same as GetElementDofs(), but with a user-allocated DofTransformation object. doftrans must be al...
Definition pfespace.cpp:561
void PRefineAndUpdate(const Array< pRefinement > &refs, bool want_transfer=true) override
const SparseMatrix * GetRestrictionMatrix() const override
Get the R matrix which restricts a local dof vector to true dof vector.
Definition pfespace.hpp:467
ParMesh * GetParMesh() const
Definition pfespace.hpp:338
void Update(bool want_transform=true) override
Class for parallel grid function.
Definition pgridfunc.hpp:50
void Save(std::ostream &out) const override
real_t ComputeL2Error(Coefficient *exsol[], const IntegrationRule *irs[]=NULL, const Array< int > *elems=NULL) const override
Returns ||u_ex - u_h||_L2 in parallel for H1 or L2 elements.
real_t GetValue(int i, const IntegrationPoint &ip, int vdim=1) const override
ParFiniteElementSpace * ParFESpace() const
void ProjectCoefficient(Coefficient &coeff) override
Project coeff Coefficient to this GridFunction. The projection computation depends on the choice of t...
std::unique_ptr< ParGridFunction > ProlongateToMaxOrder() const
Return a GridFunction with the values of this, prolongated to the maximum order of all elements in th...
Class for parallel linear form.
Class for parallel meshes.
Definition pmesh.hpp:34
void ParPrint(std::ostream &out, const std::string &comments="") const
Definition pmesh.cpp:6313
int GetNSharedFaces() const
Return the number of shared faces (3D), edges (2D), vertices (1D)
Definition pmesh.cpp:3152
FaceElementTransformations * GetFaceElementTransformations(int FaceNo, int mask=31) override
Definition pmesh.cpp:2896
bool FaceIsTrueInterior(int FaceNo) const
Definition mesh.hpp:571
int GetSharedFace(int sface) const
Return the local face index for the given shared face.
Definition pmesh.cpp:3171
FaceElementTransformations * GetSharedFaceTransformations(int sf, bool fill2=true)
Get the FaceElementTransformations for the given shared face (edge 2D) using the shared face index sf...
Definition pmesh.cpp:2922
Base class for solvers.
Definition operator.hpp:780
void Mult(const Vector &x, Vector &y) const override
Matrix vector multiplication.
A general vector function coefficient.
Vector data type.
Definition vector.hpp:82
int Size() const
Returns the size of the vector.
Definition vector.hpp:226
void SetSize(int s)
Resize the vector to size s.
Definition vector.hpp:558
real_t kappa
Definition ex24.cpp:54
int dim
Definition ex24.cpp:53
real_t freq
Definition ex24.cpp:54
int main()
HYPRE_Int HYPRE_BigInt
real_t b
Definition lissajous.cpp:42
real_t a
Definition lissajous.cpp:41
float real_t
Definition config.hpp:43
std::function< real_t(const Vector &)> f(real_t mass_coeff)
Definition lor_mms.hpp:30
IntegrationRules IntRules(0, Quadrature1D::GaussLegendre)
A global object with all integration rules (defined in intrules.cpp)
Definition intrules.hpp:492
const char vishost[]
int DetRand(int &seed)
Definition phpref.cpp:30
void f_exact(const Vector &x, Vector &f)
Definition phpref.cpp:426
real_t CheckH1Continuity(ParGridFunction &x)
Definition phpref.cpp:340