Nektar++
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
Public Member Functions | Static Public Member Functions | Static Public Attributes | Protected Member Functions | Protected Attributes | Friends | List of all members
Nektar::SolverUtils::FilterReynoldsStresses Class Reference

Append Reynolds stresses to the average fields. More...

#include <FilterReynoldsStresses.h>

Inheritance diagram for Nektar::SolverUtils::FilterReynoldsStresses:
Inheritance graph
[legend]
Collaboration diagram for Nektar::SolverUtils::FilterReynoldsStresses:
Collaboration graph
[legend]

Public Member Functions

SOLVER_UTILS_EXPORT FilterReynoldsStresses (const LibUtilities::SessionReaderSharedPtr &pSession, const std::map< std::string, std::string > &pParams)
 
SOLVER_UTILS_EXPORT ~FilterReynoldsStresses ()
 
- Public Member Functions inherited from Nektar::SolverUtils::FilterSampler
SOLVER_UTILS_EXPORT FilterSampler (const LibUtilities::SessionReaderSharedPtr &pSession, const ParamMap &pParams)
 
virtual SOLVER_UTILS_EXPORT ~FilterSampler ()
 
- Public Member Functions inherited from Nektar::SolverUtils::Filter
SOLVER_UTILS_EXPORT Filter (const LibUtilities::SessionReaderSharedPtr &pSession)
 
virtual SOLVER_UTILS_EXPORT ~Filter ()
 
SOLVER_UTILS_EXPORT void Initialise (const Array< OneD, const MultiRegions::ExpListSharedPtr > &pFields, const NekDouble &time)
 
SOLVER_UTILS_EXPORT void Update (const Array< OneD, const MultiRegions::ExpListSharedPtr > &pFields, const NekDouble &time)
 
SOLVER_UTILS_EXPORT void Finalise (const Array< OneD, const MultiRegions::ExpListSharedPtr > &pFields, const NekDouble &time)
 
SOLVER_UTILS_EXPORT bool IsTimeDependent ()
 

Static Public Member Functions

static FilterSharedPtr create (const LibUtilities::SessionReaderSharedPtr &pSession, const std::map< std::string, std::string > &pParams)
 Creates an instance of this class. More...
 

Static Public Attributes

static std::string className
 Name of the class. More...
 

Protected Member Functions

virtual void v_Initialise (const Array< OneD, const MultiRegions::ExpListSharedPtr > &pFields, const NekDouble &time)
 
virtual bool v_IsTimeDependent ()
 
virtual void v_ProcessSample (const Array< OneD, const MultiRegions::ExpListSharedPtr > &pFields, const NekDouble &time)
 
virtual void v_PrepareOutput (const Array< OneD, const MultiRegions::ExpListSharedPtr > &pFields, const NekDouble &time)
 
virtual std::string v_GetFileSuffix ()
 
- Protected Member Functions inherited from Nektar::SolverUtils::FilterSampler
virtual SOLVER_UTILS_EXPORT void v_Update (const Array< OneD, const MultiRegions::ExpListSharedPtr > &pFields, const NekDouble &time)
 
virtual SOLVER_UTILS_EXPORT void v_Finalise (const Array< OneD, const MultiRegions::ExpListSharedPtr > &pFields, const NekDouble &time)
 
void OutputField (const Array< OneD, const MultiRegions::ExpListSharedPtr > &pFields, int dump=-1)
 

Protected Attributes

std::vector< Array< OneD,
NekDouble > > 
m_fields
 
std::vector< Array< OneD,
NekDouble > > 
m_delta
 
NekDouble m_alpha
 
bool m_movAvg
 
- Protected Attributes inherited from Nektar::SolverUtils::FilterSampler
NekDouble m_scale
 
unsigned int m_numSamples
 
unsigned int m_outputFrequency
 
unsigned int m_sampleFrequency
 
unsigned int m_index
 
unsigned int m_outputIndex
 
std::string m_outputFile
 
LibUtilities::FieldIOSharedPtr m_fld
 
LibUtilities::FieldMetaDataMap m_fieldMetaData
 
std::vector< Array< OneD,
NekDouble > > 
m_outFields
 
std::vector< std::string > m_variables
 
- Protected Attributes inherited from Nektar::SolverUtils::Filter
LibUtilities::SessionReaderSharedPtr m_session
 

Friends

class MemoryManager< FilterReynoldsStresses >
 

Additional Inherited Members

- Public Types inherited from Nektar::SolverUtils::Filter
typedef std::map< std::string,
std::string > 
ParamMap
 

Detailed Description

Append Reynolds stresses to the average fields.

This class appends the average fields with the Reynolds stresses of the form $ \overline{u' v'} $.

For the default case, this is achieved by calculating $ C_{n} = \Sigma_{i=1}^{n} (u_i - \bar{u}_n)(v_i - \bar{v}_n)$ using the recursive relation:

\[ C_{n} = C_{n-1} + \frac{n}{n-1} (u_n - \bar{u}_n)(v_n - \bar{v}_n) \]

The FilterSampler base class then divides the result by n, leading to the Reynolds stress.

It is also possible to perform the averages using an exponential moving average, in which case either the moving average parameter $ \alpha $ or the time constant $ \tau $ must be prescribed.

Definition at line 45 of file FilterReynoldsStresses.h.

Constructor & Destructor Documentation

Nektar::SolverUtils::FilterReynoldsStresses::FilterReynoldsStresses ( const LibUtilities::SessionReaderSharedPtr pSession,
const std::map< std::string, std::string > &  pParams 
)

Definition at line 67 of file FilterReynoldsStresses.cpp.

References ASSERTL0, Nektar::LibUtilities::Equation::Evaluate(), m_alpha, m_movAvg, Nektar::SolverUtils::FilterSampler::m_sampleFrequency, and Nektar::SolverUtils::Filter::m_session.

70  : FilterSampler(pSession, pParams)
71 {
72  ParamMap::const_iterator it;
73 
74  // Check if should use moving average
75  it = pParams.find("MovingAverage");
76  if (it == pParams.end())
77  {
78  m_movAvg = false;
79  }
80  else
81  {
82  std::string sOption = it->second.c_str();
83  m_movAvg = (boost::iequals(sOption, "true")) ||
84  (boost::iequals(sOption, "yes"));
85  }
86 
87  if (!m_movAvg)
88  {
89  return;
90  }
91 
92  // Load alpha parameter for moving average
93  it = pParams.find("alpha");
94  if (it == pParams.end())
95  {
96  it = pParams.find("tau");
97  if (it == pParams.end())
98  {
99  ASSERTL0(false, "MovingAverage needs either alpha or tau.");
100  }
101  else
102  {
103  // Load time constant
104  LibUtilities::Equation equ(m_session, it->second);
105  NekDouble tau = equ.Evaluate();
106  // Load delta T between samples
107  NekDouble dT;
108  m_session->LoadParameter("TimeStep", dT);
109  dT = dT * m_sampleFrequency;
110  // Calculate alpha
111  m_alpha = dT / (tau + dT);
112  }
113  }
114  else
115  {
116  LibUtilities::Equation equ(m_session, it->second);
117  m_alpha = equ.Evaluate();
118  // Check if tau was also defined
119  it = pParams.find("tau");
120  if (it != pParams.end())
121  {
122  ASSERTL0(false,
123  "Cannot define both alpha and tau in MovingAverage.");
124  }
125  }
126  // Check bounds of m_alpha
127  ASSERTL0(m_alpha > 0 && m_alpha < 1, "Alpha out of bounds.");
128 }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
double NekDouble
LibUtilities::SessionReaderSharedPtr m_session
Definition: Filter.h:84
SOLVER_UTILS_EXPORT FilterSampler(const LibUtilities::SessionReaderSharedPtr &pSession, const ParamMap &pParams)
Nektar::SolverUtils::FilterReynoldsStresses::~FilterReynoldsStresses ( )

Definition at line 130 of file FilterReynoldsStresses.cpp.

131 {
132 }

Member Function Documentation

static FilterSharedPtr Nektar::SolverUtils::FilterReynoldsStresses::create ( const LibUtilities::SessionReaderSharedPtr pSession,
const std::map< std::string, std::string > &  pParams 
)
inlinestatic

Creates an instance of this class.

Definition at line 51 of file FilterReynoldsStresses.h.

References Nektar::MemoryManager< DataType >::AllocateSharedPtr().

54  {
55  FilterSharedPtr p =
57  pParams);
58  return p;
59  }
static boost::shared_ptr< DataType > AllocateSharedPtr()
Allocate a shared pointer from the memory pool.
boost::shared_ptr< Filter > FilterSharedPtr
A shared pointer to a Driver object.
Definition: Filter.h:50
virtual std::string Nektar::SolverUtils::FilterReynoldsStresses::v_GetFileSuffix ( )
inlineprotectedvirtual

Implements Nektar::SolverUtils::FilterSampler.

Definition at line 80 of file FilterReynoldsStresses.h.

81  {
82  return "_stress";
83  }
void Nektar::SolverUtils::FilterReynoldsStresses::v_Initialise ( const Array< OneD, const MultiRegions::ExpListSharedPtr > &  pFields,
const NekDouble time 
)
protectedvirtual

Reimplemented from Nektar::SolverUtils::FilterSampler.

Definition at line 134 of file FilterReynoldsStresses.cpp.

References ASSERTL0, m_delta, m_fields, Nektar::SolverUtils::FilterSampler::m_variables, and Nektar::SolverUtils::FilterSampler::v_Initialise().

137 {
138  int dim = pFields.num_elements() - 1;
139  int nExtraFields = dim == 2 ? 3 : 6;
140  int origFields = pFields.num_elements();
141 
142  // Fill name of variables
143  for (int n = 0; n < origFields; ++n)
144  {
145  m_variables.push_back(pFields[n]->GetSession()->GetVariable(n));
146  }
147  if (dim == 2)
148  {
149  m_variables.push_back("uu");
150  m_variables.push_back("uv");
151  m_variables.push_back("vv");
152  }
153  else if (dim == 3)
154  {
155  m_variables.push_back("uu");
156  m_variables.push_back("uv");
157  m_variables.push_back("uw");
158  m_variables.push_back("vv");
159  m_variables.push_back("vw");
160  m_variables.push_back("ww");
161  }
162  else
163  {
164  ASSERTL0(false, "Unsupported dimension");
165  }
166 
167  // Allocate storage
168  m_fields.resize(origFields + nExtraFields);
169  m_delta.resize(dim);
170 
171  for (int n = 0; n < m_fields.size(); ++n)
172  {
173  m_fields[n] = Array<OneD, NekDouble>(pFields[0]->GetTotPoints(), 0.0);
174  }
175  for (int n = 0; n < m_delta.size(); ++n)
176  {
177  m_delta[n] = Array<OneD, NekDouble>(pFields[0]->GetTotPoints(), 0.0);
178  }
179 
180  // Initialise output arrays
181  FilterSampler::v_Initialise(pFields, time);
182 }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
std::vector< std::string > m_variables
Definition: FilterSampler.h:91
virtual SOLVER_UTILS_EXPORT void v_Initialise(const Array< OneD, const MultiRegions::ExpListSharedPtr > &pFields, const NekDouble &time)
std::vector< Array< OneD, NekDouble > > m_delta
std::vector< Array< OneD, NekDouble > > m_fields
bool Nektar::SolverUtils::FilterReynoldsStresses::v_IsTimeDependent ( )
protectedvirtual

Implements Nektar::SolverUtils::Filter.

Definition at line 297 of file FilterReynoldsStresses.cpp.

298 {
299  return true;
300 }
void Nektar::SolverUtils::FilterReynoldsStresses::v_PrepareOutput ( const Array< OneD, const MultiRegions::ExpListSharedPtr > &  pFields,
const NekDouble time 
)
protectedvirtual

Reimplemented from Nektar::SolverUtils::FilterSampler.

Definition at line 265 of file FilterReynoldsStresses.cpp.

References m_fields, m_movAvg, Nektar::SolverUtils::FilterSampler::m_numSamples, Nektar::SolverUtils::FilterSampler::m_outFields, and Nektar::SolverUtils::FilterSampler::m_scale.

268 {
269  int dim = pFields.num_elements() - 1;
270 
271  if (m_movAvg)
272  {
273  m_scale = 1.0;
274  }
275  else
276  {
277  m_scale = 1.0 / m_numSamples;
278  }
279 
280  // Set wavespace to false, as calculations were performed in physical space
281  bool waveSpace = pFields[0]->GetWaveSpace();
282  pFields[0]->SetWaveSpace(false);
283 
284  // Forward transform and put into m_outFields (except pressure)
285  for (int i = 0; i < m_fields.size(); ++i)
286  {
287  if (i != dim)
288  {
289  pFields[0]->FwdTrans_IterPerExp(m_fields[i], m_outFields[i]);
290  }
291  }
292 
293  // Restore waveSpace
294  pFields[0]->SetWaveSpace(waveSpace);
295 }
std::vector< Array< OneD, NekDouble > > m_outFields
Definition: FilterSampler.h:90
std::vector< Array< OneD, NekDouble > > m_fields
void Nektar::SolverUtils::FilterReynoldsStresses::v_ProcessSample ( const Array< OneD, const MultiRegions::ExpListSharedPtr > &  pFields,
const NekDouble time 
)
protectedvirtual

Implements Nektar::SolverUtils::FilterSampler.

Definition at line 184 of file FilterReynoldsStresses.cpp.

References m_alpha, m_delta, m_fields, m_movAvg, Nektar::SolverUtils::FilterSampler::m_numSamples, Nektar::SolverUtils::FilterSampler::m_outFields, Vmath::Svtsvtp(), Vmath::Svtvm(), and Vmath::Vmul().

187 {
188  int i, j, n;
189  int nq = pFields[0]->GetTotPoints();
190  int dim = pFields.num_elements() - 1;
191  bool waveSpace = pFields[0]->GetWaveSpace();
192  NekDouble nSamples = (NekDouble)m_numSamples;
193 
194  // For moving average, take first sample as initial vector
195  NekDouble alpha = m_alpha;
196  if (m_numSamples == 1)
197  {
198  alpha = 1.0;
199  }
200 
201  // Define auxiliary constants for averages
202  NekDouble facOld, facAvg, facStress, facDelta;
203  if (m_movAvg)
204  {
205  facOld = 1.0 - alpha;
206  facAvg = alpha;
207  facStress = alpha;
208  facDelta = 1.0;
209  }
210  else
211  {
212  facOld = 1.0;
213  facAvg = 1.0;
214  facStress = nSamples / (nSamples - 1);
215  facDelta = 1.0 / nSamples;
216  }
217 
218  Array<OneD, NekDouble> vel(nq);
219  Array<OneD, NekDouble> tmp(nq);
220 
221  // Update original velocities in phys space and calculate (\bar{u} - u_n)
222  for (n = 0; n < dim; ++n)
223  {
224  if (waveSpace)
225  {
226  pFields[n]->HomogeneousBwdTrans(pFields[n]->GetPhys(), vel);
227  }
228  else
229  {
230  vel = pFields[n]->GetPhys();
231  }
233  nq, facAvg, vel, 1, facOld, m_fields[n], 1, m_fields[n], 1);
234  Vmath::Svtvm(nq, facDelta, m_fields[n], 1, vel, 1, m_delta[n], 1);
235  }
236  // Update pressure (directly to outFields)
237  Vmath::Svtsvtp(m_outFields[dim].num_elements(),
238  facAvg,
239  pFields[dim]->GetCoeffs(),
240  1,
241  facOld,
242  m_outFields[dim],
243  1,
244  m_outFields[dim],
245  1);
246 
247  // Ignore Reynolds stress for first sample (its contribution is zero)
248  if (m_numSamples == 1)
249  {
250  return;
251  }
252 
253  // Calculate C_{n} = facOld * C_{n-1} + facStress * deltaI * deltaJ
254  for (i = 0, n = dim + 1; i < dim; ++i)
255  {
256  for (j = i; j < dim; ++j, ++n)
257  {
258  Vmath::Vmul(nq, m_delta[i], 1, m_delta[j], 1, tmp, 1);
260  nq, facStress, tmp, 1, facOld, m_fields[n], 1, m_fields[n], 1);
261  }
262  }
263 }
std::vector< Array< OneD, NekDouble > > m_outFields
Definition: FilterSampler.h:90
void Svtvm(int n, const T alpha, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
svtvp (scalar times vector plus vector): z = alpha*x - y
Definition: Vmath.cpp:504
std::vector< Array< OneD, NekDouble > > m_delta
double NekDouble
std::vector< Array< OneD, NekDouble > > m_fields
void Svtsvtp(int n, const T alpha, const T *x, int incx, const T beta, const T *y, int incy, T *z, int incz)
vvtvvtp (scalar times vector plus scalar times vector):
Definition: Vmath.cpp:577
void Vmul(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Multiply vector z = x*y.
Definition: Vmath.cpp:169

Friends And Related Function Documentation

friend class MemoryManager< FilterReynoldsStresses >
friend

Definition at line 48 of file FilterReynoldsStresses.h.

Member Data Documentation

std::string Nektar::SolverUtils::FilterReynoldsStresses::className
static
Initial value:

Name of the class.

Definition at line 62 of file FilterReynoldsStresses.h.

NekDouble Nektar::SolverUtils::FilterReynoldsStresses::m_alpha
protected

Definition at line 87 of file FilterReynoldsStresses.h.

Referenced by FilterReynoldsStresses(), and v_ProcessSample().

std::vector<Array<OneD, NekDouble> > Nektar::SolverUtils::FilterReynoldsStresses::m_delta
protected

Definition at line 86 of file FilterReynoldsStresses.h.

Referenced by v_Initialise(), and v_ProcessSample().

std::vector<Array<OneD, NekDouble> > Nektar::SolverUtils::FilterReynoldsStresses::m_fields
protected

Definition at line 85 of file FilterReynoldsStresses.h.

Referenced by v_Initialise(), v_PrepareOutput(), and v_ProcessSample().

bool Nektar::SolverUtils::FilterReynoldsStresses::m_movAvg
protected