MFEM  v4.3.0
Finite element discretization library
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
nurbs_ex1p.cpp
Go to the documentation of this file.
1 // MFEM Example 1 - Parallel NURBS Version
2 //
3 // Compile with: make nurbs_ex1p
4 //
5 // Sample runs: mpirun -np 4 nurbs_ex1p -m ../../data/square-disc.mesh
6 // mpirun -np 4 nurbs_ex1p -m ../../data/star.mesh
7 // mpirun -np 4 nurbs_ex1p -m ../../data/escher.mesh
8 // mpirun -np 4 nurbs_ex1p -m ../../data/fichera.mesh
9 // mpirun -np 4 nurbs_ex1p -m ../../data/square-disc-p2.vtk -o 2
10 // mpirun -np 4 nurbs_ex1p -m ../../data/square-disc-p3.mesh -o 3
11 // mpirun -np 4 nurbs_ex1p -m ../../data/square-disc-nurbs.mesh -o -1
12 // mpirun -np 4 nurbs_ex1p -m ../../data/disc-nurbs.mesh -o -1
13 // mpirun -np 4 nurbs_ex1p -m ../../data/pipe-nurbs.mesh -o -1
14 // mpirun -np 4 nurbs_ex1p -m ../../data/ball-nurbs.mesh -o 2
15 // mpirun -np 4 nurbs_ex1p -m ../../data/star-surf.mesh
16 // mpirun -np 4 nurbs_ex1p -m ../../data/square-disc-surf.mesh
17 // mpirun -np 4 nurbs_ex1p -m ../../data/inline-segment.mesh
18 // mpirun -np 4 nurbs_ex1p -m ../../data/amr-quad.mesh
19 // mpirun -np 4 nurbs_ex1p -m ../../data/amr-hex.mesh
20 // mpirun -np 4 nurbs_ex1p -m ../../data/mobius-strip.mesh
21 // mpirun -np 4 nurbs_ex1p -m ../../data/mobius-strip.mesh -o -1 -sc
22 // mpirun -np 4 nurbs_ex1p -m ../../data/square-disc-nurbs.mesh -o -1
23 // mpirun -np 4 nurbs_ex1p -m ../../data/disc-nurbs.mesh -o -1
24 // mpirun -np 4 nurbs_ex1p -m ../../data/pipe-nurbs.mesh -o -1
25 // mpirun -np 4 nurbs_ex1p -m ../../data/square-nurbs.mesh -o 2 -no-ibp
26 // mpirun -np 4 nurbs_ex1p -m ../../data/cube-nurbs.mesh -o 2 -no-ibp
27 // mpirun -np 4 nurbs_ex1p -m ../../data/pipe-nurbs-2d.mesh -o 2 -no-ibp
28 
29 // Description: This example code demonstrates the use of MFEM to define a
30 // simple finite element discretization of the Laplace problem
31 // -Delta u = 1 with homogeneous Dirichlet boundary conditions.
32 // Specifically, we discretize using a FE space of the specified
33 // order, or if order < 1 using an isoparametric/isogeometric
34 // space (i.e. quadratic for quadratic curvilinear mesh, NURBS for
35 // NURBS mesh, etc.)
36 //
37 // The example highlights the use of mesh refinement, finite
38 // element grid functions, as well as linear and bilinear forms
39 // corresponding to the left-hand side and right-hand side of the
40 // discrete linear system. We also cover the explicit elimination
41 // of essential boundary conditions, static condensation, and the
42 // optional connection to the GLVis tool for visualization.
43 
44 #include "mfem.hpp"
45 #include <fstream>
46 #include <iostream>
47 
48 using namespace std;
49 using namespace mfem;
50 /** Class for integrating the bilinear form a(u,v) := (Q Laplace u, v) where Q
51  can be a scalar coefficient. */
52 class Diffusion2Integrator: public BilinearFormIntegrator
53 {
54 private:
55 #ifndef MFEM_THREAD_SAFE
56  Vector shape,laplace;
57 #endif
58  Coefficient *Q;
59 
60 public:
61  /// Construct a diffusion integrator with coefficient Q = 1
62  Diffusion2Integrator() { Q = NULL; }
63 
64  /// Construct a diffusion integrator with a scalar coefficient q
65  Diffusion2Integrator (Coefficient &q) : Q(&q) { }
66 
67  /** Given a particular Finite Element
68  computes the element stiffness matrix elmat. */
69  virtual void AssembleElementMatrix(const FiniteElement &el,
71  DenseMatrix &elmat)
72  {
73  int nd = el.GetDof();
74  int dim = el.GetDim();
75  double w;
76 
77 #ifdef MFEM_THREAD_SAFE
78  Vector shape(nd);
79  Vector laplace(nd);
80 #else
81  shape.SetSize(nd);
82  laplace.SetSize(nd);
83 #endif
84  elmat.SetSize(nd);
85 
86  const IntegrationRule *ir = IntRule;
87  if (ir == NULL)
88  {
89  int order;
90  if (el.Space() == FunctionSpace::Pk)
91  {
92  order = 2*el.GetOrder() - 2;
93  }
94  else
95  {
96  order = 2*el.GetOrder() + dim - 1;
97  }
98 
99  if (el.Space() == FunctionSpace::rQk)
100  {
101  ir = &RefinedIntRules.Get(el.GetGeomType(),order);
102  }
103  else
104  {
105  ir = &IntRules.Get(el.GetGeomType(),order);
106  }
107  }
108 
109  elmat = 0.0;
110  for (int i = 0; i < ir->GetNPoints(); i++)
111  {
112  const IntegrationPoint &ip = ir->IntPoint(i);
113  Trans.SetIntPoint(&ip);
114  w = -ip.weight * Trans.Weight();
115 
116  el.CalcShape(ip, shape);
117  el.CalcPhysLaplacian(Trans, laplace);
118 
119  if (Q)
120  {
121  w *= Q->Eval(Trans, ip);
122  }
123 
124  for (int j = 0; j < nd; j++)
125  {
126  for (int i = 0; i < nd; i++)
127  {
128  elmat(i, j) += w*shape(i)*laplace(j);
129  }
130  }
131  }
132  }
133 
134 };
135 
136 int main(int argc, char *argv[])
137 {
138  // 1. Initialize MPI.
139  int num_procs, myid;
140  MPI_Init(&argc, &argv);
141  MPI_Comm_size(MPI_COMM_WORLD, &num_procs);
142  MPI_Comm_rank(MPI_COMM_WORLD, &myid);
143 
144  // 2. Parse command-line options.
145  const char *mesh_file = "../../data/star.mesh";
146  int ref_levels = -1;
147  Array<int> order(1);
148  order[0] = 1;
149  bool static_cond = false;
150  bool visualization = 1;
151  bool ibp = 1;
152  bool strongBC = 1;
153  double kappa = -1;
154 
155  OptionsParser args(argc, argv);
156  args.AddOption(&mesh_file, "-m", "--mesh",
157  "Mesh file to use.");
158  args.AddOption(&ref_levels, "-r", "--refine",
159  "Number of times to refine the mesh uniformly, -1 for auto.");
160  args.AddOption(&order, "-o", "--order",
161  "Finite element order (polynomial degree) or -1 for"
162  " isoparametric space.");
163  args.AddOption(&ibp, "-ibp", "--ibp", "-no-ibp",
164  "--no-ibp",
165  "Selects the standard weak form (IBP) or the nonstandard (NO-IBP).");
166  args.AddOption(&strongBC, "-sbc", "--strong-bc", "-wbc",
167  "--weak-bc",
168  "Selects strong or weak enforcement of Dirichlet BCs.");
169  args.AddOption(&kappa, "-k", "--kappa",
170  "Sets the SIPG penalty parameters, should be positive."
171  " Negative values are replaced with (order+1)^2.");
172  args.AddOption(&static_cond, "-sc", "--static-condensation", "-no-sc",
173  "--no-static-condensation", "Enable static condensation.");
174  args.AddOption(&visualization, "-vis", "--visualization", "-no-vis",
175  "--no-visualization",
176  "Enable or disable GLVis visualization.");
177  args.Parse();
178  if (!args.Good())
179  {
180  if (myid == 0)
181  {
182  args.PrintUsage(cout);
183  }
184  MPI_Finalize();
185  return 1;
186  }
187  if (!strongBC & (kappa < 0))
188  {
189  kappa = 4*(order.Max()+1)*(order.Max()+1);
190  }
191  if (myid == 0)
192  {
193  args.PrintOptions(cout);
194  }
195 
196  // 3. Read the (serial) mesh from the given mesh file on all processors. We
197  // can handle triangular, quadrilateral, tetrahedral, hexahedral, surface
198  // and volume meshes with the same code.
199  Mesh *mesh = new Mesh(mesh_file, 1, 1);
200  int dim = mesh->Dimension();
201 
202  // 4. Refine the serial mesh on all processors to increase the resolution. In
203  // this example we do 'ref_levels' of uniform refinement. We choose
204  // 'ref_levels' to be the largest number that gives a final mesh with no
205  // more than 10,000 elements.
206  {
207  if (ref_levels < 0)
208  {
209  ref_levels =
210  (int)floor(log(5000./mesh->GetNE())/log(2.)/dim);
211  }
212 
213  for (int l = 0; l < ref_levels; l++)
214  {
215  mesh->UniformRefinement();
216  }
217  if (myid == 0)
218  {
219  mesh->PrintInfo();
220  }
221  }
222 
223  // 5. Define a parallel mesh by a partitioning of the serial mesh. Refine
224  // this mesh further in parallel to increase the resolution. Once the
225  // parallel mesh is defined, the serial mesh can be deleted.
226  ParMesh *pmesh = new ParMesh(MPI_COMM_WORLD, *mesh);
227  delete mesh;
228  if (!pmesh->NURBSext)
229  {
230  int par_ref_levels = 2;
231  for (int l = 0; l < par_ref_levels; l++)
232  {
233  pmesh->UniformRefinement();
234  }
235  }
236 
237  // 6. Define a parallel finite element space on the parallel mesh. Here we
238  // use continuous Lagrange finite elements of the specified order. If
239  // order < 1, we instead use an isoparametric/isogeometric space.
241  NURBSExtension *NURBSext = NULL;
242  int own_fec = 0;
243 
244  if (order[0] == -1) // Isoparametric
245  {
246  if (pmesh->GetNodes())
247  {
248  fec = pmesh->GetNodes()->OwnFEC();
249  own_fec = 0;
250  cout << "Using isoparametric FEs: " << fec->Name() << endl;
251  }
252  else
253  {
254  cout <<"Mesh does not have FEs --> Assume order 1.\n";
255  fec = new H1_FECollection(1, dim);
256  own_fec = 1;
257  }
258  }
259  else if (pmesh->NURBSext && (order[0] > 0) ) // Subparametric NURBS
260  {
261  fec = new NURBSFECollection(order[0]);
262  own_fec = 1;
263  int nkv = pmesh->NURBSext->GetNKV();
264 
265  if (order.Size() == 1)
266  {
267  int tmp = order[0];
268  order.SetSize(nkv);
269  order = tmp;
270  }
271  if (order.Size() != nkv ) { mfem_error("Wrong number of orders set."); }
272  NURBSext = new NURBSExtension(pmesh->NURBSext, order);
273  }
274  else
275  {
276  if (order.Size() > 1) { cout <<"Wrong number of orders set, needs one.\n"; }
277  fec = new H1_FECollection(abs(order[0]), dim);
278  own_fec = 1;
279  }
280  ParFiniteElementSpace *fespace = new ParFiniteElementSpace(pmesh,NURBSext,fec);
281  HYPRE_BigInt size = fespace->GlobalTrueVSize();
282  if (myid == 0)
283  {
284  cout << "Number of finite element unknowns: " << size << endl;
285  }
286 
287  if (!ibp)
288  {
289  if (!pmesh->NURBSext)
290  {
291  cout << "No integration by parts requires a NURBS mesh."<< endl;
292  return 2;
293  }
294  if (pmesh->NURBSext->GetNP()>1)
295  {
296  cout << "No integration by parts requires a NURBS mesh, with only 1 patch."<<
297  endl;
298  cout << "A C_1 discretisation is required."<< endl;
299  cout << "Currently only C_0 multipatch coupling implemented."<< endl;
300  return 3;
301  }
302  if (order[0]<2)
303  {
304  cout << "No integration by parts requires at least quadratic NURBS."<< endl;
305  cout << "A C_1 discretisation is required."<< endl;
306  return 4;
307  }
308  }
309 
310  // 7. Determine the list of true (i.e. parallel conforming) essential
311  // boundary dofs. In this example, the boundary conditions are defined
312  // by marking all the boundary attributes from the mesh as essential
313  // (Dirichlet) and converting them to a list of true dofs.
314  Array<int> ess_tdof_list;
315  if (pmesh->bdr_attributes.Size())
316  {
317  Array<int> ess_bdr(pmesh->bdr_attributes.Max());
318  if (strongBC)
319  {
320  ess_bdr = 1;
321  }
322  else
323  {
324  ess_bdr = 0;
325  }
326  fespace->GetEssentialTrueDofs(ess_bdr, ess_tdof_list);
327  }
328 
329  // 8. Set up the parallel linear form b(.) which corresponds to the
330  // right-hand side of the FEM linear system, which in this case is
331  // (1,phi_i) where phi_i are the basis functions in fespace.
332  ConstantCoefficient one(1.0);
333  ConstantCoefficient zero(0.0);
334 
335  ParLinearForm *b = new ParLinearForm(fespace);
337 
338  if (!strongBC)
340  new DGDirichletLFIntegrator(zero, one, -1.0, kappa));
341  b->Assemble();
342 
343  // 9. Define the solution vector x as a parallel finite element grid function
344  // corresponding to fespace. Initialize x with initial guess of zero,
345  // which satisfies the boundary conditions.
346  ParGridFunction x(fespace);
347  x = 0.0;
348 
349  // 10. Set up the parallel bilinear form a(.,.) on the finite element space
350  // corresponding to the Laplacian operator -Delta, by adding the Diffusion
351  // domain integrator.
352  ParBilinearForm *a = new ParBilinearForm(fespace);
353  if (ibp)
354  {
356  }
357  else
358  {
359  a->AddDomainIntegrator(new Diffusion2Integrator(one));
360  }
361  if (!strongBC)
362  {
363  a->AddBdrFaceIntegrator(new DGDiffusionIntegrator(one, -1.0, kappa));
364  }
365 
366  // 11. Assemble the parallel bilinear form and the corresponding linear
367  // system, applying any necessary transformations such as: parallel
368  // assembly, eliminating boundary conditions, applying conforming
369  // constraints for non-conforming AMR, static condensation, etc.
370  if (static_cond) { a->EnableStaticCondensation(); }
371  a->Assemble();
372 
373  HypreParMatrix A;
374  Vector B, X;
375  a->FormLinearSystem(ess_tdof_list, x, *b, A, X, B);
376 
377  if (myid == 0)
378  {
379  cout << "Size of linear system: " << A.GetGlobalNumRows() << endl;
380  }
381 
382  // 12. Define and apply a parallel PCG solver for AX=B with the BoomerAMG
383  // preconditioner from hypre.
384  HypreSolver *amg = new HypreBoomerAMG(A);
385  HyprePCG *pcg = new HyprePCG(A);
386  pcg->SetTol(1e-12);
387  pcg->SetMaxIter(200);
388  pcg->SetPrintLevel(2);
389  pcg->SetPreconditioner(*amg);
390  pcg->Mult(B, X);
391 
392  // 13. Recover the parallel grid function corresponding to X. This is the
393  // local finite element solution on each processor.
394  a->RecoverFEMSolution(X, *b, x);
395 
396  // 14. Save the refined mesh and the solution in parallel. This output can
397  // be viewed later using GLVis: "glvis -np <np> -m mesh -g sol".
398  {
399  ostringstream mesh_name, sol_name;
400  mesh_name << "mesh." << setfill('0') << setw(6) << myid;
401  sol_name << "sol." << setfill('0') << setw(6) << myid;
402 
403  ofstream mesh_ofs(mesh_name.str().c_str());
404  mesh_ofs.precision(8);
405  pmesh->Print(mesh_ofs);
406 
407  ofstream sol_ofs(sol_name.str().c_str());
408  sol_ofs.precision(8);
409  x.Save(sol_ofs);
410  }
411 
412  // 15. Send the solution by socket to a GLVis server.
413  if (visualization)
414  {
415  char vishost[] = "localhost";
416  int visport = 19916;
417  socketstream sol_sock(vishost, visport);
418  sol_sock << "parallel " << num_procs << " " << myid << "\n";
419  sol_sock.precision(8);
420  sol_sock << "solution\n" << *pmesh << x << flush;
421  }
422 
423  // 16. Save data in the VisIt format
424  VisItDataCollection visit_dc("Example1-Parallel", pmesh);
425  visit_dc.RegisterField("solution", &x);
426  visit_dc.Save();
427 
428  // 17. Free the used memory.
429  delete pcg;
430  delete amg;
431  delete a;
432  delete b;
433  delete fespace;
434  if (own_fec) { delete fec; }
435  delete pmesh;
436 
437  MPI_Finalize();
438 
439  return 0;
440 }
int GetNPoints() const
Returns the number of the points in the integration rule.
Definition: intrules.hpp:247
Abstract class for all finite elements.
Definition: fe.hpp:243
Arbitrary order non-uniform rational B-splines (NURBS) finite elements.
Definition: fe_coll.hpp:466
void SetTol(double tol)
Definition: hypre.cpp:3569
int Size() const
Return the logical size of the array.
Definition: array.hpp:134
Class for domain integration L(v) := (f, v)
Definition: lininteg.hpp:97
virtual void GetEssentialTrueDofs(const Array< int > &bdr_attr_is_ess, Array< int > &ess_tdof_list, int component=-1)
Definition: pfespace.cpp:790
int GetDim() const
Returns the reference space dimension for the finite element.
Definition: fe.hpp:317
Class for an integration rule - an Array of IntegrationPoint.
Definition: intrules.hpp:90
const IntegrationRule & Get(int GeomType, int Order)
Returns an integration rule for given GeomType and Order.
Definition: intrules.cpp:920
A coefficient that is constant across space and time.
Definition: coefficient.hpp:78
virtual void FormLinearSystem(const Array< int > &ess_tdof_list, Vector &x, Vector &b, OperatorHandle &A, Vector &X, Vector &B, int copy_interior=0)
Form the linear system A X = B, corresponding to this bilinear form and the linear form b(...
virtual void CalcPhysLaplacian(ElementTransformation &Trans, Vector &Laplacian) const
Evaluate the Laplacian of all shape functions of a scalar finite element in reference space at the gi...
Definition: fe.cpp:203
void SetIntPoint(const IntegrationPoint *ip)
Set the integration point ip that weights and Jacobians will be evaluated at.
Definition: eltrans.hpp:85
HYPRE_BigInt GlobalTrueVSize() const
Definition: pfespace.hpp:275
int GetNKV() const
Definition: nurbs.hpp:350
int GetOrder() const
Returns the order of the finite element. In the case of anisotropic orders, returns the maximum order...
Definition: fe.hpp:327
Data type dense matrix using column-major storage.
Definition: densemat.hpp:23
int GetNE() const
Returns number of elements.
Definition: mesh.hpp:846
virtual void Save()
Save the collection and a VisIt root file.
virtual void Save(std::ostream &out) const
Definition: pgridfunc.cpp:841
int Space() const
Returns the type of FunctionSpace on the element.
Definition: fe.hpp:337
double kappa
Definition: ex24.cpp:54
Abstract parallel finite element space.
Definition: pfespace.hpp:28
void SetPrintLevel(int print_lvl)
Definition: hypre.cpp:3589
Geometry::Type GetGeomType() const
Returns the Geometry::Type of the reference element.
Definition: fe.hpp:320
The BoomerAMG solver in hypre.
Definition: hypre.hpp:1387
IntegrationPoint & IntPoint(int i)
Returns a reference to the i-th integration point.
Definition: intrules.hpp:250
Class for parallel linear form.
Definition: plinearform.hpp:26
void Parse()
Parse the command-line options. Note that this function expects all the options provided through the ...
Definition: optparser.cpp:150
constexpr char vishost[]
void mfem_error(const char *msg)
Function called when an error is encountered. Used by the macros MFEM_ABORT, MFEM_ASSERT, MFEM_VERIFY.
Definition: error.cpp:153
double b
Definition: lissajous.cpp:42
void UniformRefinement(int i, const DSTable &, int *, int *, int *)
Definition: mesh.cpp:9143
constexpr int visport
Data collection with VisIt I/O routines.
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const =0
Evaluate the values of all shape functions of a scalar finite element in reference space at the given...
virtual void PrintInfo(std::ostream &out=mfem::out)
In serial, this method calls PrintCharacteristics(). In parallel, additional information about the pa...
Definition: mesh.hpp:1497
T Max() const
Find the maximal element in the array, using the comparison operator &lt; for class T.
Definition: array.cpp:68
void Assemble(int skip_zeros=1)
Assemble the local matrix.
virtual void Print(std::ostream &out=mfem::out) const
Definition: pmesh.cpp:4382
double Weight()
Return the weight of the Jacobian matrix of the transformation at the currently set IntegrationPoint...
Definition: eltrans.hpp:123
void AddBdrFaceIntegrator(LinearFormIntegrator *lfi)
Adds new Boundary Face Integrator. Assumes ownership of lfi.
Definition: linearform.cpp:83
int Dimension() const
Definition: mesh.hpp:911
void PrintUsage(std::ostream &out) const
Print the usage message.
Definition: optparser.cpp:457
void SetMaxIter(int max_iter)
Definition: hypre.cpp:3579
void AddDomainIntegrator(LinearFormIntegrator *lfi)
Adds new Domain Integrator. Assumes ownership of lfi.
Definition: linearform.cpp:39
Abstract base class BilinearFormIntegrator.
Definition: bilininteg.hpp:34
Array< int > bdr_attributes
A list of all unique boundary attributes used by the Mesh.
Definition: mesh.hpp:204
PCG solver in hypre.
Definition: hypre.hpp:1060
int GetDof() const
Returns the number of degrees of freedom in the finite element.
Definition: fe.hpp:323
Base class Coefficients that optionally depend on space and time. These are used by the BilinearFormI...
Definition: coefficient.hpp:39
Collection of finite elements from the same family in multiple dimensions. This class is used to matc...
Definition: fe_coll.hpp:26
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 &#39;var&#39; to receive the value. Enable/disable tags are used to set the bool...
Definition: optparser.hpp:82
int GetNP() const
Definition: nurbs.hpp:341
void SetSize(int nsize)
Change the logical size of the array, keep existing entries.
Definition: array.hpp:674
HYPRE_Int HYPRE_BigInt
virtual void RecoverFEMSolution(const Vector &X, const Vector &b, Vector &x)
virtual void RegisterField(const std::string &field_name, GridFunction *gf)
Add a grid function to the collection and update the root file.
double a
Definition: lissajous.cpp:41
NURBSExtension * NURBSext
Optional NURBS mesh extension.
Definition: mesh.hpp:206
virtual const char * Name() const
Definition: fe_coll.hpp:61
Class for integration point with weight.
Definition: intrules.hpp:25
HYPRE_BigInt GetGlobalNumRows() const
Return the global number of rows.
Definition: hypre.hpp:573
int dim
Definition: ex24.cpp:53
void SetPreconditioner(HypreSolver &precond)
Set the hypre solver to be used as a preconditioner.
Definition: hypre.cpp:3594
void AddDomainIntegrator(BilinearFormIntegrator *bfi)
Adds new Domain Integrator. Assumes ownership of bfi.
void PrintOptions(std::ostream &out) const
Print the options.
Definition: optparser.cpp:327
Class for parallel bilinear form.
Abstract class for hypre&#39;s solvers and preconditioners.
Definition: hypre.hpp:970
Vector data type.
Definition: vector.hpp:60
IntegrationRules RefinedIntRules(1, Quadrature1D::GaussLegendre)
A global object with all refined integration rules.
Definition: intrules.hpp:380
void GetNodes(Vector &node_coord) const
Definition: mesh.cpp:7343
Arbitrary order H1-conforming (continuous) finite elements.
Definition: fe_coll.hpp:216
Class for parallel grid function.
Definition: pgridfunc.hpp:32
void SetSize(int s)
Change the size of the DenseMatrix to s x s.
Definition: densemat.hpp:105
Wrapper for hypre&#39;s ParCSR matrix class.
Definition: hypre.hpp:277
virtual void Mult(const HypreParVector &b, HypreParVector &x) const
Solve Ax=b with hypre&#39;s PCG.
Definition: hypre.cpp:3617
Class for parallel meshes.
Definition: pmesh.hpp:32
IntegrationRules IntRules(0, Quadrature1D::GaussLegendre)
A global object with all integration rules (defined in intrules.cpp)
Definition: intrules.hpp:377
int main()
void EnableStaticCondensation()
Enable the use of static condensation. For details see the description for class StaticCondensation i...
void AddBdrFaceIntegrator(BilinearFormIntegrator *bfi)
Adds new boundary Face Integrator. Assumes ownership of bfi.
bool Good() const
Return true if the command line options were parsed successfully.
Definition: optparser.hpp:150