Nektar++
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
ContField3DHomogeneous1D.cpp
Go to the documentation of this file.
1 //////////////////////////////////////////////////////////////////////////////
2 //
3 // File ContField3DHomogeneous1D.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 and a 1D homogeneous direction
34 //
35 ///////////////////////////////////////////////////////////////////////////////
36 
39 
40 namespace Nektar
41 {
42  namespace MultiRegions
43  {
44 
47  {
48  }
49 
51  const ContField3DHomogeneous1D &In):
53  {
54 
55  bool False = false;
56  ContField2DSharedPtr zero_plane =
57  boost::dynamic_pointer_cast<ContField2D> (In.m_planes[0]);
58 
59  for(int n = 0; n < m_planes.num_elements(); ++n)
60  {
62  AllocateSharedPtr(*zero_plane,False);
63  }
64 
65  SetCoeffPhys();
66  }
67 
69  const ContField3DHomogeneous1D &In,
71  const std::string &variable):
72  DisContField3DHomogeneous1D (In, false)
73  {
74  ContField2DSharedPtr zero_plane_old =
75  boost::dynamic_pointer_cast<ContField2D> (In.m_planes[0]);
76 
77  ContField2DSharedPtr zero_plane =
79  AllocateSharedPtr(*zero_plane_old,graph2D,
80  variable);
81 
82  for(int n = 0; n < m_planes.num_elements(); ++n)
83  {
85  AllocateSharedPtr(*zero_plane,graph2D,
86  variable);
87  }
88 
89  SetCoeffPhys();
90 
91  if(variable.compare("DefaultVar") != 0)
92  {
95  m_lhom, bcs,variable);
96  }
97  }
98 
100  {
101  }
102 
104  const LibUtilities::SessionReaderSharedPtr &pSession,
105  const LibUtilities::BasisKey &HomoBasis,
106  const NekDouble lhom,
107  const bool useFFT,
108  const bool dealiasing,
109  const SpatialDomains::MeshGraphSharedPtr &graph2D,
110  const std::string &variable,
111  const bool CheckIfSingularSystem,
112  const Collections::ImplementationType ImpType):
113  DisContField3DHomogeneous1D(pSession,HomoBasis,lhom,
114  useFFT,dealiasing)
115  {
116  int i,n,nel;
117  ContField2DSharedPtr plane_zero;
118  ContField2DSharedPtr plane_two;
119 
121  m_graph = graph2D;
122 
123  // Plane zero (k=0 - cos) - singularaty check required for Poisson
124  // problems
126  pSession, graph2D, variable, false,
127  CheckIfSingularSystem, ImpType);
128 
130  pSession, graph2D, variable, false,
131  false, ImpType);
132 
135 
136  for(n = 0; n < m_planes.num_elements(); ++n)
137  {
138  // Plane zero and one (k=0 - cos and sin) - singularaty check
139  // required for Poisson problems
140  if(m_transposition->GetK(n) == 0)
141  {
143  ::AllocateSharedPtr(*plane_zero, graph2D, variable,
144  false, CheckIfSingularSystem);
145  }
146  else
147  {
148  // For k > 0 singularty check not required anymore -
149  // creating another ContField2D to avoid Assembly Map copy
150  // TODO: We may want to deal with it in a more efficient
151  // way in the future.
153  ::AllocateSharedPtr(*plane_two, graph2D, variable,
154  false, false);
155  }
156 
157  nel = m_planes[n]->GetExpSize();
158 
159  for(i = 0; i < nel; ++i)
160  {
161  (*m_exp).push_back(m_planes[n]->GetExp(i));
162  }
163  }
164 
165  nel = GetExpSize();
166 
168  AllocateSharedPtr(nel);
169 
170  SetCoeffPhys();
171 
172  // Do not set up BCs if default variable
173  if(variable.compare("DefaultVar") != 0)
174  {
175  SetupBoundaryConditions(HomoBasis,lhom,bcs,variable);
176  }
177  }
178 
179 
181  Array<OneD,NekDouble>& outarray)
182  {
184  int ncoeffs = m_planes[0]->GetNcoeffs();
185 
186  for(int n = 0; n < m_planes.num_elements(); ++n)
187  {
188  m_planes[n]->ImposeDirichletConditions(tmp = outarray +
189  n*ncoeffs);
190  }
191  }
192 
194  {
195  for(int n = 0; n < m_planes.num_elements(); ++n)
196  {
197  m_planes[n]->FillBndCondFromField();
198  }
199  }
200 
202  {
203  for(int n = 0; n < m_planes.num_elements(); ++n)
204  {
205  m_planes[n]->FillBndCondFromField(nreg);
206  }
207  }
208 
209  /**
210  *
211  */
213  {
214  for(int n = 0; n < m_planes.num_elements(); ++n)
215  {
216  m_planes[n]->LocalToGlobal(useComm);
217  }
218  };
219 
220 
221  /**
222  *
223  */
225  {
226  for(int n = 0; n < m_planes.num_elements(); ++n)
227  {
228  m_planes[n]->GlobalToLocal();
229  }
230  };
231 
232 
233  /**
234  *
235  */
237  Array<OneD,NekDouble> &field)
238  {
239  int cnt = 0;
241 
242  for(int n = 0; n < m_planes.num_elements(); ++n)
243  {
244  m_planes[n]->SmoothField(tmp = field + cnt);
245 
246  cnt += m_planes[n]->GetTotPoints();
247  }
248  }
249 
250 
252  const Array<OneD, const NekDouble> &inarray,
253  Array<OneD, NekDouble> &outarray,
254  const FlagList &flags,
255  const StdRegions::ConstFactorMap &factors,
256  const StdRegions::VarCoeffMap &varcoeff,
257  const Array<OneD, const NekDouble> &dirForcing,
258  const bool PhysSpaceForcing)
259  {
260 
261  int n;
262  int cnt = 0;
263  int cnt1 = 0;
264  NekDouble beta;
265  StdRegions::ConstFactorMap new_factors;
266 
268  Array<OneD, NekDouble> fce(inarray.num_elements());
270 
271  // Fourier transform forcing function
272  if(m_WaveSpace)
273  {
274  fce = inarray;
275  }
276  else
277  {
278  HomogeneousFwdTrans(inarray, fce,
279  (flags.isSet(eUseGlobal))?eGlobal:eLocal);
280  }
281 
282  bool smode = false;
283 
284  if (m_homogeneousBasis->GetBasisType() ==
286  m_homogeneousBasis->GetBasisType() ==
288  {
289  smode = true;
290  }
291 
292  for(n = 0; n < m_planes.num_elements(); ++n)
293  {
294  if(n != 1 || m_transposition->GetK(n) != 0 || smode)
295  {
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(wfce,
305  e_out = outarray + cnt1,
306  flags, new_factors, varcoeff,
307  dirForcing,
308  PhysSpaceForcing);
309  }
310 
311  cnt += m_planes[n]->GetTotPoints();
312  cnt1 += m_planes[n]->GetNcoeffs();
313  }
314  }
315 
316  /**
317  * Reset the GlobalLinSys Manager
318  */
320  {
321  for(int n = 0; n < m_planes.num_elements(); ++n)
322  {
323  m_planes[n]->ClearGlobalLinSysManager();
324  }
325  }
326 
327  } // end of namespace
328 } //end of namespace
virtual void v_SmoothField(Array< OneD, NekDouble > &field)
Local coefficients.
virtual void v_ImposeDirichletConditions(Array< OneD, NekDouble > &outarray)
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
const boost::shared_ptr< LocalRegions::ExpansionVector > GetExp() const
This function returns the vector of elements in the expansion.
Definition: ExpList.h:2067
boost::shared_ptr< ContField2D > ContField2DSharedPtr
Definition: ContField2D.h:289
std::map< ConstFactorType, NekDouble > ConstFactorMap
Definition: StdRegions.hpp:252
NekDouble m_lhom
Width of homogeneous direction.
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
Global coefficients.
void SetCoeffPhys(void)
Definition of the total number of degrees of freedom and quadrature points. Sets up the storage for m...
bool isSet(const FlagType &key) const
LibUtilities::BasisSharedPtr m_homogeneousBasis
Definition of the total number of degrees of freedom and quadrature points. Sets up the storage for m...
Fourier Modified expansions with just the real part of the first mode .
Definition: BasisType.h:59
boost::shared_ptr< LocalRegions::ExpansionVector > m_exp
The list of local expansions.
Definition: ExpList.h:1036
This class is the abstraction of a global continuous two- dimensional spectral/hp element expansion w...
Definition: ContField2D.h:56
std::map< StdRegions::VarCoeffType, Array< OneD, NekDouble > > VarCoeffMap
Definition: StdRegions.hpp:227
SpatialDomains::MeshGraphSharedPtr m_graph
Mesh associated with this expansion list.
Definition: ExpList.h:972
LibUtilities::SessionReaderSharedPtr m_session
Session.
Definition: ExpList.h:969
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)
Fourier Modified expansions with just the imaginary part of the first mode .
Definition: BasisType.h:60
void SetupBoundaryConditions(const LibUtilities::BasisKey &HomoBasis, const NekDouble lhom, SpatialDomains::BoundaryConditions &bcs, const std::string variable)
virtual void v_GlobalToLocal(void)
Template method virtual forwarded for GlobalToLocal()
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)
Solves the three-dimensional Helmholtz equation, subject to the boundary conditions specified...
boost::shared_ptr< MeshGraph > MeshGraphSharedPtr
Definition: MeshGraph.h:442
Describes the specification for a Basis.
Definition: Basis.h:50
virtual void v_LocalToGlobal(bool useComm)
Template method virtual forwarded for LocalToGlobal()