Nektar++
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
Functions
CalcDivergence.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 CalcDivergence.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 20 of file CalcDivergence.cpp.

References Nektar::MemoryManager< DataType >::AllocateSharedPtr(), ASSERTL0, Nektar::LibUtilities::SessionReader::CreateInstance(), Nektar::MultiRegions::DirCartesianMap, Nektar::LibUtilities::ePolyEvenlySpaced, Nektar::LibUtilities::Import(), Nektar::SpatialDomains::MeshGraph::Read(), Vmath::Vsub(), and Nektar::LibUtilities::Write().

21 {
22  int i,j;
23 
24  if(argc != 3)
25  {
26  fprintf(stderr,"Usage: ./CalcDivergence file.xml file.fld\n");
27  exit(1);
28  }
29 
30  LibUtilities::SessionReaderSharedPtr vSession
31  = LibUtilities::SessionReader::CreateInstance(argc, argv);
32 
33 
34  //----------------------------------------------
35  // Read in mesh from input file
36  string meshfile(argv[argc-2]);
37  SpatialDomains::MeshGraphSharedPtr graphShPt = SpatialDomains::MeshGraph::Read(vSession);//meshfile);
38  //----------------------------------------------
39 
40  //----------------------------------------------
41  // Import field file.
42  string fieldfile(argv[argc-1]);
43  vector<LibUtilities::FieldDefinitionsSharedPtr> fielddef;
44  vector<vector<NekDouble> > fielddata;
45  LibUtilities::Import(fieldfile,fielddef,fielddata);
46  bool useFFT = false;
47  bool dealiasing = false;
48  //----------------------------------------------
49 
50  //----------------------------------------------
51  // Define Expansion
52  int expdim = graphShPt->GetMeshDimension();
53  int nfields = fielddef[0]->m_fields.size();
54 
55  Array<OneD, MultiRegions::ExpListSharedPtr> Exp(nfields + 1);
56 
57  switch(expdim)
58  {
59  case 1:
60  {
61  ASSERTL0(fielddef[0]->m_numHomogeneousDir <= 2,"Quasi-3D approach is only set up for 1 or 2 homogeneous directions");
62 
63  if(fielddef[0]->m_numHomogeneousDir == 1)
64  {
65  MultiRegions::ExpList2DHomogeneous1DSharedPtr Exp2DH1;
66 
67  // Define Homogeneous expansion
68  //int nplanes = fielddef[0]->m_numModes[1];
69 
70  int nplanes;
71  vSession->LoadParameter("HomModesZ",nplanes,fielddef[0]->m_numModes[1]);
72 
73  // choose points to be at evenly spaced points at
74  const LibUtilities::PointsKey Pkey(nplanes+1,LibUtilities::ePolyEvenlySpaced);
75  const LibUtilities::BasisKey Bkey(fielddef[0]->m_basis[1],nplanes,Pkey);
76  NekDouble ly = fielddef[0]->m_homogeneousLengths[0];
77 
78  Exp2DH1 = MemoryManager<MultiRegions::ExpList2DHomogeneous1D>::AllocateSharedPtr(vSession,Bkey,ly,useFFT,dealiasing,graphShPt);
79  Exp[0] = Exp2DH1;
80 
81  for(i = 1; i < nfields + 1; ++i)
82  {
83  Exp[i] = MemoryManager<MultiRegions::ExpList2DHomogeneous1D>::AllocateSharedPtr(*Exp2DH1);
84  }
85  }
86  else if(fielddef[0]->m_numHomogeneousDir == 2)
87  {
88  MultiRegions::ExpList3DHomogeneous2DSharedPtr Exp3DH2;
89 
90  // Define Homogeneous expansion
91  //int nylines = fielddef[0]->m_numModes[1];
92  //int nzlines = fielddef[0]->m_numModes[2];
93 
94  int nylines;
95  int nzlines;
96  vSession->LoadParameter("HomModesY",nylines,fielddef[0]->m_numModes[1]);
97  vSession->LoadParameter("HomModesZ",nzlines,fielddef[0]->m_numModes[2]);
98 
99  // choose points to be at evenly spaced points at
100  const LibUtilities::PointsKey PkeyY(nylines+1,LibUtilities::ePolyEvenlySpaced);
101  const LibUtilities::BasisKey BkeyY(fielddef[0]->m_basis[1],nylines,PkeyY);
102 
103  const LibUtilities::PointsKey PkeyZ(nzlines+1,LibUtilities::ePolyEvenlySpaced);
104  const LibUtilities::BasisKey BkeyZ(fielddef[0]->m_basis[2],nzlines,PkeyZ);
105 
106  NekDouble ly = fielddef[0]->m_homogeneousLengths[0];
107  NekDouble lz = fielddef[0]->m_homogeneousLengths[1];
108 
109  Exp3DH2 = MemoryManager<MultiRegions::ExpList3DHomogeneous2D>::AllocateSharedPtr(vSession,BkeyY,BkeyZ,ly,lz,useFFT,dealiasing,graphShPt);
110  Exp[0] = Exp3DH2;
111 
112  for(i = 1; i < nfields + 1; ++i)
113  {
114  Exp[i] = MemoryManager<MultiRegions::ExpList3DHomogeneous2D>::AllocateSharedPtr(*Exp3DH2);
115  }
116  }
117  else
118  {
119  MultiRegions::ExpList1DSharedPtr Exp1D;
120  Exp1D = MemoryManager<MultiRegions::ExpList1D>
121  ::AllocateSharedPtr(vSession,graphShPt);
122  Exp[0] = Exp1D;
123  for(i = 1; i < nfields + 1; ++i)
124  {
125  Exp[i] = MemoryManager<MultiRegions::ExpList1D>
126  ::AllocateSharedPtr(*Exp1D);
127  }
128  }
129  }
130  break;
131  case 2:
132  {
133  ASSERTL0(fielddef[0]->m_numHomogeneousDir <= 1,"NumHomogeneousDir is only set up for 1");
134 
135  if(fielddef[0]->m_numHomogeneousDir == 1)
136  {
137  MultiRegions::ExpList3DHomogeneous1DSharedPtr Exp3DH1;
138 
139  // Define Homogeneous expansion
140  //int nplanes = fielddef[0]->m_numModes[2];
141 
142  int nplanes;
143  vSession->LoadParameter("HomModesZ",nplanes,fielddef[0]->m_numModes[2]);
144 
145  // choose points to be at evenly spaced points at
146  // nplanes + 1 points
147  const LibUtilities::PointsKey Pkey(nplanes+1,LibUtilities::ePolyEvenlySpaced);
148  const LibUtilities::BasisKey Bkey(fielddef[0]->m_basis[2],nplanes,Pkey);
149  NekDouble lz = fielddef[0]->m_homogeneousLengths[0];
150 
151  Exp3DH1 = MemoryManager<MultiRegions::ExpList3DHomogeneous1D>::AllocateSharedPtr(vSession,Bkey,lz,useFFT,dealiasing,graphShPt);
152  Exp[0] = Exp3DH1;
153 
154  for(i = 1; i < nfields + 1; ++i)
155  {
156  Exp[i] = MemoryManager<MultiRegions::ExpList3DHomogeneous1D>::AllocateSharedPtr(*Exp3DH1);
157  }
158  }
159  else
160  {
161  MultiRegions::ExpList2DSharedPtr Exp2D;
162  Exp2D = MemoryManager<MultiRegions::ExpList2D>
163  ::AllocateSharedPtr(vSession,graphShPt);
164  Exp[0] = Exp2D;
165 
166  for(i = 1; i < nfields + 1; ++i)
167  {
168  Exp[i] = MemoryManager<MultiRegions::ExpList2D>
169  ::AllocateSharedPtr(*Exp2D);
170  }
171  }
172  }
173  break;
174  case 3:
175  {
176  MultiRegions::ExpList3DSharedPtr Exp3D;
177  Exp3D = MemoryManager<MultiRegions::ExpList3D>
178  ::AllocateSharedPtr(vSession,graphShPt);
179  Exp[0] = Exp3D;
180 
181  for(i = 1; i < nfields + 1; ++i)
182  {
183  Exp[i] = MemoryManager<MultiRegions::ExpList3D>
184  ::AllocateSharedPtr(*Exp3D);
185  }
186  }
187  break;
188  default:
189  ASSERTL0(false,"Expansion dimension not recognised");
190  break;
191  }
192 
193  //----------------------------------------------
194  // Copy data from field file
195  for(j = 0; j < nfields; ++j)
196  {
197  for(unsigned int i = 0; i < fielddata.size(); ++i)
198  {
199  Exp[j]->ExtractDataToCoeffs(fielddef [i],
200  fielddata[i],
201  fielddef [i]->m_fields[j],
202  Exp[j]->UpdateCoeffs());
203  }
204  Exp[j]->BwdTrans(Exp[j]->GetCoeffs(),Exp[j]->UpdatePhys());
205  }
206  //----------------------------------------------
207 
208  int nq = Exp[0]->GetNpoints();
209 
210  Array<OneD, NekDouble> Ux(nq);
211  Array<OneD, NekDouble> Vy(nq);
212  Array<OneD, NekDouble> Div(nq);
213 
214  switch(expdim)
215  {
216  case 1:
217  {
218  if(fielddef[0]->m_numHomogeneousDir == 1)
219  {
220  ASSERTL0(false,"Not implemented yet");
221  }
222  else if(fielddef[0]->m_numHomogeneousDir == 2)
223  {
224  ASSERTL0(false,"Not implemented yet");
225  }
226  else
227  {
228  ASSERTL0(false,"Not implemented yet");
229  }
230  }
231  break;
232  case 2:
233  {
234  if(fielddef[0]->m_numHomogeneousDir == 1)
235  {
236  ASSERTL0(false,"Not implemented yet");
237  }
238  else
239  {
240  Exp[0]->PhysDeriv(MultiRegions::DirCartesianMap[0],Exp[0]->GetPhys(),Ux);
241  Exp[0]->PhysDeriv(MultiRegions::DirCartesianMap[1],Exp[1]->GetPhys(),Vy);
242 
243  Vmath::Vsub(nq,Ux,1,Vy,1,Div,1);
244 
245  Exp[3]->FwdTrans(Div,Exp[3]->UpdateCoeffs());
246  }
247  }
248  break;
249  case 3:
250  {
251  ASSERTL0(false,"Not implemented yet");
252  }
253  break;
254  default:
255  {
256  ASSERTL0(false,"Expansion dimension not recognised");
257  }
258  break;
259  }
260 
261  //-----------------------------------------------
262  // Write solution to file with additional computed fields
263  string fldfilename(argv[2]);
264  string out = fldfilename.substr(0, fldfilename.find_last_of("."));
265  string endfile("_with_divergence.fld");
266  out += endfile;
267  std::vector<LibUtilities::FieldDefinitionsSharedPtr> FieldDef
268  = Exp[0]->GetFieldDefinitions();
269  std::vector<std::vector<NekDouble> > FieldData(FieldDef.size());
270 
271  for(j = 0; j < nfields + 1; ++j)
272  {
273  for(i = 0; i < FieldDef.size(); ++i)
274  {
275  if (j >= nfields)
276  {
277  FieldDef[i]->m_fields.push_back("Div");
278  }
279  else
280  {
281  FieldDef[i]->m_fields.push_back(fielddef[i]->m_fields[j]);
282  }
283  Exp[j]->AppendFieldData(FieldDef[i], FieldData[i]);
284  }
285  }
286  LibUtilities::Write(out, FieldDef, FieldData);
287  //-----------------------------------------------
288 
289  return 0;
290 }
ASSERTL0
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
Nektar::MemoryManager
General purpose memory allocation routines with the ability to allocate from thread specific memory p...
Definition: NekMemoryManager.hpp:102
Nektar::LibUtilities::SessionReaderSharedPtr
boost::shared_ptr< SessionReader > SessionReaderSharedPtr
Definition: MeshPartition.h:51
Nektar::MultiRegions::ExpList3DHomogeneous2DSharedPtr
boost::shared_ptr< ExpList3DHomogeneous2D > ExpList3DHomogeneous2DSharedPtr
Shared pointer to an ExpList3DHomogeneous2D object.
Definition: ExpList3DHomogeneous2D.h:50
Nektar::LibUtilities::Import
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
Nektar::LibUtilities::ePolyEvenlySpaced
1D Evenly-spaced points using Lagrange polynomial
Definition: PointsType.h:63
Nektar::MultiRegions::ExpList1DSharedPtr
boost::shared_ptr< ExpList1D > ExpList1DSharedPtr
Shared pointer to an ExpList1D object.
Definition: ExpList1D.h:50
Nektar::LibUtilities::PointsKey
Defines a specification for a set of points.
Definition: Points.h:58
Nektar::NekDouble
double NekDouble
Definition: NektarUnivTypeDefs.hpp:44
Nektar::MultiRegions::ExpList2DSharedPtr
boost::shared_ptr< ExpList2D > ExpList2DSharedPtr
Shared pointer to an ExpList2D object.
Definition: ExpList2D.h:49
Nektar::MultiRegions::ExpList2DHomogeneous1DSharedPtr
boost::shared_ptr< ExpList2DHomogeneous1D > ExpList2DHomogeneous1DSharedPtr
Shared pointer to an ExpList2DHomogeneous1D object.
Definition: ExpList2DHomogeneous1D.h:50
Vmath::Vsub
void Vsub(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Subtract vector z = x-y.
Definition: Vmath.cpp:329
Nektar::MultiRegions::ExpList3DSharedPtr
boost::shared_ptr< ExpList3D > ExpList3DSharedPtr
Shared pointer to an ExpList3D object.
Definition: ExpList3D.h:114
Nektar::MultiRegions::DirCartesianMap
MultiRegions::Direction const DirCartesianMap[]
Definition: ExpList.h:86
Nektar::MultiRegions::ExpList3DHomogeneous1DSharedPtr
boost::shared_ptr< ExpList3DHomogeneous1D > ExpList3DHomogeneous1DSharedPtr
Shared pointer to an ExpList3DHomogeneous1D object.
Definition: ExpList3DHomogeneous1D.h:50
Nektar::LibUtilities::Write
void Write(const std::string &outFile, std::vector< FieldDefinitionsSharedPtr > &fielddefs, std::vector< std::vector< NekDouble > > &fielddata, const FieldMetaDataMap &fieldinfomap)
Write a field file in serial only.
Definition: FieldIO.cpp:81
Nektar::SpatialDomains::MeshGraphSharedPtr
boost::shared_ptr< MeshGraph > MeshGraphSharedPtr
Definition: MeshGraph.h:442
Nektar::LibUtilities::BasisKey
Describes the specification for a Basis.
Definition: Basis.h:50
Generated on Wed Mar 23 2016 20:14:24 for Nektar++ by   doxygen 1.8.8