MFEM  v4.5.2
Finite element discretization library
qspace.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 #include "qspace.hpp"
13 
14 namespace mfem
15 {
16 
18  const IntegrationRule &ir)
19  : mesh(mesh_), order(ir.GetOrder())
20 {
21  for (int g = 0; g < Geometry::NumGeom; g++)
22  {
23  int_rule[g] = NULL;
24  }
25  int_rule[geom] = &ir;
26 }
27 
29 {
31  mesh.GetGeometries(dim, geoms);
32  for (Geometry::Type geom : geoms)
33  {
34  int_rule[geom] = &IntRules.Get(geom, order);
35  }
36 }
37 
39 {
40  const int num_elem = mesh.GetNE();
41  offsets.SetSize(num_elem + 1);
42  int offset = 0;
43  for (int i = 0; i < num_elem; i++)
44  {
45  offsets[i] = offset;
46  int geom = mesh.GetElementBaseGeometry(i);
47  MFEM_ASSERT(int_rule[geom] != NULL, "Missing integration rule.");
48  offset += int_rule[geom]->GetNPoints();
49  }
50  offsets[num_elem] = size = offset;
51 }
52 
54 {
57 }
58 
59 QuadratureSpace::QuadratureSpace(Mesh *mesh_, std::istream &in)
60  : QuadratureSpaceBase(*mesh_)
61 {
62  const char *msg = "invalid input stream";
63  std::string ident;
64 
65  in >> ident; MFEM_VERIFY(ident == "QuadratureSpace", msg);
66  in >> ident; MFEM_VERIFY(ident == "Type:", msg);
67  in >> ident;
68  if (ident == "default_quadrature")
69  {
70  in >> ident; MFEM_VERIFY(ident == "Order:", msg);
71  in >> order;
72  }
73  else
74  {
75  MFEM_ABORT("unknown QuadratureSpace type: " << ident);
76  return;
77  }
78 
79  Construct();
80 }
81 
83  : QuadratureSpaceBase(mesh_, mesh_.GetElementGeometry(0), ir)
84 {
85  MFEM_VERIFY(mesh.GetNumGeometries(mesh.Dimension()) == 1,
86  "Constructor not valid for mixed meshes");
88 }
89 
90 void QuadratureSpace::Save(std::ostream &os) const
91 {
92  os << "QuadratureSpace\n"
93  << "Type: default_quadrature\n"
94  << "Order: " << order << '\n';
95 }
96 
98  FaceType face_type_)
99  : QuadratureSpaceBase(mesh_, order_),
100  face_type(face_type_),
101  num_faces(mesh.GetNFbyType(face_type))
102 {
103  Construct();
104 }
105 
107  FaceType face_type_)
108  : QuadratureSpaceBase(mesh_, mesh_.GetFaceGeometry(0), ir),
109  face_type(face_type_),
110  num_faces(mesh.GetNFbyType(face_type))
111 {
112  MFEM_VERIFY(mesh.GetNumGeometries(mesh.Dimension() - 1) == 1,
113  "Constructor not valid for mixed meshes");
114  ConstructOffsets();
115 }
116 
117 void FaceQuadratureSpace::ConstructOffsets()
118 {
119  face_indices.SetSize(num_faces);
120  offsets.SetSize(num_faces + 1);
121  int offset = 0;
122  int f_idx = 0;
123  for (int i = 0; i < mesh.GetNumFacesWithGhost(); i++)
124  {
126  if (face.IsNonconformingCoarse() || !face.IsOfFaceType(face_type))
127  {
128  continue;
129  }
130  face_indices[f_idx] = i;
131  offsets[f_idx] = offset;
133  MFEM_ASSERT(int_rule[geom] != NULL, "Missing integration rule");
134  offset += int_rule[geom]->GetNPoints();
135 
136  f_idx++;
137  }
138  offsets[num_faces] = size = offset;
139 }
140 
141 void FaceQuadratureSpace::Construct()
142 {
144  ConstructOffsets();
145 }
146 
147 int FaceQuadratureSpace::GetPermutedIndex(int idx, int iq) const
148 {
149  const int f_idx = face_indices[idx];
151  {
152  const int dim = mesh.Dimension();
153  const IntegrationRule &ir = GetIntRule(idx);
154  const int q1d = (int)floor(pow(ir.GetNPoints(), 1.0/(dim-1)) + 0.5);
155  const Mesh::FaceInformation face = mesh.GetFaceInformation(f_idx);
156  return ToLexOrdering(dim, face.element[0].local_face_id, q1d, iq);
157  }
158  else
159  {
160  return iq;
161  }
162 }
163 
164 void FaceQuadratureSpace::Save(std::ostream &os) const
165 {
166  os << "FaceQuadratureSpace\n"
167  << "Type: default_quadrature\n"
168  << "Order: " << order << '\n';
169 }
170 
171 } // namespace mfem
Geometry::Type GetGeometry(int idx) const override
Returns the geometry type of face idx.
Definition: qspace.hpp:169
int GetNPoints() const
Returns the number of the points in the integration rule.
Definition: intrules.hpp:247
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:923
Geometry::Type GetElementBaseGeometry(int i) const
Definition: mesh.hpp:1108
static const int NumGeom
Definition: geom.hpp:42
int Dimension() const
Definition: mesh.hpp:1047
Geometry::Type GetFaceGeometry(int i) const
Return the Geometry::Type associated with face i.
Definition: mesh.cpp:1413
QuadratureSpaceBase(Mesh &mesh_, int order_=0)
Protected constructor. Used by derived classes.
Definition: qspace.hpp:42
int size
Total number of quadrature points.
Definition: qspace.hpp:31
int GetPermutedIndex(int idx, int iq) const override
Returns the permuted index of the iq quadrature point in entity idx.
Definition: qspace.cpp:147
int GetNumGeometries(int dim) const
Return the number of geometries of the given dimension present in the mesh.
Definition: mesh.cpp:5921
IntegrationRules IntRules(0, Quadrature1D::GaussLegendre)
A global object with all integration rules (defined in intrules.cpp)
Definition: intrules.hpp:379
FaceType
Definition: mesh.hpp:45
void Save(std::ostream &out) const override
Write the QuadratureSpace to the stream out.
Definition: qspace.cpp:90
bool IsOfFaceType(FaceType type) const
Return true if the face is of the same type as type.
Definition: mesh.hpp:1465
struct mfem::Mesh::FaceInformation::@13 element[2]
Array< int > offsets
Entity quadrature point offset array, of size num_entities + 1.
Definition: qspace.hpp:37
Abstract base class for QuadratureSpace and FaceQuadratureSpace.
Definition: qspace.hpp:24
This structure is used as a human readable output format that decipheres the information contained in...
Definition: mesh.hpp:1417
int ToLexOrdering(const int dim, const int face_id, const int size1d, const int index)
Convert a dof face index from Native ordering to lexicographic ordering for quads and hexes...
FaceQuadratureSpace(Mesh &mesh_, int order_, FaceType face_type_)
Create a FaceQuadratureSpace based on the global rules from IntRules.
Definition: qspace.cpp:97
void SetSize(int nsize)
Change the logical size of the array, keep existing entries.
Definition: array.hpp:684
const IntegrationRule & GetIntRule(int idx) const
Return the IntegrationRule associated with entity idx.
Definition: qspace.hpp:72
void GetGeometries(int dim, Array< Geometry::Type > &el_geoms) const
Return all element geometries of the given dimension present in the mesh.
Definition: mesh.cpp:5932
int GetNE() const
Returns number of elements.
Definition: mesh.hpp:936
static bool IsTensorProduct(Type geom)
Definition: geom.hpp:108
Mesh & mesh
The underlying mesh.
Definition: qspace.hpp:29
const IntegrationRule * int_rule[Geometry::NumGeom]
The quadrature rules used for each geometry type.
Definition: qspace.hpp:39
int dim
Definition: ex24.cpp:53
void Save(std::ostream &out) const override
Write the FaceQuadratureSpace to the stream out.
Definition: qspace.cpp:164
QuadratureSpace(Mesh *mesh_, int order_)
Create a QuadratureSpace based on the global rules from IntRules.
Definition: qspace.hpp:99
FaceInformation GetFaceInformation(int f) const
Definition: mesh.cpp:1131
void ConstructOffsets()
Definition: qspace.cpp:38
bool IsNonconformingCoarse() const
Return true if the face is a nonconforming coarse face.
Definition: mesh.hpp:1496
int GetNumFacesWithGhost() const
Return the number of faces (3D), edges (2D) or vertices (1D) including ghost faces.
Definition: mesh.cpp:5400
void ConstructIntRules(int dim)
Fill the int_rule array for each geometry type using order.
Definition: qspace.cpp:28
int order
The order of integration rule.
Definition: qspace.hpp:30