Nektar++
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator 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 namespace Nektar
41 {
42  /**
43  * @class ShallowWaterSystem
44  *
45  * Provides the underlying timestepping framework for shallow water flow solvers
46  * including the general timestepping routines. This class is not intended
47  * to be directly instantiated, but rather is a base class on which to
48  * define shallow water solvers, e.g. SWE, Boussinesq, linear and nonlinear versions.
49  *
50  * For details on implementing unsteady solvers see
51  * \ref sectionADRSolverModuleImplementation
52  */
53 
54  /**
55  * Processes SolverInfo parameters from the session file and sets up
56  * timestepping-specific code.
57  * @param pSession Session object to read parameters from.
58  */
59 
62  "ShallowWaterSystem",
64  "Auxiliary functions for the shallow water system.");
65 
66 
69  : UnsteadySystem(pSession)
70  {
71  }
72 
74  {
76 
77  // if discontinuous Galerkin determine numerical flux to use
79  {
80  ASSERTL0(m_session->DefinesSolverInfo("UPWINDTYPE"),
81  "No UPWINDTYPE defined in session.");
82  }
83 
84  // Set up locations of velocity vector.
85  m_vecLocs = Array<OneD, Array<OneD, NekDouble> >(1);
86  m_vecLocs[0] = Array<OneD, NekDouble>(m_spacedim);
87  for (int i = 0; i < m_spacedim; ++i)
88  {
89  m_vecLocs[0][i] = 1 + i;
90  }
91 
92  // Load generic input parameters
93  m_session->LoadParameter("IO_InfoSteps", m_infosteps, 0);
94 
95  // Load acceleration of gravity
96  m_session->LoadParameter("Gravity", m_g, 9.81);
97 
98  // input/output in primitive variables
99  m_primitive = true;
100 
102 
103  m_constantDepth = true;
104  NekDouble depth = m_depth[0];
105  for (int i = 0; i < GetTotPoints(); ++i)
106  {
107  if (m_depth[i] != depth)
108  {
109  m_constantDepth = false;
110  break;
111  }
112  }
113 
114  // Compute the bottom slopes
115  int nq = GetTotPoints();
116  if (m_constantDepth != true)
117  {
118  m_bottomSlope = Array<OneD, Array<OneD, NekDouble> >(m_spacedim);
119  for (int i = 0; i < m_spacedim; ++i)
120  {
121  m_bottomSlope[i] = Array<OneD, NekDouble>(nq);
123  Vmath::Neg(nq,m_bottomSlope[i],1);
124  }
125  }
126 
128 
129  }
130 
131 
132  /**
133  *
134  */
136  {
137  }
138 
139 
141  {
143  if (m_constantDepth == true)
144  {
145  SolverUtils::AddSummaryItem(s, "Depth", "constant");
146  }
147  else
148  {
149  SolverUtils::AddSummaryItem(s, "Depth", "variable");
150  }
151 
152 
153  }
154 
156  {
157  ASSERTL0(false, "This function is not implemented for this equation.");
158  }
159 
161  {
162  ASSERTL0(false, "This function is not implemented for this equation.");
163  }
164 
166  {
167  EvaluateFunction("d",m_depth,"WaterDepth");
168  }
169 
170 
172  {
173  EvaluateFunction("f",m_coriolis,"Coriolis");
174  }
175 
176  void ShallowWaterSystem::CopyBoundaryTrace(const Array<OneD, NekDouble>&Fwd, Array<OneD, NekDouble>&Bwd)
177  {
178 
179  int cnt = 0;
180  // loop over Boundary Regions
181  for(int bcRegion = 0; bcRegion < m_fields[0]->GetBndConditions().num_elements(); ++bcRegion)
182  {
183 
184  // Copy the forward trace of the field to the backward trace
185  int e, id2, npts;
186 
187  for (e = 0; e < m_fields[0]->GetBndCondExpansions()[bcRegion]
188  ->GetExpSize(); ++e)
189  {
190  npts = m_fields[0]->GetBndCondExpansions()[bcRegion]->
191  GetExp(e)->GetTotPoints();
192  id2 = m_fields[0]->GetTrace()->GetPhys_Offset(
193  m_fields[0]->GetTraceMap()->
194  GetBndCondCoeffsToGlobalCoeffsMap(cnt+e));
195 
196  Vmath::Vcopy(npts, &Fwd[id2], 1, &Bwd[id2], 1);
197  }
198 
199  cnt +=m_fields[0]->GetBndCondExpansions()[bcRegion]->GetExpSize();
200  }
201  }
202 
203 }