Nektar++
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
Functions
ProbeFld.cpp File Reference
#include <cstdio>
#include <cstdlib>
#include <MultiRegions/ExpList.h>
#include <MultiRegions/ExpList1D.h>
#include <MultiRegions/ExpList2D.h>
#include <MultiRegions/ExpList3D.h>
#include <MultiRegions/ExpList2DHomogeneous1D.h>
#include <MultiRegions/ExpList3DHomogeneous1D.h>
#include <MultiRegions/ExpList3DHomogeneous2D.h>
Include dependency graph for ProbeFld.cpp:

Go to the source code of this file.

Functions

int main (int argc, char *argv[])
 

Function Documentation

int main ( int  argc,
char *  argv[] 
)

Definition at line 15 of file ProbeFld.cpp.

References ASSERTL0, Nektar::LibUtilities::ePolyEvenlySpaced, Nektar::LibUtilities::Import(), and Nektar::NekConstants::kGeomFactorsTol.

16 {
17  int i,j;
18  bool file = false;
19 
20  if(argc == 4)
21  {
22  file = true;
23  }
24  else if(argc != 10)
25  {
26  fprintf(stderr,
27  "Usage: ProbeFld meshfile fieldfile N x0 y0 z0 dx dy dz\n");
28  fprintf(stderr,
29  " Probes N points along the line from (x0,y0,z0) to "
30  "(x0+dx, y0+dy, z0+dz)\n");
31  fprintf(stderr,
32  "ProbeFld meshfile fieldfile points.txt\n");
33  fprintf(stderr,
34  " Probes the solution at the points in the points.txt file.\n");
35  fprintf(stderr,
36  " Points are given as space-separated x y z on each line.\n");
37  exit(1);
38  }
39 
41  = LibUtilities::SessionReader::CreateInstance(3, argv);
42 
43  //----------------------------------------------
44  // Read in mesh from input file
45  string meshfile(argv[1]);
46  SpatialDomains::MeshGraphSharedPtr graphShPt = SpatialDomains::MeshGraph::Read(vSession);//meshfile);
47  //----------------------------------------------
48 
49  //----------------------------------------------
50  // Import field file.
51  string fieldfile(argv[2]);
52  vector<LibUtilities::FieldDefinitionsSharedPtr> fielddef;
53  vector<vector<NekDouble> > fielddata;
54  LibUtilities::Import(fieldfile,fielddef,fielddata);
55  //----------------------------------------------
56 
57  //----------------------------------------------
58  // Set up Expansion information
59  vector< vector<LibUtilities::PointsType> > pointstype;
60  for(i = 0; i < fielddef.size(); ++i)
61  {
62  vector<LibUtilities::PointsType> ptype;
63  for(j = 0; j < 3; ++j)
64  {
65  ptype.push_back(LibUtilities::ePolyEvenlySpaced);
66  }
67  pointstype.push_back(ptype);
68  }
69  graphShPt->SetExpansions(fielddef,pointstype);
70  bool useFFT = false;
71  bool dealiasing = false;
72  //----------------------------------------------
73 
74 
75  //----------------------------------------------
76  // Define Expansion
77  int expdim = graphShPt->GetMeshDimension();
78  int nfields = fielddef[0]->m_fields.size();
80 
81  switch(expdim)
82  {
83  case 1:
84  {
85  ASSERTL0(fielddef[0]->m_numHomogeneousDir <= 2,"NumHomogeneousDir is only set up for 1 or 2");
86 
87  if(fielddef[0]->m_numHomogeneousDir == 1)
88  {
90 
91  // Define Homogeneous expansion
92  //int nplanes = fielddef[0]->m_numModes[1];
93  int nplanes;
94  vSession->LoadParameter("HomModesZ",nplanes,fielddef[0]->m_numModes[1]);
95 
96  // choose points to be at evenly spaced points at
98  const LibUtilities::BasisKey Bkey(fielddef[0]->m_basis[1],nplanes,Pkey);
99  NekDouble ly = fielddef[0]->m_homogeneousLengths[0];
100 
101  Exp2DH1 = MemoryManager<MultiRegions::ExpList2DHomogeneous1D>::AllocateSharedPtr(vSession,Bkey,ly,useFFT,dealiasing,graphShPt);
102  Exp[0] = Exp2DH1;
103 
104  for(i = 1; i < nfields; ++i)
105  {
107  }
108  }
109  else if(fielddef[0]->m_numHomogeneousDir == 2)
110  {
112 
113  // Define Homogeneous expansion
114  //int nylines = fielddef[0]->m_numModes[1];
115  //int nzlines = fielddef[0]->m_numModes[2];
116 
117  int nylines;
118  int nzlines;
119  vSession->LoadParameter("HomModesY",nylines,fielddef[0]->m_numModes[1]);
120  vSession->LoadParameter("HomModesZ",nzlines,fielddef[0]->m_numModes[2]);
121 
122  // choose points to be at evenly spaced points at
124  const LibUtilities::BasisKey BkeyY(fielddef[0]->m_basis[1],nylines,PkeyY);
125 
127  const LibUtilities::BasisKey BkeyZ(fielddef[0]->m_basis[2],nzlines,PkeyZ);
128 
129  NekDouble ly = fielddef[0]->m_homogeneousLengths[0];
130  NekDouble lz = fielddef[0]->m_homogeneousLengths[1];
131 
132  Exp3DH2 = MemoryManager<MultiRegions::ExpList3DHomogeneous2D>::AllocateSharedPtr(vSession,BkeyY,BkeyZ,ly,lz,useFFT,dealiasing,graphShPt);
133  Exp[0] = Exp3DH2;
134 
135  for(i = 1; i < nfields; ++i)
136  {
138  }
139  }
140  else
141  {
144  ::AllocateSharedPtr(vSession,graphShPt);
145  Exp[0] = Exp1D;
146  for(i = 1; i < nfields; ++i)
147  {
149  ::AllocateSharedPtr(*Exp1D);
150  }
151  }
152  }
153  break;
154  case 2:
155  {
156  ASSERTL0(fielddef[0]->m_numHomogeneousDir <= 1,"NumHomogeneousDir is only set up for 1");
157 
158  if(fielddef[0]->m_numHomogeneousDir == 1)
159  {
161 
162  // Define Homogeneous expansion
163  //int nplanes = fielddef[0]->m_numModes[2];
164 
165  int nplanes;
166  vSession->LoadParameter("HomModesZ",nplanes,fielddef[0]->m_numModes[2]);
167 
168  // choose points to be at evenly spaced points at
169  // nplanes + 1 points
171  const LibUtilities::BasisKey Bkey(fielddef[0]->m_basis[2],nplanes,Pkey);
172  NekDouble lz = fielddef[0]->m_homogeneousLengths[0];
173 
174  Exp3DH1 = MemoryManager<MultiRegions::ExpList3DHomogeneous1D>::AllocateSharedPtr(vSession,Bkey,lz,useFFT,dealiasing,graphShPt);
175  Exp[0] = Exp3DH1;
176 
177  for(i = 1; i < nfields; ++i)
178  {
180  }
181  }
182  else
183  {
186  ::AllocateSharedPtr(vSession,graphShPt);
187  Exp[0] = Exp2D;
188 
189  for(i = 1; i < nfields; ++i)
190  {
192  ::AllocateSharedPtr(*Exp2D);
193  }
194  }
195  }
196  break;
197  case 3:
198  {
201  ::AllocateSharedPtr(vSession,graphShPt);
202  Exp[0] = Exp3D;
203 
204  for(i = 1; i < nfields; ++i)
205  {
207  ::AllocateSharedPtr(*Exp3D);
208  }
209  }
210  break;
211 
212  default:
213  ASSERTL0(false,"Expansion dimension not recognised");
214  break;
215  }
216  //----------------------------------------------
217 
218  //----------------------------------------------
219  // Copy data from field file
220  for(j = 0; j < nfields; ++j)
221  {
222  for(int i = 0; i < fielddata.size(); ++i)
223  {
224  Exp[j]->ExtractDataToCoeffs(fielddef [i],
225  fielddata[i],
226  fielddef [i]->m_fields[j],
227  Exp[j]->UpdateCoeffs());
228  }
229  Exp[j]->BwdTrans(Exp[j]->GetCoeffs(),Exp[j]->UpdatePhys());
230  }
231  //----------------------------------------------
232 
233  Array<OneD, NekDouble> gloCoord(3,0.0);
234  if (file)
235  {
236  string line;
237  ifstream pts(argv[3]);
238  while (getline(pts, line))
239  {
240  stringstream ss(line);
241  ss >> gloCoord[0];
242  ss >> gloCoord[1];
243  ss >> gloCoord[2];
244  cout << gloCoord[0] << " " << gloCoord[1] << " " << gloCoord[2];
245  int ExpId = Exp[0]->GetExpIndex(gloCoord,NekConstants::kGeomFactorsTol);
246 
247  for (int j = 0; j < nfields; ++j)
248  {
249  if (ExpId == -1)
250  {
251  cout << " -";
252  }
253  else
254  {
255  Array<OneD, NekDouble> phys(Exp[j]->GetPhys() + Exp[j]->GetPhys_Offset(ExpId));
256  cout << " " << Exp[j]->GetExp(ExpId)->PhysEvaluate(gloCoord, phys);
257  }
258  }
259 
260  cout << endl;
261  }
262  }
263  else
264  {
265  //----------------------------------------------
266  // Probe data fields
267  NekDouble N = atoi(argv[3]);
268  NekDouble x0 = atof(argv[4]);
269  NekDouble y0 = atof(argv[5]);
270  NekDouble z0 = atof(argv[6]);
271  NekDouble dx = atof(argv[7])/(N>1 ? (N-1) : 1);
272  NekDouble dy = atof(argv[8])/(N>1 ? (N-1) : 1);
273  NekDouble dz = atof(argv[9])/(N>1 ? (N-1) : 1);
274 
275  for (int i = 0; i < N; ++i)
276  {
277  gloCoord[0] = x0 + i*dx;
278  gloCoord[1] = y0 + i*dy;
279  gloCoord[2] = z0 + i*dz;
280  cout << gloCoord[0] << " " << gloCoord[1] << " " << gloCoord[2];
281  int ExpId = Exp[0]->GetExpIndex(gloCoord,NekConstants::kGeomFactorsTol);
282 
283  for (int j = 0; j < nfields; ++j)
284  {
285  if (ExpId == -1)
286  {
287  cout << " -";
288  }
289  else
290  {
291  Array<OneD, NekDouble> phys(Exp[j]->GetPhys() + Exp[j]->GetPhys_Offset(ExpId));
292  cout << " " << Exp[j]->GetExp(ExpId)->PhysEvaluate(gloCoord, phys);
293  }
294  }
295 
296  cout << endl;
297  }
298  }
299 
300  //----------------------------------------------
301  return 0;
302 }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
General purpose memory allocation routines with the ability to allocate from thread specific memory p...
boost::shared_ptr< SessionReader > SessionReaderSharedPtr
Definition: MeshPartition.h:51
boost::shared_ptr< ExpList3DHomogeneous2D > ExpList3DHomogeneous2DSharedPtr
Shared pointer to an ExpList3DHomogeneous2D object.
void Import(const std::string &infilename, std::vector< FieldDefinitionsSharedPtr > &fielddefs, std::vector< std::vector< NekDouble > > &fielddata, FieldMetaDataMap &fieldinfomap, const Array< OneD, int > ElementiDs)
Imports an FLD file.
Definition: FieldIO.cpp:115
1D Evenly-spaced points using Lagrange polynomial
Definition: PointsType.h:63
boost::shared_ptr< ExpList1D > ExpList1DSharedPtr
Shared pointer to an ExpList1D object.
Definition: ExpList1D.h:50
Defines a specification for a set of points.
Definition: Points.h:58
double NekDouble
boost::shared_ptr< ExpList2D > ExpList2DSharedPtr
Shared pointer to an ExpList2D object.
Definition: ExpList2D.h:49
boost::shared_ptr< ExpList2DHomogeneous1D > ExpList2DHomogeneous1DSharedPtr
Shared pointer to an ExpList2DHomogeneous1D object.
boost::shared_ptr< ExpList3D > ExpList3DSharedPtr
Shared pointer to an ExpList3D object.
Definition: ExpList3D.h:114
static const NekDouble kGeomFactorsTol
boost::shared_ptr< ExpList3DHomogeneous1D > ExpList3DHomogeneous1DSharedPtr
Shared pointer to an ExpList3DHomogeneous1D object.
boost::shared_ptr< MeshGraph > MeshGraphSharedPtr
Definition: MeshGraph.h:442
Describes the specification for a Basis.
Definition: Basis.h:50