Nektar++
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  : ExpList3DHomogeneous1D(pSession, HomoBasis, lhom, useFFT,
98  dealiasing),
99  m_bndCondExpansions(),
100  m_bndConditions()
101  {
102  int i, n, nel;
103  DisContField2DSharedPtr plane_zero;
105 
106  // note that nzplanes can be larger than nzmodes
107  m_planes[0] = plane_zero = MemoryManager<DisContField2D>::
108  AllocateSharedPtr(pSession, graph2D, variable, true, false);
109 
112 
113  nel = m_planes[0]->GetExpSize();
114 
115  for (i = 0; i < nel; ++i)
116  {
117  (*m_exp).push_back(m_planes[0]->GetExp(i));
118  }
119 
120  for (n = 1; n < m_planes.num_elements(); ++n)
121  {
123  AllocateSharedPtr(*plane_zero, graph2D,
124  variable, true, false);
125  for(i = 0; i < nel; ++i)
126  {
127  (*m_exp).push_back((*m_exp)[i]);
128  }
129  }
130 
131  // Set up trace object.
132  Array<OneD, ExpListSharedPtr> trace(m_planes.num_elements());
133  for (n = 0; n < m_planes.num_elements(); ++n)
134  {
135  trace[n] = m_planes[n]->GetTrace();
136  }
137 
139  pSession, HomoBasis, lhom, useFFT, dealiasing, trace);
140 
141  // Setup default optimisation information
142  nel = GetExpSize();
143 
145  AllocateSharedPtr(nel);
146 
147  SetCoeffPhys();
148 
149  // Do not set up BCs if default variable
150  if(variable.compare("DefaultVar") != 0)
151  {
152  SetupBoundaryConditions(HomoBasis, lhom, bcs, variable);
153  }
154 
155  SetUpDG();
156  }
157 
158  /**
159  * @brief Default destructor
160  */
162  {
163  }
164 
166  const LibUtilities::BasisKey &HomoBasis,
167  const NekDouble lhom,
169  const std::string variable)
170  {
171  int n;
172 
173  // Setup an ExpList2DHomogeneous1D expansion for boundary
174  // conditions and link to class declared in m_planes
176  bcs.GetBoundaryRegions();
178  bcs.GetBoundaryConditions();
179  SpatialDomains::BoundaryRegionCollection::const_iterator it;
180 
181  // count the number of non-periodic boundary regions
182  int cnt = 0;
183  for (it = bregions.begin(); it != bregions.end(); ++it)
184  {
185  SpatialDomains::BoundaryConditionShPtr boundaryCondition =
186  GetBoundaryCondition(bconditions, it->first, variable);
187  if (boundaryCondition->GetBoundaryConditionType()
189  {
190  cnt++;
191  }
192  }
193 
195  ExpListSharedPtr>(cnt);
196  m_bndConditions = m_planes[0]->UpdateBndConditions();
197 
198  int nplanes = m_planes.num_elements();
200  PlanesBndCondExp(nplanes);
201 
202  cnt = 0;
203  for (it = bregions.begin(); it != bregions.end(); ++it)
204  {
205  SpatialDomains::BoundaryConditionShPtr boundaryCondition =
206  GetBoundaryCondition(bconditions, it->first, variable);
207  if(boundaryCondition->GetBoundaryConditionType() !=
209  {
210  for (n = 0; n < nplanes; ++n)
211  {
212  PlanesBndCondExp[n] = m_planes[n]->
214  }
215 
216  m_bndCondExpansions[cnt++] =
218  AllocateSharedPtr(m_session, HomoBasis, lhom,
219  m_useFFT, false,
220  PlanesBndCondExp);
221  }
222  }
223  EvaluateBoundaryConditions(0.0, variable);
224  }
225 
227  const NekDouble time,
228  const std::string varName)
229  {
230  int n;
232  Array<OneD, NekDouble> local_z(m_planes.num_elements());
233 
234  for (n = 0; n < m_planes.num_elements(); n++)
235  {
236  local_z[n] = z[m_transposition->GetPlaneID(n)];
237  }
238 
239  for (n = 0; n < m_planes.num_elements(); ++n)
240  {
241  m_planes[n]->EvaluateBoundaryConditions(
242  time, varName, 0.5*m_lhom*(1.0+local_z[n]));
243  }
244 
245  // Fourier transform coefficient space boundary values
246  // This will only be undertaken for time dependent
247  // boundary conditions unless time == 0.0 which is the
248  // case when the method is called from the constructor.
249  for (n = 0; n < m_bndCondExpansions.num_elements(); ++n)
250  {
251  if (time == 0.0 ||
252  m_bndConditions[n]->IsTimeDependent() ||
253  boost::iequals(m_bndConditions[n]->GetUserDefined(),
254  "MovingBody"))
255  {
256  m_bndCondExpansions[n]->HomogeneousFwdTrans(
259  }
260  }
261  }
262 
264  const Array<OneD, const NekDouble> &inarray,
265  Array<OneD, NekDouble> &outarray,
266  const FlagList &flags,
267  const StdRegions::ConstFactorMap &factors,
268  const StdRegions::VarCoeffMap &varcoeff,
269  const Array<OneD, const NekDouble> &dirForcing)
270  {
271  int n;
272  int cnt = 0;
273  int cnt1 = 0;
274  NekDouble beta;
275  StdRegions::ConstFactorMap new_factors;
276 
278  Array<OneD, NekDouble> fce(inarray.num_elements());
279 
280  // Transform forcing function in half-physical space if required
281  if (m_WaveSpace)
282  {
283  fce = inarray;
284  }
285  else
286  {
287  HomogeneousFwdTrans(inarray,fce);
288  }
289 
290  for (n = 0; n < m_planes.num_elements(); ++n)
291  {
292  if (n != 1 || m_transposition->GetK(n) != 0)
293  {
294  beta = 2*M_PI*(m_transposition->GetK(n))/m_lhom;
295  new_factors = factors;
296  // add in Homogeneous Fourier direction and SVV if turned on
297  new_factors[StdRegions::eFactorLambda] +=
298  beta*beta*(1+GetSpecVanVisc(n));
299  m_planes[n]->HelmSolve(
300  fce + cnt,
301  e_out = outarray + cnt1,
302  flags, new_factors, varcoeff, dirForcing);
303  }
304 
305  cnt += m_planes[n]->GetTotPoints();
306  cnt1 += m_planes[n]->GetNcoeffs();
307  }
308  }
309 
311  const NekDouble time,
312  const std::string varName,
313  const NekDouble x2_in,
314  const NekDouble x3_in)
315  {
316  EvaluateBoundaryConditions(time, varName);
317  }
318 
319  boost::shared_ptr<ExpList> &DisContField3DHomogeneous1D::
321  {
322  return UpdateBndCondExpansion(i);
323  }
324 
327  {
328  return UpdateBndConditions();
329  }
330 
332  Array<OneD, int> &ElmtID,
333  Array<OneD,int> &EdgeID)
334  {
335 
336  if(m_BCtoElmMap.num_elements() == 0)
337  {
338  Array<OneD, int> ElmtID_tmp;
339  Array<OneD, int> EdgeID_tmp;
340 
341  m_planes[0]->GetBoundaryToElmtMap(ElmtID_tmp, EdgeID_tmp);
342  int nel_per_plane = m_planes[0]->GetExpSize();
343  int nplanes = m_planes.num_elements();
344 
345  int MapSize = ElmtID_tmp.num_elements();
346 
347  ElmtID = Array<OneD, int>(nplanes*MapSize);
348  EdgeID = Array<OneD, int>(nplanes*MapSize);
349 
350  // If this mesh (or partition) has no BCs, skip this step
351  if (MapSize > 0)
352  {
353  for(int i = 0; i < nplanes; i++)
354  {
355  for(int j = 0; j < MapSize; j++)
356  {
357  ElmtID[j+i*MapSize] = ElmtID_tmp[j]+i*nel_per_plane;
358  EdgeID[j+i*MapSize] = EdgeID_tmp[j];
359  }
360  }
361  m_BCtoElmMap = Array<OneD, int>(nplanes*MapSize);
362  m_BCtoEdgMap = Array<OneD, int>(nplanes*MapSize);
363 
364  Vmath::Vcopy(nplanes*MapSize,ElmtID,1,m_BCtoElmMap,1);
365  Vmath::Vcopy(nplanes*MapSize,EdgeID,1,m_BCtoEdgMap,1);
366  }
367  }
368  else
369  {
370  int MapSize = m_BCtoElmMap.num_elements();
371 
372  ElmtID = Array<OneD, int>(MapSize);
373  EdgeID = Array<OneD, int>(MapSize);
374 
375  Vmath::Vcopy(MapSize, m_BCtoElmMap, 1, ElmtID, 1);
376  Vmath::Vcopy(MapSize, m_BCtoEdgMap, 1, EdgeID, 1);
377  }
378  }
379 
381  Array<OneD, NekDouble> &BndVals,
382  const Array<OneD, NekDouble> &TotField,
383  int BndID)
384  {
387 
389  Array<OneD, NekDouble> tmp_BC;
390 
391  int cnt = 0;
392  int pos = 0;
393  int exp_size, exp_size_per_plane, elmtID, boundaryID;
394  int offset, exp_dim;
395 
396  for (int k = 0; k < m_planes.num_elements(); k++)
397  {
398  for (int n = 0; n < m_bndConditions.num_elements(); ++n)
399  {
400  exp_size = m_bndCondExpansions[n]->GetExpSize();
401  exp_size_per_plane = exp_size/m_planes.num_elements();
402 
403  for (int i = 0; i < exp_size_per_plane; i++)
404  {
405  if(n == BndID)
406  {
407  elmtID = m_BCtoElmMap[cnt];
408  boundaryID = m_BCtoEdgMap[cnt];
409  exp_dim = m_bndCondExpansions[n]->
410  GetExp(i+k*exp_size_per_plane)->GetTotPoints();
411  offset = GetPhys_Offset(elmtID);
412  elmt = GetExp(elmtID);
413  temp_BC_exp = boost::dynamic_pointer_cast<
415  m_bndCondExpansions[n]->GetExp(
416  i+k*exp_size_per_plane));
417 
418  elmt->GetEdgePhysVals(boundaryID, temp_BC_exp,
419  tmp_Tot = TotField + offset,
420  tmp_BC = BndVals + pos);
421  pos += exp_dim;
422  }
423  cnt++;
424  }
425  }
426  }
427  }
428 
432  Array<OneD, NekDouble> &outarray,
433  int BndID)
434  {
437 
438  Array<OneD, NekDouble> tmp_V1;
439  Array<OneD, NekDouble> tmp_V2;
440  Array<OneD, NekDouble> tmp_outarray;
441 
442  int cnt = 0;
443  int exp_size, exp_size_per_plane, elmtID, Phys_offset, Coef_offset;
444 
445  for(int k = 0; k < m_planes.num_elements(); k++)
446  {
447  for(int n = 0; n < m_bndConditions.num_elements(); ++n)
448  {
449  exp_size = m_bndCondExpansions[n]->GetExpSize();
450  exp_size_per_plane = exp_size/m_planes.num_elements();
451 
452  for(int i = 0; i < exp_size_per_plane; i++)
453  {
454  if(n == BndID)
455  {
456  elmtID = m_BCtoElmMap[cnt];
457 
458  Phys_offset = m_bndCondExpansions[n]->
459  GetPhys_Offset(i+k*exp_size_per_plane);
460  Coef_offset = m_bndCondExpansions[n]->
461  GetCoeff_Offset(i+k*exp_size_per_plane);
462 
463  elmt = GetExp(elmtID);
464  temp_BC_exp = boost::dynamic_pointer_cast<
466  m_bndCondExpansions[n]->GetExp(
467  i+k*exp_size_per_plane));
468 
469  temp_BC_exp->NormVectorIProductWRTBase(
470  tmp_V1 = V1 + Phys_offset,
471  tmp_V2 = V2 + Phys_offset,
472  tmp_outarray = outarray + Coef_offset);
473  }
474  cnt++;
475  }
476  }
477  }
478  }
479 
481  Array<OneD, NekDouble> &outarray)
482  {
483  ASSERTL1(m_physState == true,
484  "Field must be in physical state to extract trace space.");
485 
486  v_ExtractTracePhys(m_phys, outarray);
487  }
488 
489  /**
490  * @brief This method extracts the trace (edges in 2D) for each plane
491  * from the field @a inarray and puts the values in @a outarray.
492  *
493  * It assumes the field is C0 continuous so that it can overwrite the
494  * edge data when visited by the two adjacent elements.
495  *
496  * @param inarray An array containing the 2D data from which we wish
497  * to extract the edge data.
498  * @param outarray The resulting edge information.
499  */
501  const Array<OneD, const NekDouble> &inarray,
502  Array<OneD, NekDouble> &outarray)
503  {
504  int nPoints_plane = m_planes[0]->GetTotPoints();
505  int nTracePts = m_planes[0]->GetTrace()->GetTotPoints();
506 
507  for (int i = 0; i < m_planes.num_elements(); ++i)
508  {
509  Array<OneD, NekDouble> inarray_plane(nPoints_plane, 0.0);
510  Array<OneD, NekDouble> outarray_plane(nPoints_plane, 0.0);
511 
512  Vmath::Vcopy(nPoints_plane,
513  &inarray[i*nPoints_plane], 1,
514  &inarray_plane[0], 1);
515 
516  m_planes[i]->ExtractTracePhys(inarray_plane, outarray_plane);
517 
518  Vmath::Vcopy(nTracePts,
519  &outarray_plane[0], 1,
520  &outarray[i*nTracePts], 1);
521  }
522  }
523 
524  /**
525  * @brief Set up all DG member variables and maps.
526  */
528  {
529  const int nPlanes = m_planes.num_elements();
530  const int nTracePlane = m_planes[0]->GetTrace()->GetExpSize();
531 
532  // Get trace map from first plane.
533  AssemblyMapDGSharedPtr traceMap = m_planes[0]->GetTraceMap();
534  const Array<OneD, const int> &traceBndMap
535  = traceMap->GetBndCondTraceToGlobalTraceMap();
536  int mapSize = traceBndMap.num_elements();
537 
538  // Set up trace boundary map
539  m_traceBndMap = Array<OneD, int>(nPlanes * mapSize);
540 
541  int i, n, e, cnt = 0, cnt1 = 0;
542 
543  for (i = 0; i < m_bndCondExpansions.num_elements(); ++i)
544  {
545  int nExp = m_bndCondExpansions[i]->GetExpSize();
546  int nPlaneExp = nExp / nPlanes;
547 
548  for (n = 0; n < nPlanes; ++n)
549  {
550  const int offset = n * nTracePlane;
551  for (e = 0; e < nPlaneExp; ++e)
552  {
553  m_traceBndMap[cnt++] = offset + traceBndMap[cnt1+e];
554  }
555  }
556 
557  cnt1 += nPlaneExp;
558  }
559  }
560  } // end of namespace
561 } //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:1775
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:1867
static boost::shared_ptr< DataType > AllocateSharedPtr()
Allocate a shared pointer from the memory pool.
NekOptimize::GlobalOptParamSharedPtr m_globalOptParam
Definition: ExpList.h:971
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:1875
const BoundaryConditionCollection & GetBoundaryConditions(void) const
Definition: Conditions.h:235
Array< OneD, NekDouble > & UpdateCoeffs()
This function returns (a reference to) the array (implemented as m_coeffs) containing all local expa...
Definition: ExpList.h:1894
void SetUpDG()
Set up all DG member variables and maps.
void NormVectorIProductWRTBase(const Array< OneD, const NekDouble > &Fx, const Array< OneD, const NekDouble > &Fy, Array< OneD, NekDouble > &outarray)
Definition: StdExpansion.h:711
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:1858
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)
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:248
NekDouble m_lhom
Width of homogeneous direction.
Array< OneD, NekDouble > m_phys
The global expansion evaluated at the quadrature points.
Definition: ExpList.h:926
int GetExpSize(void)
This function returns the number of elements in the expansion.
Definition: ExpList.h:1837
boost::shared_ptr< SessionReader > SessionReaderSharedPtr
Definition: MeshPartition.h:50
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.
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:947
std::map< StdRegions::VarCoeffType, Array< OneD, NekDouble > > VarCoeffMap
Definition: StdRegions.hpp:225
boost::shared_ptr< ExpList > ExpListSharedPtr
Shared pointer to an ExpList object.
Definition: ExpList.h:1340
bool m_physState
The state of the array m_phys.
Definition: ExpList.h:935
std::map< int, BoundaryRegionShPtr > BoundaryRegionCollection
Definition: Conditions.h:204
boost::shared_ptr< DisContField2D > DisContField2DSharedPtr
boost::shared_ptr< StdExpansion1D > StdExpansion1DSharedPtr
LibUtilities::SessionReaderSharedPtr m_session
Session.
Definition: ExpList.h:880
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:213
void SetupBoundaryConditions(const LibUtilities::BasisKey &HomoBasis, const NekDouble lhom, SpatialDomains::BoundaryConditions &bcs, const std::string variable)
Array< OneD, SpatialDomains::BoundaryConditionShPtr > & UpdateBndConditions()
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:843
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:2631
const BoundaryRegionCollection & GetBoundaryRegions(void) const
Definition: Conditions.h:225
boost::shared_ptr< BoundaryConditionBase > BoundaryConditionShPtr
Definition: Conditions.h:206
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:191
boost::shared_ptr< MeshGraph > MeshGraphSharedPtr
Definition: MeshGraph.h:432
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1038
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...