Nektar++
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
DisContField3DHomogeneous1D.cpp
Go to the documentation of this file.
1 //////////////////////////////////////////////////////////////////////////////
2 //
3 // File DisContField3DHomogeneous1D.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: Field definition for 3D domain with boundary
33 // conditions using LDG flux and a 1D homogeneous direction
34 //
35 ///////////////////////////////////////////////////////////////////////////////
36 
40 
41 
42 namespace Nektar
43 {
44  namespace MultiRegions
45  {
46 
49  m_bndCondExpansions(),
50  m_bndConditions()
51  {
52  }
53 
56  const LibUtilities::BasisKey &HomoBasis,
57  const NekDouble lhom,
58  const bool useFFT,
59  const bool dealiasing):
60  ExpList3DHomogeneous1D(pSession,HomoBasis,lhom,useFFT,dealiasing),
61  m_bndCondExpansions(),
62  m_bndConditions()
63  {
64  }
65 
68  const bool DeclarePlanesSetCoeffPhys)
69  : ExpList3DHomogeneous1D (In,false),
70  m_bndCondExpansions (In.m_bndCondExpansions),
71  m_bndConditions (In.m_bndConditions)
72  {
73  if (DeclarePlanesSetCoeffPhys)
74  {
75  DisContField2DSharedPtr zero_plane =
76  boost::dynamic_pointer_cast<DisContField2D> (In.m_planes[0]);
77 
78  for(int n = 0; n < m_planes.num_elements(); ++n)
79  {
80  m_planes[n] =
82  AllocateSharedPtr(*zero_plane, false);
83  }
84 
85  SetCoeffPhys();
86  }
87  }
88 
91  const LibUtilities::BasisKey &HomoBasis,
92  const NekDouble lhom,
93  const bool useFFT,
94  const bool dealiasing,
96  const std::string &variable,
97  const Collections::ImplementationType ImpType):
98  ExpList3DHomogeneous1D(pSession, HomoBasis, lhom, useFFT,
99  dealiasing),
100  m_bndCondExpansions(),
101  m_bndConditions()
102  {
103  int i, n, nel;
104  DisContField2DSharedPtr plane_zero;
106 
107  // note that nzplanes can be larger than nzmodes
108  m_planes[0] = plane_zero = MemoryManager<DisContField2D>::
109  AllocateSharedPtr(pSession, graph2D, variable, true, false,
110  ImpType);
111 
114 
115  nel = m_planes[0]->GetExpSize();
116 
117  for (i = 0; i < nel; ++i)
118  {
119  (*m_exp).push_back(m_planes[0]->GetExp(i));
120  }
121 
122  for (n = 1; n < m_planes.num_elements(); ++n)
123  {
125  AllocateSharedPtr(*plane_zero, graph2D,
126  variable, true, false);
127  for(i = 0; i < nel; ++i)
128  {
129  (*m_exp).push_back((*m_exp)[i]);
130  }
131  }
132 
133  // Set up trace object.
134  Array<OneD, ExpListSharedPtr> trace(m_planes.num_elements());
135  for (n = 0; n < m_planes.num_elements(); ++n)
136  {
137  trace[n] = m_planes[n]->GetTrace();
138  }
139 
141  pSession, HomoBasis, lhom, useFFT,
142  dealiasing, trace);
143 
144  // Setup default optimisation information
145  nel = GetExpSize();
146 
148  AllocateSharedPtr(nel);
149 
150  SetCoeffPhys();
151 
152  // Do not set up BCs if default variable
153  if(variable.compare("DefaultVar") != 0)
154  {
155  SetupBoundaryConditions(HomoBasis, lhom, bcs, variable);
156  }
157 
158  SetUpDG();
159  }
160 
161  /**
162  * @brief Default destructor
163  */
165  {
166  }
167 
169  const LibUtilities::BasisKey &HomoBasis,
170  const NekDouble lhom,
172  const std::string variable)
173  {
174  int n;
175 
176  // Setup an ExpList2DHomogeneous1D expansion for boundary
177  // conditions and link to class declared in m_planes
179  bcs.GetBoundaryRegions();
181  bcs.GetBoundaryConditions();
182  SpatialDomains::BoundaryRegionCollection::const_iterator it;
183 
184  // count the number of non-periodic boundary regions
185  int cnt = 0;
186  for (it = bregions.begin(); it != bregions.end(); ++it)
187  {
188  SpatialDomains::BoundaryConditionShPtr boundaryCondition =
189  GetBoundaryCondition(bconditions, it->first, variable);
190  if (boundaryCondition->GetBoundaryConditionType()
192  {
193  cnt++;
194  }
195  }
196 
198  ExpListSharedPtr>(cnt);
199  m_bndConditions = m_planes[0]->UpdateBndConditions();
200 
201  int nplanes = m_planes.num_elements();
203  PlanesBndCondExp(nplanes);
204 
205  cnt = 0;
206  for (it = bregions.begin(); it != bregions.end(); ++it)
207  {
208  SpatialDomains::BoundaryConditionShPtr boundaryCondition =
209  GetBoundaryCondition(bconditions, it->first, variable);
210  if(boundaryCondition->GetBoundaryConditionType() !=
212  {
213  for (n = 0; n < nplanes; ++n)
214  {
215  PlanesBndCondExp[n] = m_planes[n]->
217  }
218 
219  m_bndCondExpansions[cnt++] =
221  AllocateSharedPtr(m_session, HomoBasis, lhom,
222  m_useFFT, false,
223  PlanesBndCondExp);
224  }
225  }
226  EvaluateBoundaryConditions(0.0, variable);
227  }
228 
230  const NekDouble time,
231  const std::string varName)
232  {
233  int n;
235  Array<OneD, NekDouble> local_z(m_planes.num_elements());
236 
237  for (n = 0; n < m_planes.num_elements(); n++)
238  {
239  local_z[n] = z[m_transposition->GetPlaneID(n)];
240  }
241 
242  for (n = 0; n < m_planes.num_elements(); ++n)
243  {
244  m_planes[n]->EvaluateBoundaryConditions(
245  time, varName, 0.5*m_lhom*(1.0+local_z[n]));
246  }
247 
248  // Fourier transform coefficient space boundary values
249  // This will only be undertaken for time dependent
250  // boundary conditions unless time == 0.0 which is the
251  // case when the method is called from the constructor.
252  for (n = 0; n < m_bndCondExpansions.num_elements(); ++n)
253  {
254  if (time == 0.0 ||
255  m_bndConditions[n]->IsTimeDependent() )
256  {
257  m_bndCondExpansions[n]->HomogeneousFwdTrans(
260  }
261  }
262  }
263 
265  const Array<OneD, const NekDouble> &inarray,
266  Array<OneD, NekDouble> &outarray,
267  const FlagList &flags,
268  const StdRegions::ConstFactorMap &factors,
269  const StdRegions::VarCoeffMap &varcoeff,
270  const Array<OneD, const NekDouble> &dirForcing,
271  const bool PhysSpaceForcing)
272  {
273  int n;
274  int cnt = 0;
275  int cnt1 = 0;
276  NekDouble beta;
277  StdRegions::ConstFactorMap new_factors;
278 
280  Array<OneD, NekDouble> fce(inarray.num_elements());
282 
283  // Transform forcing function in half-physical space if required
284  if (m_WaveSpace)
285  {
286  fce = inarray;
287  }
288  else
289  {
290  HomogeneousFwdTrans(inarray,fce);
291  }
292 
293  for (n = 0; n < m_planes.num_elements(); ++n)
294  {
295  if (n != 1 || m_transposition->GetK(n) != 0)
296  {
297  beta = 2*M_PI*(m_transposition->GetK(n))/m_lhom;
298  new_factors = factors;
299  // add in Homogeneous Fourier direction and SVV if turned on
300  new_factors[StdRegions::eFactorLambda] +=
301  beta*beta*(1+GetSpecVanVisc(n));
302 
303  wfce = (PhysSpaceForcing)? fce+cnt:fce+cnt1;
304  m_planes[n]->HelmSolve(
305  wfce,
306  e_out = outarray + cnt1,
307  flags, new_factors, varcoeff, dirForcing,
308  PhysSpaceForcing);
309  }
310 
311  cnt += m_planes[n]->GetTotPoints();
312  cnt1 += m_planes[n]->GetNcoeffs();
313  }
314  }
315 
317  const NekDouble time,
318  const std::string varName,
319  const NekDouble x2_in,
320  const NekDouble x3_in)
321  {
322  EvaluateBoundaryConditions(time, varName);
323  }
324 
325  boost::shared_ptr<ExpList> &DisContField3DHomogeneous1D::
327  {
328  return UpdateBndCondExpansion(i);
329  }
330 
333  {
334  return UpdateBndConditions();
335  }
336 
338  Array<OneD, int> &ElmtID,
339  Array<OneD,int> &EdgeID)
340  {
341 
342  if(m_BCtoElmMap.num_elements() == 0)
343  {
344  Array<OneD, int> ElmtID_tmp;
345  Array<OneD, int> EdgeID_tmp;
346 
347  m_planes[0]->GetBoundaryToElmtMap(ElmtID_tmp, EdgeID_tmp);
348  int nel_per_plane = m_planes[0]->GetExpSize();
349  int nplanes = m_planes.num_elements();
350 
351  int MapSize = ElmtID_tmp.num_elements();
352 
353  m_BCtoElmMap = Array<OneD, int>(nplanes*MapSize);
354  m_BCtoEdgMap = Array<OneD, int>(nplanes*MapSize);
355 
356  // If this mesh (or partition) has no BCs, skip this step
357  if (MapSize > 0)
358  {
359  int i ,j, n, cnt;
360  int cntPlane = 0;
361  for (cnt=n=0; n < m_bndCondExpansions.num_elements(); ++n)
362  {
363  int planeExpSize = m_planes[0]
364  ->GetBndCondExpansions()[n]
365  ->GetExpSize();
366  for (i = 0; i < planeExpSize ; ++i, ++cntPlane)
367  {
368  for(j = 0; j < nplanes; j++)
369  {
370  m_BCtoElmMap[cnt+i+j*planeExpSize] =
371  ElmtID_tmp[cntPlane]+j*nel_per_plane;
372  m_BCtoEdgMap[cnt+i+j*planeExpSize] =
373  EdgeID_tmp[cntPlane];
374  }
375  }
376  cnt += m_bndCondExpansions[n]->GetExpSize();
377  }
378  }
379  }
380  ElmtID = m_BCtoElmMap;
381  EdgeID = m_BCtoEdgMap;
382  }
383 
385  boost::shared_ptr<ExpList> &result,
386  const bool DeclareCoeffPhysArrays)
387  {
388  int n, cnt, nq;
389  int offsetOld, offsetNew;
390 
391  std::vector<unsigned int> eIDs;
392  Array<OneD, int> ElmtID,EdgeID;
393  GetBoundaryToElmtMap(ElmtID,EdgeID);
394 
395  // Skip other boundary regions
396  for (cnt = n = 0; n < i; ++n)
397  {
398  cnt += m_bndCondExpansions[n]->GetExpSize();
399  }
400 
401  // Populate eIDs with information from BoundaryToElmtMap
402  for (n = 0; n < m_bndCondExpansions[i]->GetExpSize(); ++n)
403  {
404  eIDs.push_back(ElmtID[cnt+n]);
405  }
406 
407  // Create expansion list
408  result =
410  (*this, eIDs);
411 
412  // Copy phys and coeffs to new explist
413  if ( DeclareCoeffPhysArrays)
414  {
415  Array<OneD, NekDouble> tmp1, tmp2;
416  for (n = 0; n < result->GetExpSize(); ++n)
417  {
418  nq = GetExp(ElmtID[cnt+n])->GetTotPoints();
419  offsetOld = GetPhys_Offset(ElmtID[cnt+n]);
420  offsetNew = result->GetPhys_Offset(n);
421  Vmath::Vcopy(nq, tmp1 = GetPhys()+ offsetOld, 1,
422  tmp2 = result->UpdatePhys()+ offsetNew, 1);
423 
424  nq = GetExp(ElmtID[cnt+n])->GetNcoeffs();
425  offsetOld = GetCoeff_Offset(ElmtID[cnt+n]);
426  offsetNew = result->GetCoeff_Offset(n);
427  Vmath::Vcopy(nq, tmp1 = GetCoeffs()+ offsetOld, 1,
428  tmp2 = result->UpdateCoeffs()+ offsetNew, 1);
429  }
430  }
431 
432  // Set wavespace value
433  result->SetWaveSpace(GetWaveSpace());
434  }
435 
437  Array<OneD, NekDouble> &BndVals,
438  const Array<OneD, NekDouble> &TotField,
439  int BndID)
440  {
443 
445  Array<OneD, NekDouble> tmp_BC;
446 
447  int cnt = 0;
448  int pos = 0;
449  int exp_size, exp_size_per_plane, elmtID, boundaryID;
450  int offset, exp_dim;
451 
452  for (int k = 0; k < m_planes.num_elements(); k++)
453  {
454  for (int n = 0; n < m_bndConditions.num_elements(); ++n)
455  {
456  exp_size = m_bndCondExpansions[n]->GetExpSize();
457  exp_size_per_plane = exp_size/m_planes.num_elements();
458 
459  for (int i = 0; i < exp_size_per_plane; i++)
460  {
461  if(n == BndID)
462  {
463  elmtID = m_BCtoElmMap[cnt];
464  boundaryID = m_BCtoEdgMap[cnt];
465  exp_dim = m_bndCondExpansions[n]->
466  GetExp(i+k*exp_size_per_plane)->GetTotPoints();
467  offset = GetPhys_Offset(elmtID);
468  elmt = GetExp(elmtID);
469  temp_BC_exp = boost::dynamic_pointer_cast<
471  m_bndCondExpansions[n]->GetExp(
472  i+k*exp_size_per_plane));
473 
474  elmt->GetEdgePhysVals(boundaryID, temp_BC_exp,
475  tmp_Tot = TotField + offset,
476  tmp_BC = BndVals + pos);
477  pos += exp_dim;
478  }
479  cnt++;
480  }
481  }
482  }
483  }
484 
488  Array<OneD, NekDouble> &outarray,
489  int BndID)
490  {
493 
494  Array<OneD, NekDouble> tmp_V1;
495  Array<OneD, NekDouble> tmp_V2;
496  Array<OneD, NekDouble> tmp_outarray;
497 
498  int cnt = 0;
499  int exp_size, exp_size_per_plane, elmtID, Phys_offset, Coef_offset;
500 
501  for(int k = 0; k < m_planes.num_elements(); k++)
502  {
503  for(int n = 0; n < m_bndConditions.num_elements(); ++n)
504  {
505  exp_size = m_bndCondExpansions[n]->GetExpSize();
506  exp_size_per_plane = exp_size/m_planes.num_elements();
507 
508  for(int i = 0; i < exp_size_per_plane; i++)
509  {
510  if(n == BndID)
511  {
512  elmtID = m_BCtoElmMap[cnt];
513 
514  Phys_offset = m_bndCondExpansions[n]->
515  GetPhys_Offset(i+k*exp_size_per_plane);
516  Coef_offset = m_bndCondExpansions[n]->
517  GetCoeff_Offset(i+k*exp_size_per_plane);
518 
519  elmt = GetExp(elmtID);
520  temp_BC_exp = boost::dynamic_pointer_cast<
522  m_bndCondExpansions[n]->GetExp(
523  i+k*exp_size_per_plane));
524 
525  temp_BC_exp->NormVectorIProductWRTBase(
526  tmp_V1 = V1 + Phys_offset,
527  tmp_V2 = V2 + Phys_offset,
528  tmp_outarray = outarray + Coef_offset);
529  }
530  cnt++;
531  }
532  }
533  }
534  }
535 
537  Array<OneD, NekDouble> &outarray)
538  {
539  ASSERTL1(m_physState == true,
540  "Field must be in physical state to extract trace space.");
541 
542  v_ExtractTracePhys(m_phys, outarray);
543  }
544 
545  /**
546  * @brief This method extracts the trace (edges in 2D) for each plane
547  * from the field @a inarray and puts the values in @a outarray.
548  *
549  * It assumes the field is C0 continuous so that it can overwrite the
550  * edge data when visited by the two adjacent elements.
551  *
552  * @param inarray An array containing the 2D data from which we wish
553  * to extract the edge data.
554  * @param outarray The resulting edge information.
555  */
557  const Array<OneD, const NekDouble> &inarray,
558  Array<OneD, NekDouble> &outarray)
559  {
560  int nPoints_plane = m_planes[0]->GetTotPoints();
561  int nTracePts = m_planes[0]->GetTrace()->GetTotPoints();
562 
563  for (int i = 0; i < m_planes.num_elements(); ++i)
564  {
565  Array<OneD, NekDouble> inarray_plane(nPoints_plane, 0.0);
566  Array<OneD, NekDouble> outarray_plane(nPoints_plane, 0.0);
567 
568  Vmath::Vcopy(nPoints_plane,
569  &inarray[i*nPoints_plane], 1,
570  &inarray_plane[0], 1);
571 
572  m_planes[i]->ExtractTracePhys(inarray_plane, outarray_plane);
573 
574  Vmath::Vcopy(nTracePts,
575  &outarray_plane[0], 1,
576  &outarray[i*nTracePts], 1);
577  }
578  }
579 
580  /**
581  */
583  Array<OneD, Array<OneD, NekDouble> > &normals)
584  {
585  int j, n, cnt, nq;
586  int expdim = GetCoordim(0);
587  int coordim = 3;
590 
591  Array<OneD, int> ElmtID,EdgeID;
592  GetBoundaryToElmtMap(ElmtID,EdgeID);
593 
594  // Initialise result
595  normals = Array<OneD, Array<OneD, NekDouble> > (coordim);
596  for (j = 0; j < coordim; ++j)
597  {
598  normals[j] = Array<OneD, NekDouble> (
599  GetBndCondExpansions()[i]->GetTotPoints(), 0.0);
600  }
601 
602  // Skip other boundary regions
603  for (cnt = n = 0; n < i; ++n)
604  {
605  cnt += GetBndCondExpansions()[n]->GetExpSize();
606  }
607 
608  int offset;
609  for (n = 0; n < GetBndCondExpansions()[i]->GetExpSize(); ++n)
610  {
611  offset = GetBndCondExpansions()[i]->GetPhys_Offset(n);
612  nq = GetBndCondExpansions()[i]->GetExp(n)->GetTotPoints();
613 
614  elmt = GetExp(ElmtID[cnt+n]);
615  const Array<OneD, const Array<OneD, NekDouble> > normalsElmt
616  = elmt->GetSurfaceNormal(EdgeID[cnt+n]);
617  // Copy to result
618  for (j = 0; j < expdim; ++j)
619  {
620  Vmath::Vcopy(nq, normalsElmt[j], 1,
621  tmp = normals[j] + offset, 1);
622  }
623  }
624  }
625 
626  /**
627  * @brief Set up all DG member variables and maps.
628  */
630  {
631  const int nPlanes = m_planes.num_elements();
632  const int nTracePlane = m_planes[0]->GetTrace()->GetExpSize();
633 
634  // Get trace map from first plane.
635  AssemblyMapDGSharedPtr traceMap = m_planes[0]->GetTraceMap();
636  const Array<OneD, const int> &traceBndMap
637  = traceMap->GetBndCondTraceToGlobalTraceMap();
638  int mapSize = traceBndMap.num_elements();
639 
640  // Set up trace boundary map
641  m_traceBndMap = Array<OneD, int>(nPlanes * mapSize);
642 
643  int i, n, e, cnt = 0, cnt1 = 0;
644 
645  for (i = 0; i < m_bndCondExpansions.num_elements(); ++i)
646  {
647  int nExp = m_bndCondExpansions[i]->GetExpSize();
648  int nPlaneExp = nExp / nPlanes;
649 
650  for (n = 0; n < nPlanes; ++n)
651  {
652  const int offset = n * nTracePlane;
653  for (e = 0; e < nPlaneExp; ++e)
654  {
655  m_traceBndMap[cnt++] = offset + traceBndMap[cnt1+e];
656  }
657  }
658 
659  cnt1 += nPlaneExp;
660  }
661  }
662  } // end of namespace
663 } //end of namespace
boost::shared_ptr< AssemblyMapDG > AssemblyMapDGSharedPtr
Definition: AssemblyMapDG.h:49
const Array< OneD, const NekDouble > & GetCoeffs() const
This function returns (a reference to) the array (implemented as m_coeffs) containing all local expa...
Definition: ExpList.h:1938
int GetCoeff_Offset(int n) const
Get the start offset position for a global list of m_coeffs correspoinding to element n...
Definition: ExpList.h:2076
static boost::shared_ptr< DataType > AllocateSharedPtr()
Allocate a shared pointer from the memory pool.
NekOptimize::GlobalOptParamSharedPtr m_globalOptParam
Definition: ExpList.h:1052
LibUtilities::TranspositionSharedPtr m_transposition
int GetPhys_Offset(int n) const
Get the start offset position for a global list of m_phys correspoinding to element n...
Definition: ExpList.h:2084
const BoundaryConditionCollection & GetBoundaryConditions(void) const
Definition: Conditions.h:248
Array< OneD, NekDouble > & UpdateCoeffs()
This function returns (a reference to) the array (implemented as m_coeffs) containing all local expa...
Definition: ExpList.h:2095
void SetUpDG()
Set up all DG member variables and maps.
boost::shared_ptr< ExpList > & UpdateBndCondExpansion(int i)
const boost::shared_ptr< LocalRegions::ExpansionVector > GetExp() const
This function returns the vector of elements in the expansion.
Definition: ExpList.h:2067
void NormVectorIProductWRTBase(const Array< OneD, const NekDouble > &Fx, Array< OneD, NekDouble > &outarray)
Definition: StdExpansion.h:733
Array< OneD, int > m_BCtoElmMap
Storage space for the boundary to element and boundary to trace map. This member variable is really a...
std::map< ConstFactorType, NekDouble > ConstFactorMap
Definition: StdRegions.hpp:252
NekDouble m_lhom
Width of homogeneous direction.
Array< OneD, NekDouble > m_phys
The global expansion evaluated at the quadrature points.
Definition: ExpList.h:1015
int GetExpSize(void)
This function returns the number of elements in the expansion.
Definition: ExpList.h:2046
boost::shared_ptr< SessionReader > SessionReaderSharedPtr
Definition: MeshPartition.h:51
virtual void v_ExtractTracePhys(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
This method extracts the trace (edges in 2D) for each plane from the field inarray and puts the value...
Array< OneD, MultiRegions::ExpListSharedPtr > m_bndCondExpansions
An object which contains the discretised boundary conditions.
virtual void v_GetBndElmtExpansion(int i, boost::shared_ptr< ExpList > &result, const bool DeclareCoeffPhysArrays)
void SetCoeffPhys(void)
Definition of the total number of degrees of freedom and quadrature points. Sets up the storage for m...
LibUtilities::BasisSharedPtr m_homogeneousBasis
Definition of the total number of degrees of freedom and quadrature points. Sets up the storage for m...
virtual Array< OneD, SpatialDomains::BoundaryConditionShPtr > & v_UpdateBndConditions()
virtual boost::shared_ptr< ExpList > & v_UpdateBndCondExpansion(int i)
boost::shared_ptr< LocalRegions::ExpansionVector > m_exp
The list of local expansions.
Definition: ExpList.h:1036
std::map< StdRegions::VarCoeffType, Array< OneD, NekDouble > > VarCoeffMap
Definition: StdRegions.hpp:227
bool m_physState
The state of the array m_phys.
Definition: ExpList.h:1024
std::map< int, BoundaryRegionShPtr > BoundaryRegionCollection
Definition: Conditions.h:217
virtual void v_GetBoundaryNormals(int i, Array< OneD, Array< OneD, NekDouble > > &normals)
boost::shared_ptr< ExpList > ExpListSharedPtr
Shared pointer to an ExpList object.
boost::shared_ptr< DisContField2D > DisContField2DSharedPtr
boost::shared_ptr< StdExpansion1D > StdExpansion1DSharedPtr
LibUtilities::SessionReaderSharedPtr m_session
Session.
Definition: ExpList.h:969
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, const bool PhysSpaceForcing)
Array< OneD, ExpListSharedPtr > m_planes
double NekDouble
Defines a list of flags.
void HomogeneousFwdTrans(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate=eLocal, bool Shuff=true, bool UnShuff=true)
std::map< int, BoundaryConditionMapShPtr > BoundaryConditionCollection
Definition: Conditions.h:226
void SetupBoundaryConditions(const LibUtilities::BasisKey &HomoBasis, const NekDouble lhom, SpatialDomains::BoundaryConditions &bcs, const std::string variable)
const Array< OneD, const NekDouble > & GetPhys() const
This function returns (a reference to) the array (implemented as m_phys) containing the function ev...
Definition: ExpList.h:2037
Array< OneD, SpatialDomains::BoundaryConditionShPtr > & UpdateBndConditions()
const Array< OneD, const MultiRegions::ExpListSharedPtr > & GetBndCondExpansions()
void GetEdgePhysVals(const int edge, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Extract the physical values along edge edge from inarray into outarray following the local edge orien...
Definition: StdExpansion.h:876
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)
void GetBCValues(Array< OneD, NekDouble > &BndVals, const Array< OneD, NekDouble > &TotField, int BndID)
This funtion extract form a vector containing a full 3D-homogenous-1D field the value associated with...
void NormVectorIProductWRTBase(Array< OneD, const NekDouble > &V1, Array< OneD, const NekDouble > &V2, Array< OneD, NekDouble > &outarray, int BndID)
This function calculate the inner product of two vectors (V1 and V2) respect to the basis along a bou...
void GetBoundaryToElmtMap(Array< OneD, int > &ElmtID, Array< OneD, int > &EdgeID)
Set up a list of element ids and edge ids the link to the boundary conditions.
static SpatialDomains::BoundaryConditionShPtr GetBoundaryCondition(const SpatialDomains::BoundaryConditionCollection &collection, unsigned int index, const std::string &variable)
Definition: ExpList.cpp:3059
const BoundaryRegionCollection & GetBoundaryRegions(void) const
Definition: Conditions.h:238
bool GetWaveSpace(void) const
This function returns the third direction expansion condition, which can be in wave space (coefficien...
Definition: ExpList.h:1584
int GetCoordim(int eid)
This function returns the dimension of the coordinates of the element eid.
Definition: ExpList.h:1898
boost::shared_ptr< BoundaryConditionBase > BoundaryConditionShPtr
Definition: Conditions.h:219
boost::shared_ptr< StdExpansion > StdExpansionSharedPtr
void EvaluateBoundaryConditions(const NekDouble time=0.0, const std::string varName="")
This function evaluates the boundary conditions at a certaintime-level.
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode...
Definition: ErrorUtil.hpp:228
boost::shared_ptr< MeshGraph > MeshGraphSharedPtr
Definition: MeshGraph.h:442
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1061
Abstraction of a two-dimensional multi-elemental expansion which is merely a collection of local expa...
Describes the specification for a Basis.
Definition: Basis.h:50
Array< OneD, SpatialDomains::BoundaryConditionShPtr > m_bndConditions
An array which contains the information about the boundary condition on the different boundary region...