Nektar++
ExpList.cpp
Go to the documentation of this file.
1 ///////////////////////////////////////////////////////////////////////////////
2 //
3 // File ExpList.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: Expansion list definition
33 //
34 ///////////////////////////////////////////////////////////////////////////////
35 
36 #include <MultiRegions/ExpList.h>
39 
40 #include <StdRegions/StdSegExp.h>
41 #include <StdRegions/StdTriExp.h>
43 #include <StdRegions/StdQuadExp.h>
44 #include <StdRegions/StdTetExp.h>
45 #include <StdRegions/StdPyrExp.h>
46 #include <StdRegions/StdPrismExp.h>
47 #include <StdRegions/StdHexExp.h>
48 
49 #include <LocalRegions/MatrixKey.h> // for MatrixKey
50 #include <LocalRegions/Expansion.h> // for Expansion
51 
52 #include <MultiRegions/AssemblyMap/AssemblyMapCG.h> // for AssemblyMapCG, etc
53 #include <MultiRegions/AssemblyMap/AssemblyMapDG.h> // for AssemblyMapCG, etc
54 #include <MultiRegions/GlobalLinSysKey.h> // for GlobalLinSysKey
55 #include <MultiRegions/GlobalMatrix.h> // for GlobalMatrix, etc
56 #include <MultiRegions/GlobalMatrixKey.h> // for GlobalMatrixKey
57 
61 
62 #include <Collections/CollectionOptimisation.h>
63 #include <Collections/Operator.h>
64 
65 namespace Nektar
66 {
67  namespace MultiRegions
68  {
69  /**
70  * @class ExpList
71  * All multi-elemental expansions \f$u^{\delta}(\boldsymbol{x})\f$ can
72  * be considered as the assembly of the various elemental contributions.
73  * On a discrete level, this yields,
74  * \f[u^{\delta}(\boldsymbol{x}_i)=\sum_{e=1}^{{N_{\mathrm{el}}}}
75  * \sum_{n=0}^{N^{e}_m-1}\hat{u}_n^e\phi_n^e(\boldsymbol{x}_i).\f]
76  * where \f${N_{\mathrm{el}}}\f$ is the number of elements and
77  * \f$N^{e}_m\f$ is the local elemental number of expansion modes.
78  * As it is the lowest level class, it contains the definition of the
79  * common data and common routines to all multi-elemental expansions.
80  *
81  * The class stores a vector of expansions, \a m_exp, (each derived from
82  * StdRegions#StdExpansion) which define the constituent components of
83  * the domain. The coefficients from these expansions are concatenated
84  * in \a m_coeffs, while the expansion evaluated at the quadrature
85  * points is stored in \a m_phys.
86  */
87 
88  /**
89  * Creates an empty expansion list. The expansion list will typically be
90  * populated by a derived class (namely one of MultiRegions#ExpList1D,
91  * MultiRegions#ExpList2D or MultiRegions#ExpList3D).
92  */
94  m_comm(),
95  m_session(),
96  m_graph(),
97  m_ncoeffs(0),
98  m_npoints(0),
99  m_coeffs(),
100  m_phys(),
101  m_physState(false),
102  m_exp(MemoryManager<LocalRegions::ExpansionVector>
103  ::AllocateSharedPtr()),
104  m_coeff_offset(),
105  m_phys_offset(),
106  m_offset_elmt_id(),
107  m_blockMat(MemoryManager<BlockMatrixMap>::AllocateSharedPtr()),
108  m_WaveSpace(false)
109  {
111  }
112 
113 
114  /**
115  * Creates an empty expansion list. The expansion list will typically be
116  * populated by a derived class (namely one of MultiRegions#ExpList1D,
117  * MultiRegions#ExpList2D or MultiRegions#ExpList3D).
118  */
120  const LibUtilities::SessionReaderSharedPtr &pSession):
121  m_comm(pSession->GetComm()),
122  m_session(pSession),
123  m_graph(),
124  m_ncoeffs(0),
125  m_npoints(0),
126  m_coeffs(),
127  m_phys(),
128  m_physState(false),
129  m_exp(MemoryManager<LocalRegions::ExpansionVector>
130  ::AllocateSharedPtr()),
131  m_coeff_offset(),
132  m_phys_offset(),
133  m_offset_elmt_id(),
134  m_blockMat(MemoryManager<BlockMatrixMap>::AllocateSharedPtr()),
135  m_WaveSpace(false)
136  {
138  }
139 
140 
141  /**
142  * Creates an empty expansion list. The expansion list will typically be
143  * populated by a derived class (namely one of MultiRegions#ExpList1D,
144  * MultiRegions#ExpList2D or MultiRegions#ExpList3D).
145  */
147  const LibUtilities::SessionReaderSharedPtr &pSession,
148  const SpatialDomains::MeshGraphSharedPtr &pGraph):
149  m_comm(pSession->GetComm()),
150  m_session(pSession),
151  m_graph(pGraph),
152  m_ncoeffs(0),
153  m_npoints(0),
154  m_coeffs(),
155  m_phys(),
156  m_physState(false),
157  m_exp(MemoryManager<LocalRegions::ExpansionVector>
158  ::AllocateSharedPtr()),
159  m_coeff_offset(),
160  m_phys_offset(),
161  m_offset_elmt_id(),
162  m_blockMat(MemoryManager<BlockMatrixMap>::AllocateSharedPtr()),
163  m_WaveSpace(false)
164  {
166  }
167 
168 
169  /**
170  * Copies an existing expansion list.
171  * @param in Source expansion list.
172  */
173  ExpList::ExpList(const ExpList &in, const bool DeclareCoeffPhysArrays):
174  m_comm(in.m_comm),
175  m_session(in.m_session),
176  m_graph(in.m_graph),
177  m_ncoeffs(in.m_ncoeffs),
178  m_npoints(in.m_npoints),
179  m_physState(false),
180  m_exp(in.m_exp),
181  m_collections(in.m_collections),
182  m_coll_coeff_offset(in.m_coll_coeff_offset),
183  m_coll_phys_offset(in.m_coll_phys_offset),
184  m_coeff_offset(in.m_coeff_offset),
185  m_phys_offset(in.m_phys_offset),
186  m_offset_elmt_id(in.m_offset_elmt_id),
187  m_globalOptParam(in.m_globalOptParam),
188  m_blockMat(in.m_blockMat),
189  m_WaveSpace(false)
190  {
192 
193  if(DeclareCoeffPhysArrays)
194  {
197  }
198  }
199 
200  /**
201  *
202  */
204  {
205  return m_expType;
206  }
207 
208  /**
209  *
210  */
212  {
213  m_expType = Type;
214  }
215 
217  {
218  }
219 
220  /**
221  * The integration is evaluated locally, that is
222  * \f[\int
223  * f(\boldsymbol{x})d\boldsymbol{x}=\sum_{e=1}^{{N_{\mathrm{el}}}}
224  * \left\{\int_{\Omega_e}f(\boldsymbol{x})d\boldsymbol{x}\right\}, \f]
225  * where the integration over the separate elements is done by the
226  * function StdRegions#StdExpansion#Integral, which discretely
227  * evaluates the integral using Gaussian quadrature.
228  *
229  * Note that the array #m_phys should be filled with the values of the
230  * function \f$f(\boldsymbol{x})\f$ at the quadrature points
231  * \f$\boldsymbol{x}_i\f$.
232  *
233  * @return The value of the discretely evaluated integral
234  * \f$\int f(\boldsymbol{x})d\boldsymbol{x}\f$.
235  */
237  {
238  ASSERTL1(m_physState == true,
239  "local physical space is not true ");
240 
241  return PhysIntegral(m_phys);
242  }
243 
244 
245  /**
246  * The integration is evaluated locally, that is
247  * \f[\int
248  * f(\boldsymbol{x})d\boldsymbol{x}=\sum_{e=1}^{{N_{\mathrm{el}}}}
249  * \left\{\int_{\Omega_e}f(\boldsymbol{x})d\boldsymbol{x}\right\}, \f]
250  * where the integration over the separate elements is done by the
251  * function StdRegions#StdExpansion#Integral, which discretely
252  * evaluates the integral using Gaussian quadrature.
253  *
254  * @param inarray An array of size \f$Q_{\mathrm{tot}}\f$
255  * containing the values of the function
256  * \f$f(\boldsymbol{x})\f$ at the quadrature
257  * points \f$\boldsymbol{x}_i\f$.
258  * @return The value of the discretely evaluated integral
259  * \f$\int f(\boldsymbol{x})d\boldsymbol{x}\f$.
260  */
262  const Array<OneD, const NekDouble> &inarray)
263  {
264  int i;
265  NekDouble sum = 0.0;
266 
267  for(i = 0; i < (*m_exp).size(); ++i)
268  {
269  sum += (*m_exp)[i]->Integral(inarray + m_phys_offset[i]);
270  }
271 
272  return sum;
273  }
274 
275 
276  /**
277  * Retrieves the block matrix specified by \a bkey, and computes
278  * \f$ y=Mx \f$.
279  * @param gkey GlobalMatrixKey specifying the block matrix to
280  * use in the matrix-vector multiply.
281  * @param inarray Input vector \f$ x \f$.
282  * @param outarray Output vector \f$ y \f$.
283  */
285  const GlobalMatrixKey &gkey,
286  const Array<OneD,const NekDouble> &inarray,
287  Array<OneD, NekDouble> &outarray)
288  {
289  // Retrieve the block matrix using the given key.
290  const DNekScalBlkMatSharedPtr& blockmat = GetBlockMatrix(gkey);
291  int nrows = blockmat->GetRows();
292  int ncols = blockmat->GetColumns();
293 
294  // Create NekVectors from the given data arrays
295  NekVector<NekDouble> in (ncols,inarray, eWrapper);
296  NekVector< NekDouble> out(nrows,outarray,eWrapper);
297 
298  // Perform matrix-vector multiply.
299  out = (*blockmat)*in;
300  }
301 
302 
303  /**
304  * The operation is evaluated locally for every element by the function
305  * StdRegions#StdExpansion#IProductWRTBase.
306  *
307  * @param inarray An array of size \f$Q_{\mathrm{tot}}\f$
308  * containing the values of the function
309  * \f$f(\boldsymbol{x})\f$ at the quadrature
310  * points \f$\boldsymbol{x}_i\f$.
311  * @param outarray An array of size \f$N_{\mathrm{eof}}\f$
312  * used to store the result.
313  */
315  const Array<OneD, const NekDouble> &inarray,
316  Array<OneD, NekDouble> &outarray)
317  {
319  for (int i = 0; i < m_collections.size(); ++i)
320  {
321 
323  inarray + m_coll_phys_offset[i],
324  tmp = outarray + m_coll_coeff_offset[i]);
325  }
326  }
327 
328  /**
329  * The operation is evaluated locally for every element by the function
330  * StdRegions#StdExpansion#IProductWRTDerivBase.
331  *
332  * @param dir {0,1} is the direction in which the
333  * derivative of the basis should be taken
334  * @param inarray An array of size \f$Q_{\mathrm{tot}}\f$
335  * containing the values of the function
336  * \f$f(\boldsymbol{x})\f$ at the quadrature
337  * points \f$\boldsymbol{x}_i\f$.
338  * @param outarray An array of size \f$N_{\mathrm{eof}}\f$
339  * used to store the result.
340  */
341  void ExpList::IProductWRTDerivBase(const int dir,
342  const Array<OneD, const NekDouble> &inarray,
343  Array<OneD, NekDouble> &outarray)
344  {
345  int i;
346 
347  Array<OneD,NekDouble> e_outarray;
348 
349  for(i = 0; i < (*m_exp).size(); ++i)
350  {
351  (*m_exp)[i]->IProductWRTDerivBase(dir,inarray+m_phys_offset[i],
352  e_outarray = outarray+m_coeff_offset[i]);
353  }
354  }
355 
356 
357  /**
358  * The operation is evaluated locally for every element by the function
359  * StdRegions#StdExpansion#IProductWRTDerivBase.
360  *
361  * @param inarray An array of arrays of size \f$Q_{\mathrm{tot}}\f$
362  * containing the values of the function
363  * \f$f(\boldsymbol{x})\f$ at the quadrature
364  * points \f$\boldsymbol{x}_i\f$ in dir directions.
365  * @param outarray An array of size \f$N_{\mathrm{eof}}\f$
366  * used to store the result.
367  */
369  Array<OneD, NekDouble> &outarray)
370  {
371  Array<OneD, NekDouble> tmp0,tmp1,tmp2;
372  // assume coord dimension defines the size of Deriv Base
373  int dim = GetCoordim(0);
374 
375  ASSERTL1(inarray.num_elements() >= dim,"inarray is not of sufficient dimension");
376 
377  switch(dim)
378  {
379  case 1:
380  for (int i = 0; i < m_collections.size(); ++i)
381  {
382  m_collections[i].ApplyOperator(
384  inarray[0] + m_coll_phys_offset[i],
385  tmp0 = outarray + m_coll_coeff_offset[i]);
386  }
387  break;
388  case 2:
389  for (int i = 0; i < m_collections.size(); ++i)
390  {
391  m_collections[i].ApplyOperator(
393  inarray[0] + m_coll_phys_offset[i],
394  tmp0 = inarray[1] + m_coll_phys_offset[i],
395  tmp1 = outarray + m_coll_coeff_offset[i]);
396  }
397  break;
398  case 3:
399  for (int i = 0; i < m_collections.size(); ++i)
400  {
401  m_collections[i].ApplyOperator(
403  inarray[0] + m_coll_phys_offset[i],
404  tmp0 = inarray[1] + m_coll_phys_offset[i],
405  tmp1 = inarray[2] + m_coll_phys_offset[i],
406  tmp2 = outarray + m_coll_coeff_offset[i]);
407  }
408  break;
409  default:
410  ASSERTL0(false,"Dimension of inarray not correct");
411  break;
412  }
413  }
414  /**
415  * Given a function \f$f(\boldsymbol{x})\f$ evaluated at
416  * the quadrature points, this function calculates the
417  * derivatives \f$\frac{d}{dx_1}\f$, \f$\frac{d}{dx_2}\f$
418  * and \f$\frac{d}{dx_3}\f$ of the function
419  * \f$f(\boldsymbol{x})\f$ at the same quadrature
420  * points. The local distribution of the quadrature points
421  * allows an elemental evaluation of the derivative. This
422  * is done by a call to the function
423  * StdRegions#StdExpansion#PhysDeriv.
424  *
425  * @param inarray An array of size \f$Q_{\mathrm{tot}}\f$
426  * containing the values of the function
427  * \f$f(\boldsymbol{x})\f$ at the quadrature
428  * points \f$\boldsymbol{x}_i\f$.
429  * @param out_d0 The discrete evaluation of the
430  * derivative\f$\frac{d}{dx_1}\f$ will
431  * be stored in this array of size
432  * \f$Q_{\mathrm{tot}}\f$.
433  * @param out_d1 The discrete evaluation of the
434  * derivative\f$\frac{d}{dx_2}\f$ will be
435  * stored in this array of size
436  * \f$Q_{\mathrm{tot}}\f$. Note that if no
437  * memory is allocated for \a out_d1, the
438  * derivative \f$\frac{d}{dx_2}\f$ will not be
439  * calculated.
440  * @param out_d2 The discrete evaluation of the
441  * derivative\f$\frac{d}{dx_3}\f$ will be
442  * stored in this array of size
443  * \f$Q_{\mathrm{tot}}\f$. Note that if no
444  * memory is allocated for \a out_d2, the
445  * derivative \f$\frac{d}{dx_3}\f$ will not be
446  * calculated.
447  */
449  Array<OneD, NekDouble> &out_d0,
450  Array<OneD, NekDouble> &out_d1,
451  Array<OneD, NekDouble> &out_d2)
452  {
453  Array<OneD, NekDouble> e_out_d0;
454  Array<OneD, NekDouble> e_out_d1;
455  Array<OneD, NekDouble> e_out_d2;
456  for (int i = 0; i < m_collections.size(); ++i)
457  {
458  int offset = m_coll_phys_offset[i];
459  e_out_d0 = out_d0 + offset;
460  e_out_d1 = out_d1 + offset;
461  e_out_d2 = out_d2 + offset;
462 
463  m_collections[i].ApplyOperator(Collections::ePhysDeriv,
464  inarray + offset,
465  e_out_d0,e_out_d1, e_out_d2);
466 
467  }
468  }
469 
470  void ExpList::v_PhysDeriv(const int dir,
471  const Array<OneD, const NekDouble> &inarray,
472  Array<OneD, NekDouble> &out_d)
473  {
474  Direction edir = DirCartesianMap[dir];
475  v_PhysDeriv(edir, inarray,out_d);
476  }
477 
479  Array<OneD, NekDouble> &out_d)
480  {
481  int i;
482  if(edir==MultiRegions::eS)
483  {
484  Array<OneD, NekDouble> e_out_ds;
485  for(i=0; i<(*m_exp).size(); ++i)
486  {
487  e_out_ds = out_d + m_phys_offset[i];
488  (*m_exp)[i]->PhysDeriv_s(inarray+m_phys_offset[i],e_out_ds);
489  }
490  }
491  else if(edir==MultiRegions::eN)
492  {
493  Array<OneD, NekDouble > e_out_dn;
494  for(i=0; i<(*m_exp).size(); i++)
495  {
496  e_out_dn = out_d +m_phys_offset[i];
497  (*m_exp)[i]->PhysDeriv_n(inarray+m_phys_offset[i],e_out_dn);
498  }
499  }
500  else
501  {
502  // convert enum into int
503  int intdir= (int)edir;
504  Array<OneD, NekDouble> e_out_d;
505  for(i= 0; i < (*m_exp).size(); ++i)
506  {
507  e_out_d = out_d + m_phys_offset[i];
508  (*m_exp)[i]->PhysDeriv(intdir, inarray+m_phys_offset[i], e_out_d);
509  }
510 
511  }
512  }
513 
514 
515  /**
516  * The coefficients of the function to be acted upon
517  * should be contained in the \param inarray. The
518  * resulting coefficients are stored in \param outarray
519  *
520  * @param inarray An array of size \f$N_{\mathrm{eof}}\f$
521  * containing the inner product.
522  */
524  const Array<OneD, const NekDouble> &inarray,
525  Array<OneD, NekDouble> &outarray)
526  {
528  const DNekScalBlkMatSharedPtr& InvMass = GetBlockMatrix(mkey);
529 
530  // Inverse mass matrix
532  if(inarray.get() == outarray.get())
533  {
534  NekVector<NekDouble> in(m_ncoeffs,inarray); // copy data
535  out = (*InvMass)*in;
536  }
537  else
538  {
540  out = (*InvMass)*in;
541  }
542  }
543 
544  /**
545  * Given a function \f$u(\boldsymbol{x})\f$ defined at the
546  * quadrature points, this function determines the
547  * transformed elemental coefficients \f$\hat{u}_n^e\f$
548  * employing a discrete elemental Galerkin projection from
549  * physical space to coefficient space. For each element,
550  * the operation is evaluated locally by the function
551  * StdRegions#StdExpansion#IproductWRTBase followed by a
552  * call to #MultiRegions#MultiplyByElmtInvMass.
553  *
554  * @param inarray An array of size \f$Q_{\mathrm{tot}}\f$
555  * containing the values of the function
556  * \f$f(\boldsymbol{x})\f$ at the quadrature
557  * points \f$\boldsymbol{x}_i\f$.
558  * @param outarray The resulting coefficients
559  * \f$\hat{u}_n^e\f$ will be stored in this
560  * array of size \f$N_{\mathrm{eof}}\f$.
561  */
563  Array<OneD, NekDouble> &outarray)
564  {
566 
567  IProductWRTBase_IterPerExp(inarray,f);
568  MultiplyByElmtInvMass(f,outarray);
569 
570  }
571 
573  const Array<OneD, const NekDouble>& inarray,
574  Array<OneD, NekDouble> &outarray)
575  {
576  int i;
577 
578  Array<OneD,NekDouble> e_outarray;
579 
580  for(i= 0; i < (*m_exp).size(); ++i)
581  {
582  (*m_exp)[i]->FwdTrans_BndConstrained(inarray+m_phys_offset[i],
583  e_outarray = outarray+m_coeff_offset[i]);
584  }
585  }
586 
587  /**
588  * This function smooth a field after some calculaitons which have
589  * been done elementally.
590  *
591  * @param field An array containing the field in physical space
592  *
593  */
595  {
596  // Do nothing unless the method is implemented in the appropriate
597  // class, i.e. ContField1D,ContField2D, etc.
598 
599  // So far it has been implemented just for ContField2D and
600  // ContField3DHomogeneous1D
601 
602  // Block in case users try the smoothing with a modal expansion.
603  // Maybe a different techique for the smoothing require
604  // implementation for modal basis.
605 
606  ASSERTL0((*m_exp)[0]->GetBasisType(0)
608  (*m_exp)[0]->GetBasisType(0)
610  "Smoothing is currently not allowed unless you are using "
611  "a nodal base for efficiency reasons. The implemented "
612  "smoothing technique requires the mass matrix inversion "
613  "which is trivial just for GLL_LAGRANGE_SEM and "
614  "GAUSS_LAGRANGE_SEMexpansions.");
615  }
616 
617 
618  /**
619  * This function assembles the block diagonal matrix
620  * \f$\underline{\boldsymbol{M}}^e\f$, which is the
621  * concatenation of the local matrices
622  * \f$\boldsymbol{M}^e\f$ of the type \a mtype, that is
623  *
624  * \f[
625  * \underline{\boldsymbol{M}}^e = \left[
626  * \begin{array}{cccc}
627  * \boldsymbol{M}^1 & 0 & \hspace{3mm}0 \hspace{3mm}& 0 \\
628  * 0 & \boldsymbol{M}^2 & 0 & 0 \\
629  * 0 & 0 & \ddots & 0 \\
630  * 0 & 0 & 0 & \boldsymbol{M}^{N_{\mathrm{el}}} \end{array}\right].\f]
631  *
632  * @param mtype the type of matrix to be assembled
633  * @param scalar an optional parameter
634  * @param constant an optional parameter
635  */
637  const GlobalMatrixKey &gkey)
638  {
639  int i,cnt1;
640  int n_exp = 0;
641  DNekScalMatSharedPtr loc_mat;
642  DNekScalBlkMatSharedPtr BlkMatrix;
643  map<int,int> elmt_id;
645 
646  if(ShapeType != LibUtilities::eNoShapeType)
647  {
648  for(i = 0 ; i < (*m_exp).size(); ++i)
649  {
650  if((*m_exp)[m_offset_elmt_id[i]]->DetShapeType()
651  == ShapeType)
652  {
653  elmt_id[n_exp++] = m_offset_elmt_id[i];
654  }
655  }
656  }
657  else
658  {
659  n_exp = (*m_exp).size();
660  for(i = 0; i < n_exp; ++i)
661  {
662  elmt_id[i] = m_offset_elmt_id[i];
663  }
664  }
665 
666  Array<OneD,unsigned int> nrows(n_exp);
667  Array<OneD,unsigned int> ncols(n_exp);
668 
669  switch(gkey.GetMatrixType())
670  {
672  {
673  // set up an array of integers for block matrix construction
674  for(i = 0; i < n_exp; ++i)
675  {
676  nrows[i] = (*m_exp)[elmt_id.find(i)->second]->GetTotPoints();
677  ncols[i] = (*m_exp)[elmt_id.find(i)->second]->GetNcoeffs();
678  }
679  }
680  break;
682  {
683  // set up an array of integers for block matrix construction
684  for(i = 0; i < n_exp; ++i)
685  {
686  nrows[i] = (*m_exp)[elmt_id.find(i)->second]->GetNcoeffs();
687  ncols[i] = (*m_exp)[elmt_id.find(i)->second]->GetTotPoints();
688  }
689  }
690  break;
691  case StdRegions::eMass:
696  {
697  // set up an array of integers for block matrix construction
698  for(i = 0; i < n_exp; ++i)
699  {
700  nrows[i] = (*m_exp)[elmt_id.find(i)->second]->GetNcoeffs();
701  ncols[i] = (*m_exp)[elmt_id.find(i)->second]->GetNcoeffs();
702  }
703  }
704  break;
705 
707  {
708  // set up an array of integers for block matrix construction
709  for(i = 0; i < n_exp; ++i)
710  {
711  nrows[i] = (*m_exp)[elmt_id.find(i)->second]->GetNcoeffs();
712  ncols[i] = (*m_exp)[elmt_id.find(i)->second]->NumDGBndryCoeffs();
713  }
714  }
715  break;
716 
718  {
719  // set up an array of integers for block matrix construction
720  for(i = 0; i < n_exp; ++i)
721  {
722  nrows[i] = (*m_exp)[elmt_id.find(i)->second]->NumDGBndryCoeffs();
723  ncols[i] = (*m_exp)[elmt_id.find(i)->second]->NumDGBndryCoeffs();
724  }
725  }
726  break;
727 
728  default:
729  {
731  "Global Matrix creation not defined for this type "
732  "of matrix");
733  }
734  }
735 
736  MatrixStorage blkmatStorage = eDIAGONAL;
738  ::AllocateSharedPtr(nrows,ncols,blkmatStorage);
739 
740  int nvarcoeffs = gkey.GetNVarCoeffs();
741  int eid;
742  Array<OneD, NekDouble> varcoeffs_wk;
743 
744  for(i = cnt1 = 0; i < n_exp; ++i)
745  {
746  // need to be initialised with zero size for non variable coefficient case
747  StdRegions::VarCoeffMap varcoeffs;
748 
749  eid = elmt_id[i];
750  if(nvarcoeffs>0)
751  {
752  StdRegions::VarCoeffMap::const_iterator x;
753  for (x = gkey.GetVarCoeffs().begin(); x != gkey.GetVarCoeffs().end(); ++x)
754  {
755  varcoeffs[x->first] = x->second + m_phys_offset[eid];
756  }
757  }
758 
760  (*m_exp)[eid]->DetShapeType(),
761  *(*m_exp)[eid],
762  gkey.GetConstFactors(),
763  varcoeffs );
764 
765  loc_mat = boost::dynamic_pointer_cast<LocalRegions::Expansion>((*m_exp)[elmt_id.find(i)->second])->GetLocMatrix(matkey);
766  BlkMatrix->SetBlock(i,i,loc_mat);
767  }
768 
769  return BlkMatrix;
770  }
771 
773  const GlobalMatrixKey &gkey)
774  {
775  BlockMatrixMap::iterator matrixIter = m_blockMat->find(gkey);
776 
777  if(matrixIter == m_blockMat->end())
778  {
779  return ((*m_blockMat)[gkey] = GenBlockMatrix(gkey));
780  }
781  else
782  {
783  return matrixIter->second;
784  }
785  }
786 
788  const GlobalMatrixKey &gkey,
789  const Array<OneD,const NekDouble> &inarray,
790  Array<OneD, NekDouble> &outarray)
791  {
792  const Array<OneD, const bool> doBlockMatOp
793  = m_globalOptParam->DoBlockMatOp(gkey.GetMatrixType());
794  const Array<OneD, const int> num_elmts
795  = m_globalOptParam->GetShapeNumElements();
796 
797  Array<OneD,NekDouble> tmp_outarray;
798  int cnt = 0;
799  int eid;
800  for(int n = 0; n < num_elmts.num_elements(); ++n)
801  {
802  if(doBlockMatOp[n])
803  {
804  const LibUtilities::ShapeType vType
805  = m_globalOptParam->GetShapeList()[n];
806  const MultiRegions::GlobalMatrixKey vKey(gkey, vType);
807  if (cnt < m_offset_elmt_id.num_elements())
808  {
809  eid = m_offset_elmt_id[cnt];
810  MultiplyByBlockMatrix(vKey,inarray + m_coeff_offset[eid],
811  tmp_outarray = outarray + m_coeff_offset[eid]);
812  cnt += num_elmts[n];
813  }
814  }
815  else
816  {
817  int i;
818  int nvarcoeffs = gkey.GetNVarCoeffs();
819 
820  for(i= 0; i < num_elmts[n]; ++i)
821  {
822  // need to be initialised with zero size for non variable coefficient case
823  StdRegions::VarCoeffMap varcoeffs;
824 
825  eid = m_offset_elmt_id[cnt++];
826  if(nvarcoeffs>0)
827  {
828  StdRegions::VarCoeffMap::const_iterator x;
829  for (x = gkey.GetVarCoeffs().begin(); x != gkey.GetVarCoeffs().end(); ++x)
830  {
831  varcoeffs[x->first] = x->second + m_phys_offset[eid];
832  }
833  }
834 
836  (*m_exp)[eid]->DetShapeType(),
837  *((*m_exp)[eid]),
838  gkey.GetConstFactors(),varcoeffs);
839 
840  (*m_exp)[eid]->GeneralMatrixOp(inarray + m_coeff_offset[eid],
841  tmp_outarray = outarray+m_coeff_offset[eid],
842  mkey);
843  }
844  }
845  }
846  }
847 
848  /**
849  * Retrieves local matrices from each expansion in the expansion list
850  * and combines them together to generate a global matrix system.
851  * @param mkey Matrix key for the matrix to be generated.
852  * @param locToGloMap Local to global mapping.
853  * @returns Shared pointer to the generated global matrix.
854  */
856  const GlobalMatrixKey &mkey,
857  const AssemblyMapCGSharedPtr &locToGloMap)
858  {
859  int i,j,n,gid1,gid2,cntdim1,cntdim2;
860  NekDouble sign1,sign2;
861  DNekScalMatSharedPtr loc_mat;
862 
863  unsigned int glob_rows;
864  unsigned int glob_cols;
865  unsigned int loc_rows;
866  unsigned int loc_cols;
867 
868  bool assembleFirstDim;
869  bool assembleSecondDim;
870 
871  switch(mkey.GetMatrixType())
872  {
874  {
875  glob_rows = m_npoints;
876  glob_cols = locToGloMap->GetNumGlobalCoeffs();
877 
878  assembleFirstDim = false;
879  assembleSecondDim = true;
880  }
881  break;
883  {
884  glob_rows = locToGloMap->GetNumGlobalCoeffs();
885  glob_cols = m_npoints;
886 
887  assembleFirstDim = true;
888  assembleSecondDim = false;
889  }
890  break;
891  case StdRegions::eMass:
895  {
896  glob_rows = locToGloMap->GetNumGlobalCoeffs();
897  glob_cols = locToGloMap->GetNumGlobalCoeffs();
898 
899  assembleFirstDim = true;
900  assembleSecondDim = true;
901  }
902  break;
903  default:
904  {
906  "Global Matrix creation not defined for this type "
907  "of matrix");
908  }
909  }
910 
911  COOMatType spcoomat;
912  CoordType coord;
913 
914  int nvarcoeffs = mkey.GetNVarCoeffs();
915  int eid;
916 
917  // fill global matrix
918  for(n = cntdim1 = cntdim2 = 0; n < (*m_exp).size(); ++n)
919  {
920  // need to be initialised with zero size for non variable coefficient case
921  StdRegions::VarCoeffMap varcoeffs;
922 
923  eid = m_offset_elmt_id[n];
924  if(nvarcoeffs>0)
925  {
926  StdRegions::VarCoeffMap::const_iterator x;
927  for (x = mkey.GetVarCoeffs().begin(); x != mkey.GetVarCoeffs().end(); ++x)
928  {
929  varcoeffs[x->first] = x->second + m_phys_offset[eid];
930  }
931  }
932 
934  (*m_exp)[eid]->DetShapeType(),
935  *((*m_exp)[eid]),
936  mkey.GetConstFactors(),varcoeffs);
937 
938  loc_mat = boost::dynamic_pointer_cast<LocalRegions::Expansion>((*m_exp)[m_offset_elmt_id[n]])->GetLocMatrix(matkey);
939 
940  loc_rows = loc_mat->GetRows();
941  loc_cols = loc_mat->GetColumns();
942 
943  for(i = 0; i < loc_rows; ++i)
944  {
945  if(assembleFirstDim)
946  {
947  gid1 = locToGloMap->GetLocalToGlobalMap (cntdim1 + i);
948  sign1 = locToGloMap->GetLocalToGlobalSign(cntdim1 + i);
949  }
950  else
951  {
952  gid1 = cntdim1 + i;
953  sign1 = 1.0;
954  }
955 
956  for(j = 0; j < loc_cols; ++j)
957  {
958  if(assembleSecondDim)
959  {
960  gid2 = locToGloMap
961  ->GetLocalToGlobalMap(cntdim2 + j);
962  sign2 = locToGloMap
963  ->GetLocalToGlobalSign(cntdim2 + j);
964  }
965  else
966  {
967  gid2 = cntdim2 + j;
968  sign2 = 1.0;
969  }
970 
971  // sparse matrix fill
972  coord = make_pair(gid1,gid2);
973  if( spcoomat.count(coord) == 0 )
974  {
975  spcoomat[coord] = sign1*sign2*(*loc_mat)(i,j);
976  }
977  else
978  {
979  spcoomat[coord] += sign1*sign2*(*loc_mat)(i,j);
980  }
981  }
982  }
983  cntdim1 += loc_rows;
984  cntdim2 += loc_cols;
985  }
986 
988  ::AllocateSharedPtr(m_session,glob_rows,glob_cols,spcoomat);
989  }
990 
991 
993  {
994  int i,j,n,gid1,gid2,loc_lda,eid;
995  NekDouble sign1,sign2,value;
996  DNekScalMatSharedPtr loc_mat;
997 
998  int totDofs = locToGloMap->GetNumGlobalCoeffs();
999  int NumDirBCs = locToGloMap->GetNumGlobalDirBndCoeffs();
1000 
1001  unsigned int rows = totDofs - NumDirBCs;
1002  unsigned int cols = totDofs - NumDirBCs;
1003  NekDouble zero = 0.0;
1004 
1005  DNekMatSharedPtr Gmat;
1006  int bwidth = locToGloMap->GetFullSystemBandWidth();
1007 
1008  int nvarcoeffs = mkey.GetNVarCoeffs();
1009  MatrixStorage matStorage;
1010 
1011  map<int, RobinBCInfoSharedPtr> RobinBCInfo = GetRobinBCInfo();
1012 
1013  switch(mkey.GetMatrixType())
1014  {
1015  // case for all symmetric matices
1018  if( (2*(bwidth+1)) < rows)
1019  {
1021  Gmat = MemoryManager<DNekMat>::AllocateSharedPtr(rows,cols,zero,matStorage,bwidth,bwidth);
1022  }
1023  else
1024  {
1025  matStorage = ePOSITIVE_DEFINITE_SYMMETRIC;
1026  Gmat = MemoryManager<DNekMat>::AllocateSharedPtr(rows,cols,zero,matStorage);
1027  }
1028 
1029  break;
1030  default: // Assume general matrix - currently only set up for full invert
1031  {
1032  matStorage = eFULL;
1033  Gmat = MemoryManager<DNekMat>::AllocateSharedPtr(rows,cols,zero,matStorage);
1034  }
1035  }
1036 
1037  // fill global symmetric matrix
1038  for(n = 0; n < (*m_exp).size(); ++n)
1039  {
1040  // need to be initialised with zero size for non variable coefficient case
1041  StdRegions::VarCoeffMap varcoeffs;
1042 
1043  eid = m_offset_elmt_id[n];
1044  if(nvarcoeffs>0)
1045  {
1046  StdRegions::VarCoeffMap::const_iterator x;
1047  for (x = mkey.GetVarCoeffs().begin(); x != mkey.GetVarCoeffs().end(); ++x)
1048  {
1049  varcoeffs[x->first] = x->second + m_phys_offset[eid];
1050  }
1051  }
1052 
1053  LocalRegions::MatrixKey matkey(mkey.GetMatrixType(),
1054  (*m_exp)[eid]->DetShapeType(),
1055  *((*m_exp)[eid]),
1056  mkey.GetConstFactors(),varcoeffs);
1057 
1058  loc_mat = boost::dynamic_pointer_cast<LocalRegions::Expansion>((*m_exp)[n])->GetLocMatrix(matkey);
1059 
1060 
1061  if(RobinBCInfo.count(n) != 0) // add robin mass matrix
1062  {
1064 
1065  // declare local matrix from scaled matrix.
1066  int rows = loc_mat->GetRows();
1067  int cols = loc_mat->GetColumns();
1068  const NekDouble *dat = loc_mat->GetRawPtr();
1070  Blas::Dscal(rows*cols,loc_mat->Scale(),new_mat->GetRawPtr(),1);
1071 
1072  // add local matrix contribution
1073  for(rBC = RobinBCInfo.find(n)->second;rBC; rBC = rBC->next)
1074  {
1075  (*m_exp)[n]->AddRobinMassMatrix(rBC->m_robinID,rBC->m_robinPrimitiveCoeffs,new_mat);
1076  }
1077 
1078  NekDouble one = 1.0;
1079  // redeclare loc_mat to point to new_mat plus the scalar.
1080  loc_mat = MemoryManager<DNekScalMat>::AllocateSharedPtr(one,new_mat);
1081  }
1082 
1083  loc_lda = loc_mat->GetColumns();
1084 
1085  for(i = 0; i < loc_lda; ++i)
1086  {
1087  gid1 = locToGloMap->GetLocalToGlobalMap(m_coeff_offset[n] + i) - NumDirBCs;
1088  sign1 = locToGloMap->GetLocalToGlobalSign(m_coeff_offset[n] + i);
1089  if(gid1 >= 0)
1090  {
1091  for(j = 0; j < loc_lda; ++j)
1092  {
1093  gid2 = locToGloMap->GetLocalToGlobalMap(m_coeff_offset[n] + j) - NumDirBCs;
1094  sign2 = locToGloMap->GetLocalToGlobalSign(m_coeff_offset[n] + j);
1095  if(gid2 >= 0)
1096  {
1097  // When global matrix is symmetric,
1098  // only add the value for the upper
1099  // triangular part in order to avoid
1100  // entries to be entered twice
1101  if((matStorage == eFULL)||(gid2 >= gid1))
1102  {
1103  value = Gmat->GetValue(gid1,gid2) + sign1*sign2*(*loc_mat)(i,j);
1104  Gmat->SetValue(gid1,gid2,value);
1105  }
1106  }
1107  }
1108  }
1109  }
1110  }
1111 
1112  return Gmat;
1113  }
1114 
1115 
1116  /**
1117  * Consider a linear system
1118  * \f$\boldsymbol{M\hat{u}}_g=\boldsymbol{f}\f$ to be solved. Dependent
1119  * on the solution method, this function constructs
1120  * - <b>The full linear system</b><BR>
1121  * A call to the function #GenGlobalLinSysFullDirect
1122  * - <b>The statically condensed linear system</b><BR>
1123  * A call to the function #GenGlobalLinSysStaticCond
1124  *
1125  * @param mkey A key which uniquely defines the global
1126  * matrix to be constructed.
1127  * @param locToGloMap Contains the mapping array and required
1128  * information for the transformation from
1129  * local to global degrees of freedom.
1130  * @return (A shared pointer to) the global linear system in
1131  * required format.
1132  */
1134  const GlobalLinSysKey &mkey,
1135  const AssemblyMapCGSharedPtr &locToGloMap)
1136  {
1137  GlobalLinSysSharedPtr returnlinsys;
1138  boost::shared_ptr<ExpList> vExpList = GetSharedThisPtr();
1139 
1141 
1142  if (vType >= eSIZE_GlobalSysSolnType)
1143  {
1144  ASSERTL0(false,"Matrix solution type not defined");
1145  }
1146  std::string vSolnType = MultiRegions::GlobalSysSolnTypeMap[vType];
1147 
1148  return GetGlobalLinSysFactory().CreateInstance( vSolnType, mkey,
1149  vExpList, locToGloMap);
1150  }
1151 
1153  const GlobalLinSysKey &mkey,
1154  const AssemblyMapSharedPtr &locToGloMap)
1155  {
1156  boost::shared_ptr<ExpList> vExpList = GetSharedThisPtr();
1157  const map<int,RobinBCInfoSharedPtr> vRobinBCInfo = GetRobinBCInfo();
1158 
1160 
1161  if (vType >= eSIZE_GlobalSysSolnType)
1162  {
1163  ASSERTL0(false,"Matrix solution type not defined");
1164  }
1165  std::string vSolnType = MultiRegions::GlobalSysSolnTypeMap[vType];
1166 
1167  return GetGlobalLinSysFactory().CreateInstance(vSolnType,mkey,
1168  vExpList,locToGloMap);
1169  }
1170 
1171 
1172  /**
1173  * Given the elemental coefficients \f$\hat{u}_n^e\f$ of
1174  * an expansion, this function evaluates the spectral/hp
1175  * expansion \f$u^{\delta}(\boldsymbol{x})\f$ at the
1176  * quadrature points \f$\boldsymbol{x}_i\f$. The operation
1177  * is evaluated locally by the elemental function
1178  * StdRegions#StdExpansion#BwdTrans.
1179  *
1180  * @param inarray An array of size \f$N_{\mathrm{eof}}\f$
1181  * containing the local coefficients
1182  * \f$\hat{u}_n^e\f$.
1183  * @param outarray The resulting physical values at the
1184  * quadrature points
1185  * \f$u^{\delta}(\boldsymbol{x}_i)\f$
1186  * will be stored in this array of size
1187  * \f$Q_{\mathrm{tot}}\f$.
1188  */
1190  Array<OneD, NekDouble> &outarray)
1191  {
1193  for (int i = 0; i < m_collections.size(); ++i)
1194  {
1195  m_collections[i].ApplyOperator(Collections::eBwdTrans,
1196  inarray + m_coll_coeff_offset[i],
1197  tmp = outarray + m_coll_phys_offset[i]);
1198  }
1199  }
1200 
1202  const Array<OneD, const NekDouble> &gloCoord)
1203  {
1204  Array<OneD, NekDouble> stdCoord(GetCoordim(0),0.0);
1205  for (int i = 0; i < (*m_exp).size(); ++i)
1206  {
1207  if ((*m_exp)[i]->GetGeom()->ContainsPoint(gloCoord))
1208  {
1209  return (*m_exp)[i];
1210  }
1211  }
1212  ASSERTL0(false, "Cannot find element for this point.");
1213  return (*m_exp)[0]; // avoid warnings
1214  }
1215 
1216 
1217  /**
1218  * @todo need a smarter search here that first just looks at bounding
1219  * vertices - suggest first seeing if point is within 10% of
1220  * region defined by vertices. The do point search.
1221  */
1223  const Array<OneD, const NekDouble> &gloCoord,
1224  NekDouble tol,
1225  bool returnNearestElmt)
1226  {
1227  Array<OneD, NekDouble> Lcoords(gloCoord.num_elements());
1228 
1229  return GetExpIndex(gloCoord,Lcoords,tol,returnNearestElmt);
1230  }
1231 
1232 
1234  Array<OneD, NekDouble> &locCoords,
1235  NekDouble tol,
1236  bool returnNearestElmt)
1237  {
1238  NekDouble nearpt = 1e6;
1239 
1240  if (GetNumElmts() == 0)
1241  {
1242  return -1;
1243  }
1244  std::vector<std::pair<int,NekDouble> > elmtIdDist;
1245 
1246  // Manifold case (point may match multiple elements)
1247  if (GetExp(0)->GetCoordim() > GetExp(0)->GetShapeDimension())
1248  {
1251  NekDouble dist = 0.0;
1252 
1253  // Scan all elements and store those which may contain the point
1254  for (int i = 0; i < (*m_exp).size(); ++i)
1255  {
1256  if ((*m_exp)[i]->GetGeom()->ContainsPoint(gloCoords,
1257  locCoords,
1258  tol, nearpt))
1259  {
1260  w.SetX(gloCoords[0]);
1261  w.SetY(gloCoords[1]);
1262  w.SetZ(gloCoords[2]);
1263 
1264  // Find closest vertex
1265  for (int j = 0; j < (*m_exp)[i]->GetNverts(); ++j) {
1266  v = m_graph->GetVertex(
1267  (*m_exp)[i]->GetGeom()->GetVid(j));
1268  if (j == 0 || dist > v->dist(w))
1269  {
1270  dist = v->dist(w);
1271  }
1272  }
1273  elmtIdDist.push_back(
1274  std::pair<int, NekDouble>(i, dist));
1275  }
1276  }
1277 
1278  // Find nearest element
1279  if (!elmtIdDist.empty())
1280  {
1281  int min_id = elmtIdDist[0].first;
1282  NekDouble min_d = elmtIdDist[0].second;
1283 
1284  for (int i = 1; i < elmtIdDist.size(); ++i)
1285  {
1286  if (elmtIdDist[i].second < min_d) {
1287  min_id = elmtIdDist[i].first;
1288  min_d = elmtIdDist[i].second;
1289  }
1290  }
1291 
1292  // retrieve local coordinate of point
1293  (*m_exp)[min_id]->GetGeom()->GetLocCoords(gloCoords,
1294  locCoords);
1295  return min_id;
1296  }
1297  else
1298  {
1299  return -1;
1300  }
1301  }
1302  // non-embedded mesh (point can only match one element)
1303  else
1304  {
1305  static int start = 0;
1306  int min_id = 0;
1307  NekDouble nearpt_min = 1e6;
1308  Array<OneD, NekDouble> savLocCoords(locCoords.num_elements());
1309 
1310  // restart search from last found value
1311  for (int i = start; i < (*m_exp).size(); ++i)
1312  {
1313  if ((*m_exp)[i]->GetGeom()->ContainsPoint(gloCoords,
1314  locCoords,
1315  tol, nearpt))
1316  {
1317  start = i;
1318  return i;
1319  }
1320  else
1321  {
1322  if(nearpt < nearpt_min)
1323  {
1324  min_id = i;
1325  nearpt_min = nearpt;
1326  Vmath::Vcopy(locCoords.num_elements(),locCoords,1,savLocCoords,1);
1327  }
1328  }
1329  }
1330 
1331  for (int i = 0; i < start; ++i)
1332  {
1333  if ((*m_exp)[i]->GetGeom()->ContainsPoint(gloCoords,
1334  locCoords,
1335  tol, nearpt))
1336  {
1337  start = i;
1338  return i;
1339  }
1340  else
1341  {
1342  if(nearpt < nearpt_min)
1343  {
1344  min_id = i;
1345  nearpt_min = nearpt;
1346  Vmath::Vcopy(locCoords.num_elements(),
1347  locCoords,1,savLocCoords,1);
1348  }
1349  }
1350  }
1351 
1352  std::string msg = "Failed to find point within element to tolerance of "
1353  + boost::lexical_cast<std::string>(tol)
1354  + " using local point ("
1355  + boost::lexical_cast<std::string>(locCoords[0]) +","
1356  + boost::lexical_cast<std::string>(locCoords[1]) +","
1357  + boost::lexical_cast<std::string>(locCoords[1])
1358  + ") in element: "
1359  + boost::lexical_cast<std::string>(min_id);
1360  WARNINGL1(false,msg.c_str());
1361 
1362  if(returnNearestElmt)
1363  {
1364  Vmath::Vcopy(locCoords.num_elements(),savLocCoords,1,locCoords,1);
1365  return min_id;
1366  }
1367  else
1368  {
1369  return -1;
1370  }
1371 
1372  }
1373  }
1374 
1375 
1376  /**
1377  * Configures geometric info, such as tangent direction, on each
1378  * expansion.
1379  * @param graph2D Mesh
1380  */
1382  {
1383 
1384  }
1385 
1386  /**
1387  * @brief Reset geometry information, metrics, matrix managers and
1388  * geometry information.
1389  *
1390  * This routine clears all matrix managers and resets all geometry
1391  * information, which allows the geometry information to be dynamically
1392  * updated as the solver is run.
1393  */
1395  {
1396  // Reset matrix managers.
1398  DNekScalMat, LocalRegions::MatrixKey::opLess>::ClearManager();
1399  LibUtilities::NekManager<LocalRegions::MatrixKey,
1401 
1402  // Loop over all elements and reset geometry information.
1403  for (int i = 0; i < m_exp->size(); ++i)
1404  {
1405  (*m_exp)[i]->GetGeom()->Reset(m_graph->GetCurvedEdges(),
1406  m_graph->GetCurvedFaces());
1407  }
1408 
1409  // Loop over all elements and rebuild geometric factors.
1410  for (int i = 0; i < m_exp->size(); ++i)
1411  {
1412  (*m_exp)[i]->Reset();
1413  }
1414  }
1415 
1416  /**
1417  * Write Tecplot Files Header
1418  * @param outfile Output file name.
1419  * @param var variables names
1420  */
1421  void ExpList::v_WriteTecplotHeader(std::ostream &outfile,
1422  std::string var)
1423  {
1424  if (GetNumElmts() == 0)
1425  {
1426  return;
1427  }
1428 
1429  int coordim = GetExp(0)->GetCoordim();
1430  char vars[3] = { 'x', 'y', 'z' };
1431 
1432  if (m_expType == e3DH1D)
1433  {
1434  coordim += 1;
1435  }
1436  else if (m_expType == e3DH2D)
1437  {
1438  coordim += 2;
1439  }
1440 
1441  outfile << "Variables = x";
1442  for (int i = 1; i < coordim; ++i)
1443  {
1444  outfile << ", " << vars[i];
1445  }
1446 
1447  if (var.size() > 0)
1448  {
1449  outfile << ", " << var;
1450  }
1451 
1452  outfile << std::endl << std::endl;
1453  }
1454 
1455  /**
1456  * Write Tecplot Files Zone
1457  * @param outfile Output file name.
1458  * @param expansion Expansion that is considered
1459  */
1460  void ExpList::v_WriteTecplotZone(std::ostream &outfile, int expansion)
1461  {
1462  int i, j;
1463  int coordim = GetCoordim(0);
1464  int nPoints = GetTotPoints();
1465  int nBases = (*m_exp)[0]->GetNumBases();
1466  int numBlocks = 0;
1467 
1469 
1470  if (expansion == -1)
1471  {
1472  nPoints = GetTotPoints();
1473 
1474  coords[0] = Array<OneD, NekDouble>(nPoints);
1475  coords[1] = Array<OneD, NekDouble>(nPoints);
1476  coords[2] = Array<OneD, NekDouble>(nPoints);
1477 
1478  GetCoords(coords[0], coords[1], coords[2]);
1479 
1480  for (i = 0; i < m_exp->size(); ++i)
1481  {
1482  int numInt = 1;
1483 
1484  for (j = 0; j < nBases; ++j)
1485  {
1486  numInt *= (*m_exp)[i]->GetNumPoints(j)-1;
1487  }
1488 
1489  numBlocks += numInt;
1490  }
1491  }
1492  else
1493  {
1494  nPoints = (*m_exp)[expansion]->GetTotPoints();
1495 
1496  coords[0] = Array<OneD, NekDouble>(nPoints);
1497  coords[1] = Array<OneD, NekDouble>(nPoints);
1498  coords[2] = Array<OneD, NekDouble>(nPoints);
1499 
1500  (*m_exp)[expansion]->GetCoords(coords[0], coords[1], coords[2]);
1501 
1502  numBlocks = 1;
1503  for (j = 0; j < nBases; ++j)
1504  {
1505  numBlocks *= (*m_exp)[expansion]->GetNumPoints(j)-1;
1506  }
1507  }
1508 
1509  if (m_expType == e3DH1D)
1510  {
1511  nBases += 1;
1512  coordim += 1;
1513  int nPlanes = GetZIDs().num_elements();
1514  NekDouble tmp = numBlocks * (nPlanes-1.0) / nPlanes;
1515  numBlocks = (int)tmp;
1516  }
1517  else if (m_expType == e3DH2D)
1518  {
1519  nBases += 2;
1520  coordim += 1;
1521  }
1522 
1523  outfile << "Zone, N=" << nPoints << ", E="
1524  << numBlocks << ", F=FEBlock" ;
1525 
1526  switch(nBases)
1527  {
1528  case 2:
1529  outfile << ", ET=QUADRILATERAL" << std::endl;
1530  break;
1531  case 3:
1532  outfile << ", ET=BRICK" << std::endl;
1533  break;
1534  default:
1535  ASSERTL0(false,"Not set up for this type of output");
1536  break;
1537  }
1538 
1539  // Write out coordinates
1540  for (j = 0; j < coordim; ++j)
1541  {
1542  for (i = 0; i < nPoints; ++i)
1543  {
1544  outfile << coords[j][i] << " ";
1545  if (i % 1000 == 0 && i)
1546  {
1547  outfile << std::endl;
1548  }
1549  }
1550  outfile << std::endl;
1551  }
1552  }
1553 
1554  void ExpList::v_WriteTecplotConnectivity(std::ostream &outfile,
1555  int expansion)
1556  {
1557  int i,j,k,l;
1558  int nbase = (*m_exp)[0]->GetNumBases();
1559  int cnt = 0;
1560 
1561  boost::shared_ptr<LocalRegions::ExpansionVector> exp = m_exp;
1562 
1563  if (expansion != -1)
1564  {
1565  exp = boost::shared_ptr<LocalRegions::ExpansionVector>(
1567  (*exp)[0] = (*m_exp)[expansion];
1568  }
1569 
1570  if (nbase == 2)
1571  {
1572  for(i = 0; i < (*exp).size(); ++i)
1573  {
1574  const int np0 = (*exp)[i]->GetNumPoints(0);
1575  const int np1 = (*exp)[i]->GetNumPoints(1);
1576 
1577  for(j = 1; j < np1; ++j)
1578  {
1579  for(k = 1; k < np0; ++k)
1580  {
1581  outfile << cnt + (j-1)*np0 + k << " ";
1582  outfile << cnt + (j-1)*np0 + k+1 << " ";
1583  outfile << cnt + j *np0 + k+1 << " ";
1584  outfile << cnt + j *np0 + k << endl;
1585  }
1586  }
1587 
1588  cnt += np0*np1;
1589  }
1590  }
1591  else if (nbase == 3)
1592  {
1593  for(i = 0; i < (*exp).size(); ++i)
1594  {
1595  const int np0 = (*exp)[i]->GetNumPoints(0);
1596  const int np1 = (*exp)[i]->GetNumPoints(1);
1597  const int np2 = (*exp)[i]->GetNumPoints(2);
1598  const int np01 = np0*np1;
1599 
1600  for(j = 1; j < np2; ++j)
1601  {
1602  for(k = 1; k < np1; ++k)
1603  {
1604  for(l = 1; l < np0; ++l)
1605  {
1606  outfile << cnt + (j-1)*np01 + (k-1)*np0 + l << " ";
1607  outfile << cnt + (j-1)*np01 + (k-1)*np0 + l+1 << " ";
1608  outfile << cnt + (j-1)*np01 + k *np0 + l+1 << " ";
1609  outfile << cnt + (j-1)*np01 + k *np0 + l << " ";
1610  outfile << cnt + j *np01 + (k-1)*np0 + l << " ";
1611  outfile << cnt + j *np01 + (k-1)*np0 + l+1 << " ";
1612  outfile << cnt + j *np01 + k *np0 + l+1 << " ";
1613  outfile << cnt + j *np01 + k *np0 + l << endl;
1614  }
1615  }
1616  }
1617  cnt += np0*np1*np2;
1618  }
1619  }
1620  else
1621  {
1622  ASSERTL0(false,"Not set up for this dimension");
1623  }
1624  }
1625 
1626  /**
1627  * Write Tecplot Files Field
1628  * @param outfile Output file name.
1629  * @param expansion Expansion that is considered
1630  */
1631  void ExpList::v_WriteTecplotField(std::ostream &outfile, int expansion)
1632  {
1633  if (expansion == -1)
1634  {
1635  int totpoints = GetTotPoints();
1636  if(m_physState == false)
1637  {
1639  }
1640 
1641  for(int i = 0; i < totpoints; ++i)
1642  {
1643  outfile << m_phys[i] << " ";
1644  if(i % 1000 == 0 && i)
1645  {
1646  outfile << std::endl;
1647  }
1648  }
1649  outfile << std::endl;
1650 
1651  }
1652  else
1653  {
1654  int nPoints = (*m_exp)[expansion]->GetTotPoints();
1655 
1656  for (int i = 0; i < nPoints; ++i)
1657  {
1658  outfile << m_phys[i + m_phys_offset[expansion]] << " ";
1659  }
1660 
1661  outfile << std::endl;
1662  }
1663  }
1664 
1665  void ExpList::WriteVtkHeader(std::ostream &outfile)
1666  {
1667  outfile << "<?xml version=\"1.0\"?>" << endl;
1668  outfile << "<VTKFile type=\"UnstructuredGrid\" version=\"0.1\" "
1669  << "byte_order=\"LittleEndian\">" << endl;
1670  outfile << " <UnstructuredGrid>" << endl;
1671  }
1672 
1673  void ExpList::WriteVtkFooter(std::ostream &outfile)
1674  {
1675  outfile << " </UnstructuredGrid>" << endl;
1676  outfile << "</VTKFile>" << endl;
1677  }
1678 
1679  void ExpList::v_WriteVtkPieceHeader(std::ostream &outfile, int expansion, int istrip)
1680  {
1681  ASSERTL0(false, "Routine not implemented for this expansion.");
1682  }
1683 
1684  void ExpList::WriteVtkPieceFooter(std::ostream &outfile, int expansion)
1685  {
1686  outfile << " </PointData>" << endl;
1687  outfile << " </Piece>" << endl;
1688  }
1689 
1690  void ExpList::v_WriteVtkPieceData(std::ostream &outfile, int expansion,
1691  std::string var)
1692  {
1693  int i;
1694  int nq = (*m_exp)[expansion]->GetTotPoints();
1695 
1696  // printing the fields of that zone
1697  outfile << " <DataArray type=\"Float64\" Name=\""
1698  << var << "\">" << endl;
1699  outfile << " ";
1700  const Array<OneD, NekDouble> phys = m_phys + m_phys_offset[expansion];
1701  for(i = 0; i < nq; ++i)
1702  {
1703  outfile << (fabs(phys[i]) < NekConstants::kNekZeroTol ? 0 : phys[i]) << " ";
1704  }
1705  outfile << endl;
1706  outfile << " </DataArray>" << endl;
1707  }
1708 
1709  /**
1710  * Given a spectral/hp approximation
1711  * \f$u^{\delta}(\boldsymbol{x})\f$ evaluated at the quadrature points
1712  * (which should be contained in #m_phys), this function calculates the
1713  * \f$L_\infty\f$ error of this approximation with respect to an exact
1714  * solution. The local distribution of the quadrature points allows an
1715  * elemental evaluation of this operation through the functions
1716  * StdRegions#StdExpansion#Linf.
1717  *
1718  * The exact solution, also evaluated at the quadrature
1719  * points, should be contained in the variable #m_phys of
1720  * the ExpList object \a Sol.
1721  *
1722  * @param soln A 1D array, containing the discrete
1723  * evaluation of the exact solution at the
1724  * quadrature points in its array #m_phys.
1725  * @return The \f$L_\infty\f$ error of the approximation.
1726  */
1728  const Array<OneD, const NekDouble> &inarray,
1729  const Array<OneD, const NekDouble> &soln)
1730  {
1731  NekDouble err = 0.0;
1732 
1733  if (soln == NullNekDouble1DArray)
1734  {
1735  err = Vmath::Vmax(m_npoints, inarray, 1);
1736  }
1737  else
1738  {
1739  for (int i = 0; i < m_npoints; ++i)
1740  {
1741  err = max(err, abs(inarray[i] - soln[i]));
1742  }
1743  }
1744 
1745  m_comm->GetRowComm()->AllReduce(err, LibUtilities::ReduceMax);
1746 
1747  return err;
1748  }
1749 
1750  /**
1751  * Given a spectral/hp approximation \f$u^{\delta}(\boldsymbol{x})\f$
1752  * evaluated at the quadrature points (which should be contained in
1753  * #m_phys), this function calculates the \f$L_2\f$ error of this
1754  * approximation with respect to an exact solution. The local
1755  * distribution of the quadrature points allows an elemental evaluation
1756  * of this operation through the functions StdRegions#StdExpansion#L2.
1757  *
1758  * The exact solution, also evaluated at the quadrature points, should
1759  * be contained in the variable #m_phys of the ExpList object \a Sol.
1760  *
1761  * @param Sol An ExpList, containing the discrete
1762  * evaluation of the exact solution at the
1763  * quadrature points in its array #m_phys.
1764  * @return The \f$L_2\f$ error of the approximation.
1765  */
1767  const Array<OneD, const NekDouble> &inarray,
1768  const Array<OneD, const NekDouble> &soln)
1769  {
1770  NekDouble err = 0.0, errl2;
1771  int i;
1772 
1773  if (soln == NullNekDouble1DArray)
1774  {
1775  for (i = 0; i < (*m_exp).size(); ++i)
1776  {
1777  errl2 = (*m_exp)[i]->L2(inarray + m_phys_offset[i]);
1778  err += errl2*errl2;
1779  }
1780  }
1781  else
1782  {
1783  for (i = 0; i < (*m_exp).size(); ++i)
1784  {
1785  errl2 = (*m_exp)[i]->L2(inarray + m_phys_offset[i],
1786  soln + m_phys_offset[i]);
1787  err += errl2*errl2;
1788  }
1789  }
1790 
1791  m_comm->GetRowComm()->AllReduce(err, LibUtilities::ReduceSum);
1792 
1793  return sqrt(err);
1794  }
1795 
1797  {
1798  NekDouble err = 0.0;
1799  int i = 0;
1800 
1801  for (i = 0; i < (*m_exp).size(); ++i)
1802  {
1803  err += (*m_exp)[m_offset_elmt_id[i]]->Integral(inarray + m_phys_offset[i]);
1804  }
1805  m_comm->GetRowComm()->AllReduce(err, LibUtilities::ReduceSum);
1806 
1807  return err;
1808  }
1809 
1811  {
1812  ASSERTL0(false,
1813  "This method is not defined or valid for this class type");
1814  Array<OneD, NekDouble> NoEnergy(1,0.0);
1815  return NoEnergy;
1816  }
1817 
1819  {
1820  ASSERTL0(false,
1821  "This method is not defined or valid for this class type");
1823  return trans;
1824  }
1825 
1827  {
1828  ASSERTL0(false,
1829  "This method is not defined or valid for this class type");
1830  NekDouble len = 0.0;
1831  return len;
1832  }
1833 
1835  {
1836  ASSERTL0(false,
1837  "This method is not defined or valid for this class type");
1838  Array<OneD, unsigned int> NoModes(1);
1839  return NoModes;
1840  }
1841 
1843  {
1844  ASSERTL0(false,
1845  "This method is not defined or valid for this class type");
1846  Array<OneD, unsigned int> NoModes(1);
1847  return NoModes;
1848  }
1849 
1850 
1852  {
1853  ASSERTL0(false,
1854  "This method is not defined or valid for this class type");
1855  }
1856 
1858  ASSERTL0(false,
1859  "This method is not defined or valid for this class type");
1860  }
1861 
1863  const std::string &fileName,
1864  const std::string &varName,
1865  const boost::shared_ptr<ExpList> locExpList)
1866  {
1867  string varString = fileName.substr(0, fileName.find_last_of("."));
1868  int j, k, len = varString.length();
1869  varString = varString.substr(len-1, len);
1870 
1871  std::vector<LibUtilities::FieldDefinitionsSharedPtr> FieldDef;
1872  std::vector<std::vector<NekDouble> > FieldData;
1873 
1874  LibUtilities::FieldIO f(m_session->GetComm());
1875  f.Import(fileName, FieldDef, FieldData);
1876 
1877  bool found = false;
1878  for (j = 0; j < FieldDef.size(); ++j)
1879  {
1880  for (k = 0; k < FieldDef[j]->m_fields.size(); ++k)
1881  {
1882  if (FieldDef[j]->m_fields[k] == varName)
1883  {
1884  // Copy FieldData into locExpList
1885  locExpList->ExtractDataToCoeffs(
1886  FieldDef[j], FieldData[j],
1887  FieldDef[j]->m_fields[k],
1888  locExpList->UpdateCoeffs());
1889  found = true;
1890  }
1891  }
1892  }
1893 
1894  ASSERTL0(found, "Could not find variable '"+varName+
1895  "' in file boundary condition "+fileName);
1896  locExpList->BwdTrans_IterPerExp(
1897  locExpList->GetCoeffs(),
1898  locExpList->UpdatePhys());
1899  }
1900 
1901  /**
1902  * Given a spectral/hp approximation
1903  * \f$u^{\delta}(\boldsymbol{x})\f$ evaluated at the quadrature points
1904  * (which should be contained in #m_phys), this function calculates the
1905  * \f$H^1_2\f$ error of this approximation with respect to an exact
1906  * solution. The local distribution of the quadrature points allows an
1907  * elemental evaluation of this operation through the functions
1908  * StdRegions#StdExpansion#H1.
1909  *
1910  * The exact solution, also evaluated at the quadrature points, should
1911  * be contained in the variable #m_phys of the ExpList object \a Sol.
1912  *
1913  * @param soln An 1D array, containing the discrete evaluation
1914  * of the exact solution at the quadrature points.
1915  *
1916  * @return The \f$H^1_2\f$ error of the approximation.
1917  */
1919  const Array<OneD, const NekDouble> &inarray,
1920  const Array<OneD, const NekDouble> &soln)
1921  {
1922  NekDouble err = 0.0, errh1;
1923  int i;
1924 
1925  for (i = 0; i < (*m_exp).size(); ++i)
1926  {
1927  errh1 = (*m_exp)[i]->H1(inarray + m_phys_offset[i],
1928  soln + m_phys_offset[i]);
1929  err += errh1*errh1;
1930  }
1931 
1932  m_comm->GetRowComm()->AllReduce(err, LibUtilities::ReduceSum);
1933 
1934  return sqrt(err);
1935  }
1936 
1937  void ExpList::GeneralGetFieldDefinitions(std::vector<LibUtilities::FieldDefinitionsSharedPtr> &fielddef,
1938  int NumHomoDir,
1939  int NumHomoStrip,
1941  std::vector<NekDouble> &HomoLen,
1942  std::vector<unsigned int> &HomoZIDs,
1943  std::vector<unsigned int> &HomoYIDs)
1944  {
1945  int startenum = (int) LibUtilities::eSegment;
1946  int endenum = (int) LibUtilities::eHexahedron;
1947  int s = 0;
1949 
1950  ASSERTL1(NumHomoDir == HomoBasis.num_elements(),"Homogeneous basis is not the same length as NumHomoDir");
1951  ASSERTL1(NumHomoDir == HomoLen.size(),"Homogeneous length vector is not the same length as NumHomDir");
1952 
1953  // count number of shapes
1954  switch((*m_exp)[0]->GetShapeDimension())
1955  {
1956  case 1:
1957  startenum = (int) LibUtilities::eSegment;
1958  endenum = (int) LibUtilities::eSegment;
1959  break;
1960  case 2:
1961  startenum = (int) LibUtilities::eTriangle;
1962  endenum = (int) LibUtilities::eQuadrilateral;
1963  break;
1964  case 3:
1965  startenum = (int) LibUtilities::eTetrahedron;
1966  endenum = (int) LibUtilities::eHexahedron;
1967  break;
1968  }
1969 
1970  for(s = startenum; s <= endenum; ++s)
1971  {
1972  std::vector<unsigned int> elementIDs;
1973  std::vector<LibUtilities::BasisType> basis;
1974  std::vector<unsigned int> numModes;
1975  std::vector<std::string> fields;
1976 
1977  bool first = true;
1978  bool UniOrder = true;
1979  int n;
1980 
1981  shape = (LibUtilities::ShapeType) s;
1982 
1983  for(int i = 0; i < (*m_exp).size(); ++i)
1984  {
1985  if((*m_exp)[i]->GetGeom()->GetShapeType() == shape)
1986  {
1987  elementIDs.push_back((*m_exp)[i]->GetGeom()->GetGlobalID());
1988  if(first)
1989  {
1990  for(int j = 0; j < (*m_exp)[i]->GetNumBases(); ++j)
1991  {
1992  basis.push_back((*m_exp)[i]->GetBasis(j)->GetBasisType());
1993  numModes.push_back((*m_exp)[i]->GetBasis(j)->GetNumModes());
1994  }
1995 
1996  // add homogeneous direction details if defined
1997  for(n = 0 ; n < NumHomoDir; ++n)
1998  {
1999  basis.push_back(HomoBasis[n]->GetBasisType());
2000  numModes.push_back(HomoBasis[n]->GetNumModes());
2001  }
2002 
2003  first = false;
2004  }
2005  else
2006  {
2007  ASSERTL0((*m_exp)[i]->GetBasis(0)->GetBasisType() == basis[0],"Routine is not set up for multiple bases definitions");
2008 
2009  for(int j = 0; j < (*m_exp)[i]->GetNumBases(); ++j)
2010  {
2011  numModes.push_back((*m_exp)[i]->GetBasis(j)->GetNumModes());
2012  if(numModes[j] != (*m_exp)[i]->GetBasis(j)->GetNumModes())
2013  {
2014  UniOrder = false;
2015  }
2016  }
2017  // add homogeneous direction details if defined
2018  for(n = 0 ; n < NumHomoDir; ++n)
2019  {
2020  numModes.push_back(HomoBasis[n]->GetNumModes());
2021  }
2022  }
2023  }
2024  }
2025 
2026 
2027  if(elementIDs.size() > 0)
2028  {
2029  for(int i = 0; i < NumHomoStrip; ++i)
2030  {
2033  AllocateSharedPtr(shape, elementIDs, basis,
2034  UniOrder, numModes,fields,
2035  NumHomoDir, HomoLen, HomoZIDs,
2036  HomoYIDs);
2037  fielddef.push_back(fdef);
2038  }
2039  }
2040  }
2041  }
2042 
2043 
2044  //
2045  // Virtual functions
2046  //
2047  std::vector<LibUtilities::FieldDefinitionsSharedPtr> ExpList::v_GetFieldDefinitions()
2048  {
2049  std::vector<LibUtilities::FieldDefinitionsSharedPtr> returnval;
2050  v_GetFieldDefinitions(returnval);
2051  return returnval;
2052  }
2053 
2054  void ExpList::v_GetFieldDefinitions(std::vector<LibUtilities::FieldDefinitionsSharedPtr> &fielddef)
2055  {
2056  GeneralGetFieldDefinitions(fielddef);
2057  }
2058 
2059  //Append the element data listed in elements
2060  //fielddef->m_ElementIDs onto fielddata
2061  void ExpList::v_AppendFieldData(LibUtilities::FieldDefinitionsSharedPtr &fielddef, std::vector<NekDouble> &fielddata)
2062  {
2063  v_AppendFieldData(fielddef,fielddata,m_coeffs);
2064  }
2065 
2066  void ExpList::v_AppendFieldData(LibUtilities::FieldDefinitionsSharedPtr &fielddef, std::vector<NekDouble> &fielddata, Array<OneD, NekDouble> &coeffs)
2067  {
2068  int i;
2069  // Determine mapping from element ids to location in
2070  // expansion list
2071  // Determine mapping from element ids to location in
2072  // expansion list
2073  map<int, int> ElmtID_to_ExpID;
2074  for(i = 0; i < (*m_exp).size(); ++i)
2075  {
2076  ElmtID_to_ExpID[(*m_exp)[i]->GetGeom()->GetGlobalID()] = i;
2077  }
2078 
2079  for(i = 0; i < fielddef->m_elementIDs.size(); ++i)
2080  {
2081  int eid = ElmtID_to_ExpID[fielddef->m_elementIDs[i]];
2082  int datalen = (*m_exp)[eid]->GetNcoeffs();
2083  fielddata.insert(fielddata.end(),&coeffs[m_coeff_offset[eid]],&coeffs[m_coeff_offset[eid]]+datalen);
2084  }
2085 
2086  }
2087 
2088  /// Extract the data in fielddata into the coeffs
2091  std::vector<NekDouble> &fielddata,
2092  std::string &field,
2093  Array<OneD, NekDouble> &coeffs)
2094  {
2095  v_ExtractDataToCoeffs(fielddef,fielddata,field,coeffs);
2096  }
2097 
2098  void ExpList::ExtractCoeffsToCoeffs(const boost::shared_ptr<ExpList> &fromExpList, const Array<OneD, const NekDouble> &fromCoeffs, Array<OneD, NekDouble> &toCoeffs)
2099  {
2100  v_ExtractCoeffsToCoeffs(fromExpList,fromCoeffs,toCoeffs);
2101  }
2102 
2103  /**
2104  * @brief Extract data from raw field data into expansion list.
2105  *
2106  * @param fielddef Field definitions.
2107  * @param fielddata Data for associated field.
2108  * @param field Field variable name.
2109  * @param coeffs Resulting coefficient array.
2110  */
2113  std::vector<NekDouble> &fielddata,
2114  std::string &field,
2115  Array<OneD, NekDouble> &coeffs)
2116  {
2117  int i, expId;
2118  int offset = 0;
2119  int modes_offset = 0;
2120  int datalen = fielddata.size()/fielddef->m_fields.size();
2121 
2122  // Find data location according to field definition
2123  for(i = 0; i < fielddef->m_fields.size(); ++i)
2124  {
2125  if(fielddef->m_fields[i] == field)
2126  {
2127  break;
2128  }
2129  offset += datalen;
2130  }
2131 
2132  ASSERTL0(i != fielddef->m_fields.size(),
2133  "Field (" + field + ") not found in file.");
2134 
2135  // Determine mapping from element ids to location in expansion list
2136  map<int, int> elmtToExpId;
2137 
2138  // Loop in reverse order so that in case where using a Homogeneous
2139  // expansion it sets geometry ids to first part of m_exp
2140  // list. Otherwise will set to second (complex) expansion
2141  for(i = (*m_exp).size()-1; i >= 0; --i)
2142  {
2143  elmtToExpId[(*m_exp)[i]->GetGeom()->GetGlobalID()] = i;
2144  }
2145 
2146  for (i = 0; i < fielddef->m_elementIDs.size(); ++i)
2147  {
2148  // Reset modes_offset in the case where all expansions of
2149  // the same order.
2150  if (fielddef->m_uniOrder == true)
2151  {
2152  modes_offset = 0;
2153  }
2154 
2155  datalen = LibUtilities::GetNumberOfCoefficients(fielddef->m_shapeType,
2156  fielddef->m_numModes, modes_offset);
2157 
2158  const int elmtId = fielddef->m_elementIDs[i];
2159  if (elmtToExpId.count(elmtId) == 0)
2160  {
2161  offset += datalen;
2162  continue;
2163  }
2164 
2165  expId = elmtToExpId[elmtId];
2166 
2167  if (datalen == (*m_exp)[expId]->GetNcoeffs())
2168  {
2169  Vmath::Vcopy(datalen, &fielddata[offset], 1,
2170  &coeffs[m_coeff_offset[expId]], 1);
2171  }
2172  else
2173  {
2174  (*m_exp)[expId]->ExtractDataToCoeffs(
2175  &fielddata[offset], fielddef->m_numModes,
2176  modes_offset, &coeffs[m_coeff_offset[expId]]);
2177  }
2178 
2179  offset += datalen;
2180  modes_offset += (*m_exp)[0]->GetNumBases();
2181  }
2182 
2183  return;
2184  }
2185 
2186  void ExpList::v_ExtractCoeffsToCoeffs(const boost::shared_ptr<ExpList> &fromExpList, const Array<OneD, const NekDouble> &fromCoeffs, Array<OneD, NekDouble> &toCoeffs)
2187  {
2188  int i;
2189  int offset = 0;
2190 
2191  // check if the same and if so just copy over coeffs
2192  if(fromExpList->GetNcoeffs() == m_ncoeffs)
2193  {
2194  Vmath::Vcopy(m_ncoeffs,fromCoeffs,1,toCoeffs,1);
2195  }
2196  else
2197  {
2198  std::vector<unsigned int> nummodes;
2199  for(i = 0; i < (*m_exp).size(); ++i)
2200  {
2201  int eid = m_offset_elmt_id[i];
2202  for(int j= 0; j < fromExpList->GetExp(eid)->GetNumBases(); ++j)
2203  {
2204  nummodes.push_back(fromExpList->GetExp(eid)->GetBasisNumModes(j));
2205  }
2206 
2207  (*m_exp)[eid]->ExtractDataToCoeffs(&fromCoeffs[offset], nummodes,0,
2208  &toCoeffs[m_coeff_offset[eid]]);
2209 
2210  offset += fromExpList->GetExp(eid)->GetNcoeffs();
2211  }
2212  }
2213  }
2214 
2215 
2218  {
2219  ASSERTL0(false,
2220  "This method is not defined or valid for this class type");
2222  return result;
2223  }
2224 
2225  boost::shared_ptr<ExpList> &ExpList::v_UpdateBndCondExpansion(int i)
2226  {
2227  ASSERTL0(false,
2228  "This method is not defined or valid for this class type");
2229  static boost::shared_ptr<ExpList> result;
2230  return result;
2231  }
2232 
2234  const Array<OneD, const Array<OneD, NekDouble> > &Vec,
2235  const Array<OneD, const NekDouble> &Fwd,
2236  const Array<OneD, const NekDouble> &Bwd,
2237  Array<OneD, NekDouble> &Upwind)
2238  {
2239  ASSERTL0(false,
2240  "This method is not defined or valid for this class type");
2241  }
2242 
2244  const Array<OneD, const NekDouble> &Vn,
2245  const Array<OneD, const NekDouble> &Fwd,
2246  const Array<OneD, const NekDouble> &Bwd,
2247  Array<OneD, NekDouble> &Upwind)
2248  {
2249  ASSERTL0(false,
2250  "This method is not defined or valid for this class type");
2251  }
2252 
2253  boost::shared_ptr<ExpList> &ExpList::v_GetTrace()
2254  {
2255  ASSERTL0(false,
2256  "This method is not defined or valid for this class type");
2257  static boost::shared_ptr<ExpList> returnVal;
2258  return returnVal;
2259  }
2260 
2261  boost::shared_ptr<AssemblyMapDG> &ExpList::v_GetTraceMap()
2262  {
2263  ASSERTL0(false,
2264  "This method is not defined or valid for this class type");
2265  static boost::shared_ptr<AssemblyMapDG> result;
2266  return result;
2267  }
2268 
2270  {
2271  return GetTraceMap()->GetBndCondTraceToGlobalTraceMap();
2272  }
2273 
2275  Array<OneD, Array<OneD, NekDouble> > &normals)
2276  {
2277  ASSERTL0(false,
2278  "This method is not defined or valid for this class type");
2279  }
2280 
2282  const Array<OneD, const NekDouble> &Fx,
2283  const Array<OneD, const NekDouble> &Fy,
2284  Array<OneD, NekDouble> &outarray)
2285  {
2286  ASSERTL0(false,
2287  "This method is not defined or valid for this class type");
2288  }
2289 
2291  const Array<OneD, const NekDouble> &Fn,
2292  Array<OneD, NekDouble> &outarray)
2293  {
2294  ASSERTL0(false,
2295  "This method is not defined or valid for this class type");
2296  }
2297 
2299  const Array<OneD, const NekDouble> &Fwd,
2300  const Array<OneD, const NekDouble> &Bwd,
2301  Array<OneD, NekDouble> &outarray)
2302  {
2303  ASSERTL0(false,
2304  "This method is not defined or valid for this class type");
2305  }
2306 
2308  Array<OneD,NekDouble> &Bwd)
2309  {
2310  ASSERTL0(false,
2311  "This method is not defined or valid for this class type");
2312  }
2313 
2315  const Array<OneD,const NekDouble> &field,
2316  Array<OneD,NekDouble> &Fwd,
2317  Array<OneD,NekDouble> &Bwd)
2318  {
2319  ASSERTL0(false,
2320  "This method is not defined or valid for this class type");
2321  }
2322 
2324  {
2325  ASSERTL0(false,
2326  "This method is not defined or valid for this class type");
2327  }
2328 
2330  const Array<OneD, const NekDouble> &inarray,
2331  Array<OneD,NekDouble> &outarray)
2332  {
2333  ASSERTL0(false,
2334  "This method is not defined or valid for this class type");
2335  }
2336 
2338  const Array<OneD,const NekDouble> &inarray,
2339  Array<OneD, NekDouble> &outarray,
2340  CoeffState coeffstate)
2341  {
2342  ASSERTL0(false,
2343  "This method is not defined or valid for this class type");
2344  }
2345 
2347  const Array<OneD, const NekDouble> &inarray,
2348  Array<OneD, NekDouble> &outarray,
2349  const FlagList &flags,
2350  const StdRegions::ConstFactorMap &factors,
2351  const StdRegions::VarCoeffMap &varcoeff,
2352  const Array<OneD, const NekDouble> &dirForcing)
2353  {
2354  ASSERTL0(false, "HelmSolve not implemented.");
2355  }
2356 
2358  const Array<OneD, Array<OneD, NekDouble> > &velocity,
2359  const Array<OneD, const NekDouble> &inarray,
2360  Array<OneD, NekDouble> &outarray,
2361  const NekDouble lambda,
2362  CoeffState coeffstate,
2363  const Array<OneD, const NekDouble>& dirForcing)
2364  {
2365  ASSERTL0(false,
2366  "This method is not defined or valid for this class type");
2367  }
2368 
2370  const Array<OneD, Array<OneD, NekDouble> > &velocity,
2371  const Array<OneD, const NekDouble> &inarray,
2372  Array<OneD, NekDouble> &outarray,
2373  const NekDouble lambda,
2374  CoeffState coeffstate,
2375  const Array<OneD, const NekDouble>& dirForcing)
2376  {
2377  ASSERTL0(false,
2378  "This method is not defined or valid for this class type");
2379  }
2380 
2382  Array<OneD, NekDouble> &outarray,
2383  CoeffState coeffstate,
2384  bool Shuff,
2385  bool UnShuff)
2386  {
2387  ASSERTL0(false,
2388  "This method is not defined or valid for this class type");
2389  }
2390 
2392  Array<OneD, NekDouble> &outarray,
2393  CoeffState coeffstate,
2394  bool Shuff,
2395  bool UnShuff)
2396  {
2397  ASSERTL0(false,
2398  "This method is not defined or valid for this class type");
2399  }
2400 
2402  {
2403  ASSERTL0(false,
2404  "This method is not defined or valid for this class type");
2405  }
2406 
2408  const Array<OneD, NekDouble> &TotField,
2409  int BndID)
2410  {
2411  ASSERTL0(false,
2412  "This method is not defined or valid for this class type");
2413  }
2414 
2417  Array<OneD, NekDouble> &outarray,
2418  int BndID)
2419  {
2420  ASSERTL0(false,
2421  "This method is not defined or valid for this class type");
2422  }
2423 
2425  {
2426  ASSERTL0(false,
2427  "This method is not defined or valid for this class type");
2428  }
2429 
2430  /**
2431  */
2433  {
2434  ASSERTL0(false,
2435  "This method is not defined or valid for this class type");
2436  }
2437 
2439  {
2440  ASSERTL0(false,
2441  "This method is not defined or valid for this class type");
2442  }
2443 
2445  {
2446  ASSERTL0(false,
2447  "This method is not defined or valid for this class type");
2448  }
2449 
2450 
2452  Array<OneD, NekDouble> &outarray,
2453  CoeffState coeffstate)
2454  {
2455  v_BwdTrans_IterPerExp(inarray,outarray);
2456  }
2457 
2459  Array<OneD, NekDouble> &outarray,
2460  CoeffState coeffstate)
2461  {
2462  v_FwdTrans_IterPerExp(inarray,outarray);
2463  }
2464 
2466  const Array<OneD, const NekDouble> &inarray,
2467  Array<OneD, NekDouble> &outarray,
2468  CoeffState coeffstate)
2469  {
2471  for (int i = 0; i < m_collections.size(); ++i)
2472  {
2473 
2475  inarray + m_coll_phys_offset[i],
2476  tmp = outarray + m_coll_coeff_offset[i]);
2477  }
2478  }
2479 
2481  const GlobalMatrixKey &gkey,
2482  const Array<OneD,const NekDouble> &inarray,
2483  Array<OneD, NekDouble> &outarray,
2484  CoeffState coeffstate)
2485  {
2486  GeneralMatrixOp_IterPerExp(gkey,inarray,outarray);
2487  }
2488 
2489  /**
2490  * The operation is evaluated locally by the elemental
2491  * function StdRegions#StdExpansion#GetCoords.
2492  *
2493  * @param coord_0 After calculation, the \f$x_1\f$ coordinate
2494  * will be stored in this array.
2495  * @param coord_1 After calculation, the \f$x_2\f$ coordinate
2496  * will be stored in this array.
2497  * @param coord_2 After calculation, the \f$x_3\f$ coordinate
2498  * will be stored in this array.
2499  */
2501  Array<OneD, NekDouble> &coord_1,
2502  Array<OneD, NekDouble> &coord_2)
2503  {
2504  if (GetNumElmts() == 0)
2505  {
2506  return;
2507  }
2508 
2509  int i;
2510  Array<OneD, NekDouble> e_coord_0;
2511  Array<OneD, NekDouble> e_coord_1;
2512  Array<OneD, NekDouble> e_coord_2;
2513 
2514  switch(GetExp(0)->GetCoordim())
2515  {
2516  case 1:
2517  for(i= 0; i < (*m_exp).size(); ++i)
2518  {
2519  e_coord_0 = coord_0 + m_phys_offset[i];
2520  (*m_exp)[i]->GetCoords(e_coord_0);
2521  }
2522  break;
2523  case 2:
2524  ASSERTL0(coord_1.num_elements() != 0,
2525  "output coord_1 is not defined");
2526 
2527  for(i= 0; i < (*m_exp).size(); ++i)
2528  {
2529  e_coord_0 = coord_0 + m_phys_offset[i];
2530  e_coord_1 = coord_1 + m_phys_offset[i];
2531  (*m_exp)[i]->GetCoords(e_coord_0,e_coord_1);
2532  }
2533  break;
2534  case 3:
2535  ASSERTL0(coord_1.num_elements() != 0,
2536  "output coord_1 is not defined");
2537  ASSERTL0(coord_2.num_elements() != 0,
2538  "output coord_2 is not defined");
2539 
2540  for(i= 0; i < (*m_exp).size(); ++i)
2541  {
2542  e_coord_0 = coord_0 + m_phys_offset[i];
2543  e_coord_1 = coord_1 + m_phys_offset[i];
2544  e_coord_2 = coord_2 + m_phys_offset[i];
2545  (*m_exp)[i]->GetCoords(e_coord_0,e_coord_1,e_coord_2);
2546  }
2547  break;
2548  }
2549  }
2550 
2551  /**
2552  */
2554  {
2555  ASSERTL0(false,
2556  "This method is not defined or valid for this class type");
2557  }
2558 
2559  /**
2560  */
2562  Array<OneD,int> &EdgeID)
2563  {
2564  ASSERTL0(false,
2565  "This method is not defined or valid for this class type");
2566  }
2567 
2568  /**
2569  */
2571  {
2572  ASSERTL0(false,
2573  "This method is not defined or valid for this class type");
2574  }
2575 
2576  /**
2577  */
2580  {
2581  ASSERTL0(false,
2582  "This method is not defined or valid for this class type");
2584  result;
2585  return result;
2586  }
2587 
2588  /**
2589  */
2591  {
2592  ASSERTL0(false,
2593  "This method is not defined or valid for this class type");
2595  return result;
2596  }
2597 
2598  /**
2599  */
2601  const NekDouble time,
2602  const std::string varName,
2603  const NekDouble x2_in,
2604  const NekDouble x3_in)
2605  {
2606  ASSERTL0(false,
2607  "This method is not defined or valid for this class type");
2608  }
2609 
2610  /**
2611  */
2612  map<int, RobinBCInfoSharedPtr> ExpList::v_GetRobinBCInfo(void)
2613  {
2614  ASSERTL0(false,
2615  "This method is not defined or valid for this class type");
2616  static map<int,RobinBCInfoSharedPtr> result;
2617  return result;
2618  }
2619 
2620  /**
2621  */
2623  PeriodicMap &periodicVerts,
2624  PeriodicMap &periodicEdges,
2625  PeriodicMap &periodicFaces)
2626  {
2627  ASSERTL0(false,
2628  "This method is not defined or valid for this class type");
2629  }
2630 
2633  unsigned int regionId,
2634  const std::string& variable)
2635  {
2636  SpatialDomains::BoundaryConditionCollection::const_iterator collectionIter = collection.find(regionId);
2637  ASSERTL1(collectionIter != collection.end(), "Unable to locate collection "+boost::lexical_cast<string>(regionId));
2638  const SpatialDomains::BoundaryConditionMapShPtr boundaryConditionMap = (*collectionIter).second;
2639  SpatialDomains::BoundaryConditionMap::const_iterator conditionMapIter = boundaryConditionMap->find(variable);
2640  ASSERTL1(conditionMapIter != boundaryConditionMap->end(), "Unable to locate condition map.");
2641  const SpatialDomains::BoundaryConditionShPtr boundaryCondition = (*conditionMapIter).second;
2642  return boundaryCondition;
2643  }
2644 
2646  {
2647  ASSERTL0(false,
2648  "This method is not defined or valid for this class type");
2649  return NullExpListSharedPtr;
2650  }
2651 
2652 
2654  {
2655 
2657 
2658  switch(exp->DetShapeType())
2659  {
2662  ::AllocateSharedPtr(exp->GetBasis(0)->GetBasisKey());
2663  break;
2665  {
2667  if((nexp = exp->as<StdRegions::StdNodalTriExp>()))
2668  {
2670  ::AllocateSharedPtr(exp->GetBasis(0)->GetBasisKey(),
2671  exp->GetBasis(1)->GetBasisKey(),
2672  nexp->GetNodalPointsKey().GetPointsType());
2673  }
2674  else
2675  {
2677  ::AllocateSharedPtr(exp->GetBasis(0)->GetBasisKey(),
2678  exp->GetBasis(1)->GetBasisKey());
2679  }
2680  }
2681  break;
2684  ::AllocateSharedPtr(exp->GetBasis(0)->GetBasisKey(),
2685  exp->GetBasis(1)->GetBasisKey());
2686  break;
2689  ::AllocateSharedPtr(exp->GetBasis(0)->GetBasisKey(),
2690  exp->GetBasis(1)->GetBasisKey(),
2691  exp->GetBasis(2)->GetBasisKey());
2692  break;
2695  ::AllocateSharedPtr(exp->GetBasis(0)->GetBasisKey(),
2696  exp->GetBasis(1)->GetBasisKey(),
2697  exp->GetBasis(2)->GetBasisKey());
2698  break;
2699  case LibUtilities::ePrism:
2701  ::AllocateSharedPtr(exp->GetBasis(0)->GetBasisKey(),
2702  exp->GetBasis(1)->GetBasisKey(),
2703  exp->GetBasis(2)->GetBasisKey());
2704  break;
2707  ::AllocateSharedPtr(exp->GetBasis(0)->GetBasisKey(),
2708  exp->GetBasis(1)->GetBasisKey(),
2709  exp->GetBasis(2)->GetBasisKey());
2710  break;
2711  default:
2712  ASSERTL0(false,"Shape type not setup");
2713  break;
2714  }
2715 
2716  return stdExp;
2717  }
2718 
2719  /**
2720  * @brief Construct collections of elements containing a single element
2721  * type and polynomial order from the list of expansions.
2722  */
2723  void ExpList::CreateCollections(Collections::ImplementationType ImpType)
2724  {
2726  vector<std::pair<LocalRegions::ExpansionSharedPtr,int> > > collections;
2728  vector<std::pair<LocalRegions::ExpansionSharedPtr,int> > >::iterator it;
2729 
2730  // Figure out optimisation parameters if provided in
2731  // session file or default given
2732  Collections::CollectionOptimisation colOpt(m_session, ImpType);
2733  ImpType = colOpt.GetDefaultImplementationType();
2734 
2735  bool autotuning = colOpt.IsUsingAutotuning();
2736  bool verbose = (m_session->DefinesCmdLineArgument("verbose")) &&
2737  (m_comm->GetRank() == 0);
2738  int collmax = (colOpt.GetMaxCollectionSize() > 0
2739  ? colOpt.GetMaxCollectionSize()
2740  : 2*m_exp->size());
2741 
2742  // clear vectors in case previously called
2743  m_collections.clear();
2744  m_coll_coeff_offset.clear();
2745  m_coll_phys_offset.clear();
2746 
2747  // Loop over expansions, and create collections for each element type
2748  for (int i = 0; i < m_exp->size(); ++i)
2749  {
2750  collections[(*m_exp)[i]->DetShapeType()].push_back(
2751  std::pair<LocalRegions::ExpansionSharedPtr,int> ((*m_exp)[i],i));
2752  }
2753 
2754  for (it = collections.begin(); it != collections.end(); ++it)
2755  {
2756  LocalRegions::ExpansionSharedPtr exp = it->second[0].first;
2757 
2758  Collections::OperatorImpMap impTypes = colOpt.GetOperatorImpMap(exp);
2759  vector<StdRegions::StdExpansionSharedPtr> collExp;
2760 
2761  int prevCoeffOffset = m_coeff_offset[it->second[0].second];
2762  int prevPhysOffset = m_phys_offset [it->second[0].second];
2763  int collcnt;
2764 
2765  m_coll_coeff_offset.push_back(prevCoeffOffset);
2766  m_coll_phys_offset .push_back(prevPhysOffset);
2767 
2768  if(it->second.size() == 1) // single element case
2769  {
2770  collExp.push_back(it->second[0].first);
2771 
2772  // if no Imp Type provided and No settign in xml file.
2773  // reset impTypes using timings
2774  if(autotuning)
2775  {
2776  impTypes = colOpt.SetWithTimings(collExp,
2777  impTypes, verbose);
2778  }
2779 
2780  Collections::Collection tmp(collExp, impTypes);
2781  m_collections.push_back(tmp);
2782  }
2783  else
2784  {
2785  // set up first geometry
2786  collExp.push_back(it->second[0].first);
2787  int prevnCoeff = it->second[0].first->GetNcoeffs();
2788  int prevnPhys = it->second[0].first->GetTotPoints();
2789  collcnt = 1;
2790 
2791  for (int i = 1; i < it->second.size(); ++i)
2792  {
2793  int nCoeffs = it->second[i].first->GetNcoeffs();
2794  int nPhys = it->second[i].first->GetTotPoints();
2795  int coeffOffset = m_coeff_offset[it->second[i].second];
2796  int physOffset = m_phys_offset [it->second[i].second];
2797 
2798  // check to see if next elmt is different or
2799  // collmax reached and if so end collection
2800  // and start new one
2801  if(prevCoeffOffset + nCoeffs != coeffOffset ||
2802  prevnCoeff != nCoeffs ||
2803  prevPhysOffset + nPhys != physOffset ||
2804  prevnPhys != nPhys || collcnt >= collmax)
2805  {
2806 
2807  // if no Imp Type provided and No
2808  // settign in xml file. reset
2809  // impTypes using timings
2810  if(autotuning)
2811  {
2812  impTypes = colOpt.SetWithTimings(collExp,
2813  impTypes,
2814  verbose);
2815  }
2816 
2817  Collections::Collection tmp(collExp, impTypes);
2818  m_collections.push_back(tmp);
2819 
2820 
2821  // start new geom list
2822  collExp.clear();
2823 
2824  m_coll_coeff_offset.push_back(coeffOffset);
2825  m_coll_phys_offset .push_back(physOffset);
2826  collExp.push_back(it->second[i].first);
2827  collcnt = 1;
2828  }
2829  else // add to list of collections
2830  {
2831  collExp.push_back(it->second[i].first);
2832  collcnt++;
2833  }
2834 
2835  // if end of list finish up collection
2836  if (i == it->second.size() - 1)
2837  {
2838  // if no Imp Type provided and No
2839  // settign in xml file.
2840  if(autotuning)
2841  {
2842  impTypes = colOpt.SetWithTimings(collExp,
2843  impTypes,verbose);
2844  }
2845 
2846  Collections::Collection tmp(collExp, impTypes);
2847  m_collections.push_back(tmp);
2848  collExp.clear();
2849  collcnt = 0;
2850 
2851  }
2852 
2853  prevCoeffOffset = coeffOffset;
2854  prevPhysOffset = physOffset;
2855  prevnCoeff = nCoeffs;
2856  prevnPhys = nPhys;
2857  }
2858  }
2859  }
2860  }
2861  } //end of namespace
2862 } //end of namespace
2863 
virtual NekDouble v_L2(const Array< OneD, const NekDouble > &phys, const Array< OneD, const NekDouble > &soln=NullNekDouble1DArray)
Definition: ExpList.cpp:1766
virtual void v_GlobalToLocal(void)
Definition: ExpList.cpp:2444
virtual void v_AddFwdBwdTraceIntegral(const Array< OneD, const NekDouble > &Fwd, const Array< OneD, const NekDouble > &Bwd, Array< OneD, NekDouble > &outarray)
Definition: ExpList.cpp:2298
virtual void v_SmoothField(Array< OneD, NekDouble > &field)
Definition: ExpList.cpp:594
virtual void v_PhysDeriv(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &out_d0, Array< OneD, NekDouble > &out_d1, Array< OneD, NekDouble > &out_d2)
Definition: ExpList.cpp:448
virtual void v_GeneralMatrixOp(const GlobalMatrixKey &gkey, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate)
Definition: ExpList.cpp:2480
const StdRegions::VarCoeffMap & GetVarCoeffs() const
GlobalSysSolnType GetGlobalSysSolnType() const
Return the associated solution type.
const DNekScalBlkMatSharedPtr & GetBlockMatrix(const GlobalMatrixKey &gkey)
Definition: ExpList.cpp:772
virtual void v_ExtractDataToCoeffs(LibUtilities::FieldDefinitionsSharedPtr &fielddef, std::vector< NekDouble > &fielddata, std::string &field, Array< OneD, NekDouble > &coeffs)
Extract data from raw field data into expansion list.
Definition: ExpList.cpp:2111
void SetY(typename boost::call_traits< DataType >::const_reference val)
Definition: NekPoint.hpp:224
int GetExpIndex(const Array< OneD, const NekDouble > &gloCoord, NekDouble tol=0.0, bool returnNearestElmt=false)
Definition: ExpList.cpp:1222
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
virtual void v_Reset()
Reset geometry information, metrics, matrix managers and geometry information.
Definition: ExpList.cpp:1394
boost::shared_ptr< GlobalMatrix > GenGlobalMatrix(const GlobalMatrixKey &mkey, const boost::shared_ptr< AssemblyMapCG > &locToGloMap)
Generates a global matrix from the given key and map.
Definition: ExpList.cpp:855
#define NEKERROR(type, msg)
Assert Level 0 – Fundamental assert which is used whether in FULLDEBUG, DEBUG or OPT compilation mod...
Definition: ErrorUtil.hpp:158
void GeneralMatrixOp_IterPerExp(const GlobalMatrixKey &gkey, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Definition: ExpList.cpp:787
virtual const Array< OneD, const boost::shared_ptr< ExpList > > & v_GetBndCondExpansions(void)
Definition: ExpList.cpp:2217
boost::shared_ptr< Transposition > TranspositionSharedPtr
ExpList()
The default constructor.
Definition: ExpList.cpp:93
static ExpListSharedPtr NullExpListSharedPtr
Definition: ExpList.h:1343
tBaseSharedPtr CreateInstance(tKey idKey BOOST_PP_COMMA_IF(MAX_PARAM) BOOST_PP_ENUM_BINARY_PARAMS(MAX_PARAM, tParam, x))
Create an instance of the class referred to by idKey.
Definition: NekFactory.hpp:162
static Array< OneD, NekDouble > NullNekDouble1DArray
static boost::shared_ptr< DataType > AllocateSharedPtr()
Allocate a shared pointer from the memory pool.
std::map< int, vector< PeriodicEntity > > PeriodicMap
NekOptimize::GlobalOptParamSharedPtr m_globalOptParam
Definition: ExpList.h:971
boost::shared_ptr< AssemblyMap > AssemblyMapSharedPtr
Definition: AssemblyMap.h:53
std::map< CoordType, NekDouble > COOMatType
virtual boost::shared_ptr< AssemblyMapDG > & v_GetTraceMap()
Definition: ExpList.cpp:2261
T Vmax(int n, const T *x, const int incx)
Return the maximum element in x – called vmax to avoid conflict with max.
Definition: Vmath.cpp:756
General purpose memory allocation routines with the ability to allocate from thread specific memory p...
virtual void v_FillBndCondFromField()
Definition: ExpList.cpp:2432
void ExtractFileBCs(const std::string &fileName, const std::string &varName, const boost::shared_ptr< ExpList > locExpList)
Definition: ExpList.cpp:1862
boost::shared_ptr< RobinBCInfo > RobinBCInfoSharedPtr
boost::shared_ptr< DNekMat > GenGlobalMatrixFull(const GlobalLinSysKey &mkey, const boost::shared_ptr< AssemblyMapCG > &locToGloMap)
Definition: ExpList.cpp:992
NekDouble PhysIntegral(void)
This function integrates a function over the domain consisting of all the elements of the expansion...
Definition: ExpList.cpp:236
virtual void v_WriteTecplotZone(std::ostream &outfile, int expansion)
Definition: ExpList.cpp:1460
void Import(const std::string &infilename, std::vector< FieldDefinitionsSharedPtr > &fielddefs, std::vector< std::vector< NekDouble > > &fielddata=NullVectorNekDoubleVector, FieldMetaDataMap &fieldinfomap=NullFieldMetaDataMap, const Array< OneD, int > ElementiDs=NullInt1DArray)
Imports an FLD file.
Definition: FieldIO.cpp:383
ExpansionType GetExpType(void)
Returns the type of the expansion.
Definition: ExpList.cpp:203
boost::shared_ptr< FieldDefinitions > FieldDefinitionsSharedPtr
Definition: FieldIO.h:131
virtual void v_IProductWRTBase(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate)
Definition: ExpList.cpp:2465
const boost::shared_ptr< LocalRegions::ExpansionVector > GetExp() const
This function returns the vector of elements in the expansion.
Definition: ExpList.h:1858
virtual void v_BwdTrans_IterPerExp(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Definition: ExpList.cpp:1189
StdRegions::StdExpansionSharedPtr GetStdExp(StdRegions::StdExpansionSharedPtr exp)
Definition: ExpList.cpp:2653
Lagrange Polynomials using the Gauss points .
Definition: BasisType.h:54
void BwdTrans(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate=eLocal)
Definition: ExpList.h:1557
virtual void v_GetNormals(Array< OneD, Array< OneD, NekDouble > > &normals)
Definition: ExpList.cpp:2274
Array< OneD, const unsigned int > GetZIDs(void)
This function returns a vector containing the wave numbers in z-direction associated with the 3D homo...
Definition: ExpList.h:503
map< int, RobinBCInfoSharedPtr > GetRobinBCInfo()
Definition: ExpList.h:752
std::map< ConstFactorType, NekDouble > ConstFactorMap
Definition: StdRegions.hpp:248
void GeneralGetFieldDefinitions(std::vector< LibUtilities::FieldDefinitionsSharedPtr > &fielddef, int NumHomoDir=0, int NumHomoStrip=1, Array< OneD, LibUtilities::BasisSharedPtr > &HomoBasis=LibUtilities::NullBasisSharedPtr1DArray, std::vector< NekDouble > &HomoLen=LibUtilities::NullNekDoubleVector, std::vector< unsigned int > &HomoZIDs=LibUtilities::NullUnsignedIntVector, std::vector< unsigned int > &HomoYIDs=LibUtilities::NullUnsignedIntVector)
Definition: ExpList.cpp:1937
BlockMatrixMapShPtr m_blockMat
Definition: ExpList.h:973
Array< OneD, NekDouble > m_phys
The global expansion evaluated at the quadrature points.
Definition: ExpList.h:926
void ApplyGeomInfo()
Apply geometry information to each expansion.
Definition: ExpList.cpp:1381
const StdRegions::ConstFactorMap & GetConstFactors() const
Returns all the constants.
Array< OneD, NekDouble > m_coeffs
Concatenation of all local expansion coefficients.
Definition: ExpList.h:909
virtual void v_WriteVtkPieceHeader(std::ostream &outfile, int expansion, int istrip)
Definition: ExpList.cpp:1679
virtual void v_NormVectorIProductWRTBase(Array< OneD, const NekDouble > &V1, Array< OneD, const NekDouble > &V2, Array< OneD, NekDouble > &outarray, int BndID)
Definition: ExpList.cpp:2415
virtual LibUtilities::TranspositionSharedPtr v_GetTransposition(void)
Definition: ExpList.cpp:1818
boost::shared_ptr< DNekMat > DNekMatSharedPtr
Definition: NekTypeDefs.hpp:70
boost::shared_ptr< GlobalMatrix > GlobalMatrixSharedPtr
Shared pointer to a GlobalMatrix object.
Definition: GlobalMatrix.h:89
boost::shared_ptr< SessionReader > SessionReaderSharedPtr
Definition: MeshPartition.h:50
boost::shared_ptr< StdNodalTriExp > StdNodalTriExpSharedPtr
virtual boost::shared_ptr< ExpList > & v_UpdateBndCondExpansion(int i)
Definition: ExpList.cpp:2225
boost::shared_ptr< GlobalLinSys > GenGlobalBndLinSys(const GlobalLinSysKey &mkey, const AssemblyMapSharedPtr &locToGloMap)
Generate a GlobalLinSys from information provided by the key "mkey" and the mapping provided in LocTo...
Definition: ExpList.cpp:1152
std::vector< ExpansionSharedPtr > ExpansionVector
Definition: Expansion.h:70
boost::shared_ptr< DNekScalMat > DNekScalMatSharedPtr
void MultiplyByElmtInvMass(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
This function elementally mulplies the coefficient space of Sin my the elemental inverse of the mass ...
Definition: ExpList.cpp:523
virtual void v_WriteTecplotHeader(std::ostream &outfile, std::string var="")
Definition: ExpList.cpp:1421
virtual NekDouble v_GetHomoLen(void)
Definition: ExpList.cpp:1826
void WriteVtkFooter(std::ostream &outfile)
Definition: ExpList.cpp:1673
virtual void v_GetCoords(Array< OneD, NekDouble > &coord_0, Array< OneD, NekDouble > &coord_1, Array< OneD, NekDouble > &coord_2=NullNekDouble1DArray)
Definition: ExpList.cpp:2500
int GetNumberOfCoefficients(ShapeType shape, std::vector< unsigned int > &modes, int offset)
Definition: ShapeType.hpp:312
map< GlobalMatrixKey, DNekScalBlkMatSharedPtr > BlockMatrixMap
A map between global matrix keys and their associated block matrices.
Definition: ExpList.h:95
boost::shared_ptr< ExpList > GetSharedThisPtr()
Returns a shared pointer to the current object.
Definition: ExpList.h:831
virtual Array< OneD, const NekDouble > v_HomogeneousEnergy(void)
Definition: ExpList.cpp:1810
virtual const Array< OneD, const int > & v_GetTraceBndMap()
Definition: ExpList.cpp:2269
virtual void v_WriteVtkPieceData(std::ostream &outfile, int expansion, std::string var)
Definition: ExpList.cpp:1690
const char *const GlobalSysSolnTypeMap[]
virtual void v_GetBCValues(Array< OneD, NekDouble > &BndVals, const Array< OneD, NekDouble > &TotField, int BndID)
Definition: ExpList.cpp:2407
Array< OneD, int > m_coeff_offset
Offset of elemental data into the array m_coeffs.
Definition: ExpList.h:958
int GetNumElmts(void)
This function returns the number of elements in the expansion which may be different for a homogeoeno...
Definition: ExpList.h:562
void SetZ(typename boost::call_traits< DataType >::const_reference val)
Definition: NekPoint.hpp:230
Base class for all multi-elemental spectral/hp expansions.
Definition: ExpList.h:101
static const NekDouble kNekZeroTol
void GetCoords(Array< OneD, NekDouble > &coord_0, Array< OneD, NekDouble > &coord_1=NullNekDouble1DArray, Array< OneD, NekDouble > &coord_2=NullNekDouble1DArray)
This function calculates the coordinates of all the elemental quadrature points . ...
Definition: ExpList.h:1632
virtual void v_ImposeDirichletConditions(Array< OneD, NekDouble > &outarray)
Definition: ExpList.cpp:2424
virtual void v_Upwind(const Array< OneD, const Array< OneD, NekDouble > > &Vec, const Array< OneD, const NekDouble > &Fwd, const Array< OneD, const NekDouble > &Bwd, Array< OneD, NekDouble > &Upwind)
Definition: ExpList.cpp:2233
NekDouble Linf(const Array< OneD, const NekDouble > &inarray, const Array< OneD, const NekDouble > &soln=NullNekDouble1DArray)
This function calculates the error of the global spectral/hp element approximation.
Definition: ExpList.cpp:1727
boost::shared_ptr< LocalRegions::ExpansionVector > m_exp
The list of local expansions.
Definition: ExpList.h:947
int GetTotPoints(void) const
Returns the total number of quadrature points m_npoints .
Definition: ExpList.h:1401
NekMatrix< NekMatrix< NekMatrix< NekDouble, StandardMatrixTag >, ScaledMatrixTag >, BlockMatrixTag > DNekScalBlkMat
Definition: NekTypeDefs.hpp:66
virtual void v_FwdTrans_IterPerExp(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Definition: ExpList.cpp:562
virtual ~ExpList()
The default destructor.
Definition: ExpList.cpp:216
void IProductWRTDerivBase(const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
This function calculates the inner product of a function with respect to the derivative (in directio...
Definition: ExpList.cpp:341
std::map< StdRegions::VarCoeffType, Array< OneD, NekDouble > > VarCoeffMap
Definition: StdRegions.hpp:225
boost::shared_ptr< GlobalLinSys > GenGlobalLinSys(const GlobalLinSysKey &mkey, const boost::shared_ptr< AssemblyMapCG > &locToGloMap)
This operation constructs the global linear system of type mkey.
Definition: ExpList.cpp:1133
Array< OneD, int > m_phys_offset
Offset of elemental data into the array m_phys.
Definition: ExpList.h:961
boost::shared_ptr< ExpList > ExpListSharedPtr
Shared pointer to an ExpList object.
Definition: ExpList.h:1340
void IProductWRTBase_IterPerExp(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
This function calculates the inner product of a function with respect to all {local} expansion modes...
Definition: ExpList.h:1518
virtual NekDouble v_Integral(const Array< OneD, const NekDouble > &inarray)
Definition: ExpList.cpp:1796
virtual std::vector< LibUtilities::FieldDefinitionsSharedPtr > v_GetFieldDefinitions(void)
Definition: ExpList.cpp:2047
bool m_physState
The state of the array m_phys.
Definition: ExpList.h:935
virtual void v_WriteTecplotConnectivity(std::ostream &outfile, int expansion)
Definition: ExpList.cpp:1554
int m_ncoeffs
The total number of local degrees of freedom. m_ncoeffs .
Definition: ExpList.h:887
Array< OneD, int > m_offset_elmt_id
Array containing the element id m_offset_elmt_id[n] that the n^th consecutive block of data in m_coef...
Definition: ExpList.h:969
boost::shared_ptr< DNekScalBlkMat > DNekScalBlkMatSharedPtr
Definition: NekTypeDefs.hpp:74
virtual void v_HomogeneousFwdTrans(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate=eLocal, bool Shuff=true, bool UnShuff=true)
Definition: ExpList.cpp:2381
virtual void v_LocalToGlobal(void)
Definition: ExpList.cpp:2438
LibUtilities::ShapeType GetShapeType() const
Return the expansion type associated with key.
SpatialDomains::MeshGraphSharedPtr m_graph
Mesh associated with this expansion list.
Definition: ExpList.h:883
virtual void v_EvaluateBoundaryConditions(const NekDouble time=0.0, const std::string varName="", const NekDouble x2_in=NekConstants::kNekUnsetDouble, const NekDouble x3_in=NekConstants::kNekUnsetDouble)
Definition: ExpList.cpp:2600
LibUtilities::SessionReaderSharedPtr m_session
Session.
Definition: ExpList.h:880
double NekDouble
const DNekScalBlkMatSharedPtr GenBlockMatrix(const GlobalMatrixKey &gkey)
This function assembles the block diagonal matrix of local matrices of the type mtype.
Definition: ExpList.cpp:636
Defines a list of flags.
std::map< int, BoundaryConditionMapShPtr > BoundaryConditionCollection
Definition: Conditions.h:213
virtual Array< OneD, SpatialDomains::BoundaryConditionShPtr > & v_UpdateBndConditions()
Definition: ExpList.cpp:2590
Describe a linear system.
void FwdTrans_BndConstrained(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Definition: ExpList.cpp:572
virtual boost::shared_ptr< ExpList > & v_GetTrace()
Definition: ExpList.cpp:2253
virtual void v_DealiasedProd(const Array< OneD, NekDouble > &inarray1, const Array< OneD, NekDouble > &inarray2, Array< OneD, NekDouble > &outarray, CoeffState coeffstate=eLocal)
Definition: ExpList.cpp:2401
virtual void v_ExtractCoeffsToCoeffs(const boost::shared_ptr< ExpList > &fromExpList, const Array< OneD, const NekDouble > &fromCoeffs, Array< OneD, NekDouble > &toCoeffs)
Definition: ExpList.cpp:2186
virtual boost::shared_ptr< ExpList > & v_GetPlane(int n)
Definition: ExpList.cpp:2645
StdRegions::MatrixType GetMatrixType() const
Return the matrix type.
Describes a matrix with ordering defined by a local to global map.
virtual void v_GetPeriodicEntities(PeriodicMap &periodicVerts, PeriodicMap &periodicEdges, PeriodicMap &periodicFaces)
Definition: ExpList.cpp:2622
virtual void v_PhysInterp1DScaled(const NekDouble scale, const Array< OneD, NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Definition: ExpList.cpp:1851
MultiRegions::Direction const DirCartesianMap[]
Definition: ExpList.h:86
StandardMatrixTag boost::call_traits< LhsDataType >::const_reference rhs typedef NekMatrix< LhsDataType, StandardMatrixTag >::iterator iterator
LibUtilities::CommSharedPtr m_comm
Communicator.
Definition: ExpList.h:877
virtual void v_AddTraceIntegral(const Array< OneD, const NekDouble > &Fx, const Array< OneD, const NekDouble > &Fy, Array< OneD, NekDouble > &outarray)
Definition: ExpList.cpp:2281
void ExtractCoeffsToCoeffs(const boost::shared_ptr< ExpList > &fromExpList, const Array< OneD, const NekDouble > &fromCoeffs, Array< OneD, NekDouble > &toCoeffs)
Extract the data from fromField using fromExpList the coeffs using the basic ExpList Elemental expans...
Definition: ExpList.cpp:2098
boost::shared_ptr< BoundaryConditionMap > BoundaryConditionMapShPtr
Definition: Conditions.h:212
virtual void v_LinearAdvectionReactionSolve(const Array< OneD, Array< OneD, NekDouble > > &velocity, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const NekDouble lambda, CoeffState coeffstate=eLocal, const Array< OneD, const NekDouble > &dirForcing=NullNekDouble1DArray)
Definition: ExpList.cpp:2369
virtual Array< OneD, const unsigned int > v_GetZIDs(void)
Definition: ExpList.cpp:1834
virtual map< int, RobinBCInfoSharedPtr > v_GetRobinBCInfo(void)
Definition: ExpList.cpp:2612
#define WARNINGL1(condition, msg)
Definition: ErrorUtil.hpp:192
void SetX(typename boost::call_traits< DataType >::const_reference val)
Definition: NekPoint.hpp:218
NekDouble H1(const Array< OneD, const NekDouble > &inarray, const Array< OneD, const NekDouble > &soln=NullNekDouble1DArray)
Calculates the error of the global spectral/hp element approximation.
Definition: ExpList.cpp:1918
virtual void v_IProductWRTBase_IterPerExp(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Definition: ExpList.cpp:314
boost::shared_ptr< Expansion > ExpansionSharedPtr
Definition: Expansion.h:68
void WriteVtkPieceFooter(std::ostream &outfile, int expansion)
Definition: ExpList.cpp:1684
virtual void v_LinearAdvectionDiffusionReactionSolve(const Array< OneD, Array< OneD, NekDouble > > &velocity, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const NekDouble lambda, CoeffState coeffstate=eLocal, const Array< OneD, const NekDouble > &dirForcing=NullNekDouble1DArray)
Definition: ExpList.cpp:2357
virtual void v_HelmSolve(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const FlagList &flags, const StdRegions::ConstFactorMap &factors, const StdRegions::VarCoeffMap &varcoeff, const Array< OneD, const NekDouble > &dirForcing)
Definition: ExpList.cpp:2346
std::vector< int > m_coll_phys_offset
Offset of elemental data into the array m_phys.
Definition: ExpList.h:955
virtual void v_AppendFieldData(LibUtilities::FieldDefinitionsSharedPtr &fielddef, std::vector< NekDouble > &fielddata)
Definition: ExpList.cpp:2061
Class for operating on FLD files.
Definition: FieldIO.h:151
static SpatialDomains::BoundaryConditionShPtr GetBoundaryCondition(const SpatialDomains::BoundaryConditionCollection &collection, unsigned int index, const std::string &variable)
Definition: ExpList.cpp:2631
boost::shared_ptr< GlobalLinSys > GlobalLinSysSharedPtr
Pointer to a GlobalLinSys object.
Definition: GlobalLinSys.h:51
boost::shared_ptr< AssemblyMapDG > & GetTraceMap(void)
Definition: ExpList.h:1947
void CreateCollections(Collections::ImplementationType ImpType=Collections::eNoImpType)
Construct collections of elements containing a single element type and polynomial order from the list...
Definition: ExpList.cpp:2723
void MultiplyByBlockMatrix(const GlobalMatrixKey &gkey, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Definition: ExpList.cpp:284
virtual void v_MultiplyByInvMassMatrix(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate)
Definition: ExpList.cpp:2337
virtual void v_GetBoundaryToElmtMap(Array< OneD, int > &ElmtID, Array< OneD, int > &EdgeID)
Definition: ExpList.cpp:2561
virtual Array< OneD, const unsigned int > v_GetYIDs(void)
Definition: ExpList.cpp:1842
Used to lookup the create function in NekManager.
Definition: MatrixKey.h:68
virtual void v_SetUpPhysNormals()
Definition: ExpList.cpp:2553
int GetCoordim(int eid)
This function returns the dimension of the coordinates of the element eid.
Definition: ExpList.h:1735
boost::shared_ptr< BoundaryConditionBase > BoundaryConditionShPtr
Definition: Conditions.h:206
virtual const Array< OneD, const SpatialDomains::BoundaryConditionShPtr > & v_GetBndConditions()
Definition: ExpList.cpp:2579
Lagrange for SEM basis .
Definition: BasisType.h:53
void SetExpType(ExpansionType Type)
Returns the type of the expansion.
Definition: ExpList.cpp:211
boost::shared_ptr< AssemblyMapCG > AssemblyMapCGSharedPtr
Definition: AssemblyMapCG.h:52
GlobalLinSysFactory & GetGlobalLinSysFactory()
boost::shared_ptr< StdExpansion > StdExpansionSharedPtr
std::vector< int > m_coll_coeff_offset
Offset of elemental data into the array m_coeffs.
Definition: ExpList.h:952
void WriteVtkHeader(std::ostream &outfile)
Definition: ExpList.cpp:1665
std::pair< IndexType, IndexType > CoordType
virtual void v_ExtractTracePhys(Array< OneD, NekDouble > &outarray)
Definition: ExpList.cpp:2323
virtual void v_PhysGalerkinProjection1DScaled(const NekDouble scale, const Array< OneD, NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Definition: ExpList.cpp:1857
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode...
Definition: ErrorUtil.hpp:191
boost::shared_ptr< MeshGraph > MeshGraphSharedPtr
Definition: MeshGraph.h:432
virtual void v_BwdTrans(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate)
Definition: ExpList.cpp:2451
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1038
NekMatrix< NekMatrix< NekDouble, StandardMatrixTag >, ScaledMatrixTag > DNekScalMat
virtual void v_GetFwdBwdTracePhys(Array< OneD, NekDouble > &Fwd, Array< OneD, NekDouble > &Bwd)
Definition: ExpList.cpp:2307
void ExtractDataToCoeffs(LibUtilities::FieldDefinitionsSharedPtr &fielddef, std::vector< NekDouble > &fielddata, std::string &field, Array< OneD, NekDouble > &coeffs)
Extract the data in fielddata into the coeffs.
Definition: ExpList.cpp:2089
virtual void v_WriteTecplotField(std::ostream &outfile, int expansion)
Definition: ExpList.cpp:1631
boost::shared_ptr< PointGeom > PointGeomSharedPtr
Definition: Geometry.h:60
Collections::CollectionVector m_collections
Definition: ExpList.h:949
virtual void v_ReadGlobalOptimizationParameters()
Definition: ExpList.cpp:2570
virtual void v_HomogeneousBwdTrans(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate=eLocal, bool Shuff=true, bool UnShuff=true)
Definition: ExpList.cpp:2391
virtual void v_FwdTrans(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate)
Definition: ExpList.cpp:2458