MFEM  v4.6.0
Finite element discretization library
findpts.cpp
Go to the documentation of this file.
1 // Copyright (c) 2010-2023, Lawrence Livermore National Security, LLC. Produced
2 // at the Lawrence Livermore National Laboratory. All Rights reserved. See files
3 // LICENSE and NOTICE for details. LLNL-CODE-806117.
4 //
5 // This file is part of the MFEM library. For more information and source code
6 // availability visit https://mfem.org.
7 //
8 // MFEM is free software; you can redistribute it and/or modify it under the
9 // terms of the BSD-3 license. We welcome feedback and contributions, see file
10 // CONTRIBUTING.md for details.
11 //
12 // -------------------------------------------------------------
13 // Find Points Miniapp: Evaluate grid function in physical space
14 // -------------------------------------------------------------
15 //
16 // This miniapp demonstrates the interpolation of a high-order grid function on
17 // a set of points in physical-space. The miniapp is based on GSLIB-FindPoints,
18 // which provides two key functionalities. First, for a given set of points in
19 // the physical-space, it determines the computational coordinates (element
20 // number, reference-space coordinates inside the element, and processor number
21 // [in parallel]) for each point. Second, based on computational coordinates, it
22 // interpolates a grid function in the given points. Inside GSLIB, computation
23 // of the coordinates requires use of a Hash Table to identify the candidate
24 // processor and element for each point, followed by the Newton's method to
25 // determine the reference-space coordinates inside the candidate element.
26 //
27 // Compile with: make findpts
28 //
29 // Sample runs:
30 // findpts -m ../../data/rt-2d-p4-tri.mesh -o 4
31 // findpts -m ../../data/inline-tri.mesh -o 3
32 // findpts -m ../../data/inline-quad.mesh -o 3
33 // findpts -m ../../data/inline-quad.mesh -o 3 -po 1
34 // findpts -m ../../data/inline-quad.mesh -o 3 -po 1 -fo 1 -nc 2
35 // findpts -m ../../data/inline-quad.mesh -o 3 -hr -pr -mpr -mo 2
36 // findpts -m ../../data/inline-quad.mesh -o 3 -hr -pr -mpr -mo 3
37 // findpts -m ../../data/inline-tet.mesh -o 3
38 // findpts -m ../../data/inline-hex.mesh -o 3
39 // findpts -m ../../data/inline-wedge.mesh -o 3
40 // findpts -m ../../data/amr-quad.mesh -o 2
41 // findpts -m ../../data/rt-2d-q3.mesh -o 3 -mo 4 -ft 2
42 // findpts -m ../../data/square-mixed.mesh -o 2 -mo 2
43 // findpts -m ../../data/square-mixed.mesh -o 2 -mo 2 -hr -pr -mpr
44 // findpts -m ../../data/square-mixed.mesh -o 2 -mo 3 -ft 2
45 // findpts -m ../../data/fichera-mixed.mesh -o 3 -mo 2
46 // findpts -m ../../data/inline-pyramid.mesh -o 1 -mo 1
47 // findpts -m ../../data/tinyzoo-3d.mesh -o 1 -mo 1
48 
49 #include "mfem.hpp"
50 #include "../common/mfem-common.hpp"
51 
52 using namespace mfem;
53 using namespace std;
54 
55 // Experimental - class used to update GridFunction post p-refinement.
56 // Can change as support for p-refinement evolves.
57 class PRefinementGFUpdate
58 {
59 private:
60  FiniteElementSpace *src;
61 
62 public:
63  /// @brief Used to Update GridFunction post p-refinement.
64  /// Initialize with FESpace prior to p-refinement.
65  PRefinementGFUpdate(const FiniteElementSpace& src_);
66 
67  /// Destructor
68  ~PRefinementGFUpdate();
69 
70  /// Update source FiniteElementSpace used to construct the
71  /// PRefinementTransfer operator.
72  void SetSourceFESpace(const FiniteElementSpace& src_);
73 
74  /// @brief Update GridFunction using PRefinementTransferOperator.
75  /// Do not use GridFunction->Update() prior to this method as it is
76  /// handled internally.
77  void GridFunctionUpdate(GridFunction &targf);
78 };
79 
80 PRefinementGFUpdate::PRefinementGFUpdate(const FiniteElementSpace &src_)
81 {
82  src = new FiniteElementSpace(src_);
83 }
84 
85 PRefinementGFUpdate::~PRefinementGFUpdate()
86 {
87  delete src;
88 }
89 
90 void PRefinementGFUpdate::SetSourceFESpace(const FiniteElementSpace &src_)
91 {
92  if (src) { delete src; }
93  src = new FiniteElementSpace(src_);
94 }
95 
96 void PRefinementGFUpdate::GridFunctionUpdate(GridFunction &targf)
97 {
98  MFEM_VERIFY(targf.GetSequence() != targf.FESpace()->GetSequence(),
99  "GridFunction should not be updated prior to UpdateGF.");
100  Vector srcgf = targf;
101  targf.Update();
103  PRefinementTransferOperator(*src, *(targf.FESpace()));
104  preft.Mult(srcgf, targf);
105 }
106 
107 // Experimental - required for visualizing functions on p-refined spaces.
108 GridFunction* ProlongToMaxOrder(const GridFunction *x, const int fieldtype)
109 {
110  const FiniteElementSpace *fespace = x->FESpace();
111  Mesh *mesh = fespace->GetMesh();
112  const FiniteElementCollection *fec = fespace->FEColl();
113  const int vdim = fespace->GetVDim();
114 
115  // find the max order in the space
116  int max_order = fespace->GetMaxElementOrder();
117 
118  // create a visualization space of max order for all elements
119  FiniteElementCollection *fecInt = NULL;
120  if (fieldtype == 0)
121  {
122  fecInt = new H1_FECollection(max_order, mesh->Dimension());
123  }
124  else if (fieldtype == 1)
125  {
126  fecInt = new L2_FECollection(max_order, mesh->Dimension());
127  }
128  FiniteElementSpace *spaceInt = new FiniteElementSpace(mesh, fecInt,
129  fespace->GetVDim(),
130  fespace->GetOrdering());
131 
133  DenseMatrix I;
134 
135  GridFunction *xInt = new GridFunction(spaceInt);
136 
137  // interpolate solution vector in the larger space
138  for (int i = 0; i < mesh->GetNE(); i++)
139  {
140  Geometry::Type geom = mesh->GetElementGeometry(i);
141  T.SetIdentityTransformation(geom);
142 
143  Array<int> dofs;
144  fespace->GetElementVDofs(i, dofs);
145  Vector elemvect(0), vectInt(0);
146  x->GetSubVector(dofs, elemvect);
147  DenseMatrix elemvecMat(elemvect.GetData(), dofs.Size()/vdim, vdim);
148 
149  const auto *fe = fec->GetFE(geom, fespace->GetElementOrder(i));
150  const auto *feInt = fecInt->GetFE(geom, max_order);
151 
152  feInt->GetTransferMatrix(*fe, T, I);
153 
154  spaceInt->GetElementVDofs(i, dofs);
155  vectInt.SetSize(dofs.Size());
156  DenseMatrix vectIntMat(vectInt.GetData(), dofs.Size()/vdim, vdim);
157 
158  Mult(I, elemvecMat, vectIntMat);
159  xInt->SetSubVector(dofs, vectInt);
160  }
161 
162  xInt->MakeOwner(fecInt);
163  return xInt;
164 }
165 
167  const char *title)
168 {
169  Mesh *mesh = fespace.GetMesh();
170  L2_FECollection order_coll = L2_FECollection(0, mesh->Dimension());
171  FiniteElementSpace order_space = FiniteElementSpace(mesh, &order_coll);
172  GridFunction order_gf = GridFunction(&order_space);
173 
174  for (int e = 0; e < mesh->GetNE(); e++)
175  {
176  order_gf(e) = fespace.GetElementOrder(e);
177  }
178 
179  socketstream vis1;
180  common::VisualizeField(vis1, "localhost", 19916, order_gf, title,
181  400, 0, 400, 400, "RjmAcp");
182 }
183 
184 double func_order;
185 
186 // Scalar function to project
187 double field_func(const Vector &x)
188 {
189  const int dim = x.Size();
190  double res = 0.0;
191  for (int d = 0; d < dim; d++) { res += std::pow(x(d), func_order); }
192  return res;
193 }
194 
195 void F_exact(const Vector &p, Vector &F)
196 {
197  F(0) = field_func(p);
198  for (int i = 1; i < F.Size(); i++) { F(i) = (i+1)*F(0); }
199 }
200 
201 int main (int argc, char *argv[])
202 {
203  // Set the method's default parameters.
204  const char *mesh_file = "../../data/rt-2d-q3.mesh";
205  int order = 3;
206  int mesh_poly_deg = 3;
207  int rs_levels = 0;
208  bool visualization = true;
209  int fieldtype = 0;
210  int ncomp = 1;
211  bool hrefinement = false;
212  bool prefinement = false;
213  int point_ordering = 0;
214  int gf_ordering = 0;
215  bool mesh_prefinement = false;
216 
217  // Parse command-line options.
218  OptionsParser args(argc, argv);
219  args.AddOption(&mesh_file, "-m", "--mesh",
220  "Mesh file to use.");
221  args.AddOption(&order, "-o", "--order",
222  "Finite element order (polynomial degree).");
223  args.AddOption(&mesh_poly_deg, "-mo", "--mesh-order",
224  "Polynomial degree of mesh finite element space.");
225  args.AddOption(&rs_levels, "-rs", "--refine-serial",
226  "Number of times to refine the mesh uniformly in serial.");
227  args.AddOption(&fieldtype, "-ft", "--field-type",
228  "Field type: 0 - H1, 1 - L2, 2 - H(div), 3 - H(curl).");
229  args.AddOption(&ncomp, "-nc", "--ncomp",
230  "Number of components for H1 or L2 GridFunctions");
231  args.AddOption(&visualization, "-vis", "--visualization", "-no-vis",
232  "--no-visualization",
233  "Enable or disable GLVis visualization.");
234  args.AddOption(&hrefinement, "-hr", "--h-refinement", "-no-hr",
235  "--no-h-refinement",
236  "Do random h refinements to mesh (does not work for pyramids).");
237  args.AddOption(&prefinement, "-pr", "--p-refinement", "-no-pr",
238  "--no-p-refinement",
239  "Do random p refinements to solution field (does not work for pyramids).");
240  args.AddOption(&point_ordering, "-po", "--point-ordering",
241  "Ordering of points to be found."
242  "0 (default): byNodes, 1: byVDIM");
243  args.AddOption(&gf_ordering, "-fo", "--fespace-ordering",
244  "Ordering of fespace that will be used for gridfunction to be interpolated."
245  "0 (default): byNodes, 1: byVDIM");
246  args.AddOption(&mesh_prefinement, "-mpr", "--mesh-p-refinement", "-no-mpr",
247  "--no-mesh-p-refinement",
248  "Do random p refinements to mesh Nodes.");
249 
250  args.Parse();
251  if (!args.Good())
252  {
253  args.PrintUsage(cout);
254  return 1;
255  }
256  args.PrintOptions(cout);
257 
258  func_order = std::min(order, 2);
259 
260  // Initialize and refine the starting mesh.
261  Mesh mesh(mesh_file, 1, 1, false);
262  for (int lev = 0; lev < rs_levels; lev++) { mesh.UniformRefinement(); }
263  const int dim = mesh.Dimension();
264  cout << "Mesh curvature of the original mesh: ";
265  if (mesh.GetNodes()) { cout << mesh.GetNodes()->OwnFEC()->Name(); }
266  else { cout << "(NONE)"; }
267  cout << endl;
268 
269  // Mesh bounding box.
270  Vector pos_min, pos_max;
271  MFEM_VERIFY(mesh_poly_deg > 0, "The order of the mesh must be positive.");
272  mesh.GetBoundingBox(pos_min, pos_max, mesh_poly_deg);
273  if (hrefinement || prefinement || mesh_prefinement) { mesh.EnsureNCMesh(true); }
274  cout << "--- Generating equidistant point for:\n"
275  << "x in [" << pos_min(0) << ", " << pos_max(0) << "]\n"
276  << "y in [" << pos_min(1) << ", " << pos_max(1) << "]\n";
277  if (dim == 3)
278  {
279  cout << "z in [" << pos_min(2) << ", " << pos_max(2) << "]\n";
280  }
281 
282  // Random h-refinements to mesh
283  if (hrefinement) { mesh.RandomRefinement(0.5); }
284 
285  // Curve the mesh based on the chosen polynomial degree.
286  H1_FECollection fecm(mesh_poly_deg, dim);
287  FiniteElementSpace fespace(&mesh, &fecm, dim);
288  mesh.SetNodalFESpace(&fespace);
289  GridFunction Nodes(&fespace);
290  mesh.SetNodalGridFunction(&Nodes);
291  cout << "Mesh curvature of the curved mesh: " << fecm.Name() << endl;
292 
293  PRefinementGFUpdate preft_fespace = PRefinementGFUpdate(fespace);
294 
295  if (mesh_prefinement)
296  {
297  for (int e = 0; e < mesh.GetNE(); e++)
298  {
299  if ((double) rand() / RAND_MAX < 0.2)
300  {
301  int element_order = fespace.GetElementOrder(e);
302  fespace.SetElementOrder(e, element_order + 1);
303  }
304  }
305  fespace.Update(false);
306  preft_fespace.GridFunctionUpdate(Nodes);
307  }
308 
309  MFEM_VERIFY(ncomp > 0, "Invalid number of components.");
310  int vec_dim = ncomp;
311  FiniteElementCollection *fec = NULL;
312  if (fieldtype == 0)
313  {
314  fec = new H1_FECollection(order, dim);
315  cout << "H1-GridFunction\n";
316  }
317  else if (fieldtype == 1)
318  {
319  fec = new L2_FECollection(order, dim);
320  cout << "L2-GridFunction\n";
321  }
322  else if (fieldtype == 2)
323  {
324  fec = new RT_FECollection(order, dim);
325  ncomp = 1;
326  vec_dim = dim;
327  cout << "H(div)-GridFunction\n";
328  }
329  else if (fieldtype == 3)
330  {
331  fec = new ND_FECollection(order, dim);
332  ncomp = 1;
333  vec_dim = dim;
334  cout << "H(curl)-GridFunction\n";
335  }
336  else
337  {
338  MFEM_ABORT("Invalid field type.");
339  }
340  FiniteElementSpace sc_fes(&mesh, fec, ncomp, gf_ordering);
341  GridFunction field_vals(&sc_fes);
342 
343  // Random p-refinements to the solution field
344  if (prefinement)
345  {
346  for (int e = 0; e < mesh.GetNE(); e++)
347  {
348  if ((double) rand() / RAND_MAX < 0.5)
349  {
350  int element_order = sc_fes.GetElementOrder(e);
351  sc_fes.SetElementOrder(e, element_order + 1);
352  }
353  }
354  sc_fes.Update(false);
355  field_vals.Update();
356  }
357 
358  GridFunction *mesh_nodes_pref = mesh_prefinement ?
359  ProlongToMaxOrder(&Nodes, 0) :
360  &Nodes;
361  if (mesh_prefinement && visualization)
362  {
363  mesh.SetNodalGridFunction(mesh_nodes_pref);
364  VisualizeFESpacePolynomialOrder(fespace, "Mesh Polynomial Order");
365  mesh.SetNodalGridFunction(&Nodes);
366  }
367 
368  if (prefinement && visualization)
369  {
370  mesh.SetNodalGridFunction(mesh_nodes_pref);
371  VisualizeFESpacePolynomialOrder(sc_fes, "Solution Polynomial Order");
372  mesh.SetNodalGridFunction(&Nodes);
373  }
374 
375  // Project the GridFunction using VectorFunctionCoefficient.
377  field_vals.ProjectCoefficient(F);
378 
379  GridFunction *field_vals_pref = prefinement ?
380  ProlongToMaxOrder(&field_vals, fieldtype) :
381  &field_vals;
382 
383  // Display the mesh and the field through glvis.
384  if (visualization)
385  {
386  if (mesh_prefinement) { mesh.SetNodalGridFunction(mesh_nodes_pref); }
387  socketstream vis1;
388  common::VisualizeField(vis1, "localhost", 19916, *field_vals_pref,
389  "Solution",
390  0, 0, 400, 400, "RmjA*****");
391  if (mesh_prefinement) { mesh.SetNodalGridFunction(&Nodes); }
392  }
393 
394  // Generate equidistant points in physical coordinates over the whole mesh.
395  // Note that some points might be outside, if the mesh is not a box. Note
396  // also that all tasks search the same points (not mandatory).
397  const int pts_cnt_1D = 25;
398  int pts_cnt = pow(pts_cnt_1D, dim);
399  Vector vxyz(pts_cnt * dim);
400  if (dim == 2)
401  {
403  const IntegrationRule &ir = el.GetNodes();
404  for (int i = 0; i < ir.GetNPoints(); i++)
405  {
406  const IntegrationPoint &ip = ir.IntPoint(i);
407  if (point_ordering == Ordering::byNODES)
408  {
409  vxyz(i) = pos_min(0) + ip.x * (pos_max(0)-pos_min(0));
410  vxyz(pts_cnt + i) = pos_min(1) + ip.y * (pos_max(1)-pos_min(1));
411  }
412  else
413  {
414  vxyz(i*dim + 0) = pos_min(0) + ip.x * (pos_max(0)-pos_min(0));
415  vxyz(i*dim + 1) = pos_min(1) + ip.y * (pos_max(1)-pos_min(1));
416  }
417  }
418  }
419  else
420  {
422  const IntegrationRule &ir = el.GetNodes();
423  for (int i = 0; i < ir.GetNPoints(); i++)
424  {
425  const IntegrationPoint &ip = ir.IntPoint(i);
426  if (point_ordering == Ordering::byNODES)
427  {
428  vxyz(i) = pos_min(0) + ip.x * (pos_max(0)-pos_min(0));
429  vxyz(pts_cnt + i) = pos_min(1) + ip.y * (pos_max(1)-pos_min(1));
430  vxyz(2*pts_cnt + i) = pos_min(2) + ip.z * (pos_max(2)-pos_min(2));
431  }
432  else
433  {
434  vxyz(i*dim + 0) = pos_min(0) + ip.x * (pos_max(0)-pos_min(0));
435  vxyz(i*dim + 1) = pos_min(1) + ip.y * (pos_max(1)-pos_min(1));
436  vxyz(i*dim + 2) = pos_min(2) + ip.z * (pos_max(2)-pos_min(2));
437  }
438  }
439  }
440 
441  // Find and Interpolate FE function values on the desired points.
442  Vector interp_vals(pts_cnt*vec_dim);
443  FindPointsGSLIB finder;
444  finder.Setup(mesh);
446  finder.Interpolate(vxyz, field_vals, interp_vals, point_ordering);
447  Array<unsigned int> code_out = finder.GetCode();
448  Vector dist_p_out = finder.GetDist();
449 
450  int face_pts = 0, not_found = 0, found = 0;
451  double error = 0.0, max_err = 0.0, max_dist = 0.0;
452  Vector pos(dim);
453  for (int j = 0; j < vec_dim; j++)
454  {
455  for (int i = 0; i < pts_cnt; i++)
456  {
457  if (code_out[i] < 2)
458  {
459  if (j == 0) { found++; }
460  for (int d = 0; d < dim; d++)
461  {
462  pos(d) = point_ordering == Ordering::byNODES ?
463  vxyz(d*pts_cnt + i) :
464  vxyz(i*dim + d);
465  }
466  Vector exact_val(vec_dim);
467  F_exact(pos, exact_val);
468  error = gf_ordering == Ordering::byNODES ?
469  fabs(exact_val(j) - interp_vals[i + j*pts_cnt]) :
470  fabs(exact_val(j) - interp_vals[i*vec_dim + j]);
471  max_err = std::max(max_err, error);
472  max_dist = std::max(max_dist, dist_p_out(i));
473  if (code_out[i] == 1 && j == 0) { face_pts++; }
474  }
475  else { if (j == 0) { not_found++; } }
476  }
477  }
478 
479  cout << setprecision(16)
480  << "Searched points: " << pts_cnt
481  << "\nFound points: " << found
482  << "\nMax interp error: " << max_err
483  << "\nMax dist (of found): " << max_dist
484  << "\nPoints not found: " << not_found
485  << "\nPoints on faces: " << face_pts << endl;
486 
487  // Free the internal gslib data.
488  finder.FreeData();
489 
490  if (mesh_prefinement) { delete mesh_nodes_pref; }
491  if (prefinement) { delete field_vals_pref; }
492  delete fec;
493 
494  return 0;
495 }
virtual void Interpolate(const GridFunction &field_in, Vector &field_out)
Definition: gslib.cpp:820
void SetSubVector(const Array< int > &dofs, const double value)
Set the entries listed in dofs to the given value.
Definition: vector.cpp:605
int GetNPoints() const
Returns the number of the points in the integration rule.
Definition: intrules.hpp:253
Class for an integration rule - an Array of IntegrationPoint.
Definition: intrules.hpp:96
Class for grid function - Vector with associated FE space.
Definition: gridfunc.hpp:30
virtual void Update(bool want_transform=true)
Reflect changes in the mesh: update number of DOFs, etc. Also, calculate GridFunction transformation ...
Definition: fespace.cpp:3402
void PrintOptions(std::ostream &out) const
Print the options.
Definition: optparser.cpp:331
Arbitrary order L2 elements in 3D on a cube.
Definition: fe_l2.hpp:78
int Dimension() const
Dimension of the reference space used within the elements.
Definition: mesh.hpp:1020
Matrix-free transfer operator between finite element spaces on the same mesh.
Definition: transfer.hpp:428
const FiniteElement * GetFE(Geometry::Type geom, int p) const
Variable order version of FiniteElementForGeometry().
Definition: fe_coll.hpp:191
void SetSize(int s)
Resize the vector to size s.
Definition: vector.hpp:517
void PrintUsage(std::ostream &out) const
Print the usage message.
Definition: optparser.cpp:462
void SetElementOrder(int i, int p)
Sets the order of the i&#39;th finite element.
Definition: fespace.cpp:151
void Mult(const Table &A, const Table &B, Table &C)
C = A * B (as boolean matrices)
Definition: table.cpp:475
void GetBoundingBox(Vector &min, Vector &max, int ref=2)
Returns the minimum and maximum corners of the mesh bounding box.
Definition: mesh.cpp:136
int Size() const
Returns the size of the vector.
Definition: vector.hpp:197
Data type dense matrix using column-major storage.
Definition: densemat.hpp:23
bool Good() const
Return true if the command line options were parsed successfully.
Definition: optparser.hpp:159
long GetSequence() const
Definition: gridfunc.hpp:689
void MakeOwner(FiniteElementCollection *fec_)
Make the GridFunction the owner of fec and fes.
Definition: gridfunc.hpp:122
STL namespace.
void GetSubVector(const Array< int > &dofs, Vector &elemvect) const
Extract entries listed in dofs to the output Vector elemvect.
Definition: vector.cpp:579
double field_func(const Vector &x)
Definition: findpts.cpp:187
virtual void SetL2AvgType(AvgType avgtype_)
Definition: gslib.hpp:184
IntegrationPoint & IntPoint(int i)
Returns a reference to the i-th integration point.
Definition: intrules.hpp:256
void VisualizeFESpacePolynomialOrder(FiniteElementSpace &fespace, const char *title)
Definition: findpts.cpp:166
int GetMaxElementOrder() const
Return the maximum polynomial order.
Definition: fespace.hpp:573
void Parse()
Parse the command-line options. Note that this function expects all the options provided through the ...
Definition: optparser.cpp:151
virtual void FreeData()
Definition: gslib.cpp:272
void EnsureNCMesh(bool simplices_nonconforming=false)
Definition: mesh.cpp:9888
const FiniteElementCollection * FEColl() const
Definition: fespace.hpp:723
void UniformRefinement(int i, const DSTable &, int *, int *, int *)
Definition: mesh.cpp:10232
GridFunction * ProlongToMaxOrder(const GridFunction *x, const int fieldtype)
Definition: findpts.cpp:108
void Setup(Mesh &m, const double bb_t=0.1, const double newt_tol=1.0e-12, const int npt_max=256)
Definition: gslib.cpp:107
long GetSequence() const
Definition: fespace.hpp:1271
void RandomRefinement(double prob, bool aniso=false, int nonconforming=-1, int nc_limit=0)
Refine each element with given probability. Uses GeneralRefinement.
Definition: mesh.cpp:9912
virtual void GetTransferMatrix(const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &I) const
Return interpolation matrix, I, which maps dofs from a coarse element, fe, to the fine dofs on this f...
Definition: fe_base.cpp:119
virtual void SetNodalFESpace(FiniteElementSpace *nfes)
Definition: mesh.cpp:5577
FiniteElementSpace * FESpace()
Definition: gridfunc.hpp:691
double * GetData() const
Return a pointer to the beginning of the Vector data.
Definition: vector.hpp:206
Arbitrary order H(div)-conforming Raviart-Thomas finite elements.
Definition: fe_coll.hpp:380
int GetVDim() const
Returns vector dimension.
Definition: fespace.hpp:702
A general vector function coefficient.
DofTransformation * GetElementVDofs(int i, Array< int > &vdofs) const
Returns indices of degrees of freedom for the i&#39;th element. The returned indices are offsets into an ...
Definition: fespace.cpp:281
Geometry::Type GetElementGeometry(int i) const
Definition: mesh.hpp:1228
Mesh * GetMesh() const
Returns the mesh.
Definition: fespace.hpp:555
Class FiniteElementSpace - responsible for providing FEM view of the mesh, mainly managing the set of...
Definition: fespace.hpp:219
virtual const char * Name() const
Definition: fe_coll.hpp:280
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
virtual void Update()
Transform by the Space UpdateMatrix (e.g., on Mesh change).
Definition: gridfunc.cpp:164
int GetNE() const
Returns number of elements.
Definition: mesh.hpp:1086
Class for integration point with weight.
Definition: intrules.hpp:31
int GetElementOrder(int i) const
Returns the order of the i&#39;th finite element.
Definition: fespace.cpp:178
A standard isoparametric element transformation.
Definition: eltrans.hpp:361
virtual void Mult(const Vector &x, Vector &y) const override
Interpolation or prolongation of a vector x corresponding to the coarse space to the vector y corresp...
Definition: transfer.cpp:1216
Ordering::Type GetOrdering() const
Return the ordering method.
Definition: fespace.hpp:721
int dim
Definition: ex24.cpp:53
virtual const Vector & GetDist() const
Definition: gslib.hpp:217
Arbitrary order L2 elements in 2D on a square.
Definition: fe_l2.hpp:44
void F_exact(const Vector &p, Vector &F)
Definition: findpts.cpp:195
virtual void ProjectCoefficient(Coefficient &coeff)
Project coeff Coefficient to this GridFunction. The projection computation depends on the choice of t...
Definition: gridfunc.cpp:2415
int Size() const
Return the logical size of the array.
Definition: array.hpp:141
Arbitrary order H(curl)-conforming Nedelec finite elements.
Definition: fe_coll.hpp:454
Nodes: x_i = i/(n-1), i=0,...,n-1.
Definition: fe_base.hpp:35
FindPointsGSLIB can robustly evaluate a GridFunction on an arbitrary collection of points...
Definition: gslib.hpp:52
Vector data type.
Definition: vector.hpp:58
void VisualizeField(socketstream &sock, const char *vishost, int visport, GridFunction &gf, const char *title, int x, int y, int w, int h, const char *keys, bool vec)
Definition: fem_extras.cpp:92
Arbitrary order H1-conforming (continuous) finite elements.
Definition: fe_coll.hpp:259
void GetNodes(Vector &node_coord) const
Definition: mesh.cpp:8302
void SetNodalGridFunction(GridFunction *nodes, bool make_owner=false)
Definition: mesh.cpp:5624
virtual const Array< unsigned int > & GetCode() const
Definition: gslib.hpp:208
const IntegrationRule & GetNodes() const
Get a const reference to the nodes of the element.
Definition: fe_base.hpp:389
int main(int argc, char *argv[])
Definition: findpts.cpp:201
double func_order
Definition: findpts.cpp:184
Arbitrary order "L2-conforming" discontinuous finite elements.
Definition: fe_coll.hpp:327