Nektar++
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
InputSwan.cpp
Go to the documentation of this file.
1 ////////////////////////////////////////////////////////////////////////////////
2 //
3 // File: InputSwan.cpp
4 //
5 // For more information, please see: http://www.nektar.info/
6 //
7 // The MIT License
8 //
9 // Copyright (c) 2006 Division of Applied Mathematics, Brown University (USA),
10 // Department of Aeronautics, Imperial College London (UK), and Scientific
11 // Computing and Imaging Institute, University of Utah (USA).
12 //
13 // License for the specific language governing rights and limitations under
14 // Permission is hereby granted, free of charge, to any person obtaining a
15 // copy of this software and associated documentation files (the "Software"),
16 // to deal in the Software without restriction, including without limitation
17 // the rights to use, copy, modify, merge, publish, distribute, sublicense,
18 // and/or sell copies of the Software, and to permit persons to whom the
19 // Software is furnished to do so, subject to the following conditions:
20 //
21 // The above copyright notice and this permission notice shall be included
22 // in all copies or substantial portions of the Software.
23 //
24 // THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
25 // OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
26 // FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
27 // THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
28 // LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
29 // FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
30 // DEALINGS IN THE SOFTWARE.
31 //
32 // Description: Swansea session converter.
33 //
34 ////////////////////////////////////////////////////////////////////////////////
35 
36 #include <iostream>
37 using namespace std;
38 
39 #include "MeshElements.h"
40 #include "InputSwan.h"
41 
42 namespace Nektar
43 {
44  namespace Utilities
45  {
46  ModuleKey InputSwan::className =
48  ModuleKey(eInputModule, "plt"), InputSwan::create,
49  "Reads Swansea plt format for third-order tetrahedra.");
50 
51  InputSwan::InputSwan(MeshSharedPtr m) : InputModule(m)
52  {
53 
54  }
55 
57  {
58 
59  }
60 
62  {
63  // Open the file stream.
64  OpenStream();
65 
66  vector<vector<NodeSharedPtr> > elementList;
67  vector<int> tmp, tets;
68  vector<double> pts;
69 
70  if (m_mesh->m_verbose)
71  {
72  cout << "InputSwan: Start reading file..." << endl;
73  }
74 
75  m_mesh->m_expDim = 3;
76  m_mesh->m_spaceDim = 3;
77 
78  // First read in header; 4 integers containing number of tets,
79  // number of points, nas2 (unknown) and order of the grid
80  // (i.e. GridOrder = 3 => cubic mesh).
81  tmp.resize(6);
82  m_mshFile.read(reinterpret_cast<char*>(&tmp[0]),
83  static_cast<int> (6*sizeof(int)));
84 
85  if (tmp[0] != tmp[5] || tmp[0] != 4*sizeof(int))
86  {
87  cout << "Header data broken" << endl;
88  m_mshFile.close();
89  return;
90  }
91 
92  int NB_Tet = tmp[1];
93  int NB_Points = tmp[2];
94  int GridOrder = tmp[4];
95  int ND = (GridOrder+1)*(GridOrder+2)*(GridOrder+3)/6;
96 
97  cout << "NB_Tet = " << NB_Tet << endl;
98  cout << "NB_Points = " << NB_Points << endl;
99  cout << "GridOrder = " << GridOrder << endl;
100  cout << "ND = " << ND << endl;
101 
102  // Now read in list of tetrahedra. Each tet has ND points, and
103  // data are ordered with memory traversed fastest with point
104  // number.
105  tets.resize(ND*NB_Tet+2);
106  m_mshFile.read(reinterpret_cast<char*>(&tets[0]),
107  static_cast<int> ((ND*NB_Tet+2)*sizeof(int)));
108 
109  if (tets[0] != tets[ND*NB_Tet+1] || tets[0] != ND*NB_Tet*sizeof(int))
110  {
111  cout << "ERROR [InputSwan]: Tetrahedron data broken." << endl;
112  m_mshFile.close();
113  return;
114  }
115 
116  // Finally, read point data: NB_Points tuples (x,y,z).
117  tmp.resize(2);
118  pts.resize(3*NB_Points);
119  m_mshFile.read(reinterpret_cast<char*>(&tmp[0]),
120  static_cast<int> (sizeof(int)));
121  m_mshFile.read(reinterpret_cast<char*>(&pts[0]),
122  static_cast<int> (3*NB_Points*sizeof(double)));
123  m_mshFile.read(reinterpret_cast<char*>(&tmp[1]),
124  static_cast<int> (sizeof(int)));
125 
126  if (tmp[0] != tmp[1] || tmp[0] != 3*NB_Points*sizeof(double))
127  {
128  cout << "ERROR [InputSwan]: Point data broken." << endl;
129  m_mshFile.close();
130  return;
131  }
132 
133  int vid = 0, i, j;
135 
136  // Read in list of vertices.
137  for (i = 0; i < NB_Points; ++i)
138  {
139  double x = pts [ i];
140  double y = pts [1*NB_Points+i];
141  double z = pts [2*NB_Points+i];
142  m_mesh->m_node.push_back(boost::shared_ptr<Node>(new Node(vid, x, y, z)));
143  vid++;
144  }
145 
146  // Iterate over list of tetrahedra: for each, create nodes. At the
147  // moment discard high order data and create linear mesh.
148  for (i = 0; i < NB_Tet; ++i)
149  {
150  vector<NodeSharedPtr> nodeList;
151  for (j = 0; j < 20; ++j)
152  {
153  int vid = tets[j*NB_Tet+i+1];
154  nodeList.push_back(m_mesh->m_node[vid-1]);
155  }
156 
157  vector<int> tags;
158  tags.push_back(0); // composite
159  tags.push_back(elType); // element type
160 
161  ElmtConfig conf(elType,3,true,true);
163  CreateInstance(elType,conf,nodeList,tags);
164  m_mesh->m_element[3].push_back(E);
165  }
166 
167  // Attempt to read in composites. Need to determine number of
168  // triangles from data.
169  tmp.resize(2);
170  m_mshFile.read(reinterpret_cast<char*>(&tmp[0]),
171  static_cast<int> (sizeof(int)));
172  int n_tri = tmp[0]/sizeof(int)/5;
173  tets.resize(n_tri*5);
174  m_mshFile.read(reinterpret_cast<char*>(&tets[0]),
175  static_cast<int> (tmp[0]));
176  m_mshFile.read(reinterpret_cast<char*>(&tmp[1]),
177  static_cast<int> (sizeof(int)));
178 
179  if (tmp[0] != tmp[1])
180  {
181  cout << "ERROR [InputSwan]: Surface data broken." << endl;
182  m_mshFile.close();
183  return;
184  }
185 
186  elType = LibUtilities::eTriangle;
187 
188  // Process list of triangles forming surfaces.
189  for (i = 0; i < n_tri; ++i)
190  {
191  vector<NodeSharedPtr> nodeList;
192 
193  for (j = 0; j < 3; ++j)
194  {
195  nodeList.push_back(m_mesh->m_node[tets[i+j*n_tri]-1]);
196  }
197 
198  vector<int> tags;
199  tags.push_back(1); // composite
200  tags.push_back(elType); // element type
201 
202  ElmtConfig conf(elType,1,false,false);
204  CreateInstance(elType,conf,nodeList,tags);
205  m_mesh->m_element[2].push_back(E);
206  }
207 
208  m_mshFile.close();
209 
210  // Process the rest of the mesh.
211  ProcessVertices ();
212  ProcessEdges ();
213  ProcessFaces ();
214  ProcessElements ();
216  }
217  }
218 }