Nektar++
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
Classes | Typedefs | Functions | Variables
Nektar::GlobalMapping Namespace Reference

Classes

class  Mapping
 Base class for mapping to be applied to the coordinate system. More...
 
class  MappingGeneral
 
class  MappingTranslation
 
class  MappingXofXZ
 
class  MappingXofZ
 
class  MappingXYofXY
 
class  MappingXYofZ
 

Typedefs

typedef LibUtilities::NekFactory< std::string, Mapping, const LibUtilities::SessionReaderSharedPtr &, const Array< OneD, MultiRegions::ExpListSharedPtr > &, const TiXmlElement * > MappingFactory
 Declaration of the mapping factory. More...
 

Functions

void UpdateGeometry (SpatialDomains::MeshGraphSharedPtr graph, Array< OneD, MultiRegions::ExpListSharedPtr > &fields, bool modal)
 Update geometry according to displacement that is in current fields. More...
 
MappingFactoryGetMappingFactory ()
 Declaration of the mapping factory singleton. More...
 

Variables

GLOBAL_MAPPING_EXPORT typedef std::shared_ptr< MappingMappingSharedPtr
 A shared pointer to a Mapping object. More...
 

Typedef Documentation

◆ MappingFactory

Declaration of the mapping factory.

Definition at line 59 of file Mapping.h.

Function Documentation

◆ GetMappingFactory()

GLOBAL_MAPPING_EXPORT MappingFactory & Nektar::GlobalMapping::GetMappingFactory ( )

Declaration of the mapping factory singleton.

Definition at line 52 of file Mapping.cpp.

53 {
54  static MappingFactory instance;
55  return instance;
56 }
LibUtilities::NekFactory< std::string, Mapping, const LibUtilities::SessionReaderSharedPtr &, const Array< OneD, MultiRegions::ExpListSharedPtr > &, const TiXmlElement * > MappingFactory
Declaration of the mapping factory.
Definition: Mapping.h:59

Referenced by Nektar::GlobalMapping::Mapping::Load().

◆ UpdateGeometry()

GLOBAL_MAPPING_EXPORT void Nektar::GlobalMapping::UpdateGeometry ( SpatialDomains::MeshGraphSharedPtr  graph,
Array< OneD, MultiRegions::ExpListSharedPtr > &  fields,
bool  modal 
)

Update geometry according to displacement that is in current fields.

Adds a summary item to the summary info list.

Parameters
graphThe MeshGraph of the current geometry.
fieldsThe fields containing the displacement.

Definition at line 59 of file Deform.cpp.

63  {
64  // Clear existing curvature.
65  SpatialDomains::CurveMap &curvedEdges = graph->GetCurvedEdges();
66  SpatialDomains::CurveMap &curvedFaces = graph->GetCurvedFaces();
67  curvedEdges.clear();
68  curvedFaces.clear();
69 
70  int i, j, k, l, dim;
71 
72  // Sets to hold IDs of updated vertices to avoid duplicating effort.
73  set<int> updatedVerts, updatedEdges, updatedFaces;
74 
75  dim = graph->GetSpaceDimension();
76  Array<OneD, Array<OneD, NekDouble> > phys (dim);
77  Array<OneD, Array<OneD, NekDouble> > coord(dim);
78 
79  for (i = 0; i < fields[0]->GetExpSize(); ++i)
80  {
81  LocalRegions::ExpansionSharedPtr exp = fields[0]->GetExp(i);
82  int offset = fields[0]->GetPhys_Offset(i);
83  int nquad = exp->GetTotPoints();
84 
85  // Extract displacement for this element, allocate storage for
86  // elemental coordinates.
87  for (j = 0; j < dim; ++j)
88  {
89  phys[j] = Array<OneD, NekDouble>(
90  nquad, fields[j]->UpdatePhys() + offset);
91  coord[j] = Array<OneD, NekDouble>(nquad);
92  }
93 
94  // In 2D loop over edges.
95  if (dim == 2)
96  {
97  exp->GetCoords(coord[0], coord[1]);
98 
100  std::dynamic_pointer_cast<SpatialDomains::Geometry2D>(
101  exp->GetGeom());
102 
103  for (j = 0; j < exp->GetGeom()->GetNumEdges(); ++j)
104  {
105  SpatialDomains::Geometry1DSharedPtr edge = geom->GetEdge(j);
106 
107  // This edge has already been processed.
108  if (updatedEdges.find(edge->GetGlobalID()) !=
109  updatedEdges.end())
110  {
111  continue;
112  }
113 
114  // Extract edge displacement.
115  int nEdgePts = exp->GetTraceNumPoints(j);
116  Array<OneD, Array<OneD, NekDouble> > edgePhys (dim);
117  Array<OneD, Array<OneD, NekDouble> > edgeCoord(dim);
118 
119  const LibUtilities::BasisKey B(
120  LibUtilities::eModified_A, nEdgePts,
121  LibUtilities::PointsKey(
123  StdRegions::StdExpansion1DSharedPtr seg = MemoryManager<
124  StdRegions::StdSegExp>::AllocateSharedPtr(B);
125 
126  for (k = 0; k < dim; ++k)
127  {
128  edgePhys [k] = Array<OneD, NekDouble>(nEdgePts);
129  edgeCoord[k] = Array<OneD, NekDouble>(nEdgePts);
130  exp->GetTracePhysVals(j, seg, phys [k], edgePhys [k]);
131  exp->GetTracePhysVals(j, seg, coord[k], edgeCoord[k]);
132  }
133 
134  // Update verts
135  for (k = 0; k < 2; ++k)
136  {
137  int id = edge->GetVid(k);
138  if (updatedVerts.find(id) != updatedVerts.end())
139  {
140  continue;
141  }
142 
144  edge->GetVertex(k);
145 
146  pt->UpdatePosition(
147  (*pt)(0) + edgePhys[0][k*(nEdgePts-1)],
148  (*pt)(1) + edgePhys[1][k*(nEdgePts-1)],
149  (*pt)(2));
150 
151  updatedVerts.insert(id);
152  }
153 
154  // Update curve
155  SpatialDomains::CurveSharedPtr curve = MemoryManager<
156  SpatialDomains::Curve>::AllocateSharedPtr(
157  edge->GetGlobalID(),
159 
160  for (k = 0; k < nEdgePts; ++k)
161  {
163  MemoryManager<SpatialDomains::PointGeom>
164  ::AllocateSharedPtr(
165  dim, edge->GetGlobalID(),
166  edgeCoord[0][k] + edgePhys[0][k],
167  edgeCoord[1][k] + edgePhys[1][k], 0.0);
168 
169  curve->m_points.push_back(vert);
170  }
171 
172  curvedEdges[edge->GetGlobalID()] = curve;
173  updatedEdges.insert(edge->GetGlobalID());
174  }
175  }
176  else if (dim == 3)
177  {
178  exp->GetCoords(coord[0], coord[1], coord[2]);
179 
181  std::dynamic_pointer_cast<SpatialDomains::Geometry3D>(
182  exp->GetGeom());
183 
184  for (j = 0; j < exp->GetNtraces(); ++j)
185  {
186  SpatialDomains::Geometry2DSharedPtr face = geom->GetFace(j);
187 
189  exp->as<LocalRegions::Expansion3D>();
190 
191  // This edge has already been processed.
192  if (updatedFaces.find(face->GetGlobalID()) !=
193  updatedFaces.end())
194  {
195  continue;
196  }
197 
198  // Extract face displacement.
199  LibUtilities::BasisKey B0 = exp->GetTraceBasisKey(j,0);
200  LibUtilities::BasisKey B1 = exp->GetTraceBasisKey(j,1);
201  int nq0 = B0.GetNumPoints();
202  int nq1 = B1.GetNumPoints();
203 
204  ASSERTL1(B0.GetPointsType()
206  B1.GetPointsType()
208  "Deformation requires GLL points in both "
209  "directions on a face.");
210 
211  Array<OneD, Array<OneD, NekDouble> > newPos(dim);
212 
214  StdRegions::Orientation orient = exp->GetTraceOrient(j);
215 
216  if (face->GetShapeType() == LibUtilities::eTriangle)
217  {
218  faceexp = MemoryManager<StdRegions::StdTriExp>::
219  AllocateSharedPtr(B0, B1);
220  }
221  else
222  {
223  faceexp = MemoryManager<StdRegions::StdQuadExp>::
224  AllocateSharedPtr(B0, B1);
225  }
226 
227  for (k = 0; k < dim; ++k)
228  {
229  Array<OneD, NekDouble> tmp(nq0*nq1);
230  newPos[k] = Array<OneD, NekDouble>(nq0*nq1);
231  exp3d->GetTracePhysVals(
232  j, faceexp, phys [k], tmp, orient);
233  exp3d->GetTracePhysVals(
234  j, faceexp, coord[k], newPos[k], orient);
235  Vmath::Vadd(
236  nq0*nq1, tmp, 1, newPos[k], 1, newPos[k], 1);
237  }
238 
239  // Now interpolate face onto a more reasonable set of
240  // points.
241  int nq = max(nq0, nq1);
242  if(!modal)
243  nq--;
244 
245  LibUtilities::PointsKey edgePts(
247  LibUtilities::PointsKey triPts(
249 
250  Array<OneD, Array<OneD, NekDouble> > intPos(dim);
251 
252  for (k = 0; k < dim; ++k)
253  {
254  intPos[k] = Array<OneD, NekDouble>(nq*nq);
256  faceexp->GetPointsKeys()[0],
257  faceexp->GetPointsKeys()[1],
258  newPos[k], edgePts, edgePts, intPos[k]);
259  }
260 
261  int edgeOff[2][4][2] = {
262  {
263  {0, 1},
264  {nq-1, nq},
265  {nq*(nq-1), -nq},
266  {-1,-1}
267  },
268  {
269  {0, 1},
270  {nq-1, nq},
271  {nq*nq-1, -1},
272  {nq*(nq-1), -nq}
273  }
274  };
275 
276  for (k = 0; k < face->GetNumVerts(); ++k)
277  {
278  // Update verts
279  int id = face->GetVid(k);
280  const int o =
281  face->GetShapeType() - LibUtilities::eTriangle;
282 
283  if (updatedVerts.find(id) == updatedVerts.end())
284  {
286  face->GetVertex(k);
287  pt->UpdatePosition(
288  intPos[0][edgeOff[o][k][0]],
289  intPos[1][edgeOff[o][k][0]],
290  intPos[2][edgeOff[o][k][0]]);
291  updatedVerts.insert(id);
292  }
293 
294  // Update edges
295  id = face->GetEid(k);
296  if (updatedEdges.find(id) == updatedEdges.end())
297  {
299  = face->GetEdge(k);
301  MemoryManager<SpatialDomains::Curve>
302  ::AllocateSharedPtr(
303  edge->GetGlobalID(),
305 
306  const int offset = edgeOff[o][k][0];
307  const int pos = edgeOff[o][k][1];
308 
309  if (face->GetEorient(k) == StdRegions::eBackwards)
310  {
311  for (l = nq-1; l >= 0; --l)
312  {
313  int m = offset + pos*l;
315  MemoryManager<SpatialDomains::PointGeom>
316  ::AllocateSharedPtr(
317  dim, edge->GetGlobalID(),
318  intPos[0][m], intPos[1][m],
319  intPos[2][m]);
320  curve->m_points.push_back(vert);
321  }
322  }
323  else
324  {
325  for (l = 0; l < nq; ++l)
326  {
327  int m = offset + pos*l;
329  MemoryManager<SpatialDomains::PointGeom>
330  ::AllocateSharedPtr(
331  dim, edge->GetGlobalID(),
332  intPos[0][m], intPos[1][m],
333  intPos[2][m]);
334  curve->m_points.push_back(vert);
335  }
336  }
337 
338  curvedEdges[edge->GetGlobalID()] = curve;
339  updatedEdges.insert(edge->GetGlobalID());
340  }
341  }
342 
343  // Update face-interior curvature
345  face->GetShapeType() == LibUtilities::eTriangle ?
348 
349  SpatialDomains::CurveSharedPtr curve = MemoryManager<
350  SpatialDomains::Curve>::AllocateSharedPtr(
351  face->GetGlobalID(),
352  pType);
353 
354  if (face->GetShapeType() == LibUtilities::eTriangle)
355  {
356  // This code is probably pretty crappy. Have to go from
357  // GLL-GLL points -> GLL-Gauss-Radau -> nodal triangle
358  // points.
359  const LibUtilities::BasisKey B0(
361  LibUtilities::PointsKey(
363  const LibUtilities::BasisKey B1(
365  LibUtilities::PointsKey(
367  StdRegions::StdNodalTriExp nodalTri(B0, B1, pType);
368  StdRegions::StdTriExp tri (B0, B1);
369 
370  for (k = 0; k < dim; ++k)
371  {
372  Array<OneD, NekDouble> nodal(nq*nq);
373 
375  faceexp->GetBasis(0)->GetBasisKey(),
376  faceexp->GetBasis(1)->GetBasisKey(),
377  newPos[k], B0, B1, nodal);
378 
379  Array<OneD, NekDouble> tmp1(nq*(nq+1)/2);
380  Array<OneD, NekDouble> tmp2(nq*(nq+1)/2);
381 
382  tri.FwdTrans(nodal, tmp1);
383  nodalTri.ModalToNodal(tmp1, tmp2);
384  newPos[k] = tmp2;
385  }
386 
387  for (l = 0; l < nq*(nq+1)/2; ++l)
388  {
390  MemoryManager<SpatialDomains::PointGeom>
391  ::AllocateSharedPtr(
392  dim, face->GetGlobalID(),
393  newPos[0][l], newPos[1][l], newPos[2][l]);
394  curve->m_points.push_back(vert);
395  }
396  }
397  else
398  {
399  for (l = 0; l < nq*nq; ++l)
400  {
402  MemoryManager<SpatialDomains::PointGeom>
403  ::AllocateSharedPtr(
404  dim, face->GetGlobalID(),
405  intPos[0][l], intPos[1][l], intPos[2][l]);
406  curve->m_points.push_back(vert);
407  }
408  }
409 
410  curvedFaces[face->GetGlobalID()] = curve;
411  updatedFaces.insert(face->GetGlobalID());
412  }
413  }
414  }
415 
416  // Reset geometry information
417  for (i = 0; i < fields.size(); ++i)
418  {
419  fields[i]->Reset();
420  }
421  }
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode....
Definition: ErrorUtil.hpp:250
void Interp2D(const BasisKey &fbasis0, const BasisKey &fbasis1, const Array< OneD, const NekDouble > &from, const BasisKey &tbasis0, const BasisKey &tbasis1, Array< OneD, NekDouble > &to)
this function interpolates a 2D function evaluated at the quadrature points of the 2D basis,...
Definition: Interp.cpp:115
@ eNodalTriElec
2D Nodal Electrostatic Points on a Triangle
Definition: PointsType.h:69
@ eGaussLobattoLegendre
1D Gauss-Lobatto-Legendre quadrature points
Definition: PointsType.h:51
@ eGaussRadauMAlpha1Beta0
Gauss Radau pinned at x=-1, .
Definition: PointsType.h:58
@ eOrtho_A
Principle Orthogonal Functions .
Definition: BasisType.h:45
@ eOrtho_B
Principle Orthogonal Functions .
Definition: BasisType.h:46
@ eModified_A
Principle Modified Functions .
Definition: BasisType.h:48
std::shared_ptr< Expansion > ExpansionSharedPtr
Definition: Expansion.h:68
std::shared_ptr< Expansion3D > Expansion3DSharedPtr
Definition: Expansion2D.h:49
std::shared_ptr< Curve > CurveSharedPtr
Definition: Curve.hpp:61
std::unordered_map< int, CurveSharedPtr > CurveMap
Definition: Curve.hpp:62
std::shared_ptr< PointGeom > PointGeomSharedPtr
Definition: Geometry.h:59
std::shared_ptr< Geometry2D > Geometry2DSharedPtr
Definition: Geometry.h:65
std::shared_ptr< Geometry1D > Geometry1DSharedPtr
Definition: Geometry.h:63
std::shared_ptr< Geometry3D > Geometry3DSharedPtr
Definition: Geometry3D.h:52
std::shared_ptr< StdExpansion2D > StdExpansion2DSharedPtr
std::shared_ptr< StdExpansion1D > StdExpansion1DSharedPtr
void Vadd(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Add vector z = x+y.
Definition: Vmath.cpp:322

References ASSERTL1, Nektar::StdRegions::eBackwards, Nektar::LibUtilities::eGaussLobattoLegendre, Nektar::LibUtilities::eGaussRadauMAlpha1Beta0, Nektar::LibUtilities::eModified_A, Nektar::LibUtilities::eNodalTriElec, Nektar::LibUtilities::eOrtho_A, Nektar::LibUtilities::eOrtho_B, Nektar::LibUtilities::eTriangle, Nektar::StdRegions::StdExpansion::FwdTrans(), Nektar::LibUtilities::BasisKey::GetNumPoints(), Nektar::LibUtilities::BasisKey::GetPointsType(), Nektar::LibUtilities::Interp2D(), Nektar::StdRegions::StdNodalTriExp::ModalToNodal(), and Vmath::Vadd().

Referenced by Nektar::FieldUtils::ProcessDeform::Process(), and Nektar::IterativeElasticSystem::v_DoSolve().

Variable Documentation

◆ MappingSharedPtr

GLOBAL_MAPPING_EXPORT typedef std::shared_ptr<Mapping> Nektar::GlobalMapping::MappingSharedPtr