Nektar++
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
SplitModes.cpp
Go to the documentation of this file.
1 #include <cstdio>
2 #include <cstdlib>
3 
4 #include <StdRegions/StdTriExp.h>
5 #include <SpatialDomains/MeshGraph.h> // for FieldDefinitions, etc
6 
7 using namespace Nektar;
8 
9 int main(int argc, char *argv[])
10 {
11  if(argc != 2)
12  {
13  fprintf(stderr,"Usage: Splitmodes fieldfile \n");
14  exit(1);
15  }
16 
17  //default meshgraph
19 
20 
21  //----------------------------------------------
22  // Import fieldfile2.
23  string fieldfile(argv[argc-1]);
24  vector<LibUtilities::FieldDefinitionsSharedPtr> fielddef;
25  vector<vector<NekDouble> > fielddata;
26  LibUtilities::Import(fieldfile,fielddef,fielddata);
27  //----------------------------------------------
28 
29 
30  ASSERTL0(fielddef[0]->m_numModes[2] > 1,"Expected Fourier field to have at least 2 modes");
31 
32  ASSERTL0(fielddef[0]->m_numHomogeneousDir == 1,"Expected second fld to have one homogeneous direction");
33 
34  int nmodes = fielddef[0]->m_numModes[2];
35 
36  // take off modes from fielddef
37  vector<unsigned int> newNumModes;
38  newNumModes.push_back(fielddef[0]->m_numModes[0]);
39  newNumModes.push_back(fielddef[0]->m_numModes[1]);
40  vector<LibUtilities::BasisType> newBasis;
41  newBasis.push_back(fielddef[0]->m_basis[0]);
42  newBasis.push_back(fielddef[0]->m_basis[1]);
43  for(int i = 0; i < fielddata.size(); ++i)
44  {
45  fielddef[i]->m_numModes = newNumModes;
46  fielddef[i]->m_basis = newBasis;
47  fielddef[i]->m_numHomogeneousDir = 0;
48  }
49 
50 
51  for(int m = 0; m < nmodes; ++m)
52  {
53  vector<vector<NekDouble> > writedata;
54 
55  // outfile
56  string outfile(argv[argc-1]);
57  string out = outfile.substr(0, outfile.find_last_of("."));
58  char num[16] ="";
59  sprintf(num,"%d",(m/2));
60 
61  if(m%2 == 0)
62  {
63  out = out + "_mode_" + num + "_real.fld";
64  }
65  else
66  {
67  out = out + "_mode_" + num + "_imag.fld";
68  }
69 
70  for(int i = 0; i < fielddata.size(); ++i)
71  {
72  vector<NekDouble> newdata;
73 
74  // Determine the number of coefficients per element
75  int ncoeffs;
76  switch(fielddef[i]->m_shapeType)
77  {
79  ncoeffs = LibUtilities::StdTriData::getNumberOfCoefficients(fielddef[i]->m_numModes[0], fielddef[i]->m_numModes[1]);
80  break;
82  ncoeffs = fielddef[i]->m_numModes[0]*fielddef[i]->m_numModes[1];
83  break;
84  default:
85  ASSERTL0(false,"Shape not recognised");
86  break;
87  }
88 
90  vec_iter = fielddata[i].begin();
91 
92  for(int k = 0; k < fielddef[i]->m_fields.size(); ++k)
93  {
94  for(int n = 0; n < fielddef[i]->m_elementIDs.size(); ++n)
95  {
96  // Put orginal mode in here.
97  vec_iter += m*ncoeffs;
98  newdata.insert(newdata.end(),vec_iter,vec_iter+ncoeffs);
99  vec_iter += (nmodes-m)*ncoeffs;
100  }
101  }
102  writedata.push_back(newdata);
103  }
104 
105  //-----------------------------------------------
106  // Write out datafile.
107  LibUtilities::Write(out.c_str(), fielddef, writedata);
108  //-----------------------------------------------
109 
110  }
111  //----------------------------------------------
112 
113 
114  return 0;
115 }
116