Nektar++
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
RectangularMesh.cpp
Go to the documentation of this file.
1 #include <cstdlib>
2 #include <cstdio>
3 #include <iomanip>
4 #include <iostream>
5 #include <vector>
6 
8 
9 // Use the stl version, primarily for string.
10 #ifndef TIXML_USE_STL
11 #define TIXML_USE_STL
12 #endif
13 
14 #include <tinyxml.h>
15 
16 
17 void Header(FILE *, int nel);
18 void Middle(FILE *);
19 void End(FILE *);
20 
21 using namespace std;
22 
23 void PrintConditions(void);
24 
25 int main(int argc, char *argv[])
26 {
27  vector<double> xc,yc;
28  int nx = 0, ny = 0;
29  int i,j,k;
30 
31  if(argc != 2)
32  {
33  fprintf(stderr,"Usage RectanuglarMesh file\n");
34  exit(1);
35  }
36 
37  try{
38 
39  TiXmlDocument doc(argv[argc-1]);
40  bool loadOkay = doc.LoadFile();
41 
42  std::stringstream errstr;
43  errstr << "Unable to load file: " << argv[argc-1] << " (";
44  errstr << doc.ErrorDesc() << ", line " << doc.ErrorRow()
45  << ", column " << doc.ErrorCol() << ")";
46  ASSERTL0(loadOkay, errstr.str());
47 
48  TiXmlHandle docHandle(&doc);
49  TiXmlElement* master = NULL;
50  TiXmlElement* block = NULL;
51 
52  master = doc.FirstChildElement("NEKBLOCK");
53  ASSERTL0(master, "Unable to find NEKBLOCK tag in file.");
54 
55  // Find the Mesh tag and same the dim and space attributes
56  block = master->FirstChildElement("XBLOCK");
57 
58  ASSERTL0(block, "Unable to find XBLOCK tag in file.");
59  TiXmlElement *val = block->FirstChildElement("X");
60  while (val)
61  {
62  TiXmlNode *xval = val->FirstChild();
63 
64  std::istringstream valDataStrm(xval->ToText()->Value());
65 
66  try
67  {
68  while(!valDataStrm.fail())
69  {
70  double x_val;
71  valDataStrm >> x_val;
72 
73  if (!valDataStrm.fail())
74  {
75  xc.push_back(x_val);
76  }
77  }
78  }
79  catch(...)
80  {
81  ASSERTL0(false, "Unable to read Xval data.");
82  }
83 
84  val= val->NextSiblingElement("X");
85  }
86 
87  block = master->FirstChildElement("YBLOCK");
88 
89  val = block->FirstChildElement("Y");
90  while (val)
91  {
92  TiXmlNode *yval = val->FirstChild();
93 
94  std::istringstream valDataStrm(yval->ToText()->Value());
95 
96  try
97  {
98  while(!valDataStrm.fail())
99  {
100  double y_val;
101  valDataStrm >> y_val;
102 
103  if (!valDataStrm.fail())
104  {
105  yc.push_back(y_val);
106  }
107  }
108  }
109  catch(...)
110  {
111  ASSERTL0(false, "Unable to read Yval data.");
112  }
113 
114  val= val->NextSiblingElement("Y");
115  }
116 
117  cout << "<?xml version=\"1.0\" encoding=\"utf-8\" ?>" << endl;
118  cout << "<NEKTAR>" << endl;
119 
120 
121  cout << "<EXPANSIONS>" << endl;
122  cout << "<E COMPOSITE=\"C[0]\" NUMMODES=\"7\" FIELDS=\"u\" TYPE=\"MODIFIED\" />" <<endl;
123  cout << "</EXPANSIONS>\n" << endl;
124 
125  PrintConditions();
126 
127  //Vertices
128  cout << "<GEOMETRY DIM=\"2\" SPACE=\"2\">" << endl;
129  cout << " <VERTEX>" << endl;
130 
131  nx = xc.size();
132  ny = yc.size();
133  for(j = 0; j < ny; ++j)
134  {
135  for(k = 0; k < nx; ++k)
136  {
137  cout << " <V ID=\"" << j*nx+k << "\">\t";
138  cout << std::setprecision(8)<< xc[k] << " " << yc[j] << " 0.0";
139  cout << " </V>" << endl;
140  }
141  }
142 
143  cout << " </VERTEX>\n" << endl;
144 
145  cout << " <EDGE>" << endl;
146  int cnt = 0;
147  for(j = 0; j < ny-1; ++j)
148  {
149  for(i = 0; i < nx-1; ++i)
150  {
151  cout << " <E ID=\"" << cnt++ << "\">\t";
152  cout << j*nx+i <<" " << j*nx + i+1 ;
153  cout << " </E>" << endl;
154  }
155 
156  for(i = 0; i < nx; ++i)
157  {
158  cout << " <E ID=\"" << cnt++ << "\">\t";
159  cout << j*nx+i <<" " << (j+1)*nx+i;
160  cout << " </E>" << endl;
161  }
162  }
163 
164  for(i = 0; i < nx-1; ++i)
165  {
166  cout << " <E ID=\"" << cnt++ << "\">\t";
167  cout << j*nx+i <<" " << j*nx + i+1 ;
168  cout << " </E>" << endl;
169  }
170  cout << " </EDGE>\n" << endl;
171 
172  cout << " <ELEMENT>" << endl;
173  cnt = 0;
174  for(j = 0; j < ny-1; ++j)
175  {
176  for(i = 0; i < nx-1; ++i)
177  {
178  cout << " <Q ID=\"" << cnt++ << "\">\t";
179  cout << j*((nx-1)+nx)+i <<" " << j*((nx-1)+nx) + (nx-1)+i+1 << " " ;
180  cout << (j+1)*((nx-1)+nx)+i <<" " << j*((nx-1)+nx) + (nx-1)+i ;
181  cout << " </Q>" << endl;
182  }
183  }
184  cout << " </ELEMENT>\n" << endl;
185 
186 
187  cout << "<COMPOSITE>" << endl;
188  cout << "<C ID=\"0\"> Q[0-" << (nx-1)*(ny-1)-1 << "] </C>" << endl;
189 
190  cout << "<C ID=\"1\"> E[";
191  for(i = 0; i < nx-1; ++i)
192  {
193  cout << i;
194  if(i != nx-2)
195  {
196  cout << ",";
197  }
198  }
199  cout << "] </C> // south border" << endl;
200 
201  cout << "<C ID=\"2\"> E[";
202  for(i = 0; i < ny-1; ++i)
203  {
204  cout << (nx-1)*(i+1) + nx*i;
205  if(i != ny-2)
206  {
207  cout << ",";
208  }
209  }
210  cout << "] </C> // west border" << endl;
211 
212  cout << "<C ID=\"3\"> E[";
213  for(i = 0; i < nx-1; ++i)
214  {
215  cout << (nx-1)*(ny-1)+ nx*(ny-1)+ i;
216  if(i != nx-2)
217  {
218  cout << ",";
219  }
220  }
221  cout << "] </C> // north border" << endl;
222 
223  cout << "<C ID=\"4\"> E[";
224  for(i = 0; i < ny-1; ++i)
225  {
226  cout << (nx-1)*(i+1) + nx*i + nx-1;
227  if(i != ny-2)
228  {
229  cout << ",";
230  }
231  }
232  cout << "] </C> // East border" << endl;
233 
234 
235  cout << "</COMPOSITE>\n" << endl;
236 
237 
238  cout << "<DOMAIN> C[0] </DOMAIN>\n" << endl;
239  cout << "</GEOMETRY>\n" << endl;
240 
241  cout << "</NEKTAR>" << endl;
242 
243  }
244  catch(...)
245  {
246  return 1;
247  }
248 
249  return 0;
250 
251 }
252 
253 
254 
255 void PrintConditions(void)
256 {
257  cout << "<CONDITIONS>" << endl;
258 
259  cout << "<SOLVERINFO>" << endl;
260  cout << "<I PROPERTY=\"SolverType\" VALUE=\" \"/>" << endl;
261  cout << "</SOLVERINFO>\n" << endl;
262 
263  cout << "<PARAMETERS>" << endl;
264  cout << "<P> TimeStep = 0.002 </P>" << endl;
265  cout << "</PARAMETERS>\n" << endl;
266 
267  cout << "<VARIABLES>" << endl;
268  cout << " <V ID=\"0\"> u </V>" << endl;
269  cout << "</VARIABLES>\n" << endl;
270 
271  cout << "<BOUNDARYREGIONS>" << endl;
272  cout << "<B ID=\"0\"> C[1] </B>" << endl;
273  cout << "<B ID=\"1\"> C[2] </B>" << endl;
274  cout << "<B ID=\"2\"> C[3] </B>" << endl;
275  cout << "<B ID=\"3\"> C[4] </B>" << endl;
276  cout << "</BOUNDARYREGIONS>\n" << endl;
277 
278  cout << "<BOUNDARYCONDITIONS>" << endl;
279  cout << " <REGION REF=\"0\"> // South border " << endl;
280  cout << " <D VAR=\"u\" VALUE=\"0\" />" << endl;
281  cout << " </REGION>" << endl;
282 
283  cout << " <REGION REF=\"1\"> // West border " << endl;
284  cout << " <D VAR=\"u\" VALUE=\"0\" />" << endl;
285  cout << " </REGION>" << endl;
286 
287  cout << " <REGION REF=\"2\"> // North border " << endl;
288  cout << " <D VAR=\"u\" VALUE=\"0\" />" << endl;
289  cout << " </REGION>" << endl;
290 
291  cout << " <REGION REF=\"3\"> // East border " << endl;
292  cout << " <D VAR=\"u\" VALUE=\"0\" />" << endl;
293  cout << " </REGION>" << endl;
294  cout << "</BOUNDARYCONDITIONS>" << endl;
295 
296  cout << "</CONDITIONS>" << endl;
297 }
298 
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:188
void PrintConditions(void)
void End(FILE *)
STL namespace.
int main(int argc, char *argv[])
void Middle(FILE *)
void Header(FILE *, int nel)