Nektar++
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
Expansion.cpp
Go to the documentation of this file.
1 ///////////////////////////////////////////////////////////////////////////////
2 //
3 // File Expansion.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: File for Expansion routines
33 //
34 ///////////////////////////////////////////////////////////////////////////////
35 
36 #include <LibUtilities/Foundations/Interp.h>
37 #include <LocalRegions/Expansion.h>
38 #include <LocalRegions/MatrixKey.h>
39 
40 #include <SpatialDomains/MeshComponents.h>
41 
42 namespace Nektar
43 {
44  namespace LocalRegions
45  {
46  Expansion::Expansion(SpatialDomains::GeometrySharedPtr pGeom) :
47  m_geom(pGeom),
48  m_metricinfo(m_geom->GetGeomFactors())
49  {
50  if (!m_metricinfo)
51  {
52  return;
53  }
54 
55  if (!m_metricinfo->IsValid())
56  {
57  int nDim = m_base.num_elements();
58  string type = "regular";
59  if (m_metricinfo->GetGtype() == SpatialDomains::eDeformed)
60  {
61  type = "deformed";
62  }
63 
64  stringstream err;
65  err << nDim << "D " << type << " Jacobian not positive "
66  << "(element ID = " << m_geom->GetGlobalID() << ") "
67  << "(first vertex ID = " << m_geom->GetVid(0) << ")";
68  NEKERROR(ErrorUtil::ewarning, err.str());
69  }
70  }
71 
72  Expansion::Expansion(const Expansion &pSrc) :
73  m_geom(pSrc.m_geom),
74  m_metricinfo(pSrc.m_metricinfo)
75  {
76 
77  }
78 
79  Expansion::~Expansion()
80  {
81  }
82 
83  DNekScalMatSharedPtr Expansion::GetLocMatrix(const LocalRegions::MatrixKey &mkey)
84  {
85  return v_GetLocMatrix(mkey);
86  }
87 
88  DNekMatSharedPtr Expansion::BuildTransformationMatrix(
89  const DNekScalMatSharedPtr &r_bnd,
90  const StdRegions::MatrixType matrixType)
91  {
92  return v_BuildTransformationMatrix(r_bnd,matrixType);
93  }
94 
95 
96  DNekMatSharedPtr Expansion::BuildVertexMatrix(
97  const DNekScalMatSharedPtr &r_bnd)
98  {
99  return v_BuildVertexMatrix(r_bnd);
100  }
101 
102  void Expansion::AddEdgeNormBoundaryInt(
103  const int edge,
104  const boost::shared_ptr<Expansion> &EdgeExp,
105  const Array<OneD, const NekDouble> &Fx,
106  const Array<OneD, const NekDouble> &Fy,
107  Array<OneD, NekDouble> &outarray)
108  {
109  v_AddEdgeNormBoundaryInt(edge, EdgeExp, Fx, Fy, outarray);
110  }
111 
112  void Expansion::AddEdgeNormBoundaryInt(
113  const int edge,
114  const boost::shared_ptr<Expansion> &EdgeExp,
115  const Array<OneD, const NekDouble> &Fn,
116  Array<OneD, NekDouble> &outarray)
117  {
118  v_AddEdgeNormBoundaryInt(edge, EdgeExp, Fn, outarray);
119  }
120 
121  void Expansion::AddFaceNormBoundaryInt(
122  const int face,
123  const boost::shared_ptr<Expansion> &FaceExp,
124  const Array<OneD, const NekDouble> &Fn,
125  Array<OneD, NekDouble> &outarray)
126  {
127  v_AddFaceNormBoundaryInt(face, FaceExp, Fn, outarray);
128  }
129 
130  void Expansion::DGDeriv(
131  const int dir,
132  const Array<OneD, const NekDouble>& inarray,
133  Array<OneD, ExpansionSharedPtr> &EdgeExp,
134  Array<OneD, Array<OneD, NekDouble> > &coeffs,
135  Array<OneD, NekDouble> &outarray)
136  {
137  v_DGDeriv(dir, inarray, EdgeExp, coeffs, outarray);
138  }
139 
140  DNekScalMatSharedPtr Expansion::GetLocMatrix(const StdRegions::MatrixType mtype,
141  const StdRegions::ConstFactorMap &factors,
142  const StdRegions::VarCoeffMap &varcoeffs)
143  {
144  MatrixKey mkey(mtype, DetShapeType(), *this, factors, varcoeffs);
145  return GetLocMatrix(mkey);
146  }
147 
148  SpatialDomains::GeometrySharedPtr Expansion::GetGeom() const
149  {
150  return m_geom;
151  }
152 
153  const SpatialDomains::GeomFactorsSharedPtr& Expansion::v_GetMetricInfo() const
154  {
155  return m_metricinfo;
156  }
157 
158 
159  DNekScalMatSharedPtr Expansion::v_GetLocMatrix(const LocalRegions::MatrixKey &mkey)
160  {
161  NEKERROR(ErrorUtil::efatal, "This function is only valid for LocalRegions");
162  return NullDNekScalMatSharedPtr;
163  }
164 
165  void Expansion::v_MultiplyByQuadratureMetric(const Array<OneD, const NekDouble>& inarray,
166  Array<OneD, NekDouble> &outarray)
167  {
168  const int nqtot = GetTotPoints();
169 
170  if (m_metrics.count(MetricQuadrature) == 0)
171  {
172  ComputeQuadratureMetric();
173  }
174 
175  Vmath::Vmul(nqtot, m_metrics[MetricQuadrature], 1, inarray, 1, outarray, 1);
176  }
177 
178  void Expansion::ComputeLaplacianMetric()
179  {
180  v_ComputeLaplacianMetric();
181  }
182 
183  void Expansion::ComputeQuadratureMetric()
184  {
185  unsigned int nqtot = GetTotPoints();
186  SpatialDomains::GeomType type = m_metricinfo->GetGtype();
187  LibUtilities::PointsKeyVector p = GetPointsKeys();
188  if (type == SpatialDomains::eRegular ||
189  type == SpatialDomains::eMovingRegular)
190  {
191  m_metrics[MetricQuadrature] = Array<OneD, NekDouble>(nqtot, m_metricinfo->GetJac(p)[0]);
192  }
193  else
194  {
195  m_metrics[MetricQuadrature] = m_metricinfo->GetJac(p);
196  }
197 
198  MultiplyByStdQuadratureMetric(m_metrics[MetricQuadrature],
199  m_metrics[MetricQuadrature]);
200  }
201 
202  void Expansion::v_GetCoords(
203  Array<OneD, NekDouble> &coords_0,
204  Array<OneD, NekDouble> &coords_1,
205  Array<OneD, NekDouble> &coords_2)
206  {
207  ASSERTL1(m_geom, "m_geom not defined");
208 
209  // get physical points defined in Geom
210  m_geom->FillGeom();
211 
212  const int expDim = m_base.num_elements();
213  int nqGeom = 1;
214  bool doCopy = true;
215 
216  Array<OneD, LibUtilities::BasisSharedPtr> CBasis(expDim);
217  Array<OneD, Array<OneD, NekDouble> > tmp(3);
218 
219  for (int i = 0; i < expDim; ++i)
220  {
221  CBasis[i] = m_geom->GetBasis(i);
222  nqGeom *= CBasis[i]->GetNumPoints();
223  doCopy = doCopy && m_base[i]->GetBasisKey().SamePoints(
224  CBasis[i]->GetBasisKey());
225  }
226 
227  tmp[0] = coords_0;
228  tmp[1] = coords_1;
229  tmp[2] = coords_2;
230 
231  if (doCopy)
232  {
233  for (int i = 0; i < m_geom->GetCoordim(); ++i)
234  {
235  m_geom->GetXmap()->BwdTrans(m_geom->GetCoeffs(i), tmp[i]);
236  }
237  }
238  else
239  {
240  for (int i = 0; i < m_geom->GetCoordim(); ++i)
241  {
242  Array<OneD, NekDouble> tmpGeom(nqGeom);
243  m_geom->GetXmap()->BwdTrans(m_geom->GetCoeffs(i), tmpGeom);
244 
245  switch (expDim)
246  {
247  case 1:
248  {
249  LibUtilities::Interp1D(
250  CBasis[0]->GetPointsKey(), &tmpGeom[0],
251  m_base[0]->GetPointsKey(), &tmp[i][0]);
252  break;
253  }
254  case 2:
255  {
256  LibUtilities::Interp2D(
257  CBasis[0]->GetPointsKey(),
258  CBasis[1]->GetPointsKey(),
259  &tmpGeom[0],
260  m_base[0]->GetPointsKey(),
261  m_base[1]->GetPointsKey(),
262  &tmp[i][0]);
263  break;
264  }
265  case 3:
266  {
267  LibUtilities::Interp3D(
268  CBasis[0]->GetPointsKey(),
269  CBasis[1]->GetPointsKey(),
270  CBasis[2]->GetPointsKey(),
271  &tmpGeom[0],
272  m_base[0]->GetPointsKey(),
273  m_base[1]->GetPointsKey(),
274  m_base[2]->GetPointsKey(),
275  &tmp[i][0]);
276  break;
277  }
278  }
279  }
280  }
281  }
282 
283  DNekMatSharedPtr Expansion::v_BuildTransformationMatrix(
284  const DNekScalMatSharedPtr &r_bnd,
285  const StdRegions::MatrixType matrixType)
286  {
287  NEKERROR(ErrorUtil::efatal, "This function is only valid for LocalRegions");
288  return NullDNekMatSharedPtr;
289  }
290 
291  DNekMatSharedPtr Expansion::v_BuildVertexMatrix(
292  const DNekScalMatSharedPtr &r_bnd)
293  {
294  NEKERROR(ErrorUtil::efatal, "This function is only valid for LocalRegions");
295  return NullDNekMatSharedPtr;
296  }
297 
298  void Expansion::v_AddEdgeNormBoundaryInt(
299  const int edge,
300  const boost::shared_ptr<Expansion> &EdgeExp,
301  const Array<OneD, const NekDouble> &Fx,
302  const Array<OneD, const NekDouble> &Fy,
303  Array<OneD, NekDouble> &outarray)
304  {
305  NEKERROR(ErrorUtil::efatal, "This function is only valid for LocalRegions");
306  }
307 
308  void Expansion::v_AddEdgeNormBoundaryInt(
309  const int edge,
310  const boost::shared_ptr<Expansion> &EdgeExp,
311  const Array<OneD, const NekDouble> &Fn,
312  Array<OneD, NekDouble> &outarray)
313  {
314  NEKERROR(ErrorUtil::efatal, "This function is only valid for LocalRegions");
315  }
316 
317  void Expansion::v_AddFaceNormBoundaryInt(
318  const int face,
319  const boost::shared_ptr<Expansion> &FaceExp,
320  const Array<OneD, const NekDouble> &Fn,
321  Array<OneD, NekDouble> &outarray)
322  {
323  NEKERROR(ErrorUtil::efatal, "This function is only valid for LocalRegions");
324  }
325 
326  void Expansion::v_DGDeriv(
327  const int dir,
328  const Array<OneD, const NekDouble>& inarray,
329  Array<OneD, ExpansionSharedPtr> &EdgeExp,
330  Array<OneD, Array<OneD, NekDouble> > &coeffs,
331  Array<OneD, NekDouble> &outarray)
332  {
333  NEKERROR(ErrorUtil::efatal, "This function is only valid for LocalRegions");
334  }
335  } //end of namespace
336 } //end of namespace
337 
Generated on Wed Mar 4 2015 16:33:13 for Nektar++ by   doxygen 1.8.1.2