Nektar++
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
StdSegExp.cpp
Go to the documentation of this file.
1 ///////////////////////////////////////////////////////////////////////////////
2 //
3 // File StdSegExp.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: Routines within Standard Segment Expansions
33 //
34 ///////////////////////////////////////////////////////////////////////////////
35 
36 #include <StdRegions/StdSegExp.h>
37 
38 using namespace std;
39 
40 namespace Nektar
41 {
42  namespace StdRegions
43  {
44 
45  /** \brief Default constructor */
46 
47  StdSegExp::StdSegExp()
48  {
49  }
50 
51 
52  /** \brief Constructor using BasisKey class for quadrature points and
53  * order definition
54  *
55  * \param Ba BasisKey class definition containing order and quadrature
56  * points.
57  */
58 
59  StdSegExp::StdSegExp(const LibUtilities::BasisKey &Ba):
60  StdExpansion(Ba.GetNumModes(), 1, Ba),
61  StdExpansion1D(Ba.GetNumModes(),Ba)
62  {
63  }
64 
65 
66  /** \brief Copy Constructor */
67 
69  StdExpansion(T),
71  {
72  }
73 
74 
76  {
77  }
78 
79  /** \brief Return Shape of region, using ShapeType enum list.
80  * i.e. Segment
81  */
83  {
85  }
86 
88  {
89 
90  bool returnval = false;
91 
93  {
94  returnval = true;
95  }
96 
98  {
99  returnval = true;
100  }
101 
102  return returnval;
103  }
104 
105 
106 
107 
108  //---------------------------------------------------------------------
109  // Integration Methods
110  //---------------------------------------------------------------------
111 
112  /** \brief Integrate the physical point list \a inarray over region
113  * and return the value
114  *
115  * \param inarray definition of function to be integrated evauluated at
116  * quadrature point of expansion.
117  * \return returns \f$\int^1_{-1} u(\xi_1)d \xi_1 \f$ where \f$inarray[i]
118  * = u(\xi_{1i}) \f$
119  */
121  {
122  NekDouble Int = 0.0;
123  int nquad0 = m_base[0]->GetNumPoints();
124  Array<OneD, NekDouble> tmp(nquad0);
125  Array<OneD, const NekDouble> z = m_base[0]->GetZ();
126  Array<OneD, const NekDouble> w0 = m_base[0]->GetW();
127 
128  // multiply by integration constants
129  Vmath::Vmul(nquad0, inarray, 1, w0, 1, tmp, 1);
130 
131  Int = Vmath::Vsum(nquad0, tmp, 1);
132 
133  return Int;
134  }
135 
136 
137 
138 
139  //---------------------------------------------------------------------
140  // Differentiation Methods
141  //---------------------------------------------------------------------
142 
143 
144  /** \brief Evaluate the derivative \f$ d/d{\xi_1} \f$ at the physical
145  * quadrature points given by \a inarray and return in \a outarray.
146  *
147  * This is a wrapper around StdExpansion1D::Tensor_Deriv
148  * \param inarray array of a function evaluated at the quadrature points
149  * \param outarray the resulting array of the derivative \f$
150  * du/d_{\xi_1}|_{\xi_{1i}} \f$ will be stored in the array \a outarra
151  */
152 
154  Array<OneD, NekDouble> &out_d0,
155  Array<OneD, NekDouble> &out_d1,
156  Array<OneD, NekDouble> &out_d2)
157  {
158  PhysTensorDeriv(inarray,out_d0);
159  }
160 
161 
162  void StdSegExp::v_PhysDeriv(const int dir,
163  const Array<OneD, const NekDouble>& inarray,
164  Array<OneD, NekDouble> &outarray)
165  {
166  ASSERTL1(dir==0,"input dir is out of range");
167  PhysTensorDeriv(inarray,outarray);
168  // PhysDeriv(inarray, outarray);
169  }
170 
172  const Array<OneD, const NekDouble>& inarray,
173  Array<OneD, NekDouble> &out_d0,
174  Array<OneD, NekDouble> &out_d1,
175  Array<OneD, NekDouble> &out_d2)
176  {
177  PhysTensorDeriv(inarray,out_d0);
178  // PhysDeriv(inarray, out_d0);
179  }
180 
182  const int dir,
183  const Array<OneD, const NekDouble>& inarray,
184  Array<OneD, NekDouble> &outarray)
185  {
186  ASSERTL1(dir==0,"input dir is out of range");
187  PhysTensorDeriv(inarray,outarray);
188  // PhysDeriv(inarray, outarray);
189  }
190 
191 
192 
193  //---------------------------------------------------------------------
194  // Transforms
195  //---------------------------------------------------------------------
196 
197  /** \brief Backward transform from coefficient space given
198  * in \a inarray and evaluate at the physical quadrature
199  * points \a outarray
200  *
201  * Operation can be evaluated as \f$ u(\xi_{1i}) =
202  * \sum_{p=0}^{order-1} \hat{u}_p \phi_p(\xi_{1i}) \f$ or equivalently
203  * \f$ {\bf u} = {\bf B}^T {\bf \hat{u}} \f$ where
204  * \f${\bf B}[i][j] = \phi_i(\xi_{1j}), \mbox{\_coeffs}[p] = {\bf
205  * \hat{u}}[p] \f$
206  *
207  * The function takes the coefficient array \a inarray as
208  * input for the transformation
209  *
210  * \param inarray: the coeffficients of the expansion
211  *
212  * \param outarray: the resulting array of the values of the function at
213  * the physical quadrature points will be stored in the array \a outarray
214  */
215 
217  const Array<OneD, const NekDouble>& inarray,
218  Array<OneD, NekDouble> &outarray)
219  {
220  int nquad = m_base[0]->GetNumPoints();
221 
222  if(m_base[0]->Collocation())
223  {
224  Vmath::Vcopy(nquad, inarray, 1, outarray, 1);
225  }
226  else
227  {
228 
229 #ifdef NEKTAR_USING_DIRECT_BLAS_CALLS
230 
231  Blas::Dgemv('N',nquad,m_base[0]->GetNumModes(),1.0, (m_base[0]->GetBdata()).get(),
232  nquad,&inarray[0],1,0.0,&outarray[0],1);
233 
234 #else //NEKTAR_USING_DIRECT_BLAS_CALLS
235 
237  NekVector<NekDouble> out(nquad,outarray,eWrapper);
238  NekMatrix<NekDouble> B(nquad,m_ncoeffs,m_base[0]->GetBdata(),eWrapper);
239  out = B * in;
240 
241 #endif //NEKTAR_USING_DIRECT_BLAS_CALLS
242 
243  }
244  }
245 
246  /** \brief Forward transform from physical quadrature space stored in
247  * \a inarray and evaluate the expansion coefficients and store in
248  * \a outarray
249  *
250  * Perform a forward transform using a Galerkin projection by taking the
251  * inner product of the physical points and multiplying by the inverse of
252  * the mass matrix using the Solve method of the standard matrix
253  * container holding the local mass matrix, i.e. \f$ {\bf \hat{u}} =
254  * {\bf M}^{-1} {\bf I} \f$ where \f$ {\bf I}[p] = \int^1_{-1}
255  * \phi_p(\xi_1) u(\xi_1) d\xi_1 \f$
256  *
257  * This function stores the expansion coefficients calculated by the
258  * transformation in the coefficient space array \a outarray
259  *
260  * \param inarray: array of physical quadrature points to be transformed
261  *
262  * \param outarray: the coeffficients of the expansion
263  */
264 
266  Array<OneD, NekDouble> &outarray)
267  {
268  if(m_base[0]->Collocation())
269  {
270  Vmath::Vcopy(m_ncoeffs, inarray, 1, outarray, 1);
271  }
272  else
273  {
274  v_IProductWRTBase(inarray,outarray);
275 
276  // get Mass matrix inverse
277  StdMatrixKey masskey(eInvMass,v_DetShapeType(),*this);
278  DNekMatSharedPtr matsys = GetStdMatrix(masskey);
279 
280  NekVector<NekDouble> in(m_ncoeffs,outarray,eCopy);
282 
283  out = (*matsys)*in;
284  }
285  }
286 
288  const Array<OneD, const NekDouble>& inarray,
289  Array<OneD, NekDouble> &outarray)
290  {
291  if(m_base[0]->Collocation())
292  {
293  Vmath::Vcopy(m_ncoeffs, inarray, 1, outarray, 1);
294  }
295  else
296  {
297  int nInteriorDofs = m_ncoeffs-2;
298  int offset;
299 
300  switch(m_base[0]->GetBasisType())
301  {
303  {
304  offset = 1;
305  }
306  break;
308  {
309  nInteriorDofs = m_ncoeffs;
310  offset = 0;
311  }
312  break;
315  {
316  offset = 2;
317  }
318  break;
319  default:
320  ASSERTL0(false,"This type of FwdTrans is not defined for this expansion type");
321  }
322 
323  fill(outarray.get(), outarray.get()+m_ncoeffs, 0.0 );
324 
326  {
327  outarray[GetVertexMap(0)] = inarray[0];
328  outarray[GetVertexMap(1)] = inarray[m_base[0]->GetNumPoints()-1];
329 
330  if(m_ncoeffs>2)
331  {
332  // ideally, we would like to have tmp0 to be replaced by
333  // outarray (currently MassMatrixOp does not allow aliasing)
336 
337  StdMatrixKey masskey(eMass,v_DetShapeType(),*this);
338  MassMatrixOp(outarray,tmp0,masskey);
339  v_IProductWRTBase(inarray,tmp1);
340 
341  Vmath::Vsub(m_ncoeffs, tmp1, 1, tmp0, 1, tmp1, 1);
342 
343  // get Mass matrix inverse (only of interior DOF)
344  DNekMatSharedPtr matsys =
345  (m_stdStaticCondMatrixManager[masskey])-> GetBlock(1,1);
346 
347  Blas::Dgemv('N',nInteriorDofs,nInteriorDofs,1.0,
348  &(matsys->GetPtr())[0],nInteriorDofs,tmp1.get()
349  +offset,1,0.0,outarray.get()+offset,1);
350  }
351  }
352  else
353  {
354  StdSegExp::v_FwdTrans(inarray, outarray);
355  }
356  }
357 
358  }
359 
360 
362  Array<OneD, NekDouble> &outarray)
363  {
364  v_BwdTrans(inarray, outarray);
365  }
366 
367 
368 
369  //---------------------------------------------------------------------
370  // Inner product functions
371  //---------------------------------------------------------------------
372 
373 
374 
375  /** \brief Inner product of \a inarray over region with respect to
376  * expansion basis \a base and return in \a outarray
377  *
378  * Calculate \f$ I[p] = \int^{1}_{-1} \phi_p(\xi_1) u(\xi_1) d\xi_1
379  * = \sum_{i=0}^{nq-1} \phi_p(\xi_{1i}) u(\xi_{1i}) w_i \f$ where
380  * \f$ outarray[p] = I[p], inarray[i] = u(\xi_{1i}), base[p*nq+i] =
381  * \phi_p(\xi_{1i}) \f$.
382  *
383  * \param base an array defining the local basis for the inner product
384  * usually passed from Basis->GetBdata() or Basis->GetDbdata()
385  * \param inarray: physical point array of function to be integrated
386  * \f$ u(\xi_1) \f$
387  * \param coll_check flag to identify when a Basis->Collocation() call
388  * should be performed to see if this is a GLL_Lagrange basis with a
389  * collocation property. (should be set to 0 if taking the inner
390  * product with respect to the derivative of basis)
391  * \param outarray the values of the inner product with respect to
392  * each basis over region will be stored in the array \a outarray as
393  * output of the function
394  */
395 
397  const Array<OneD, const NekDouble> &base,
398  const Array<OneD, const NekDouble> &inarray,
399  Array<OneD, NekDouble> &outarray,
400  int coll_check)
401  {
402  int nquad = m_base[0]->GetNumPoints();
403  Array<OneD, NekDouble> tmp(nquad);
404  Array<OneD, const NekDouble> w = m_base[0]->GetW();
405 
406  Vmath::Vmul(nquad, inarray, 1, w, 1, tmp, 1);
407 
408  /* Comment below was a bug for collocated basis
409  if(coll_check&&m_base[0]->Collocation())
410  {
411  Vmath::Vcopy(nquad, tmp, 1, outarray, 1);
412  }
413  else
414  {
415  Blas::Dgemv('T',nquad,m_ncoeffs,1.0,base.get(),nquad,
416  &tmp[0],1,0.0,outarray.get(),1);
417  }*/
418 
419  // Correct implementation
420  Blas::Dgemv('T',nquad,m_ncoeffs,1.0,base.get(),nquad,
421  &tmp[0],1,0.0,outarray.get(),1);
422  }
423 
424  /** \brief Inner product of \a inarray over region with respect to the
425  * expansion basis (this)->m_base[0] and return in \a outarray
426  *
427  * Wrapper call to StdSegExp::IProductWRTBase
428  * \param inarray array of function values evaluated at the physical
429  * collocation points
430  * \param outarray the values of the inner product with respect to
431  * each basis over region will be stored in the array \a outarray as
432  * output of the function
433  */
435  Array<OneD, NekDouble> &outarray)
436  {
437  v_IProductWRTBase(m_base[0]->GetBdata(),inarray,outarray,1);
438  }
439 
441  const int dir,
442  const Array<OneD, const NekDouble>& inarray,
443  Array<OneD, NekDouble> & outarray)
444  {
445  ASSERTL1(dir >= 0 && dir < 1,"input dir is out of range");
446  v_IProductWRTBase(m_base[0]->GetDbdata(),inarray,outarray,1);
447  }
448 
450  const Array<OneD, const NekDouble>& inarray,
451  Array<OneD, NekDouble> &outarray,
452  bool multiplybyweights)
453  {
454  int nquad = m_base[0]->GetNumPoints();
455  Array<OneD, NekDouble> tmp(nquad);
456  Array<OneD, const NekDouble> w = m_base[0]->GetW();
457  Array<OneD, const NekDouble> base = m_base[0]->GetBdata();
458 
459  if(multiplybyweights)
460  {
461  Vmath::Vmul(nquad, inarray, 1, w, 1, tmp, 1);
462 
463  Blas::Dgemv('T',nquad,m_ncoeffs,1.0,base.get(),nquad,
464  &tmp[0],1,0.0,outarray.get(),1);
465  }
466  else
467  {
468  Blas::Dgemv('T',nquad,m_ncoeffs,1.0,base.get(),nquad,
469  &inarray[0],1,0.0,outarray.get(),1);
470  }
471  }
472 
473 
474 
475 
476 
477  //----------------------------
478  // Evaluation
479  //----------------------------
480 
481  void StdSegExp::v_FillMode(const int mode, Array<OneD, NekDouble> &outarray)
482  {
483  int nquad = m_base[0]->GetNumPoints();
484  const NekDouble * base = m_base[0]->GetBdata().get();
485 
486  ASSERTL2(mode <= m_ncoeffs,
487  "calling argument mode is larger than total expansion order");
488 
489  Vmath::Vcopy(nquad,(NekDouble *)base+mode*nquad,1, &outarray[0],1);
490  }
491 
492 
494  const Array<OneD, const NekDouble>& coords,
495  const Array<OneD, const NekDouble>& physvals)
496  {
497  return StdExpansion1D::v_PhysEvaluate(coords, physvals);
498  }
499 
501  const Array<OneD, const NekDouble> &inarray,
502  Array<OneD, NekDouble> &outarray,
503  const StdMatrixKey &mkey)
504  {
505  int nquad = m_base[0]->GetNumPoints();
506 
507  Array<OneD,NekDouble> physValues(nquad);
508  Array<OneD,NekDouble> dPhysValuesdx(nquad);
509 
510  v_BwdTrans(inarray,physValues);
511 
512  // Laplacian matrix operation
513  v_PhysDeriv(physValues,dPhysValuesdx);
514  v_IProductWRTBase(m_base[0]->GetDbdata(),dPhysValuesdx,outarray,1);
515  }
516 
517 
519  const Array<OneD, const NekDouble> &inarray,
520  Array<OneD, NekDouble> &outarray,
521  const StdMatrixKey &mkey)
522  {
523  int nquad = m_base[0]->GetNumPoints();
524 
525  Array<OneD,NekDouble> physValues(nquad);
526  Array<OneD,NekDouble> dPhysValuesdx(nquad);
528 
529  v_BwdTrans(inarray,physValues);
530 
531  // mass matrix operation
532  v_IProductWRTBase((m_base[0]->GetBdata()),physValues,wsp,1);
533 
534  // Laplacian matrix operation
535  v_PhysDeriv(physValues,dPhysValuesdx);
536  v_IProductWRTBase(m_base[0]->GetDbdata(),dPhysValuesdx,outarray,1);
537  Blas::Daxpy(m_ncoeffs, mkey.GetConstFactor(eFactorLambda), wsp.get(), 1, outarray.get(), 1);
538  }
539 
540  //up to here
542  const Array<OneD, const NekDouble> &inarray,
543  Array<OneD, NekDouble> &outarray)
544  {
545  int nquad0 = m_base[0]->GetNumPoints();
546 
547  const Array<OneD, const NekDouble>& w0 = m_base[0]->GetW();
548 
549  Vmath::Vmul(nquad0, inarray.get(),1,
550  w0.get(),1,outarray.get(),1);
551  }
552 
554  Array<OneD, NekDouble> &coords_0,
555  Array<OneD, NekDouble> &coords_1,
556  Array<OneD, NekDouble> &coords_2)
557  {
558  Blas::Dcopy(GetNumPoints(0),(m_base[0]->GetZ()).get(),
559  1,&coords_0[0],1);
560  }
561 
562 
563 
564 
565 
566  //---------------------------------------------------------------------
567  // Helper functions
568  //---------------------------------------------------------------------
569 
571  {
572  return 2;
573  }
574 
576  {
577  return 2;
578  }
579 
581  {
582  return 2;
583  }
584 
586  const std::vector<unsigned int> &nummodes,
587  int &modes_offset)
588  {
589  int nmodes = nummodes[modes_offset];
590  modes_offset += 1;
591 
592  return nmodes;
593  }
594 
595  //---------------------------------------------------------------------
596  // Wrapper functions
597  //---------------------------------------------------------------------
598 
600  {
601  DNekMatSharedPtr Mat;
602  MatrixType mattype;
603 
604  switch(mattype = mkey.GetMatrixType())
605  {
606  case eFwdTrans:
607  {
609  StdMatrixKey iprodkey(eIProductWRTBase,v_DetShapeType(),*this);
610  DNekMat &Iprod = *GetStdMatrix(iprodkey);
611  StdMatrixKey imasskey(eInvMass,v_DetShapeType(),*this);
612  DNekMat &Imass = *GetStdMatrix(imasskey);
613 
614  (*Mat) = Imass*Iprod;
615  }
616  break;
617  default:
618  {
620 
621  if(mattype == eMass)
622  {
623  // For Fourier basis set the imaginary component
624  // of mean mode to have a unit diagonal component
625  // in mass matrix
627  {
628  (*Mat)(1,1) = 1.0;
629  }
630  }
631  }
632  break;
633  }
634 
635  return Mat;
636  }
637 
638 
640  {
641  return v_GenMatrix(mkey);
642  }
643 
644  //---------------------------------------------------------------------
645  // Mappings
646  //---------------------------------------------------------------------
647 
648 
650  {
651  if(outarray.num_elements() != NumBndryCoeffs())
652  {
654  }
655  const LibUtilities::BasisType Btype = GetBasisType(0);
656  int nummodes = m_base[0]->GetNumModes();
657 
658  outarray[0] = 0;
659 
660  switch(Btype)
661  {
666  outarray[1]= nummodes-1;
667  break;
670  outarray[1] = 1;
671  break;
672  default:
673  ASSERTL0(0,"Mapping array is not defined for this expansion");
674  break;
675  }
676  }
677 
679  {
680  int i;
681  if(outarray.num_elements()!=GetNcoeffs()-NumBndryCoeffs())
682  {
684  }
685  const LibUtilities::BasisType Btype = GetBasisType(0);
686 
687  switch(Btype)
688  {
693  for(i = 0 ; i < GetNcoeffs()-2;i++)
694  {
695  outarray[i] = i+1;
696  }
697  break;
700  for(i = 0 ; i < GetNcoeffs()-2;i++)
701  {
702  outarray[i] = i+2;
703  }
704  break;
705  default:
706  ASSERTL0(0,"Mapping array is not defined for this expansion");
707  break;
708  }
709  }
710 
711  int StdSegExp::v_GetVertexMap(int localVertexId,bool useCoeffPacking)
712  {
713  ASSERTL0((localVertexId==0)||(localVertexId==1),"local vertex id"
714  "must be between 0 or 1");
715 
716  int localDOF = localVertexId;
717 
719  (localVertexId==1) )
720  {
721  localDOF = m_base[0]->GetNumModes()-1;
722  }
723  return localDOF;
724  }
725 
726 
727  }//end namespace
728 }//end namespace
729 
virtual NekDouble v_Integral(const Array< OneD, const NekDouble > &inarray)
Integrate the physical point list inarray over region and return the value.
Definition: StdSegExp.cpp:120
virtual void v_IProductWRTDerivBase(const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Definition: StdSegExp.cpp:440
NekDouble GetConstFactor(const ConstFactorType &factor) const
Definition: StdMatrixKey.h:122
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:188
virtual void v_GetBoundaryMap(Array< OneD, unsigned int > &outarray)
Definition: StdSegExp.cpp:649
virtual int v_NumDGBndryCoeffs() const
Definition: StdSegExp.cpp:580
virtual void v_StdPhysDeriv(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &out_d0, Array< OneD, NekDouble > &out_d1=NullNekDouble1DArray, Array< OneD, NekDouble > &out_d2=NullNekDouble1DArray)
Definition: StdSegExp.cpp:171
virtual void v_GetCoords(Array< OneD, NekDouble > &coords_0, Array< OneD, NekDouble > &coords_1, Array< OneD, NekDouble > &coords_2)
Definition: StdSegExp.cpp:553
MatrixType GetMatrixType() const
Definition: StdMatrixKey.h:82
void MassMatrixOp(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
Definition: StdExpansion.h:971
static boost::shared_ptr< DataType > AllocateSharedPtr()
Allocate a shared pointer from the memory pool.
virtual void v_BwdTrans(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Backward transform from coefficient space given in inarray and evaluate at the physical quadrature po...
Definition: StdSegExp.cpp:216
virtual void v_BwdTrans_SumFac(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Definition: StdSegExp.cpp:361
virtual int v_NumBndryCoeffs() const
Definition: StdSegExp.cpp:575
int GetNumPoints(const int dir) const
This function returns the number of quadrature points in the dir direction.
Definition: StdExpansion.h:229
virtual void v_FwdTrans(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Forward transform from physical quadrature space stored in inarray and evaluate the expansion coeffic...
Definition: StdSegExp.cpp:265
Principle Modified Functions .
Definition: BasisType.h:49
Lagrange Polynomials using the Gauss points .
Definition: BasisType.h:54
STL namespace.
virtual void v_FillMode(const int mode, Array< OneD, NekDouble > &outarray)
Definition: StdSegExp.cpp:481
StdSegExp()
Default constructor.
Definition: StdSegExp.cpp:47
Fourier Expansion .
Definition: BasisType.h:52
Chebyshev Polynomials .
Definition: BasisType.h:56
LibUtilities::NekManager< StdMatrixKey, DNekBlkMat, StdMatrixKey::opLess > m_stdStaticCondMatrixManager
boost::shared_ptr< DNekMat > DNekMatSharedPtr
Definition: NekTypeDefs.hpp:70
DNekMatSharedPtr GetStdMatrix(const StdMatrixKey &mkey)
Definition: StdExpansion.h:700
virtual bool v_IsBoundaryInteriorExpansion()
Definition: StdSegExp.cpp:87
virtual int v_GetVertexMap(int localVertexId, bool useCoeffPacking=false)
Definition: StdSegExp.cpp:711
virtual int v_CalcNumberOfCoefficients(const std::vector< unsigned int > &nummodes, int &modes_offset)
Definition: StdSegExp.cpp:585
Class representing a segment element in reference space.
Definition: StdSegExp.h:54
The base class for all shapes.
Definition: StdExpansion.h:69
int GetVertexMap(const int localVertexId, bool useCoeffPacking=false)
Definition: StdExpansion.h:826
virtual int v_GetNverts() const
Definition: StdSegExp.cpp:570
virtual DNekMatSharedPtr v_GenMatrix(const StdMatrixKey &mkey)
Definition: StdSegExp.cpp:599
Principle Modified Functions .
Definition: BasisType.h:50
virtual DNekMatSharedPtr v_CreateStdMatrix(const StdMatrixKey &mkey)
Definition: StdSegExp.cpp:639
DNekMatSharedPtr CreateGeneralMatrix(const StdMatrixKey &mkey)
this function generates the mass matrix
double NekDouble
virtual void v_IProductWRTBase(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Inner product of inarray over region with respect to the expansion basis (this)->m_base[0] and return...
Definition: StdSegExp.cpp:434
virtual NekDouble v_PhysEvaluate(const Array< OneD, const NekDouble > &Lcoords, const Array< OneD, const NekDouble > &physvals)
Definition: StdSegExp.cpp:493
LibUtilities::BasisType GetBasisType(const int dir) const
This function returns the type of basis used in the dir direction.
Definition: StdExpansion.h:165
void Vsub(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Subtract vector z = x-y.
Definition: Vmath.cpp:329
void PhysTensorDeriv(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Evaluate the derivative at the physical quadrature points given by inarray and return in outarray...
#define ASSERTL2(condition, msg)
Assert Level 2 – Debugging which is used FULLDEBUG compilation mode. This level assert is designed t...
Definition: ErrorUtil.hpp:240
virtual void v_HelmholtzMatrixOp(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
Definition: StdSegExp.cpp:518
virtual void v_FwdTrans_BndConstrained(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Definition: StdSegExp.cpp:287
int GetNcoeffs(void) const
This function returns the total number of coefficients used in the expansion.
Definition: StdExpansion.h:131
virtual void v_MultiplyByStdQuadratureMetric(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Definition: StdSegExp.cpp:541
virtual LibUtilities::ShapeType v_DetShapeType() const
Return Shape of region, using ShapeType enum list. i.e. Segment.
Definition: StdSegExp.cpp:82
Lagrange for SEM basis .
Definition: BasisType.h:53
T Vsum(int n, const T *x, const int incx)
Subtract return sum(x)
Definition: Vmath.cpp:723
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode...
Definition: ErrorUtil.hpp:218
Array< OneD, LibUtilities::BasisSharedPtr > m_base
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1047
virtual void v_LaplacianMatrixOp(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
Definition: StdSegExp.cpp:500
Describes the specification for a Basis.
Definition: Basis.h:50
virtual NekDouble v_PhysEvaluate(const Array< OneD, const NekDouble > &coords, const Array< OneD, const NekDouble > &physvals)
virtual void v_PhysDeriv(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &out_d0, Array< OneD, NekDouble > &out_d1=NullNekDouble1DArray, Array< OneD, NekDouble > &out_d2=NullNekDouble1DArray)
Evaluate the derivative at the physical quadrature points given by inarray and return in outarray...
Definition: StdSegExp.cpp:153
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
virtual void v_IProductWRTBase_SumFac(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, bool multiplybyweights=true)
Definition: StdSegExp.cpp:449
virtual void v_GetInteriorMap(Array< OneD, unsigned int > &outarray)
Definition: StdSegExp.cpp:678