Nektar++
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 points0 = max(m_edges[0]->GetBasis(0)->GetNumPoints(),
223  m_edges[2]->GetBasis(0)->GetNumPoints());
224  int order1 = max(m_edges[1]->GetBasis(0)->GetNumModes(),
225  m_edges[3]->GetBasis(0)->GetNumModes());
226  int points1 = max(m_edges[1]->GetBasis(0)->GetNumPoints(),
227  m_edges[3]->GetBasis(0)->GetNumPoints());
228 
229  const LibUtilities::BasisKey B0(
233  const LibUtilities::BasisKey B1(
237 
239  ::AllocateSharedPtr(B0,B1);
240  }
241 
242  /**
243  *
244  */
246  const Array<OneD, const NekDouble> &Lcoord)
247  {
249  "Geometry is not in physical space");
250 
251  Array<OneD, NekDouble> tmp(m_xmap->GetTotPoints());
252  m_xmap->BwdTrans(m_coeffs[i], tmp);
253 
254  return m_xmap->PhysEvaluate(Lcoord, tmp);
255  }
256 
258  const QuadGeom &face1,
259  const QuadGeom &face2)
260  {
261  return GetFaceOrientation(face1.m_verts, face2.m_verts);
262  }
263 
264  /**
265  * Calculate the orientation of face2 to face1 (note this is
266  * not face1 to face2!).
267  */
269  const PointGeomVector &face1,
270  const PointGeomVector &face2)
271  {
272  int i, j, vmap[4] = {-1,-1,-1,-1};
273  NekDouble x, y, z, x1, y1, z1, cx = 0.0, cy = 0.0, cz = 0.0;
274 
275  // For periodic faces, we calculate the vector between the centre
276  // points of the two faces. (For connected faces this will be
277  // zero). We can then use this to determine alignment later in the
278  // algorithm.
279  for (i = 0; i < 4; ++i)
280  {
281  cx += (*face2[i])(0) - (*face1[i])(0);
282  cy += (*face2[i])(1) - (*face1[i])(1);
283  cz += (*face2[i])(2) - (*face1[i])(2);
284  }
285  cx /= 4;
286  cy /= 4;
287  cz /= 4;
288 
289  // Now construct a mapping which takes us from the vertices of one
290  // face to the other. That is, vertex j of face2 corresponds to
291  // vertex vmap[j] of face1.
292  for (i = 0; i < 4; ++i)
293  {
294  x = (*face1[i])(0);
295  y = (*face1[i])(1);
296  z = (*face1[i])(2);
297  for (j = 0; j < 4; ++j)
298  {
299  x1 = (*face2[j])(0)-cx;
300  y1 = (*face2[j])(1)-cy;
301  z1 = (*face2[j])(2)-cz;
302  if (sqrt((x1-x)*(x1-x)+(y1-y)*(y1-y)+(z1-z)*(z1-z)) < 1e-8)
303  {
304  vmap[j] = i;
305  break;
306  }
307  }
308  }
309 
310  // Use the mapping to determine the eight alignment options between
311  // faces.
312  if (vmap[1] == (vmap[0]+1) % 4)
313  {
314  switch (vmap[0])
315  {
316  case 0:
318  break;
319  case 1:
321  break;
322  case 2:
324  break;
325  case 3:
327  break;
328  }
329  }
330  else
331  {
332  switch (vmap[0])
333  {
334  case 0:
336  break;
337  case 1:
339  break;
340  case 2:
342  break;
343  case 3:
345  break;
346  }
347  }
348 
349  ASSERTL0(false, "unable to determine face orientation");
351  }
352 
353 
354  /**
355  *
356  */
357  void QuadGeom::v_AddElmtConnected(int gvo_id, int locid)
358  {
359  CompToElmt ee(gvo_id,locid);
360  m_elmtMap.push_back(ee);
361  }
362 
363 
364  /**
365  *
366  */
368  {
369  return int(m_elmtMap.size());
370  }
371 
372 
373  /**
374  *
375  */
376  bool QuadGeom::v_IsElmtConnected(int gvo_id, int locid) const
377  {
378  std::list<CompToElmt>::const_iterator def;
379  CompToElmt ee(gvo_id,locid);
380 
381  def = find(m_elmtMap.begin(),m_elmtMap.end(),ee);
382 
383  // Found the element connectivity object in the list
384  if(def != m_elmtMap.end())
385  {
386  return(true);
387  }
388 
389  return(false);
390  }
391 
392 
393  /**
394  *
395  */
396  int QuadGeom::v_GetFid() const
397  {
398  return m_fid;
399  }
400 
401 
402  /**
403  *
404  */
406  {
407  return m_coordim;
408  }
409 
410 
411  /**
412  *
413  */
415  {
416  return m_xmap->GetBasis(i);
417  }
418 
419 
420  /**
421  *
422  */
424  {
425  ASSERTL1(i <= 3,"edge is out of range");
426 
427  if (i == 0 || i == 2)
428  {
429  return m_xmap->GetBasis(0);
430  }
431  else
432  {
433  return m_xmap->GetBasis(1);
434  }
435  }
436 
437  /**
438  *
439  */
441  {
442  return GetCoord(i,Lcoord);
443  }
444 
445 
446  /**
447  * Set up GeoFac for this geometry using Coord quadrature distribution
448  */
450  {
452  {
453  int i;
454  GeomType Gtype = eRegular;
455 
457 
458  // We will first check whether we have a regular or deformed
459  // geometry. We will define regular as those cases where the
460  // Jacobian and the metric terms of the derivative are constants
461  // (i.e. not coordinate dependent)
462 
463  // Check to see if expansions are linear
464  // If not linear => deformed geometry
465  for(i = 0; i < m_coordim; ++i)
466  {
467  if((m_xmap->GetBasisNumModes(0) != 2)||
468  (m_xmap->GetBasisNumModes(1) != 2))
469  {
470  Gtype = eDeformed;
471  }
472  }
473 
474  // For linear expansions, the mapping from standard to local
475  // element is given by the relation:
476  // x_i = 0.25 * [ ( x_i^A + x_i^B + x_i^C + x_i^D) +
477  // (-x_i^A + x_i^B + x_i^C - x_i^D)*xi_1 +
478  // (-x_i^A - x_i^B + x_i^C + x_i^D)*xi_2 +
479  // ( x_i^A - x_i^B + x_i^C - x_i^D)*xi_1*xi_2 ]
480  //
481  // The jacobian of the transformation and the metric terms
482  // dxi_i/dx_j, involve only terms of the form dx_i/dxi_j (both
483  // for coordim == 2 or 3). Inspecting the formula above, it can
484  // be appreciated that the derivatives dx_i/dxi_j will be
485  // constant, if the coefficient of the non-linear term is zero.
486  //
487  // That is why for regular geometry, we require
488  //
489  // x_i^A - x_i^B + x_i^C - x_i^D = 0
490  //
491  // or equivalently
492  //
493  // x_i^A - x_i^B = x_i^D - x_i^C
494  //
495  // This corresponds to quadrilaterals which are paralellograms.
496  if(Gtype == eRegular)
497  {
498  for(i = 0; i < m_coordim; i++)
499  {
500  if( fabs( (*m_verts[0])(i) - (*m_verts[1])(i) +
501  (*m_verts[2])(i) - (*m_verts[3])(i) )
503  {
504  Gtype = eDeformed;
505  break;
506  }
507  }
508  }
509 
511  Gtype, m_coordim, m_xmap, m_coeffs);
513  }
514  }
515 
516 
517  /**
518  *
519  */
521  {
522  m_ownData = true;
523  }
524 
525 
526  /**
527  * Note verts and edges are listed according to anticlockwise
528  * convention but points in _coeffs have to be in array format from
529  * left to right.
530  */
532  {
533  // check to see if geometry structure is already filled
534  if(m_state != ePtsFilled)
535  {
536  int i,j,k;
537  int nEdgeCoeffs;
538 
539  if (m_curve)
540  {
541  int npts = m_curve->m_points.size();
542  int nEdgePts = (int)sqrt(static_cast<NekDouble>(npts));
543  Array<OneD, NekDouble> tmp(npts);
544  Array<OneD, NekDouble> tmp2(m_xmap->GetTotPoints());
545  LibUtilities::PointsKey curveKey(
546  nEdgePts, m_curve->m_ptype);
547 
548  // Sanity checks:
549  // - Curved faces should have square number of points;
550  // - Each edge should have sqrt(npts) points.
551  ASSERTL0(nEdgePts * nEdgePts == npts,
552  "NUMPOINTS should be a square number in"
553  " quadrilteral "
554  + boost::lexical_cast<string>(m_globalID));
555 
556  for (i = 0; i < kNedges; ++i)
557  {
558  ASSERTL0(
559  m_edges[i]->GetXmap()->GetNcoeffs() == nEdgePts,
560  "Number of edge points does not correspond to "
561  "number of face points in quadrilateral "
562  + boost::lexical_cast<string>(m_globalID));
563  }
564 
565  for (i = 0; i < m_coordim; ++i)
566  {
567  for (j = 0; j < npts; ++j)
568  {
569  tmp[j] = (m_curve->m_points[j]->GetPtr())[i];
570  }
571 
572  // Interpolate m_curve points to GLL points
574  curveKey, curveKey, tmp,
575  m_xmap->GetBasis(0)->GetPointsKey(),
576  m_xmap->GetBasis(1)->GetPointsKey(),
577  tmp2);
578 
579  // Forwards transform to get coefficient space.
580  m_xmap->FwdTrans(tmp2, m_coeffs[i]);
581  }
582  }
583 
584  // Now fill in edges.
585  Array<OneD, unsigned int> mapArray;
586  Array<OneD, int> signArray;
587 
588  for(i = 0; i < kNedges; i++)
589  {
590  m_edges[i]->FillGeom();
591  m_xmap->GetEdgeToElementMap(i,m_eorient[i],
592  mapArray,signArray);
593 
594  nEdgeCoeffs = m_edges[i]->GetXmap()->GetNcoeffs();
595 
596  for(j = 0; j < m_coordim; j++)
597  {
598  for(k = 0; k < nEdgeCoeffs; k++)
599  {
600  m_coeffs[j][mapArray[k]]
601  = signArray[k]*(m_edges[i]->GetCoeffs(j))[k];
602  }
603  }
604  }
605 
607  }
608  }
609 
610  /**
611  *
612  */
614  Array<OneD,NekDouble> &Lcoords)
615  {
616  NekDouble resid = 0.0;
617  if(GetMetricInfo()->GetGtype() == eRegular)
618  {
619  NekDouble coords2 = (m_coordim == 3)? coords[2]: 0.0;
620  PointGeom dv1, dv2, norm, orth1, orth2;
621  PointGeom xin(m_coordim,0,coords[0],coords[1],coords2);
622 
623  // Calculate edge vectors from 0-1 and 0-3 edges.
624  dv1.Sub(*m_verts[1],*m_verts[0]);
625  dv2.Sub(*m_verts[3],*m_verts[0]);
626 
627  // Obtain normal to plane in which dv1 and dv2 lie
628  norm.Mult(dv1,dv2);
629 
630  // Obtain vector which are normal to dv1 and dv2.
631  orth1.Mult(norm,dv1);
632  orth2.Mult(norm,dv2);
633 
634  // Start with vector of desired points minus vertex_0
635  xin -= *m_verts[0];
636 
637  // Calculate length using L/|dv1| = (x-v0).n1/(dv1.n1) for coordiante 1
638  // Then rescale to [-1,1].
639  Lcoords[0] = xin.dot(orth2)/dv1.dot(orth2);
640  Lcoords[0] = 2*Lcoords[0]-1;
641  Lcoords[1] = xin.dot(orth1)/dv2.dot(orth1);
642  Lcoords[1] = 2*Lcoords[1]-1;
643  }
644  else
645  {
647 
648  // Determine nearest point of coords to values in m_xmap
649  int npts = m_xmap->GetTotPoints();
650  Array<OneD, NekDouble> ptsx(npts), ptsy(npts);
651  Array<OneD, NekDouble> tmpx(npts), tmpy(npts);
652 
653  m_xmap->BwdTrans(m_coeffs[0], ptsx);
654  m_xmap->BwdTrans(m_coeffs[1], ptsy);
655 
656  const Array<OneD, const NekDouble> za = m_xmap->GetPoints(0);
657  const Array<OneD, const NekDouble> zb = m_xmap->GetPoints(1);
658 
659  //guess the first local coords based on nearest point
660  Vmath::Sadd(npts, -coords[0], ptsx,1,tmpx,1);
661  Vmath::Sadd(npts, -coords[1], ptsy,1,tmpy,1);
662  Vmath::Vmul (npts, tmpx,1,tmpx,1,tmpx,1);
663  Vmath::Vvtvp(npts, tmpy,1,tmpy,1,tmpx,1,tmpx,1);
664 
665  int min_i = Vmath::Imin(npts,tmpx,1);
666 
667  Lcoords[0] = za[min_i%za.num_elements()];
668  Lcoords[1] = zb[min_i/za.num_elements()];
669 
670  // Perform newton iteration to find local coordinates
671  NewtonIterationForLocCoord(coords, ptsx, ptsy, Lcoords,resid);
672  }
673  return resid;
674  }
675 
676 
677  /**
678  *
679  */
680  int QuadGeom::v_GetEid(int i) const
681  {
682  ASSERTL2((i >=0) && (i <= 3),"Edge id must be between 0 and 3");
683  return m_edges[i]->GetEid();
684  }
685 
686 
687  /**
688  *
689  */
690  int QuadGeom::v_GetVid(int i) const
691  {
692  ASSERTL2((i >=0) && (i <= 3),"Verted id must be between 0 and 3");
693  return m_verts[i]->GetVid();
694  }
695 
696 
697  /**
698  *
699  */
701  {
702  ASSERTL2((i >=0) && (i <= 3),"Vertex id must be between 0 and 3");
703  return m_verts[i];
704  }
705 
706 
707  /**
708  *
709  */
711  {
712  ASSERTL2((i >=0) && (i <= 3),"Edge id must be between 0 and 3");
713  return m_edges[i];
714  }
715 
716 
717  /**
718  *
719  */
721  {
722  ASSERTL2((i >=0) && (i <= 3),"Edge id must be between 0 and 3");
723  return m_eorient[i];
724  }
725 
726 
727  /**
728  *
729  */
731  {
732  ASSERTL2((i >=0) && (i <= 3),"Edge id must be between 0 and 3");
733  if(i < 2)
734  {
735  return m_eorient[i];
736  }
737  else
738  {
740  {
741  return StdRegions::eBackwards;
742  }
743  else
744  {
745  return StdRegions::eForwards;
746  }
747  }
748  }
749 
750 
751  /**
752  *
753  */
755  {
756  int returnval = -1;
757 
758  SegGeomVector::iterator edgeIter;
759  int i;
760 
761  for (i=0,edgeIter = m_edges.begin(); edgeIter != m_edges.end(); ++edgeIter,++i)
762  {
763  if (*edgeIter == edge)
764  {
765  returnval = i;
766  break;
767  }
768  }
769 
770  return returnval;
771  }
772 
773 
774  /**
775  *
776  */
778  {
779  return kNverts;
780  }
781 
782 
783  /**
784  *
785  */
787  {
788  return kNedges;
789  }
790 
791 
793  const Array<OneD, const NekDouble> &gloCoord, NekDouble tol)
794  {
795  Array<OneD,NekDouble> locCoord(GetCoordim(),0.0);
796  return v_ContainsPoint(gloCoord,locCoord,tol);
797 
798  }
799  /**
800  *
801  */
803  Array<OneD, NekDouble> &stdCoord,
804  NekDouble tol)
805  {
806  NekDouble resid;
807  return v_ContainsPoint(gloCoord,stdCoord,tol,resid);
808  }
809 
811  Array<OneD, NekDouble> &stdCoord,
812  NekDouble tol,
813  NekDouble &resid)
814  {
815  ASSERTL1(gloCoord.num_elements() >= 2,
816  "Two dimensional geometry expects at least two coordinates.");
817 
818  resid = GetLocCoords(gloCoord, stdCoord);
819  if (stdCoord[0] >= -(1+tol) && stdCoord[1] >= -(1+tol)
820  && stdCoord[0] <= (1+tol) && stdCoord[1] <= (1+tol))
821  {
822  return true;
823  }
824  return false;
825  }
826 
827  void QuadGeom::v_Reset(CurveMap &curvedEdges,
828  CurveMap &curvedFaces)
829  {
830  Geometry::v_Reset(curvedEdges, curvedFaces);
831  CurveMap::iterator it = curvedFaces.find(m_globalID);
832 
833  if (it != curvedFaces.end())
834  {
835  m_curve = it->second;
836  }
837 
838  for (int i = 0; i < 4; ++i)
839  {
840  m_edges[i]->Reset(curvedEdges, curvedFaces);
841  }
842 
843  SetUpXmap();
844  SetUpCoeffs(m_xmap->GetNcoeffs());
845  }
846  }; //end of namespace
847 }; //end of namespace
StdRegions::StdExpansionSharedPtr m_xmap
Definition: Geometry.h:172
virtual const Geometry1DSharedPtr v_GetEdge(int i) const
Definition: QuadGeom.cpp:710
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
virtual int v_GetCoordim() const
Definition: QuadGeom.cpp:405
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:792
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:824
virtual int v_NumElmtConnected() const
Definition: QuadGeom.cpp:367
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:700
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:376
virtual int v_GetVid(int i) const
Definition: QuadGeom.cpp:690
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:786
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:245
static StdRegions::Orientation GetFaceOrientation(const QuadGeom &face1, const QuadGeom &face2)
Get the orientation of face1.
Definition: QuadGeom.cpp:257
static const NekDouble kNekZeroTol
virtual NekDouble v_GetCoord(const int i, const Array< OneD, const NekDouble > &Lcoord)
Definition: QuadGeom.cpp:440
virtual void v_AddElmtConnected(int gvo_id, int locid)
Definition: QuadGeom.cpp:357
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:827
virtual StdRegions::Orientation v_GetCartesianEorient(const int i) const
Definition: QuadGeom.cpp:730
virtual const LibUtilities::BasisSharedPtr v_GetEdgeBasis(const int i)
Definition: QuadGeom.cpp:423
static std::string npts
Definition: InputFld.cpp:43
virtual int v_GetNumVerts() const
Definition: QuadGeom.cpp:777
virtual int v_GetFid() const
Definition: QuadGeom.cpp:396
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:754
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:531
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:312
virtual int v_GetEid(int i) const
Definition: QuadGeom.cpp:680
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:613
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:720
virtual const LibUtilities::BasisSharedPtr v_GetBasis(const int i)
Definition: QuadGeom.cpp:414
#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