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