Nektar++
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
ShallowWaterSystem.cpp
Go to the documentation of this file.
1 ///////////////////////////////////////////////////////////////////////////////
2 //
3 // File ShallowWaterSystem.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: Generic timestepping for shallow water solvers
33 //
34 ///////////////////////////////////////////////////////////////////////////////
35 
36 #include <iostream>
37 
39 
40 using namespace std;
41 
42 namespace Nektar
43 {
44  /**
45  * @class ShallowWaterSystem
46  *
47  * Provides the underlying timestepping framework for shallow water flow solvers
48  * including the general timestepping routines. This class is not intended
49  * to be directly instantiated, but rather is a base class on which to
50  * define shallow water solvers, e.g. SWE, Boussinesq, linear and nonlinear versions.
51  *
52  * For details on implementing unsteady solvers see
53  * \ref sectionADRSolverModuleImplementation
54  */
55 
56  /**
57  * Processes SolverInfo parameters from the session file and sets up
58  * timestepping-specific code.
59  * @param pSession Session object to read parameters from.
60  */
61 
62  string ShallowWaterSystem::className =
64  "ShallowWaterSystem",
65  ShallowWaterSystem::create,
66  "Auxiliary functions for the shallow water system.");
67 
68 
69  ShallowWaterSystem::ShallowWaterSystem(
71  : UnsteadySystem(pSession)
72  {
73  }
74 
76  {
78 
79  // if discontinuous Galerkin determine numerical flux to use
81  {
82  ASSERTL0(m_session->DefinesSolverInfo("UPWINDTYPE"),
83  "No UPWINDTYPE defined in session.");
84  }
85 
86  // Set up locations of velocity vector.
89  for (int i = 0; i < m_spacedim; ++i)
90  {
91  m_vecLocs[0][i] = 1 + i;
92  }
93 
94  // Load generic input parameters
95  m_session->LoadParameter("IO_InfoSteps", m_infosteps, 0);
96 
97  // Load acceleration of gravity
98  m_session->LoadParameter("Gravity", m_g, 9.81);
99 
100  // input/output in primitive variables
101  m_primitive = true;
102 
104 
105  m_constantDepth = true;
106  NekDouble depth = m_depth[0];
107  for (int i = 0; i < GetTotPoints(); ++i)
108  {
109  if (m_depth[i] != depth)
110  {
111  m_constantDepth = false;
112  break;
113  }
114  }
115 
116  // Compute the bottom slopes
117  int nq = GetTotPoints();
118  if (m_constantDepth != true)
119  {
121  for (int i = 0; i < m_spacedim; ++i)
122  {
125  Vmath::Neg(nq,m_bottomSlope[i],1);
126  }
127  }
128 
130 
131  }
132 
133 
134  /**
135  *
136  */
138  {
139  }
140 
141 
143  {
145  if (m_constantDepth == true)
146  {
147  SolverUtils::AddSummaryItem(s, "Depth", "constant");
148  }
149  else
150  {
151  SolverUtils::AddSummaryItem(s, "Depth", "variable");
152  }
153 
154 
155  }
156 
158  {
159  ASSERTL0(false, "This function is not implemented for this equation.");
160  }
161 
163  {
164  ASSERTL0(false, "This function is not implemented for this equation.");
165  }
166 
168  {
169  EvaluateFunction("d",m_depth,"WaterDepth");
170  }
171 
172 
174  {
175  EvaluateFunction("f",m_coriolis,"Coriolis");
176  }
177 
179  {
180 
181  int cnt = 0;
182  // loop over Boundary Regions
183  for(int bcRegion = 0; bcRegion < m_fields[0]->GetBndConditions().num_elements(); ++bcRegion)
184  {
185 
186  // Copy the forward trace of the field to the backward trace
187  int e, id2, npts;
188 
189  for (e = 0; e < m_fields[0]->GetBndCondExpansions()[bcRegion]
190  ->GetExpSize(); ++e)
191  {
192  npts = m_fields[0]->GetBndCondExpansions()[bcRegion]->
193  GetExp(e)->GetTotPoints();
194  id2 = m_fields[0]->GetTrace()->GetPhys_Offset(
195  m_fields[0]->GetTraceMap()->
196  GetBndCondCoeffsToGlobalCoeffsMap(cnt+e));
197 
198  Vmath::Vcopy(npts, &Fwd[id2], 1, &Bwd[id2], 1);
199  }
200 
201  cnt +=m_fields[0]->GetBndCondExpansions()[bcRegion]->GetExpSize();
202  }
203  }
204 
205 }
Array< OneD, NekDouble > m_coriolis
Coriolis force.
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:198
bool m_primitive
Indicates if variables are primitive or conservative.
Array< OneD, NekDouble > m_depth
Still water depth.
std::vector< std::pair< std::string, std::string > > SummaryList
Definition: Misc.h:47
enum MultiRegions::ProjectionType m_projectionType
Type of projection; e.g continuous or discontinuous.
STL namespace.
boost::shared_ptr< SessionReader > SessionReaderSharedPtr
Definition: MeshPartition.h:51
SOLVER_UTILS_EXPORT int GetTotPoints()
virtual SOLVER_UTILS_EXPORT void v_GenerateSummary(SummaryList &s)
Print a summary of time stepping parameters.
Array< OneD, Array< OneD, NekDouble > > m_bottomSlope
virtual void v_InitObject()
Init object for UnsteadySystem class.
Base class for unsteady solvers.
void AddSummaryItem(SummaryList &l, const std::string &name, const std::string &value)
Adds a summary item to the summary info list.
Definition: Misc.cpp:50
static std::string npts
Definition: InputFld.cpp:43
int m_spacedim
Spatial dimension (>= expansion dim).
void CopyBoundaryTrace(const Array< OneD, NekDouble > &Fwd, Array< OneD, NekDouble > &Bwd)
virtual void v_GenerateSummary(SolverUtils::SummaryList &s)
Print a summary of time stepping parameters.
void Neg(int n, T *x, const int incx)
Negate x = -x.
Definition: Vmath.cpp:396
double NekDouble
virtual SOLVER_UTILS_EXPORT void v_InitObject()
Init object for UnsteadySystem class.
SOLVER_UTILS_EXPORT void EvaluateFunction(Array< OneD, Array< OneD, NekDouble > > &pArray, std::string pFunctionName, const NekDouble pTime=0.0, const int domain=0)
Evaluates a function as specified in the session file.
EquationSystemFactory & GetEquationSystemFactory()
MultiRegions::Direction const DirCartesianMap[]
Definition: ExpList.h:86
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Array holding all dependent variables.
LibUtilities::SessionReaderSharedPtr m_session
The session reader.
virtual ~ShallowWaterSystem()
Destructor.
bool m_constantDepth
Indicates if constant depth case.
int m_infosteps
Number of time steps between outputting status information.
Array< OneD, Array< OneD, NekDouble > > m_vecLocs
NekDouble m_g
Acceleration of gravity.
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1061
tKey RegisterCreatorFunction(tKey idKey, CreatorFunction classCreator, tDescription pDesc="")
Register a class with the factory.
Definition: NekFactory.hpp:215