Nektar++
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
QuadGeom.cpp
Go to the documentation of this file.
1 ////////////////////////////////////////////////////////////////////////////////
2 //
3 // File: QuadGeom.cpp
4 //
5 // For more information, please see: http://www.nektar.info/
6 //
7 // The MIT License
8 //
9 // Copyright (c) 2006 Division of Applied Mathematics, Brown University (USA),
10 // Department of Aeronautics, Imperial College London (UK), and Scientific
11 // Computing and Imaging Institute, University of Utah (USA).
12 //
13 // License for the specific language governing rights and limitations under
14 // Permission is hereby granted, free of charge, to any person obtaining a
15 // copy of this software and associated documentation files (the "Software"),
16 // to deal in the Software without restriction, including without limitation
17 // the rights to use, copy, modify, merge, publish, distribute, sublicense,
18 // and/or sell copies of the Software, and to permit persons to whom the
19 // Software is furnished to do so, subject to the following conditions:
20 //
21 // The above copyright notice and this permission notice shall be included
22 // in all copies or substantial portions of the Software.
23 //
24 // THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
25 // OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
26 // FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
27 // THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
28 // LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
29 // FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
30 // DEALINGS IN THE SOFTWARE.
31 //
32 // Description:
33 //
34 //
35 ////////////////////////////////////////////////////////////////////////////////
36 
39 
40 #include <StdRegions/StdQuadExp.h>
41 #include <SpatialDomains/SegGeom.h>
42 #include <SpatialDomains/Curve.hpp>
44 
45 
46 namespace Nektar
47 {
48  namespace SpatialDomains
49  {
50  /**
51  *
52  */
54  {
56  }
57 
58 
59  /**
60  *
61  */
62  QuadGeom::QuadGeom(const int id,
63  const PointGeomSharedPtr verts[],
64  const SegGeomSharedPtr edges[],
65  const StdRegions::Orientation eorient[]):
66  Geometry2D(verts[0]->GetCoordim()),
67  m_fid(id)
68  {
69  m_globalID = id;
71 
72  /// Copy the vert shared pointers.
73  m_verts.insert(m_verts.begin(), verts, verts+QuadGeom::kNverts);
74 
75  /// Copy the edge shared pointers.
76  m_edges.insert(m_edges.begin(), edges, edges+QuadGeom::kNedges);
77 
78 
79  for (int j=0; j<kNedges; ++j)
80  {
81  m_eorient[j] = eorient[j];
82  }
83 
84  m_coordim = verts[0]->GetCoordim();
85  ASSERTL0(m_coordim > 1,
86  "Cannot call function with dim == 1");
87 
88  SetUpXmap();
89  SetUpCoeffs(m_xmap->GetNcoeffs());
90  }
91 
92 
93  /**
94  *
95  */
96  QuadGeom::QuadGeom(const int id,
97  const SegGeomSharedPtr edges[],
98  const StdRegions::Orientation eorient[],
99  const CurveSharedPtr &curve) :
100  Geometry2D(edges[0]->GetVertex(0)->GetCoordim()),
101  m_fid(id),
102  m_curve(curve)
103  {
104  int j;
105 
106  m_globalID = m_fid;
107 
109 
110  /// Copy the edge shared pointers.
111  m_edges.insert(m_edges.begin(), edges, edges+QuadGeom::kNedges);
112 
113  for(j=0; j <kNedges; ++j)
114  {
115  if(eorient[j] == StdRegions::eForwards)
116  {
117  m_verts.push_back(edges[j]->GetVertex(0));
118  }
119  else
120  {
121  m_verts.push_back(edges[j]->GetVertex(1));
122  }
123  }
124 
125  for (j=0; j<kNedges; ++j)
126  {
127  m_eorient[j] = eorient[j];
128  }
129 
130  m_coordim = edges[0]->GetVertex(0)->GetCoordim();
131  ASSERTL0(m_coordim > 1,
132  "Cannot call function with dim == 1");
133 
134  SetUpXmap();
135  SetUpCoeffs(m_xmap->GetNcoeffs());
136  }
137 
138 
139  /**
140  *
141  */
142  QuadGeom::QuadGeom(const int id,
143  const SegGeomSharedPtr edges[],
144  const StdRegions::Orientation eorient[]):
145  Geometry2D(edges[0]->GetVertex(0)->GetCoordim()),
146  m_fid(id)
147  {
148  int j;
149 
150  m_globalID = m_fid;
152 
153  /// Copy the edge shared pointers.
154  m_edges.insert(m_edges.begin(), edges, edges+QuadGeom::kNedges);
155 
156  for(j=0; j <kNedges; ++j)
157  {
158  if(eorient[j] == StdRegions::eForwards)
159  {
160  m_verts.push_back(edges[j]->GetVertex(0));
161  }
162  else
163  {
164  m_verts.push_back(edges[j]->GetVertex(1));
165  }
166  }
167 
168  for (j=0; j<kNedges; ++j)
169  {
170  m_eorient[j] = eorient[j];
171  }
172 
173  m_coordim = edges[0]->GetVertex(0)->GetCoordim();
174  ASSERTL0(m_coordim > 1,
175  "Cannot call function with dim == 1");
176 
177  SetUpXmap();
178  SetUpCoeffs(m_xmap->GetNcoeffs());
179  }
180 
181 
182  /**
183  *
184  */
186  {
187  // From Geometry
189 
190  // From QuadFaceComponent
191  m_fid = in.m_fid;
192  m_globalID = m_fid;
193  m_ownVerts = in.m_ownVerts;
194  std::list<CompToElmt>::const_iterator def;
195  for(def = in.m_elmtMap.begin(); def != in.m_elmtMap.end(); def++)
196  {
197  m_elmtMap.push_back(*def);
198  }
199 
200  // From QuadGeom
201  m_verts = in.m_verts;
202  m_edges = in.m_edges;
203  for (int i = 0; i < kNedges; i++)
204  {
205  m_eorient[i] = in.m_eorient[i];
206  }
207  m_ownData = in.m_ownData;
208  }
209 
210 
211  /**
212  *
213  */
215  {
216  }
217 
219  {
220  int order0 = max(m_edges[0]->GetBasis(0)->GetNumModes(),
221  m_edges[2]->GetBasis(0)->GetNumModes());
222  int order1 = max(m_edges[1]->GetBasis(0)->GetNumModes(),
223  m_edges[3]->GetBasis(0)->GetNumModes());
224 
225  const LibUtilities::BasisKey B0(
229  const LibUtilities::BasisKey B1(
233 
235  ::AllocateSharedPtr(B0,B1);
236  }
237 
238  /**
239  *
240  */
242  const Array<OneD, const NekDouble> &Lcoord)
243  {
245  "Geometry is not in physical space");
246 
247  Array<OneD, NekDouble> tmp(m_xmap->GetTotPoints());
248  m_xmap->BwdTrans(m_coeffs[i], tmp);
249 
250  return m_xmap->PhysEvaluate(Lcoord, tmp);
251  }
252 
254  const QuadGeom &face1,
255  const QuadGeom &face2)
256  {
257  return GetFaceOrientation(face1.m_verts, face2.m_verts);
258  }
259 
260  /**
261  * Calculate the orientation of face2 to face1 (note this is
262  * not face1 to face2!).
263  */
265  const PointGeomVector &face1,
266  const PointGeomVector &face2)
267  {
268  int i, j, vmap[4] = {-1,-1,-1,-1};
269  NekDouble x, y, z, x1, y1, z1, cx = 0.0, cy = 0.0, cz = 0.0;
270 
271  // For periodic faces, we calculate the vector between the centre
272  // points of the two faces. (For connected faces this will be
273  // zero). We can then use this to determine alignment later in the
274  // algorithm.
275  for (i = 0; i < 4; ++i)
276  {
277  cx += (*face2[i])(0) - (*face1[i])(0);
278  cy += (*face2[i])(1) - (*face1[i])(1);
279  cz += (*face2[i])(2) - (*face1[i])(2);
280  }
281  cx /= 4;
282  cy /= 4;
283  cz /= 4;
284 
285  // Now construct a mapping which takes us from the vertices of one
286  // face to the other. That is, vertex j of face2 corresponds to
287  // vertex vmap[j] of face1.
288  for (i = 0; i < 4; ++i)
289  {
290  x = (*face1[i])(0);
291  y = (*face1[i])(1);
292  z = (*face1[i])(2);
293  for (j = 0; j < 4; ++j)
294  {
295  x1 = (*face2[j])(0)-cx;
296  y1 = (*face2[j])(1)-cy;
297  z1 = (*face2[j])(2)-cz;
298  if (sqrt((x1-x)*(x1-x)+(y1-y)*(y1-y)+(z1-z)*(z1-z)) < 1e-8)
299  {
300  vmap[j] = i;
301  break;
302  }
303  }
304  }
305 
306  // Use the mapping to determine the eight alignment options between
307  // faces.
308  if (vmap[1] == (vmap[0]+1) % 4)
309  {
310  switch (vmap[0])
311  {
312  case 0:
314  break;
315  case 1:
317  break;
318  case 2:
320  break;
321  case 3:
323  break;
324  }
325  }
326  else
327  {
328  switch (vmap[0])
329  {
330  case 0:
332  break;
333  case 1:
335  break;
336  case 2:
338  break;
339  case 3:
341  break;
342  }
343  }
344 
345  ASSERTL0(false, "unable to determine face orientation");
347  }
348 
349 
350  /**
351  *
352  */
353  void QuadGeom::v_AddElmtConnected(int gvo_id, int locid)
354  {
355  CompToElmt ee(gvo_id,locid);
356  m_elmtMap.push_back(ee);
357  }
358 
359 
360  /**
361  *
362  */
364  {
365  return int(m_elmtMap.size());
366  }
367 
368 
369  /**
370  *
371  */
372  bool QuadGeom::v_IsElmtConnected(int gvo_id, int locid) const
373  {
374  std::list<CompToElmt>::const_iterator def;
375  CompToElmt ee(gvo_id,locid);
376 
377  def = find(m_elmtMap.begin(),m_elmtMap.end(),ee);
378 
379  // Found the element connectivity object in the list
380  if(def != m_elmtMap.end())
381  {
382  return(true);
383  }
384 
385  return(false);
386  }
387 
388 
389  /**
390  *
391  */
392  int QuadGeom::v_GetFid() const
393  {
394  return m_fid;
395  }
396 
397 
398  /**
399  *
400  */
402  {
403  return m_coordim;
404  }
405 
406 
407  /**
408  *
409  */
411  {
412  return m_xmap->GetBasis(i);
413  }
414 
415 
416  /**
417  *
418  */
420  {
421  ASSERTL1(i <= 3,"edge is out of range");
422 
423  if (i == 0 || i == 2)
424  {
425  return m_xmap->GetBasis(0);
426  }
427  else
428  {
429  return m_xmap->GetBasis(1);
430  }
431  }
432 
433  /**
434  *
435  */
437  {
438  return GetCoord(i,Lcoord);
439  }
440 
441 
442  /**
443  * Set up GeoFac for this geometry using Coord quadrature distribution
444  */
446  {
448  {
449  int i;
450  GeomType Gtype = eRegular;
451 
453 
454  // We will first check whether we have a regular or deformed
455  // geometry. We will define regular as those cases where the
456  // Jacobian and the metric terms of the derivative are constants
457  // (i.e. not coordinate dependent)
458 
459  // Check to see if expansions are linear
460  // If not linear => deformed geometry
461  for(i = 0; i < m_coordim; ++i)
462  {
463  if((m_xmap->GetBasisNumModes(0) != 2)||
464  (m_xmap->GetBasisNumModes(1) != 2))
465  {
466  Gtype = eDeformed;
467  }
468  }
469 
470  // For linear expansions, the mapping from standard to local
471  // element is given by the relation:
472  // x_i = 0.25 * [ ( x_i^A + x_i^B + x_i^C + x_i^D) +
473  // (-x_i^A + x_i^B + x_i^C - x_i^D)*xi_1 +
474  // (-x_i^A - x_i^B + x_i^C + x_i^D)*xi_2 +
475  // ( x_i^A - x_i^B + x_i^C - x_i^D)*xi_1*xi_2 ]
476  //
477  // The jacobian of the transformation and the metric terms
478  // dxi_i/dx_j, involve only terms of the form dx_i/dxi_j (both
479  // for coordim == 2 or 3). Inspecting the formula above, it can
480  // be appreciated that the derivatives dx_i/dxi_j will be
481  // constant, if the coefficient of the non-linear term is zero.
482  //
483  // That is why for regular geometry, we require
484  //
485  // x_i^A - x_i^B + x_i^C - x_i^D = 0
486  //
487  // or equivalently
488  //
489  // x_i^A - x_i^B = x_i^D - x_i^C
490  //
491  // This corresponds to quadrilaterals which are paralellograms.
492  if(Gtype == eRegular)
493  {
494  for(i = 0; i < m_coordim; i++)
495  {
496  if( fabs( (*m_verts[0])(i) - (*m_verts[1])(i) +
497  (*m_verts[2])(i) - (*m_verts[3])(i) )
499  {
500  Gtype = eDeformed;
501  break;
502  }
503  }
504  }
505 
507  Gtype, m_coordim, m_xmap, m_coeffs);
509  }
510  }
511 
512 
513  /**
514  *
515  */
517  {
518  m_ownData = true;
519  }
520 
521 
522  /**
523  * Note verts and edges are listed according to anticlockwise
524  * convention but points in _coeffs have to be in array format from
525  * left to right.
526  */
528  {
529  // check to see if geometry structure is already filled
530  if(m_state != ePtsFilled)
531  {
532  int i,j,k;
533  int nEdgeCoeffs;
534 
535  if (m_curve)
536  {
537  int npts = m_curve->m_points.size();
538  int nEdgePts = (int)sqrt(static_cast<NekDouble>(npts));
539  Array<OneD, NekDouble> tmp(npts);
540  Array<OneD, NekDouble> tmp2(m_xmap->GetTotPoints());
541  LibUtilities::PointsKey curveKey(
542  nEdgePts, m_curve->m_ptype);
543 
544  // Sanity checks:
545  // - Curved faces should have square number of points;
546  // - Each edge should have sqrt(npts) points.
547  ASSERTL0(nEdgePts * nEdgePts == npts,
548  "NUMPOINTS should be a square number in"
549  " quadrilteral "
550  + boost::lexical_cast<string>(m_globalID));
551 
552  for (i = 0; i < kNedges; ++i)
553  {
554  ASSERTL0(
555  m_edges[i]->GetXmap()->GetNcoeffs() == nEdgePts,
556  "Number of edge points does not correspond to "
557  "number of face points in quadrilateral "
558  + boost::lexical_cast<string>(m_globalID));
559  }
560 
561  for (i = 0; i < m_coordim; ++i)
562  {
563  for (j = 0; j < npts; ++j)
564  {
565  tmp[j] = (m_curve->m_points[j]->GetPtr())[i];
566  }
567 
568  // Interpolate m_curve points to GLL points
570  curveKey, curveKey, tmp,
571  m_xmap->GetBasis(0)->GetPointsKey(),
572  m_xmap->GetBasis(1)->GetPointsKey(),
573  tmp2);
574 
575  // Forwards transform to get coefficient space.
576  m_xmap->FwdTrans(tmp2, m_coeffs[i]);
577  }
578  }
579 
580  // Now fill in edges.
581  Array<OneD, unsigned int> mapArray;
582  Array<OneD, int> signArray;
583 
584  for(i = 0; i < kNedges; i++)
585  {
586  m_edges[i]->FillGeom();
587  m_xmap->GetEdgeToElementMap(i,m_eorient[i],
588  mapArray,signArray);
589 
590  nEdgeCoeffs = m_edges[i]->GetXmap()->GetNcoeffs();
591 
592  for(j = 0; j < m_coordim; j++)
593  {
594  for(k = 0; k < nEdgeCoeffs; k++)
595  {
596  m_coeffs[j][mapArray[k]]
597  = signArray[k]*(m_edges[i]->GetCoeffs(j))[k];
598  }
599  }
600  }
601 
603  }
604  }
605 
606  /**
607  *
608  */
610  Array<OneD,NekDouble> &Lcoords)
611  {
612  NekDouble resid = 0.0;
613  if(GetMetricInfo()->GetGtype() == eRegular)
614  {
615  NekDouble coords2 = (m_coordim == 3)? coords[2]: 0.0;
616  PointGeom dv1, dv2, norm, orth1, orth2;
617  PointGeom xin(m_coordim,0,coords[0],coords[1],coords2);
618 
619  // Calculate edge vectors from 0-1 and 0-3 edges.
620  dv1.Sub(*m_verts[1],*m_verts[0]);
621  dv2.Sub(*m_verts[3],*m_verts[0]);
622 
623  // Obtain normal to plane in which dv1 and dv2 lie
624  norm.Mult(dv1,dv2);
625 
626  // Obtain vector which are normal to dv1 and dv2.
627  orth1.Mult(norm,dv1);
628  orth2.Mult(norm,dv2);
629 
630  // Start with vector of desired points minus vertex_0
631  xin -= *m_verts[0];
632 
633  // Calculate length using L/|dv1| = (x-v0).n1/(dv1.n1) for coordiante 1
634  // Then rescale to [-1,1].
635  Lcoords[0] = xin.dot(orth2)/dv1.dot(orth2);
636  Lcoords[0] = 2*Lcoords[0]-1;
637  Lcoords[1] = xin.dot(orth1)/dv2.dot(orth1);
638  Lcoords[1] = 2*Lcoords[1]-1;
639  }
640  else
641  {
643 
644  // Determine nearest point of coords to values in m_xmap
645  int npts = m_xmap->GetTotPoints();
646  Array<OneD, NekDouble> ptsx(npts), ptsy(npts);
647  Array<OneD, NekDouble> tmpx(npts), tmpy(npts);
648 
649  m_xmap->BwdTrans(m_coeffs[0], ptsx);
650  m_xmap->BwdTrans(m_coeffs[1], ptsy);
651 
652  const Array<OneD, const NekDouble> za = m_xmap->GetPoints(0);
653  const Array<OneD, const NekDouble> zb = m_xmap->GetPoints(1);
654 
655  //guess the first local coords based on nearest point
656  Vmath::Sadd(npts, -coords[0], ptsx,1,tmpx,1);
657  Vmath::Sadd(npts, -coords[1], ptsy,1,tmpy,1);
658  Vmath::Vmul (npts, tmpx,1,tmpx,1,tmpx,1);
659  Vmath::Vvtvp(npts, tmpy,1,tmpy,1,tmpx,1,tmpx,1);
660 
661  int min_i = Vmath::Imin(npts,tmpx,1);
662 
663  Lcoords[0] = za[min_i%za.num_elements()];
664  Lcoords[1] = zb[min_i/za.num_elements()];
665 
666  // Perform newton iteration to find local coordinates
667  NewtonIterationForLocCoord(coords, ptsx, ptsy, Lcoords,resid);
668  }
669  return resid;
670  }
671 
672 
673  /**
674  *
675  */
676  int QuadGeom::v_GetEid(int i) const
677  {
678  ASSERTL2((i >=0) && (i <= 3),"Edge id must be between 0 and 3");
679  return m_edges[i]->GetEid();
680  }
681 
682 
683  /**
684  *
685  */
686  int QuadGeom::v_GetVid(int i) const
687  {
688  ASSERTL2((i >=0) && (i <= 3),"Verted id must be between 0 and 3");
689  return m_verts[i]->GetVid();
690  }
691 
692 
693  /**
694  *
695  */
697  {
698  ASSERTL2((i >=0) && (i <= 3),"Vertex id must be between 0 and 3");
699  return m_verts[i];
700  }
701 
702 
703  /**
704  *
705  */
707  {
708  ASSERTL2((i >=0) && (i <= 3),"Edge id must be between 0 and 3");
709  return m_edges[i];
710  }
711 
712 
713  /**
714  *
715  */
717  {
718  ASSERTL2((i >=0) && (i <= 3),"Edge id must be between 0 and 3");
719  return m_eorient[i];
720  }
721 
722 
723  /**
724  *
725  */
727  {
728  ASSERTL2((i >=0) && (i <= 3),"Edge id must be between 0 and 3");
729  if(i < 2)
730  {
731  return m_eorient[i];
732  }
733  else
734  {
736  {
737  return StdRegions::eBackwards;
738  }
739  else
740  {
741  return StdRegions::eForwards;
742  }
743  }
744  }
745 
746 
747  /**
748  *
749  */
751  {
752  int returnval = -1;
753 
754  SegGeomVector::iterator edgeIter;
755  int i;
756 
757  for (i=0,edgeIter = m_edges.begin(); edgeIter != m_edges.end(); ++edgeIter,++i)
758  {
759  if (*edgeIter == edge)
760  {
761  returnval = i;
762  break;
763  }
764  }
765 
766  return returnval;
767  }
768 
769 
770  /**
771  *
772  */
774  {
775  return kNverts;
776  }
777 
778 
779  /**
780  *
781  */
783  {
784  return kNedges;
785  }
786 
787 
789  const Array<OneD, const NekDouble> &gloCoord, NekDouble tol)
790  {
791  Array<OneD,NekDouble> locCoord(GetCoordim(),0.0);
792  return v_ContainsPoint(gloCoord,locCoord,tol);
793 
794  }
795  /**
796  *
797  */
799  Array<OneD, NekDouble> &stdCoord,
800  NekDouble tol)
801  {
802  NekDouble resid;
803  return v_ContainsPoint(gloCoord,stdCoord,tol,resid);
804  }
805 
807  Array<OneD, NekDouble> &stdCoord,
808  NekDouble tol,
809  NekDouble &resid)
810  {
811  ASSERTL1(gloCoord.num_elements() >= 2,
812  "Two dimensional geometry expects at least two coordinates.");
813 
814  resid = GetLocCoords(gloCoord, stdCoord);
815  if (stdCoord[0] >= -(1+tol) && stdCoord[1] >= -(1+tol)
816  && stdCoord[0] <= (1+tol) && stdCoord[1] <= (1+tol))
817  {
818  return true;
819  }
820  return false;
821  }
822 
823  void QuadGeom::v_Reset(CurveMap &curvedEdges,
824  CurveMap &curvedFaces)
825  {
826  Geometry::v_Reset(curvedEdges, curvedFaces);
827  CurveMap::iterator it = curvedFaces.find(m_globalID);
828 
829  if (it != curvedFaces.end())
830  {
831  m_curve = it->second;
832  }
833 
834  for (int i = 0; i < 4; ++i)
835  {
836  m_edges[i]->Reset(curvedEdges, curvedFaces);
837  }
838 
839  SetUpXmap();
840  SetUpCoeffs(m_xmap->GetNcoeffs());
841  }
842  }; //end of namespace
843 }; //end of namespace
StdRegions::StdExpansionSharedPtr m_xmap
Definition: Geometry.h:172
virtual const Geometry1DSharedPtr v_GetEdge(int i) const
Definition: QuadGeom.cpp:706
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
virtual int v_GetCoordim() const
Definition: QuadGeom.cpp:401
static boost::shared_ptr< DataType > AllocateSharedPtr()
Allocate a shared pointer from the memory pool.
virtual bool v_ContainsPoint(const Array< OneD, const NekDouble > &gloCoord, NekDouble tol=0.0)
Definition: QuadGeom.cpp:788
std::vector< PointGeomSharedPtr > PointGeomVector
Definition: Geometry3D.h:61
GeomFactorsSharedPtr m_geomFactors
Definition: Geometry.h:170
Structure holding graphvertexobject id and local element facet id.
NekDouble dot(PointGeom &a)
Definition: PointGeom.cpp:192
int Imin(int n, const T *x, const int incx)
Return the index of the minimum element in x.
Definition: Vmath.cpp:833
virtual int v_NumElmtConnected() const
Definition: QuadGeom.cpp:363
void Vvtvp(int n, const T *w, const int incw, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
vvtvp (vector times vector plus vector): z = w*x + y
Definition: Vmath.cpp:428
virtual PointGeomSharedPtr v_GetVertex(int i) const
Definition: QuadGeom.cpp:696
Principle Modified Functions .
Definition: BasisType.h:49
void Sub(PointGeom &a, PointGeom &b)
Definition: PointGeom.cpp:168
virtual void v_Reset(CurveMap &curvedEdges, CurveMap &curvedFaces)
Reset this geometry object: unset the current state and remove allocated GeomFactors.
Definition: Geometry.cpp:307
StdRegions::StdExpansionSharedPtr GetXmap() const
Definition: Geometry.h:383
boost::shared_ptr< Curve > CurveSharedPtr
Definition: Curve.hpp:62
virtual bool v_IsElmtConnected(int gvo_id, int locid) const
Definition: QuadGeom.cpp:372
virtual int v_GetVid(int i) const
Definition: QuadGeom.cpp:686
const LibUtilities::BasisSharedPtr GetBasis(const int i)
Return the j-th basis of the i-th co-ordinate dimension.
Definition: Geometry.h:475
virtual int v_GetNumEdges() const
Definition: QuadGeom.cpp:782
std::list< CompToElmt > m_elmtMap
Definition: QuadGeom.h:110
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:116
NekDouble GetCoord(const int i, const Array< OneD, const NekDouble > &Lcoord)
Definition: QuadGeom.cpp:241
static StdRegions::Orientation GetFaceOrientation(const QuadGeom &face1, const QuadGeom &face2)
Get the orientation of face1.
Definition: QuadGeom.cpp:253
static const NekDouble kNekZeroTol
virtual NekDouble v_GetCoord(const int i, const Array< OneD, const NekDouble > &Lcoord)
Definition: QuadGeom.cpp:436
virtual void v_AddElmtConnected(int gvo_id, int locid)
Definition: QuadGeom.cpp:353
boost::shared_ptr< SegGeom > SegGeomSharedPtr
Definition: Geometry2D.h:60
void NewtonIterationForLocCoord(const Array< OneD, const NekDouble > &coords, const Array< OneD, const NekDouble > &ptsx, const Array< OneD, const NekDouble > &ptsy, Array< OneD, NekDouble > &Lcoords, NekDouble &resid)
Definition: Geometry2D.cpp:100
virtual void v_Reset(CurveMap &curvedEdges, CurveMap &curvedFaces)
Reset this geometry object: unset the current state and remove allocated GeomFactors.
Definition: QuadGeom.cpp:823
virtual StdRegions::Orientation v_GetCartesianEorient(const int i) const
Definition: QuadGeom.cpp:726
virtual const LibUtilities::BasisSharedPtr v_GetEdgeBasis(const int i)
Definition: QuadGeom.cpp:419
static std::string npts
Definition: InputFld.cpp:43
virtual int v_GetNumVerts() const
Definition: QuadGeom.cpp:773
virtual int v_GetFid() const
Definition: QuadGeom.cpp:392
Defines a specification for a set of points.
Definition: Points.h:58
double NekDouble
void Sadd(int n, const T alpha, const T *x, const int incx, T *y, const int incy)
Add vector y = alpha + x.
Definition: Vmath.cpp:301
PointGeomSharedPtr GetVertex(int i) const
Definition: Geometry.h:348
virtual int v_WhichEdge(SegGeomSharedPtr edge)
Return the edge number of the given edge.
Definition: QuadGeom.cpp:750
void Mult(PointGeom &a, PointGeom &b)
Definition: PointGeom.cpp:177
StandardMatrixTag boost::call_traits< LhsDataType >::const_reference rhs typedef NekMatrix< LhsDataType, StandardMatrixTag >::iterator iterator
virtual void v_FillGeom()
Put all quadrature information into edge structure.
Definition: QuadGeom.cpp:527
Array< OneD, Array< OneD, NekDouble > > m_coeffs
Definition: Geometry.h:180
2D geometry information
Definition: Geometry2D.h:65
#define ASSERTL2(condition, msg)
Assert Level 2 – Debugging which is used FULLDEBUG compilation mode. This level assert is designed t...
Definition: ErrorUtil.hpp:213
Geometry is straight-sided with constant geometric factors.
boost::shared_ptr< Geometry1D > Geometry1DSharedPtr
Definition: Geometry1D.h:48
NekDouble GetLocCoords(const Array< OneD, const NekDouble > &coords, Array< OneD, NekDouble > &Lcoords)
Definition: Geometry.h:450
LibUtilities::ShapeType m_shapeType
Definition: Geometry.h:177
InputIterator find(InputIterator first, InputIterator last, InputIterator startingpoint, const EqualityComparable &value)
Definition: StdRegions.hpp:315
virtual int v_GetEid(int i) const
Definition: QuadGeom.cpp:676
boost::unordered_map< int, CurveSharedPtr > CurveMap
Definition: Curve.hpp:63
Geometric information has been generated.
void SetUpCoeffs(const int nCoeffs)
Initialise the m_coeffs array.
Definition: Geometry.h:484
StdRegions::Orientation m_eorient[kNedges]
Definition: QuadGeom.h:107
GeomFactorsSharedPtr GetMetricInfo()
Definition: Geometry.h:299
GeomType
Indicates the type of element geometry.
bool m_ownData
Boolean indicating whether object owns the data.
Definition: QuadGeom.h:193
boost::shared_ptr< Basis > BasisSharedPtr
virtual NekDouble v_GetLocCoords(const Array< OneD, const NekDouble > &coords, Array< OneD, NekDouble > &Lcoords)
Definition: QuadGeom.cpp:609
GeomState m_state
enum identifier to determine if quad points are filled
Definition: Geometry.h:175
virtual StdRegions::Orientation v_GetEorient(const int i) const
Definition: QuadGeom.cpp:716
virtual const LibUtilities::BasisSharedPtr v_GetBasis(const int i)
Definition: QuadGeom.cpp:410
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode...
Definition: ErrorUtil.hpp:191
Geometry is curved or has non-constant factors.
int m_coordim
coordinate dimension
Definition: Geometry.h:169
Describes the specification for a Basis.
Definition: Basis.h:50
boost::shared_ptr< PointGeom > PointGeomSharedPtr
Definition: Geometry.h:60
1D Gauss-Lobatto-Legendre quadrature points
Definition: PointsType.h:50
void Vmul(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Multiply vector z = x*y.
Definition: Vmath.cpp:169