Nektar++
Public Member Functions | Static Public Member Functions | Static Public Attributes | Protected Member Functions | List of all members
Nektar::FieldUtils::ProcessMultiShear Class Reference

This processing module calculates the shear stress metrics and writes it to a surface output file. More...

#include <ProcessMultiShear.h>

Inheritance diagram for Nektar::FieldUtils::ProcessMultiShear:
[legend]

Public Member Functions

 ProcessMultiShear (FieldSharedPtr f)
 
virtual ~ProcessMultiShear ()
 
- 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)
 
virtual ~Module ()=default
 
void Process (po::variables_map &vm)
 
std::string GetModuleName ()
 
std::string GetModuleDescription ()
 
const ConfigOptionGetConfigOption (const std::string &key) const
 
ModulePriority GetModulePriority ()
 
FIELD_UTILS_EXPORT void RegisterConfig (std::string key, std::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 void AddFile (std::string fileType, std::string fileName)
 
FIELD_UTILS_EXPORT void EvaluateTriFieldAtEquiSpacedPts (LocalRegions::ExpansionSharedPtr &exp, const Array< OneD, const NekDouble > &infield, Array< OneD, NekDouble > &outfield)
 

Static Public Member Functions

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

Static Public Attributes

static ModuleKey className
 

Protected Member Functions

virtual void v_Process (po::variables_map &vm) override
 Write mesh to output file. More...
 
virtual std::string v_GetModuleName () override
 
virtual std::string v_GetModuleDescription () override
 
virtual ModulePriority v_GetModulePriority () override
 
- Protected Member Functions inherited from Nektar::FieldUtils::Module
 Module ()
 

Additional Inherited Members

- Public Attributes inherited from Nektar::FieldUtils::Module
FieldSharedPtr m_f
 Field object. More...
 
- Protected Attributes inherited from Nektar::FieldUtils::Module
std::map< std::string, ConfigOptionm_config
 List of configuration values. More...
 
std::set< std::string > m_allowedFiles
 List of allowed file formats. More...
 

Detailed Description

This processing module calculates the shear stress metrics and writes it to a surface output file.

Definition at line 49 of file ProcessMultiShear.h.

Constructor & Destructor Documentation

◆ ProcessMultiShear()

Nektar::FieldUtils::ProcessMultiShear::ProcessMultiShear ( FieldSharedPtr  f)

Definition at line 57 of file ProcessMultiShear.cpp.

57  : ProcessModule(f)
58 {
59  m_config["N"] = ConfigOption(false, "1", "Number of chk or fld files");
60  m_config["fromfld"] = ConfigOption(
61  false, "NotSet",
62  "First fld file. First underscore flags position of id in name.");
63 }
std::map< std::string, ConfigOption > m_config
List of configuration values.
Definition: Module.h:263

References Nektar::FieldUtils::Module::m_config.

◆ ~ProcessMultiShear()

Nektar::FieldUtils::ProcessMultiShear::~ProcessMultiShear ( )
virtual

Definition at line 65 of file ProcessMultiShear.cpp.

66 {
67 }

Member Function Documentation

◆ create()

static std::shared_ptr<Module> Nektar::FieldUtils::ProcessMultiShear::create ( FieldSharedPtr  f)
inlinestatic

Creates an instance of this class.

Definition at line 53 of file ProcessMultiShear.h.

54  {
56  }
static std::shared_ptr< DataType > AllocateSharedPtr(const Args &...args)
Allocate a shared pointer from the memory pool.

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

◆ v_GetModuleDescription()

virtual std::string Nektar::FieldUtils::ProcessMultiShear::v_GetModuleDescription ( )
inlineoverrideprotectedvirtual

Reimplemented from Nektar::FieldUtils::Module.

Definition at line 71 of file ProcessMultiShear.h.

72  {
73  return "Calculating shear stress metrics";
74  }

◆ v_GetModuleName()

virtual std::string Nektar::FieldUtils::ProcessMultiShear::v_GetModuleName ( )
inlineoverrideprotectedvirtual

Reimplemented from Nektar::FieldUtils::Module.

Definition at line 66 of file ProcessMultiShear.h.

67  {
68  return "ProcessMultiShear";
69  }

◆ v_GetModulePriority()

virtual ModulePriority Nektar::FieldUtils::ProcessMultiShear::v_GetModulePriority ( )
inlineoverrideprotectedvirtual

Reimplemented from Nektar::FieldUtils::Module.

Definition at line 76 of file ProcessMultiShear.h.

77  {
78  return eModifyExp;
79  }

References Nektar::FieldUtils::eModifyExp.

◆ v_Process()

void Nektar::FieldUtils::ProcessMultiShear::v_Process ( po::variables_map &  vm)
overrideprotectedvirtual

Write mesh to output file.

Reimplemented from Nektar::FieldUtils::Module.

Definition at line 69 of file ProcessMultiShear.cpp.

70 {
71  m_f->SetUpExp(vm);
72 
73  // Skip in case of empty partition
74  if (m_f->m_exp[0]->GetNumElmts() == 0)
75  {
76  return;
77  }
78 
79  ASSERTL0(m_config["fromfld"].as<string>().compare("NotSet") != 0,
80  "Need to specify fromfld=file.fld ");
81 
82  int nstart, i, j, nfields = 0;
83  bool wssg = false;
84  NekDouble nfld = m_config["N"].as<NekDouble>();
85  string fromfld, basename, endname, nstartStr;
86  stringstream filename;
87  vector<string> infiles(nfld);
88  vector<std::shared_ptr<Field>> fromField(nfld);
89 
90  // Set up list of input fld files.
91  fromfld = m_config["fromfld"].as<string>();
92  basename = fromfld.substr(0, fromfld.find_first_of("_") + 1);
93  filename << fromfld.substr(fromfld.find_first_of("_") + 1, fromfld.size());
94  filename >> nstart;
95  filename.str("");
96  filename << nstart;
97  filename >> nstartStr;
98  filename.str("");
99  endname = fromfld.substr(fromfld.find(nstartStr) + nstartStr.size(),
100  fromfld.size());
101 
102  for (i = 0; i < nfld; ++i)
103  {
104  stringstream filename;
105  filename << basename << i + nstart << endname;
106  filename >> infiles[i];
107  cout << infiles[i] << endl;
108  }
109 
110  for (i = 0; i < nfld; ++i)
111  {
112  fromField[i] = std::shared_ptr<Field>(new Field());
113  fromField[i]->m_session = m_f->m_session;
114  fromField[i]->m_graph = m_f->m_graph;
115  }
116 
117  // Import all fld files.
118  for (i = 0; i < nfld; ++i)
119  {
120  if (m_f->m_exp.size())
121  {
122  // Set up ElementGIDs in case of parallel processing
123  Array<OneD, int> ElementGIDs(m_f->m_exp[0]->GetExpSize());
124  for (j = 0; j < m_f->m_exp[0]->GetExpSize(); ++j)
125  {
126  ElementGIDs[j] =
127  m_f->m_exp[0]->GetExp(j)->GetGeom()->GetGlobalID();
128  }
129  m_f->FieldIOForFile(infiles[i])
130  ->Import(infiles[i], fromField[i]->m_fielddef,
131  fromField[i]->m_data,
133  }
134  else
135  {
136  m_f->FieldIOForFile(infiles[i])
137  ->Import(infiles[i], fromField[i]->m_fielddef,
138  fromField[i]->m_data,
140  }
141 
142  nfields = fromField[i]->m_fielddef[0]->m_fields.size();
143  int NumHomogeneousDir =
144  fromField[i]->m_fielddef[0]->m_numHomogeneousDir;
145 
146  if (nfields == 5 || nfields == 7)
147  {
148  wssg = true;
149  }
150 
151  // Set up Expansion information to use mode order from field
152  fromField[i]->m_graph->SetExpansionInfo(fromField[i]->m_fielddef);
153 
154  // Set up expansions, and extract data.
155  fromField[i]->m_exp.resize(nfields);
156  fromField[i]->m_exp[0] =
157  fromField[i]->SetUpFirstExpList(NumHomogeneousDir, true);
158 
159  for (j = 1; j < nfields; ++j)
160  {
161  fromField[i]->m_exp[j] = m_f->AppendExpList(NumHomogeneousDir);
162  }
163 
164  for (j = 0; j < nfields; ++j)
165  {
166  for (int k = 0; k < fromField[i]->m_data.size(); ++k)
167  {
168  fromField[i]->m_exp[j]->ExtractDataToCoeffs(
169  fromField[i]->m_fielddef[k], fromField[i]->m_data[k],
170  fromField[i]->m_fielddef[k]->m_fields[j],
171  fromField[i]->m_exp[j]->UpdateCoeffs());
172  }
173  fromField[i]->m_exp[j]->BwdTrans(
174  fromField[i]->m_exp[j]->GetCoeffs(),
175  fromField[i]->m_exp[j]->UpdatePhys());
176  }
177  }
178 
179  int spacedim = m_f->m_graph->GetSpaceDimension();
180  if ((fromField[0]->m_fielddef[0]->m_numHomogeneousDir) == 1 ||
181  (fromField[0]->m_fielddef[0]->m_numHomogeneousDir) == 2)
182  {
183  spacedim = 3;
184  }
185 
186  int nout = 5; // TAWSS, OSI, transWSS, AFI, CFI (optional WSSG)
187  if (wssg)
188  {
189  nout = 6;
190  }
191 
192  int npoints = fromField[0]->m_exp[0]->GetNpoints();
193  Array<OneD, Array<OneD, NekDouble>> normTemporalMeanVec(spacedim),
194  normCrossDir(spacedim), outfield(nout), dtemp(spacedim);
195  Array<OneD, NekDouble> TemporalMeanMag(npoints, 0.0),
196  DotProduct(npoints, 0.0), temp(npoints, 0.0);
197 
198  for (i = 0; i < spacedim; ++i)
199  {
200  normTemporalMeanVec[i] = Array<OneD, NekDouble>(npoints);
201  normCrossDir[i] = Array<OneD, NekDouble>(npoints);
202  dtemp[i] = Array<OneD, NekDouble>(npoints);
203  Vmath::Zero(npoints, normTemporalMeanVec[i], 1);
204  Vmath::Zero(npoints, normCrossDir[i], 1);
205  }
206 
207  for (i = 0; i < nout; ++i)
208  {
209  outfield[i] = Array<OneD, NekDouble>(npoints, 0.0);
210  }
211 
212  // -----------------------------------------------------
213  // Compute temporal average wall shear stress vector,
214  // it's spatial average, and normalise it.
215  for (i = 0; i < nfld; ++i)
216  {
217  for (j = 0; j < spacedim; ++j)
218  {
219  Vmath::Vadd(npoints, fromField[i]->m_exp[j]->GetPhys(), 1,
220  normTemporalMeanVec[j], 1, normTemporalMeanVec[j], 1);
221  }
222  }
223 
224  for (i = 0; i < spacedim; ++i)
225  {
226  Vmath::Smul(npoints, 1.0 / nfld, normTemporalMeanVec[i], 1,
227  normTemporalMeanVec[i], 1);
228  Vmath::Vvtvp(npoints, normTemporalMeanVec[i], 1, normTemporalMeanVec[i],
229  1, TemporalMeanMag, 1, TemporalMeanMag, 1);
230  }
231 
232  Vmath::Vsqrt(npoints, TemporalMeanMag, 1, TemporalMeanMag, 1);
233 
234  for (i = 0; i < spacedim; ++i)
235  {
236  Vmath::Vdiv(npoints, normTemporalMeanVec[i], 1, TemporalMeanMag, 1,
237  normTemporalMeanVec[i], 1);
238  }
239  //---------------------------------------------------
240 
241  if (wssg) // cross product with normals to obtain direction parallel to
242  // temporal mean vector.
243  {
244  Vmath::Vmul(npoints, fromField[0]->m_exp[nfields - 1]->GetPhys(), 1,
245  normTemporalMeanVec[1], 1, normCrossDir[0], 1);
246  Vmath::Vvtvm(npoints, fromField[0]->m_exp[nfields - 2]->GetPhys(), 1,
247  normTemporalMeanVec[2], 1, normCrossDir[0], 1,
248  normCrossDir[0], 1);
249  Vmath::Vmul(npoints, fromField[0]->m_exp[nfields - 3]->GetPhys(), 1,
250  normTemporalMeanVec[2], 1, normCrossDir[1], 1);
251  Vmath::Vvtvm(npoints, fromField[0]->m_exp[nfields - 1]->GetPhys(), 1,
252  normTemporalMeanVec[0], 1, normCrossDir[1], 1,
253  normCrossDir[1], 1);
254  Vmath::Vmul(npoints, fromField[0]->m_exp[nfields - 2]->GetPhys(), 1,
255  normTemporalMeanVec[0], 1, normCrossDir[2], 1);
256  Vmath::Vvtvm(npoints, fromField[0]->m_exp[nfields - 3]->GetPhys(), 1,
257  normTemporalMeanVec[1], 1, normCrossDir[2], 1,
258  normCrossDir[2], 1);
259  }
260 
261  // Compute tawss, trs, osi, taafi, tacfi, WSSG.
262  for (i = 0; i < nfld; ++i)
263  {
264  for (j = 0; j < spacedim; ++j)
265  {
266  Vmath::Vvtvp(npoints, fromField[i]->m_exp[j]->GetPhys(), 1,
267  normTemporalMeanVec[j], 1, DotProduct, 1, DotProduct,
268  1);
269  }
270  // TAWSS
271  Vmath::Vadd(npoints, fromField[i]->m_exp[spacedim]->GetPhys(), 1,
272  outfield[0], 1, outfield[0], 1);
273  // transWSS
274  Vmath::Vmul(npoints, DotProduct, 1, DotProduct, 1, temp, 1);
275  Vmath::Vvtvm(npoints, fromField[i]->m_exp[spacedim]->GetPhys(), 1,
276  fromField[i]->m_exp[spacedim]->GetPhys(), 1, temp, 1, temp,
277  1);
278 
279  for (j = 0; j < npoints; ++j)
280  {
281  if (temp[j] > 0.0)
282  {
283  outfield[1][j] = outfield[1][j] + sqrt(temp[j]);
284  }
285  }
286 
287  // TAAFI
288  Vmath::Vdiv(npoints, DotProduct, 1,
289  fromField[i]->m_exp[spacedim]->GetPhys(), 1, temp, 1);
290  Vmath::Vadd(npoints, temp, 1, outfield[3], 1, outfield[3], 1);
291 
292  // TACFI
293  for (j = 0; j < npoints; ++j)
294  {
295  temp[j] = 1 - temp[j] * temp[j];
296  if (temp[j] > 0.0)
297  {
298  outfield[4][j] = outfield[4][j] + sqrt(temp[j]);
299  }
300  }
301 
302  // WSSG
303  if (wssg)
304  {
305  // parallel component:
306  Vmath::Zero(npoints, temp, 1);
307 
308  fromField[i]->m_exp[0]->PhysDeriv(DotProduct, dtemp[0], dtemp[1],
309  dtemp[2]);
310  for (j = 0; j < spacedim; j++)
311  {
312  Vmath::Vvtvp(npoints, dtemp[j], 1, normTemporalMeanVec[j], 1,
313  temp, 1, temp, 1);
314  }
315  Vmath::Vmul(npoints, temp, 1, temp, 1, temp, 1);
316 
317  // perpendicular component:
318  Vmath::Zero(npoints, DotProduct, 1);
319 
320  for (j = 0; j < spacedim; ++j)
321  {
322  Vmath::Vvtvp(npoints, fromField[i]->m_exp[j]->GetPhys(), 1,
323  normCrossDir[j], 1, DotProduct, 1, DotProduct, 1);
324  }
325  fromField[i]->m_exp[0]->PhysDeriv(DotProduct, dtemp[0], dtemp[1],
326  dtemp[2]);
327  Vmath::Zero(npoints, DotProduct, 1);
328 
329  for (j = 0; j < spacedim; j++)
330  {
331  Vmath::Vvtvp(npoints, dtemp[j], 1, normCrossDir[j], 1,
332  DotProduct, 1, DotProduct, 1);
333  }
334  Vmath::Vvtvp(npoints, DotProduct, 1, DotProduct, 1, temp, 1, temp,
335  1);
336 
337  // Final wssg computation
338  Vmath::Vsqrt(npoints, temp, 1, temp, 1);
339  Vmath::Vadd(npoints, temp, 1, outfield[5], 1, outfield[5], 1);
340  }
341 
342  Vmath::Zero(npoints, DotProduct, 1);
343  }
344 
345  // Divide by nfld
346  Vmath::Smul(npoints, 1.0 / nfld, outfield[0], 1, outfield[0], 1);
347  Vmath::Smul(npoints, 1.0 / nfld, outfield[1], 1, outfield[1], 1);
348  Vmath::Smul(npoints, 1.0 / nfld, outfield[3], 1, outfield[3], 1);
349  Vmath::Smul(npoints, 1.0 / nfld, outfield[4], 1, outfield[4], 1);
350 
351  // OSI
352  for (i = 0; i < npoints; ++i)
353  {
354  outfield[2][i] = 0.5 * (1 - TemporalMeanMag[i] / outfield[0][i]);
355  }
356 
357  /* TAWSS = sum(wss)/nfld
358  * transWSS = sum( sqrt( wss^2 - (wss . normTempMean)^2) )/nfld.
359  * OSI = 0.5*(1-TemporalMeanMag/TAWSS)
360  * TAAFI = sum(cos)/nfld
361  * TACFI = sum(sin)/nfld = sum( sqrt(1-cos^2) )/nfld.
362  */
363 
364  m_f->m_exp.resize(nout);
365  m_f->m_fielddef = fromField[0]->m_fielddef;
366  m_f->m_exp[0] =
367  m_f->SetUpFirstExpList(m_f->m_fielddef[0]->m_numHomogeneousDir, true);
368 
369  for (i = 1; i < nout; ++i)
370  {
371  m_f->m_exp[i] =
372  m_f->AppendExpList(m_f->m_fielddef[0]->m_numHomogeneousDir);
373  }
374 
375  m_f->m_fielddef[0]->m_fields.resize(nout);
376  m_f->m_fielddef[0]->m_fields[0] = "TAWSS";
377  m_f->m_fielddef[0]->m_fields[1] = "transWSS";
378  m_f->m_fielddef[0]->m_fields[2] = "OSI";
379  m_f->m_fielddef[0]->m_fields[3] = "TAAFI";
380  m_f->m_fielddef[0]->m_fields[4] = "TACFI";
381 
382  if (wssg)
383  {
384  m_f->m_fielddef[0]->m_fields[5] = "|WSSG|";
385  Vmath::Smul(npoints, 1.0 / nfld, outfield[5], 1, outfield[5], 1);
386  }
387 
388  m_f->m_variables = m_f->m_fielddef[0]->m_fields;
389 
390  for (i = 0; i < nout; ++i)
391  {
392  m_f->m_exp[i]->FwdTrans(outfield[i], m_f->m_exp[i]->UpdateCoeffs());
393  m_f->m_exp[i]->BwdTrans(m_f->m_exp[i]->GetCoeffs(),
394  m_f->m_exp[i]->UpdatePhys());
395  }
396 }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:215
FieldSharedPtr m_f
Field object.
Definition: Module.h:234
static FieldMetaDataMap NullFieldMetaDataMap
Definition: FieldIO.h:53
double NekDouble
void Vsqrt(int n, const T *x, const int incx, T *y, const int incy)
sqrt y = sqrt(x)
Definition: Vmath.cpp:534
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:209
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:574
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:359
void Vvtvm(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)
vvtvm (vector times vector minus vector): z = w*x - y
Definition: Vmath.cpp:598
void Smul(int n, const T alpha, const T *x, const int incx, T *y, const int incy)
Scalar multiply y = alpha*x.
Definition: Vmath.cpp:248
void Vdiv(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:284
void Zero(int n, T *x, const int incx)
Zero vector.
Definition: Vmath.cpp:492
scalarT< T > sqrt(scalarT< T > in)
Definition: scalar.hpp:294

References ASSERTL0, Nektar::FieldUtils::Module::m_config, Nektar::FieldUtils::Module::m_f, Nektar::LibUtilities::NullFieldMetaDataMap, Vmath::Smul(), tinysimd::sqrt(), Vmath::Vadd(), Vmath::Vdiv(), Vmath::Vmul(), Vmath::Vsqrt(), Vmath::Vvtvm(), Vmath::Vvtvp(), and Vmath::Zero().

Member Data Documentation

◆ className

ModuleKey Nektar::FieldUtils::ProcessMultiShear::className
static
Initial value:
=
"Computes shear stress metrics.")
static std::shared_ptr< Module > create(FieldSharedPtr f)
Creates an instance of this class.
tKey RegisterCreatorFunction(tKey idKey, CreatorFunction classCreator, std::string pDesc="")
Register a class with the factory.
Definition: NekFactory.hpp:198
std::pair< ModuleType, std::string > ModuleKey
Definition: Module.h:317
ModuleFactory & GetModuleFactory()
Definition: Module.cpp:49

Definition at line 57 of file ProcessMultiShear.h.