MFEM  v4.5.1
Finite element discretization library
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
ex15p.cpp
Go to the documentation of this file.
1 // MFEM Example 15 - Parallel Version
2 //
3 // Compile with: make ex15p
4 //
5 // Sample runs: mpirun -np 4 ex15p
6 // mpirun -np 4 ex15p -o 1 -y 0.2
7 // mpirun -np 4 ex15p -o 4 -y 0.1
8 // mpirun -np 4 ex15p -n 5
9 // mpirun -np 4 ex15p -p 1 -n 3
10 //
11 // Other meshes:
12 //
13 // mpirun -np 4 ex15p -m ../data/square-disc-nurbs.mesh
14 // mpirun -np 4 ex15p -m ../data/disc-nurbs.mesh
15 // mpirun -np 4 ex15p -m ../data/fichera.mesh -tf 0.5
16 // mpirun -np 4 ex15p -m ../data/ball-nurbs.mesh -tf 0.5
17 // mpirun -np 4 ex15p -m ../data/mobius-strip.mesh
18 // mpirun -np 4 ex15p -m ../data/amr-quad.mesh
19 // mpirun -np 4 ex15p -m ../data/square-disc.mesh
20 // mpirun -np 4 ex15p -m ../data/escher.mesh -r 2 -tf 0.3
21 //
22 // Different estimators:
23 //
24 // mpirun -np 4 ex15p -est 0 -e 1e-4
25 // mpirun -np 4 ex15p -est 1 -e 1e-6
26 // mpirun -np 4 ex15p -est 1 -o 3 -tf 0.3
27 // mpirun -np 4 ex15p -est 2 -o 2
28 //
29 // Description: Building on Example 6, this example demonstrates dynamic AMR.
30 // The mesh is adapted to a time-dependent solution by refinement
31 // as well as by derefinement. For simplicity, the solution is
32 // prescribed and no time integration is done. However, the error
33 // estimation and refinement/derefinement decisions are realistic.
34 //
35 // At each outer iteration the right hand side function is changed
36 // to mimic a time dependent problem. Within each inner iteration
37 // the problem is solved on a sequence of meshes which are locally
38 // refined according to a chosen error estimator. Currently there
39 // are three error estimators supported: A L2 formulation of the
40 // Zienkiewicz-Zhu error estimator (0), a Kelly error indicator (1)
41 // and a traditional Zienkiewicz-Zhu error estimator (2). At the
42 // end of the inner iteration the error estimates are also used to
43 // identify any elements which may be over-refined and a single
44 // derefinement step is performed. After each refinement or
45 // derefinement step a rebalance operation is performed to keep
46 // the mesh evenly distributed among the available processors.
47 //
48 // The example demonstrates MFEM's capability to refine, derefine
49 // and load balance nonconforming meshes, in 2D and 3D, and on
50 // linear, curved and surface meshes. Interpolation of functions
51 // between coarse and fine meshes, persistent GLVis visualization,
52 // and saving of time-dependent fields for external visualization
53 // with VisIt (visit.llnl.gov) are also illustrated.
54 //
55 // We recommend viewing Examples 1, 6 and 9 before viewing this
56 // example.
57 
58 #include "mfem.hpp"
59 #include <fstream>
60 #include <iostream>
61 
62 using namespace std;
63 using namespace mfem;
64 
65 // Choices for the problem setup. Affect bdr_func and rhs_func.
66 int problem;
68 
69 // Prescribed time-dependent boundary and right-hand side functions.
70 double bdr_func(const Vector &pt, double t);
71 double rhs_func(const Vector &pt, double t);
72 
73 // Update the finite element space, interpolate the solution and perform
74 // parallel load balancing.
77  ParLinearForm &b);
78 
79 
80 int main(int argc, char *argv[])
81 {
82  // 1. Initialize MPI and HYPRE.
83  Mpi::Init(argc, argv);
84  int num_procs = Mpi::WorldSize();
85  int myid = Mpi::WorldRank();
86  Hypre::Init();
87 
88  // 2. Parse command-line options.
89  problem = 0;
90  nfeatures = 1;
91  const char *mesh_file = "../data/star-hilbert.mesh";
92  int order = 2;
93  double t_final = 1.0;
94  double max_elem_error = 1.0e-4;
95  double hysteresis = 0.25; // derefinement safety coefficient
96  int ref_levels = 0;
97  int nc_limit = 3; // maximum level of hanging nodes
98  bool visualization = true;
99  bool visit = false;
100  int which_estimator = 0;
101 
102  OptionsParser args(argc, argv);
103  args.AddOption(&mesh_file, "-m", "--mesh",
104  "Mesh file to use.");
105  args.AddOption(&problem, "-p", "--problem",
106  "Problem setup to use: 0 = spherical front, 1 = ball.");
107  args.AddOption(&nfeatures, "-n", "--nfeatures",
108  "Number of solution features (fronts/balls).");
109  args.AddOption(&order, "-o", "--order",
110  "Finite element order (polynomial degree).");
111  args.AddOption(&max_elem_error, "-e", "--max-err",
112  "Maximum element error");
113  args.AddOption(&hysteresis, "-y", "--hysteresis",
114  "Derefinement safety coefficient.");
115  args.AddOption(&ref_levels, "-r", "--ref-levels",
116  "Number of initial uniform refinement levels.");
117  args.AddOption(&nc_limit, "-l", "--nc-limit",
118  "Maximum level of hanging nodes.");
119  args.AddOption(&t_final, "-tf", "--t-final",
120  "Final time; start time is 0.");
121  args.AddOption(&which_estimator, "-est", "--estimator",
122  "Which estimator to use: "
123  "0 = L2ZZ, 1 = Kelly, 2 = ZZ. Defaults to L2ZZ.");
124  args.AddOption(&visualization, "-vis", "--visualization", "-no-vis",
125  "--no-visualization",
126  "Enable or disable GLVis visualization.");
127  args.AddOption(&visit, "-visit", "--visit-datafiles", "-no-visit",
128  "--no-visit-datafiles",
129  "Save data files for VisIt (visit.llnl.gov) visualization.");
130  args.Parse();
131  if (!args.Good())
132  {
133  if (myid == 0)
134  {
135  args.PrintUsage(cout);
136  }
137  return 1;
138  }
139  if (myid == 0)
140  {
141  args.PrintOptions(cout);
142  }
143 
144  // 3. Read the (serial) mesh from the given mesh file on all processors. We
145  // can handle triangular, quadrilateral, tetrahedral, hexahedral, surface
146  // and volume meshes with the same code.
147  Mesh *mesh = new Mesh(mesh_file, 1, 1);
148  int dim = mesh->Dimension();
149  int sdim = mesh->SpaceDimension();
150 
151  // 4. Project a NURBS mesh to a piecewise-quadratic curved mesh. Make sure
152  // that the mesh is non-conforming if it has quads or hexes and refine it.
153  if (mesh->NURBSext)
154  {
155  mesh->UniformRefinement();
156  if (ref_levels > 0) { ref_levels--; }
157  mesh->SetCurvature(2);
158  }
159  mesh->EnsureNCMesh(true);
160  for (int l = 0; l < ref_levels; l++)
161  {
162  mesh->UniformRefinement();
163  }
164  // Make sure tet-only meshes are marked for local refinement.
165  mesh->Finalize(true);
166 
167  // 5. Define a parallel mesh by partitioning the serial mesh. Once the
168  // parallel mesh is defined, the serial mesh can be deleted.
169  ParMesh pmesh(MPI_COMM_WORLD, *mesh);
170  delete mesh;
171 
172  MFEM_VERIFY(pmesh.bdr_attributes.Size() > 0,
173  "Boundary attributes required in the mesh.");
174  Array<int> ess_bdr(pmesh.bdr_attributes.Max());
175  ess_bdr = 1;
176 
177  // 6. Define a finite element space on the mesh. The polynomial order is one
178  // (linear) by default, but this can be changed on the command line.
179  H1_FECollection fec(order, dim);
180  ParFiniteElementSpace fespace(&pmesh, &fec);
181 
182  // 7. As in Example 1p, we set up bilinear and linear forms corresponding to
183  // the Laplace problem -\Delta u = 1. We don't assemble the discrete
184  // problem yet, this will be done in the inner loop.
185  ParBilinearForm a(&fespace);
186  ParLinearForm b(&fespace);
187 
188  ConstantCoefficient one(1.0);
191 
193  a.AddDomainIntegrator(integ);
195 
196  // 8. The solution vector x and the associated finite element grid function
197  // will be maintained over the AMR iterations.
198  ParGridFunction x(&fespace);
199 
200  // 9. Connect to GLVis. Prepare for VisIt output.
201  char vishost[] = "localhost";
202  int visport = 19916;
203 
204  socketstream sout;
205  if (visualization)
206  {
207  sout.open(vishost, visport);
208  if (!sout)
209  {
210  if (myid == 0)
211  {
212  cout << "Unable to connect to GLVis server at "
213  << vishost << ':' << visport << endl;
214  cout << "GLVis visualization disabled.\n";
215  }
216  visualization = false;
217  }
218  sout.precision(8);
219  }
220 
221  VisItDataCollection visit_dc("Example15-Parallel", &pmesh);
222  visit_dc.RegisterField("solution", &x);
223  int vis_cycle = 0;
224 
225  // 10. As in Example 6p, we set up an estimator that will be used to obtain
226  // element error indicators. The integrator needs to provide the method
227  // ComputeElementFlux. We supply an L2 space for the discontinuous flux
228  // and an H(div) space for the smoothed flux.
229  L2_FECollection flux_fec(order, dim);
230  RT_FECollection smooth_flux_fec(order-1, dim);
231  ErrorEstimator* estimator{nullptr};
232 
233  switch (which_estimator)
234  {
235  case 1:
236  {
237  auto flux_fes = new ParFiniteElementSpace(&pmesh, &flux_fec, sdim);
238  estimator = new KellyErrorEstimator(*integ, x, flux_fes);
239  break;
240  }
241  case 2:
242  {
243  auto flux_fes = new ParFiniteElementSpace(&pmesh, &fec, sdim);
244  estimator = new ZienkiewiczZhuEstimator(*integ, x, flux_fes);
245  break;
246  }
247 
248  default:
249  if (myid == 0)
250  {
251  std::cout << "Unknown estimator. Falling back to L2ZZ." << std::endl;
252  }
253  case 0:
254  {
255  auto flux_fes = new ParFiniteElementSpace(&pmesh, &flux_fec, sdim);
256  auto smooth_flux_fes = new ParFiniteElementSpace(&pmesh, &smooth_flux_fec);
257  estimator = new L2ZienkiewiczZhuEstimator(*integ, x, flux_fes, smooth_flux_fes);
258  break;
259  }
260  }
261 
262  // 11. As in Example 6p, we also need a refiner. This time the refinement
263  // strategy is based on a fixed threshold that is applied locally to each
264  // element. The global threshold is turned off by setting the total error
265  // fraction to zero. We also enforce a maximum refinement ratio between
266  // adjacent elements.
267  ThresholdRefiner refiner(*estimator);
268  refiner.SetTotalErrorFraction(0.0); // use purely local threshold
269  refiner.SetLocalErrorGoal(max_elem_error);
270  refiner.PreferConformingRefinement();
271  refiner.SetNCLimit(nc_limit);
272 
273  // 12. A derefiner selects groups of elements that can be coarsened to form
274  // a larger element. A conservative enough threshold needs to be set to
275  // prevent derefining elements that would immediately be refined again.
276  ThresholdDerefiner derefiner(*estimator);
277  derefiner.SetThreshold(hysteresis * max_elem_error);
278  derefiner.SetNCLimit(nc_limit);
279 
280  // 13. The outer time loop. In each iteration we update the right hand side,
281  // solve the problem on the current mesh, visualize the solution and
282  // refine the mesh as many times as necessary. Then we derefine any
283  // elements which have very small errors.
284  for (double time = 0.0; time < t_final + 1e-10; time += 0.01)
285  {
286  if (myid == 0)
287  {
288  cout << "\nTime " << time << "\n\nRefinement:" << endl;
289  }
290 
291  // Set the current time in the coefficients
292  bdr.SetTime(time);
293  rhs.SetTime(time);
294 
295  // Make sure errors will be recomputed in the following.
296  refiner.Reset();
297  derefiner.Reset();
298 
299  // 14. The inner refinement loop. At the end we want to have the current
300  // time step resolved to the prescribed tolerance in each element.
301  for (int ref_it = 1; ; ref_it++)
302  {
303  HYPRE_BigInt global_dofs = fespace.GlobalTrueVSize();
304  if (myid == 0)
305  {
306  cout << "Iteration: " << ref_it << ", number of unknowns: "
307  << global_dofs << flush;
308  }
309 
310  // 15. Recompute the field on the current mesh: assemble the stiffness
311  // matrix and the right-hand side.
312  a.Assemble();
313  b.Assemble();
314 
315  // 16. Project the exact solution to the essential DOFs.
316  x.ProjectBdrCoefficient(bdr, ess_bdr);
317 
318  // 17. Create and solve the parallel linear system.
319  Array<int> ess_tdof_list;
320  fespace.GetEssentialTrueDofs(ess_bdr, ess_tdof_list);
321 
322  HypreParMatrix A;
323  Vector B, X;
324  a.FormLinearSystem(ess_tdof_list, x, b, A, X, B);
325 
326  HypreBoomerAMG amg(A);
327  amg.SetPrintLevel(0);
328  HyprePCG pcg(A);
329  pcg.SetTol(1e-12);
330  pcg.SetMaxIter(200);
331  pcg.SetPrintLevel(0);
332  pcg.SetPreconditioner(amg);
333  pcg.Mult(B, X);
334 
335  // 18. Extract the local solution on each processor.
336  a.RecoverFEMSolution(X, b, x);
337 
338  // 19. Send the solution by socket to a GLVis server and optionally
339  // save it in VisIt format.
340  if (visualization)
341  {
342  sout << "parallel " << num_procs << " " << myid << "\n";
343  sout << "solution\n" << pmesh << x << flush;
344  }
345  if (visit)
346  {
347  visit_dc.SetCycle(vis_cycle++);
348  visit_dc.SetTime(time);
349  visit_dc.Save();
350  }
351 
352  // 20. Apply the refiner on the mesh. The refiner calls the error
353  // estimator to obtain element errors, then it selects elements to
354  // be refined and finally it modifies the mesh. The Stop() method
355  // determines if all elements satisfy the local threshold.
356  refiner.Apply(pmesh);
357  if (myid == 0)
358  {
359  cout << ", total error: " << estimator->GetTotalError() << endl;
360  }
361 
362  // 21. Quit the AMR loop if the termination criterion has been met
363  if (refiner.Stop())
364  {
365  a.Update(); // Free the assembled data
366  break;
367  }
368 
369  // 22. Update the space, interpolate the solution, rebalance the mesh.
370  UpdateAndRebalance(pmesh, fespace, x, a, b);
371  }
372 
373  // 23. Use error estimates from the last inner iteration to check for
374  // possible derefinements. The derefiner works similarly as the
375  // refiner. The errors are not recomputed because the mesh did not
376  // change (and also the estimator was not Reset() at this time).
377  if (derefiner.Apply(pmesh))
378  {
379  if (myid == 0)
380  {
381  cout << "\nDerefined elements." << endl;
382  }
383 
384  // 24. Update the space and the solution, rebalance the mesh.
385  UpdateAndRebalance(pmesh, fespace, x, a, b);
386  }
387  }
388 
389  delete estimator;
390 
391  // 25. Exit
392  return 0;
393 }
394 
395 
398  ParLinearForm &b)
399 {
400  // Update the space: recalculate the number of DOFs and construct a matrix
401  // that will adjust any GridFunctions to the new mesh state.
402  fespace.Update();
403 
404  // Interpolate the solution on the new mesh by applying the transformation
405  // matrix computed in the finite element space. Multiple GridFunctions could
406  // be updated here.
407  x.Update();
408 
409  if (pmesh.Nonconforming())
410  {
411  // Load balance the mesh.
412  pmesh.Rebalance();
413 
414  // Update the space again, this time a GridFunction redistribution matrix
415  // is created. Apply it to the solution.
416  fespace.Update();
417  x.Update();
418  }
419 
420  // Inform the linear and bilinear forms that the space has changed.
421  a.Update();
422  b.Update();
423 
424  // Free any transformation matrices to save memory.
425  fespace.UpdatesFinished();
426 }
427 
428 
429 const double alpha = 0.02;
430 
431 // Spherical front with a Gaussian cross section and radius t
432 double front(double x, double y, double z, double t, int)
433 {
434  double r = sqrt(x*x + y*y + z*z);
435  return exp(-0.5*pow((r - t)/alpha, 2));
436 }
437 
438 double front_laplace(double x, double y, double z, double t, int dim)
439 {
440  double x2 = x*x, y2 = y*y, z2 = z*z, t2 = t*t;
441  double r = sqrt(x2 + y2 + z2);
442  double a2 = alpha*alpha, a4 = a2*a2;
443  return -exp(-0.5*pow((r - t)/alpha, 2)) / a4 *
444  (-2*t*(x2 + y2 + z2 - (dim-1)*a2/2)/r + x2 + y2 + z2 + t2 - dim*a2);
445 }
446 
447 // Smooth spherical step function with radius t
448 double ball(double x, double y, double z, double t, int)
449 {
450  double r = sqrt(x*x + y*y + z*z);
451  return -atan(2*(r - t)/alpha);
452 }
453 
454 double ball_laplace(double x, double y, double z, double t, int dim)
455 {
456  double x2 = x*x, y2 = y*y, z2 = z*z, t2 = 4*t*t;
457  double r = sqrt(x2 + y2 + z2);
458  double a2 = alpha*alpha;
459  double den = pow(-a2 - 4*(x2 + y2 + z2 - 2*r*t) - t2, 2.0);
460  return (dim == 2) ? 2*alpha*(a2 + t2 - 4*x2 - 4*y2)/r/den
461  /* */ : 4*alpha*(a2 + t2 - 4*r*t)/r/den;
462 }
463 
464 // Composes several features into one function
465 template<typename F0, typename F1>
466 double composite_func(const Vector &pt, double t, F0 f0, F1 f1)
467 {
468  int dim = pt.Size();
469  double x = pt(0), y = pt(1), z = 0.0;
470  if (dim == 3) { z = pt(2); }
471 
472  if (problem == 0)
473  {
474  if (nfeatures <= 1)
475  {
476  return f0(x, y, z, t, dim);
477  }
478  else
479  {
480  double sum = 0.0;
481  for (int i = 0; i < nfeatures; i++)
482  {
483  double x0 = 0.5*cos(2*M_PI * i / nfeatures);
484  double y0 = 0.5*sin(2*M_PI * i / nfeatures);
485  sum += f0(x - x0, y - y0, z, t, dim);
486  }
487  return sum;
488  }
489  }
490  else
491  {
492  double sum = 0.0;
493  for (int i = 0; i < nfeatures; i++)
494  {
495  double x0 = 0.5*cos(2*M_PI * i / nfeatures + M_PI*t);
496  double y0 = 0.5*sin(2*M_PI * i / nfeatures + M_PI*t);
497  sum += f1(x - x0, y - y0, z, 0.25, dim);
498  }
499  return sum;
500  }
501 }
502 
503 // Exact solution, used for the Dirichlet BC.
504 double bdr_func(const Vector &pt, double t)
505 {
506  return composite_func(pt, t, front, ball);
507 }
508 
509 // Laplace of the exact solution, used for the right hand side.
510 double rhs_func(const Vector &pt, double t)
511 {
512  return composite_func(pt, t, front_laplace, ball_laplace);
513 }
void SetTol(double tol)
Definition: hypre.cpp:3995
double front_laplace(double x, double y, double z, double t, int dim)
Definition: ex15.cpp:378
int Size() const
Return the logical size of the array.
Definition: array.hpp:138
Class for domain integration L(v) := (f, v)
Definition: lininteg.hpp:108
virtual void GetEssentialTrueDofs(const Array< int > &bdr_attr_is_ess, Array< int > &ess_tdof_list, int component=-1)
Definition: pfespace.cpp:1032
double bdr_func(const Vector &pt, double t)
Definition: ex15.cpp:444
void SetCycle(int c)
Set time cycle (for time-dependent simulations)
A coefficient that is constant across space and time.
Definition: coefficient.hpp:84
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(...
The KellyErrorEstimator class provides a fast error indication strategy for smooth scalar parallel pr...
Definition: estimators.hpp:555
void SetNCLimit(int nc_limit_)
Set the maximum ratio of refinement levels of adjacent elements (0 = unlimited).
virtual void Update(bool want_transform=true)
Definition: pfespace.cpp:3287
HYPRE_BigInt GlobalTrueVSize() const
Definition: pfespace.hpp:285
int Size() const
Returns the size of the vector.
Definition: vector.hpp:200
bool Stop() const
Check if STOP action is requested, e.g. stopping criterion is satisfied.
Abstract parallel finite element space.
Definition: pfespace.hpp:28
virtual void SetTime(double t)
Set the time for time dependent coefficients.
Definition: coefficient.hpp:50
double rhs_func(const Vector &pt, double t)
Definition: ex15.cpp:450
The L2ZienkiewiczZhuEstimator class implements the Zienkiewicz-Zhu error estimation procedure where t...
Definition: estimators.hpp:328
void SetPrintLevel(int print_lvl)
Definition: hypre.cpp:4015
void Update(ParFiniteElementSpace *pf=NULL)
Update the object according to the given new FE space *pf.
Definition: plinearform.cpp:21
virtual void Update(FiniteElementSpace *nfes=NULL)
Update the FiniteElementSpace and delete all data associated with the old one.
The BoomerAMG solver in hypre.
Definition: hypre.hpp:1579
bool Apply(Mesh &mesh)
Perform the mesh operation.
void Rebalance()
Definition: pmesh.cpp:3989
bool Nonconforming() const
Definition: mesh.hpp:1651
Class for parallel linear form.
Definition: plinearform.hpp:26
virtual void Reset()
Reset the associated estimator.
void Parse()
Parse the command-line options. Note that this function expects all the options provided through the ...
Definition: optparser.cpp:151
constexpr char vishost[]
void EnsureNCMesh(bool simplices_nonconforming=false)
Definition: mesh.cpp:9473
virtual void Update()
Transform by the Space UpdateMatrix (e.g., on Mesh change).
Definition: pgridfunc.cpp:81
void Assemble()
Assembles the ParLinearForm i.e. sums over all domain/bdr integrators.
Definition: plinearform.cpp:46
double b
Definition: lissajous.cpp:42
void UniformRefinement(int i, const DSTable &, int *, int *, int *)
Definition: mesh.cpp:9816
virtual void Reset()
Reset the associated estimator.
void SetPrintLevel(int print_level)
Definition: hypre.hpp:1659
constexpr int visport
Mesh refinement operator using an error threshold.
Data collection with VisIt I/O routines.
virtual void SetCurvature(int order, bool discont=false, int space_dim=-1, int ordering=1)
Definition: mesh.cpp:5343
T Max() const
Find the maximal element in the array, using the comparison operator &lt; for class T.
Definition: array.cpp:68
int nfeatures
Definition: ex15.cpp:63
void Assemble(int skip_zeros=1)
Assemble the local matrix.
Base class for all element based error estimators.
Definition: estimators.hpp:41
int Dimension() const
Definition: mesh.hpp:1006
void PrintUsage(std::ostream &out) const
Print the usage message.
Definition: optparser.cpp:454
void SetTime(double t)
Set physical time (for time-dependent simulations)
void SetMaxIter(int max_iter)
Definition: hypre.cpp:4005
Arbitrary order H(div)-conforming Raviart-Thomas finite elements.
Definition: fe_coll.hpp:341
void UpdateAndRebalance(ParMesh &pmesh, ParFiniteElementSpace &fespace, ParGridFunction &x, ParBilinearForm &a, ParLinearForm &b)
Definition: ex15p.cpp:396
int SpaceDimension() const
Definition: mesh.hpp:1007
void AddDomainIntegrator(LinearFormIntegrator *lfi)
Adds new Domain Integrator. Assumes ownership of lfi.
Definition: linearform.cpp:41
Abstract base class BilinearFormIntegrator.
Definition: bilininteg.hpp:35
Array< int > bdr_attributes
A list of all unique boundary attributes used by the Mesh.
Definition: mesh.hpp:270
double ball_laplace(double x, double y, double z, double t, int dim)
Definition: ex15.cpp:394
PCG solver in hypre.
Definition: hypre.hpp:1204
void ProjectBdrCoefficient(Coefficient *coeff[], VectorCoefficient *vcoeff, Array< int > &attr)
Definition: pgridfunc.cpp:633
virtual void Save() override
Save the collection and a VisIt root file.
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 problem
Definition: ex15.cpp:62
HYPRE_Int HYPRE_BigInt
void PreferConformingRefinement()
Use conforming refinement, if possible (triangles, tetrahedra) – this is the default.
virtual void RecoverFEMSolution(const Vector &X, const Vector &b, Vector &x)
double a
Definition: lissajous.cpp:41
NURBSExtension * NURBSext
Optional NURBS mesh extension.
Definition: mesh.hpp:272
virtual void Finalize(bool refine=false, bool fix_orientation=false)
Finalize the construction of a general Mesh.
Definition: mesh.cpp:2983
int dim
Definition: ex24.cpp:53
void SetPreconditioner(HypreSolver &precond)
Set the hypre solver to be used as a preconditioner.
Definition: hypre.cpp:4020
void AddDomainIntegrator(BilinearFormIntegrator *bfi)
Adds new Domain Integrator. Assumes ownership of bfi.
void PrintOptions(std::ostream &out) const
Print the options.
Definition: optparser.cpp:324
Class for parallel bilinear form.
double ball(double x, double y, double z, double t, int)
Definition: ex15.cpp:388
int open(const char hostname[], int port)
Open the socket stream on &#39;port&#39; at &#39;hostname&#39;.
void SetThreshold(double thresh)
Set the de-refinement threshold. The default value is zero.
RefCoord t[3]
The ZienkiewiczZhuEstimator class implements the Zienkiewicz-Zhu error estimation procedure...
Definition: estimators.hpp:88
void SetLocalErrorGoal(double err_goal)
Set the local stopping criterion: stop when local_err_i &lt;= local_err_goal. The default value is zero...
const double alpha
Definition: ex15.cpp:369
A general function coefficient.
Vector data type.
Definition: vector.hpp:60
Arbitrary order H1-conforming (continuous) finite elements.
Definition: fe_coll.hpp:220
Class for parallel grid function.
Definition: pgridfunc.hpp:32
void SetTotalErrorFraction(double fraction)
Set the total fraction used in the computation of the threshold. The default value is 1/2...
Wrapper for hypre&#39;s ParCSR matrix class.
Definition: hypre.hpp:343
virtual void UpdatesFinished()
Free ParGridFunction transformation matrix (if any), to save memory.
Definition: pfespace.hpp:423
double front(double x, double y, double z, double t, int)
Definition: ex15.cpp:372
virtual void Mult(const HypreParVector &b, HypreParVector &x) const
Solve Ax=b with hypre&#39;s PCG.
Definition: hypre.cpp:4043
Class for parallel meshes.
Definition: pmesh.hpp:32
double composite_func(const Vector &pt, double t, F0 f0, F1 f1)
Definition: ex15.cpp:406
int main()
void SetNCLimit(int nc_limit_)
Set the maximum ratio of refinement levels of adjacent elements (0 = unlimited).
Arbitrary order &quot;L2-conforming&quot; discontinuous finite elements.
Definition: fe_coll.hpp:288
De-refinement operator using an error threshold.
virtual void RegisterField(const std::string &field_name, GridFunction *gf) override
Add a grid function to the collection and update the root file.
bool Good() const
Return true if the command line options were parsed successfully.
Definition: optparser.hpp:150