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

This processing module calculates the wall shear stress and adds it as an extra-field to the output file, and writes it to a surface output file. More...

#include <ProcessWSS.h>

Inheritance diagram for Nektar::FieldUtils::ProcessWSS:
Inheritance graph
[legend]
Collaboration diagram for Nektar::FieldUtils::ProcessWSS:
Collaboration graph
[legend]

Public Member Functions

 ProcessWSS (FieldSharedPtr f)
 
virtual ~ProcessWSS ()
 
virtual void Process (po::variables_map &vm)
 Write mesh to output file. More...
 
virtual std::string GetModuleName ()
 
- Public Member Functions inherited from Nektar::FieldUtils::ProcessModule
 ProcessModule ()
 
 ProcessModule (FieldSharedPtr p_f)
 
- Public Member Functions inherited from Nektar::FieldUtils::Module
FIELD_UTILS_EXPORT Module (FieldSharedPtr p_f)
 
FIELD_UTILS_EXPORT void RegisterConfig (string key, string value)
 Register a configuration option with a module. More...
 
FIELD_UTILS_EXPORT void PrintConfig ()
 Print out all configuration options for a module. More...
 
FIELD_UTILS_EXPORT void SetDefaults ()
 Sets default configuration options for those which have not been set. More...
 
FIELD_UTILS_EXPORT bool GetRequireEquiSpaced (void)
 
FIELD_UTILS_EXPORT void SetRequireEquiSpaced (bool pVal)
 
FIELD_UTILS_EXPORT void EvaluateTriFieldAtEquiSpacedPts (LocalRegions::ExpansionSharedPtr &exp, const Array< OneD, const NekDouble > &infield, Array< OneD, NekDouble > &outfield)
 

Static Public Member Functions

static boost::shared_ptr< Modulecreate (FieldSharedPtr f)
 Creates an instance of this class. More...
 

Static Public Attributes

static ModuleKey className
 

Additional Inherited Members

- Protected Member Functions inherited from Nektar::FieldUtils::Module
 Module ()
 
- Protected Attributes inherited from Nektar::FieldUtils::Module
FieldSharedPtr m_f
 Field object. More...
 
map< string, ConfigOptionm_config
 List of configuration values. More...
 
bool m_requireEquiSpaced
 

Detailed Description

This processing module calculates the wall shear stress and adds it as an extra-field to the output file, and writes it to a surface output file.

Definition at line 50 of file ProcessWSS.h.

Constructor & Destructor Documentation

Nektar::FieldUtils::ProcessWSS::ProcessWSS ( FieldSharedPtr  f)

Definition at line 56 of file ProcessWSS.cpp.

References Nektar::FieldUtils::Module::m_config, and Nektar::FieldUtils::Module::m_f.

56  : ProcessModule(f)
57 {
58  m_config["bnd"] = ConfigOption(false, "All", "Boundary to be extracted");
59  m_config["addnormals"] =
60  ConfigOption(true, "NotSet", "Add normals to output");
61  f->m_writeBndFld = true;
62  f->m_declareExpansionAsContField = true;
63  f->m_requireBoundaryExpansion = true;
64  m_f->m_fldToBnd = false;
65 }
map< string, ConfigOption > m_config
List of configuration values.
FieldSharedPtr m_f
Field object.
Nektar::FieldUtils::ProcessWSS::~ProcessWSS ( )
virtual

Definition at line 67 of file ProcessWSS.cpp.

68 {
69 }

Member Function Documentation

static boost::shared_ptr<Module> Nektar::FieldUtils::ProcessWSS::create ( FieldSharedPtr  f)
inlinestatic

Creates an instance of this class.

Definition at line 54 of file ProcessWSS.h.

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

55  {
57  }
static boost::shared_ptr< DataType > AllocateSharedPtr()
Allocate a shared pointer from the memory pool.
virtual std::string Nektar::FieldUtils::ProcessWSS::GetModuleName ( )
inlinevirtual

Implements Nektar::FieldUtils::Module.

Definition at line 66 of file ProcessWSS.h.

67  {
68  return "ProcessWSS";
69  }
void Nektar::FieldUtils::ProcessWSS::Process ( po::variables_map &  vm)
virtual

Write mesh to output file.

Implements Nektar::FieldUtils::Module.

Definition at line 71 of file ProcessWSS.cpp.

References ASSERTL0, Nektar::ParseUtils::GenerateOrderedVector(), Nektar::SpatialDomains::BoundaryConditions::GetBoundaryRegions(), Nektar::FieldUtils::Module::m_config, Nektar::FieldUtils::Module::m_f, Vmath::Neg(), Vmath::Smul(), Vmath::Vadd(), Vmath::Vsqrt(), Vmath::Vvtvp(), and Vmath::Zero().

72 {
73  if (m_f->m_verbose)
74  {
75  if (m_f->m_comm->TreatAsRankZero())
76  {
77  cout << "ProcessWSS: Calculating wall shear stress..." << endl;
78  }
79  }
80 
81  m_f->m_addNormals = m_config["addnormals"].m_beenSet;
82 
83  // Set up Field options to output boundary fld
84  string bvalues = m_config["bnd"].as<string>();
85 
86  if (bvalues.compare("All") == 0)
87  {
88  Array<OneD, const MultiRegions::ExpListSharedPtr> BndExp =
89  m_f->m_exp[0]->GetBndCondExpansions();
90 
91  for (int i = 0; i < BndExp.num_elements(); ++i)
92  {
93  m_f->m_bndRegionsToWrite.push_back(i);
94  }
95  }
96  else
97  {
99  m_f->m_bndRegionsToWrite),
100  "Failed to interpret range string");
101  }
102 
103  NekDouble kinvis = m_f->m_session->GetParameter("Kinvis");
104 
105  int i, j;
106  int spacedim = m_f->m_graph->GetSpaceDimension();
107  if ((m_f->m_fielddef[0]->m_numHomogeneousDir) == 1 ||
108  (m_f->m_fielddef[0]->m_numHomogeneousDir) == 2)
109  {
110  spacedim += m_f->m_fielddef[0]->m_numHomogeneousDir;
111  }
112 
113  int nfields = m_f->m_fielddef[0]->m_fields.size();
114  ASSERTL0(m_f->m_fielddef[0]->m_fields[0] == "u",
115  "Implicit assumption that input is in incompressible format of "
116  "(u,v,p) or (u,v,w,p)");
117 
118  if (spacedim == 1)
119  {
120  ASSERTL0(false, "Error: wss for a 1D problem cannot "
121  "be computed");
122  }
123 
124  int newfields = spacedim + 1;
125  int nshear = spacedim + 1;
126  int nstress = spacedim * spacedim;
127  int ngrad = spacedim * spacedim;
128 
129  Array<OneD, Array<OneD, NekDouble> > velocity(nfields), grad(ngrad),
130  fgrad(ngrad);
131  Array<OneD, Array<OneD, NekDouble> > stress(nstress), fstress(nstress);
132  Array<OneD, Array<OneD, NekDouble> > fshear(nshear);
133 
134  Array<OneD, MultiRegions::ExpListSharedPtr> BndExp(newfields);
135  Array<OneD, MultiRegions::ExpListSharedPtr> BndElmtExp(spacedim);
136 
137  // Extract original fields to boundary (for output)
138  for (int i = 0; i < m_f->m_exp.size(); ++i)
139  {
140  m_f->m_exp[i]->FillBndCondFromField();
141  }
142 
143  m_f->m_exp.resize(nfields + newfields);
144  string var = "u";
145  for (i = 0; i < newfields; ++i)
146  {
147  m_f->m_exp[nfields + i] =
148  m_f->AppendExpList(m_f->m_fielddef[0]->m_numHomogeneousDir, var);
149  }
150 
151  if (spacedim == 2)
152  {
153  m_f->m_fielddef[0]->m_fields.push_back("Shear_x");
154  m_f->m_fielddef[0]->m_fields.push_back("Shear_y");
155  m_f->m_fielddef[0]->m_fields.push_back("Shear_mag");
156  }
157  else
158  {
159  m_f->m_fielddef[0]->m_fields.push_back("Shear_x");
160  m_f->m_fielddef[0]->m_fields.push_back("Shear_y");
161  m_f->m_fielddef[0]->m_fields.push_back("Shear_z");
162  m_f->m_fielddef[0]->m_fields.push_back("Shear_mag");
163  }
164 
165  // Create map of boundary ids for partitioned domains
166  SpatialDomains::BoundaryConditions bcs(m_f->m_session,
167  m_f->m_exp[0]->GetGraph());
169  bcs.GetBoundaryRegions();
170  SpatialDomains::BoundaryRegionCollection::const_iterator breg_it;
171  map<int, int> BndRegionMap;
172  int cnt = 0;
173  for (breg_it = bregions.begin(); breg_it != bregions.end();
174  ++breg_it, ++cnt)
175  {
176  BndRegionMap[breg_it->first] = cnt;
177  }
178  // Loop over boundaries to Write
179  for (int b = 0; b < m_f->m_bndRegionsToWrite.size(); ++b)
180  {
181  if (BndRegionMap.count(m_f->m_bndRegionsToWrite[b]) == 1)
182  {
183  int bnd = BndRegionMap[m_f->m_bndRegionsToWrite[b]];
184  // Get expansion list for boundary and for elements containing this
185  // bnd
186  for (i = 0; i < newfields; i++)
187  {
188  BndExp[i] =
189  m_f->m_exp[nfields + i]->UpdateBndCondExpansion(bnd);
190  }
191  for (i = 0; i < spacedim; i++)
192  {
193  m_f->m_exp[i]->GetBndElmtExpansion(bnd, BndElmtExp[i]);
194  }
195 
196  // Get number of points in expansions
197  int nqb = BndExp[0]->GetTotPoints();
198  int nqe = BndElmtExp[0]->GetTotPoints();
199 
200  // Initialise local arrays for the velocity gradients, and stress
201  // components
202  // size of total number of quadrature points for elements in this
203  // bnd
204  for (i = 0; i < ngrad; ++i)
205  {
206  grad[i] = Array<OneD, NekDouble>(nqe);
207  }
208 
209  for (i = 0; i < nstress; ++i)
210  {
211  stress[i] = Array<OneD, NekDouble>(nqe);
212  }
213 
214  // initialise arrays in the boundary
215  for (i = 0; i < nstress; ++i)
216  {
217  fstress[i] = Array<OneD, NekDouble>(nqb);
218  }
219 
220  for (i = 0; i < ngrad; ++i)
221  {
222  fgrad[i] = Array<OneD, NekDouble>(nqb);
223  }
224 
225  for (i = 0; i < nshear; ++i)
226  {
227  fshear[i] = Array<OneD, NekDouble>(nqb, 0.0);
228  }
229 
230  // Extract Velocities
231  for (i = 0; i < spacedim; ++i)
232  {
233  velocity[i] = BndElmtExp[i]->GetPhys();
234  }
235 
236  // Compute gradients (velocity correction scheme method)
237  for (i = 0; i < spacedim; ++i)
238  {
239  if (spacedim == 2)
240  {
241  BndElmtExp[i]->PhysDeriv(velocity[i],
242  grad[i * spacedim + 0],
243  grad[i * spacedim + 1]);
244  }
245  else
246  {
247  BndElmtExp[i]->PhysDeriv(
248  velocity[i], grad[i * spacedim + 0],
249  grad[i * spacedim + 1], grad[i * spacedim + 2]);
250  }
251  }
252 
253  // Compute stress component terms tau_ij = mu*(u_i,j + u_j,i)
254  for (i = 0; i < spacedim; ++i)
255  {
256  for (j = 0; j < spacedim; ++j)
257  {
258  Vmath::Vadd(nqe, grad[i * spacedim + j], 1,
259  grad[j * spacedim + i], 1,
260  stress[i * spacedim + j], 1);
261 
262  Vmath::Smul(nqe, kinvis, stress[i * spacedim + j], 1,
263  stress[i * spacedim + j], 1);
264  }
265  }
266 
267  // Get boundary stress values.
268  for (j = 0; j < nstress; ++j)
269  {
270  m_f->m_exp[0]->ExtractElmtToBndPhys(bnd, stress[j], fstress[j]);
271  }
272 
273  // Get normals
274  Array<OneD, Array<OneD, NekDouble> > normals;
275  m_f->m_exp[0]->GetBoundaryNormals(bnd, normals);
276  // Reverse normals, to get correct orientation for the body
277  for (i = 0; i < spacedim; ++i)
278  {
279  Vmath::Neg(nqb, normals[i], 1);
280  }
281 
282  // calculate wss, and update coeffs in the boundary expansion
283  // S = tau_ij * n_j
284  for (i = 0; i < spacedim; ++i)
285  {
286  for (j = 0; j < spacedim; ++j)
287  {
288  Vmath::Vvtvp(nqb, normals[j], 1, fstress[i * spacedim + j],
289  1, fshear[i], 1, fshear[i], 1);
290  }
291  }
292 
293  // T = S - (S.n)n
294  for (i = 0; i < spacedim; ++i)
295  {
296  Vmath::Vvtvp(nqb, normals[i], 1, fshear[i], 1,
297  fshear[nshear - 1], 1, fshear[nshear - 1], 1);
298  }
299  Vmath::Smul(nqb, -1.0, fshear[nshear - 1], 1, fshear[nshear - 1],
300  1);
301 
302  for (i = 0; i < spacedim; i++)
303  {
304  Vmath::Vvtvp(nqb, normals[i], 1, fshear[nshear - 1], 1,
305  fshear[i], 1, fshear[i], 1);
306  BndExp[i]->FwdTrans(fshear[i], BndExp[i]->UpdateCoeffs());
307  }
308 
309  // Tw
310  Vmath::Zero(nqb, fshear[nshear - 1], 1);
311  for (i = 0; i < spacedim; ++i)
312  {
313  Vmath::Vvtvp(nqb, fshear[i], 1, fshear[i], 1,
314  fshear[nshear - 1], 1, fshear[nshear - 1], 1);
315  }
316  Vmath::Vsqrt(nqb, fshear[nshear - 1], 1, fshear[nshear - 1], 1);
317  BndExp[nshear - 1]->FwdTrans(fshear[nshear - 1],
318  BndExp[nshear - 1]->UpdateCoeffs());
319  }
320  }
321 }
map< string, ConfigOption > m_config
List of configuration values.
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:198
static bool GenerateOrderedVector(const char *const str, std::vector< unsigned int > &vec)
Definition: ParseUtils.hpp:97
void Vsqrt(int n, const T *x, const int incx, T *y, const int incy)
sqrt y = sqrt(x)
Definition: Vmath.cpp:408
void Vvtvp(int n, const T *w, const int incw, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
vvtvp (vector times vector plus vector): z = w*x + y
Definition: Vmath.cpp:442
void Smul(int n, const T alpha, const T *x, const int incx, T *y, const int incy)
Scalar multiply y = alpha*y.
Definition: Vmath.cpp:213
std::map< int, BoundaryRegionShPtr > BoundaryRegionCollection
Definition: Conditions.h:217
void Neg(int n, T *x, const int incx)
Negate x = -x.
Definition: Vmath.cpp:396
double NekDouble
void Zero(int n, T *x, const int incx)
Zero vector.
Definition: Vmath.cpp:373
void Vadd(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Add vector z = x+y.
Definition: Vmath.cpp:299
FieldSharedPtr m_f
Field object.

Member Data Documentation

ModuleKey Nektar::FieldUtils::ProcessWSS::className
static
Initial value:
"Computes wall shear stress field.")

Definition at line 58 of file ProcessWSS.h.