#include <mcnpmit.hpp>
Public Member Functions | |
McnpData () | |
~McnpData () | |
MCNPError | set_coord_system (int) |
int | get_coord_system () |
MCNPError | set_rotation_matrix (double[16]) |
double * | get_rotation_matrix () |
MCNPError | set_filename (std::string) |
std::string | get_filename () |
MCNPError | read_mcnpfile (bool) |
MCNPError | read_coord_system (std::string) |
MCNPError | read_rotation_matrix (std::string, int) |
MCNPError | make_elements (std::vector< double >[3], int *) |
MCNPError | make_adjacencies (int *) |
MCNPError | initialize_tags () |
MCNPError | extract_tally_data (std::string, moab::EntityHandle) |
MCNPError | transform_point (double *, double *, int, double *) |
Public Attributes | |
int | coord_system |
double | rotation_matrix [16] |
std::vector< moab::EntityHandle > | MCNP_vertices |
std::vector< moab::EntityHandle > | MCNP_elems |
moab::Range | vert_handles |
moab::Range | elem_handles |
moab::Tag | box_min_tag |
moab::Tag | box_max_tag |
moab::Tag | tally_tag |
moab::Tag | relerr_tag |
std::string | MCNP_filename |
Definition at line 26 of file mcnpmit.hpp.
McnpData::McnpData | ( | ) |
Definition at line 20 of file mcnpmit.cpp.
21 {
22
23 // Default value for coordinate system
24 coord_system = 0;
25
26 // Default rotation matrix is identity matrix
27 for( int i = 0; i < 4; i++ )
28 {
29 for( int j = 0; j < 4; j++ )
30 {
31 if( i == j )
32 rotation_matrix[4 * i + j] = 1;
33 else
34 rotation_matrix[4 * i + j] = 0;
35 }
36 }
37 }
References coord_system, and rotation_matrix.
McnpData::~McnpData | ( | ) |
Definition at line 40 of file mcnpmit.cpp.
41 {
42
43 // Vertices and elements
44 MCNP_vertices.clear();
45 }
References MCNP_vertices.
MCNPError McnpData::extract_tally_data | ( | std::string | s, |
moab::EntityHandle | handle | ||
) |
Definition at line 335 of file mcnpmit.cpp.
336 {
337
338 int fpos = 0;
339 double d = 0;
340
341 MCNPError result;
342 moab::ErrorCode MBresult;
343
344 // Discard first three lines
345 for( int i = 0; i < 3; i++ )
346 {
347 result = next_number( s, d, fpos );
348 if( result == MCNP_FAILURE ) return MCNP_FAILURE;
349 }
350 // Need to read in tally entry and tag ...
351 result = next_number( s, d, fpos );
352 if( result == MCNP_FAILURE ) return MCNP_FAILURE;
353 MBresult = MBI->tag_set_data( tally_tag, &handle, 1, &d );
354 if( MBresult != moab::MB_SUCCESS ) return MCNP_FAILURE;
355
356 // Need to read in relative error entry and tag ...
357 result = next_number( s, d, fpos );
358 if( result == MCNP_FAILURE ) return MCNP_FAILURE;
359 MBresult = MBI->tag_set_data( relerr_tag, &handle, 1, &d );
360 if( MBresult != moab::MB_SUCCESS ) return MCNP_FAILURE;
361
362 return MCNP_SUCCESS;
363 }
References ErrorCode, MB_SUCCESS, MBI, MCNP_FAILURE, MCNP_SUCCESS, next_number(), relerr_tag, and tally_tag.
Referenced by read_mcnpfile().
int McnpData::get_coord_system | ( | ) |
Definition at line 53 of file mcnpmit.cpp.
54 {
55 return coord_system;
56 }
References coord_system.
std::string McnpData::get_filename | ( | ) |
Definition at line 78 of file mcnpmit.cpp.
79 {
80 return MCNP_filename;
81 }
References MCNP_filename.
double * McnpData::get_rotation_matrix | ( | ) |
Definition at line 67 of file mcnpmit.cpp.
68 {
69 return rotation_matrix;
70 }
References rotation_matrix.
MCNPError McnpData::initialize_tags | ( | ) |
Definition at line 326 of file mcnpmit.cpp.
327 {
328
329 MBI->tag_get_handle( TALLY_TAG, 1, moab::MB_TYPE_DOUBLE, tally_tag, moab::MB_TAG_DENSE | moab::MB_TAG_CREAT );
330 MBI->tag_get_handle( ERROR_TAG, 1, moab::MB_TYPE_DOUBLE, relerr_tag, moab::MB_TAG_DENSE | moab::MB_TAG_CREAT );
331
332 return MCNP_SUCCESS;
333 }
References ERROR_TAG, MB_TAG_CREAT, MB_TAG_DENSE, MB_TYPE_DOUBLE, MBI, MCNP_SUCCESS, relerr_tag, TALLY_TAG, and tally_tag.
MCNPError McnpData::make_adjacencies | ( | int * | ) |
MCNPError McnpData::make_elements | ( | std::vector< double > | x[3], |
int * | n | ||
) |
Definition at line 273 of file mcnpmit.cpp.
274 {
275
276 // double v[3];
277 // MBEntityHandle dumhandle;
278 // MBEntityHandle vstart, vijk;
279 unsigned int num_verts = n[0] * n[1] * n[2];
280 double* coords;
281 coords = new double[3 * num_verts];
282
283 /*
284 // Enter the vertices ...
285 for (int k=0; k < n[2]; k++) {
286 v[2] = x[2].at(k);
287 for (int j=0; j < n[1]; j++) {
288 v[1] = x[1].at(j);
289 for (int i=0; i < n[0]; i++) {
290 v[0] = x[0].at(i);
291 MBresult = MBI->create_vertex(v, dumhandle);
292 if (MBresult != MB_SUCCESS) return MCNP_FAILURE;
293 MCNP_vertices.push_back(dumhandle);
294
295 }
296 }
297 }
298 */
299
300 // Enter the vertices ...
301 for( int k = 0; k < n[2]; k++ )
302 {
303 for( int j = 0; j < n[1]; j++ )
304 {
305 for( int i = 0; i < n[0]; i++ )
306 {
307 unsigned int ijk = 3 * ( k * n[0] * n[1] + j * n[0] + i );
308 coords[ijk] = x[0][i];
309 coords[ijk + 1] = x[1][j];
310 coords[ijk + 2] = x[2][k];
311
312 // std::cout << coords[ijk] << " " << coords[ijk+1] << " "
313 // << coords[ijk+2] << std::endl;
314
315 // MCNP_vertices.push_back(dumhandle);
316 }
317 }
318 }
319
320 MBI->create_vertices( coords, num_verts, vert_handles );
321
322 delete[] coords;
323 return MCNP_SUCCESS;
324 }
References MBI, MCNP_SUCCESS, and vert_handles.
Referenced by read_mcnpfile().
MCNPError McnpData::read_coord_system | ( | std::string | s | ) |
Definition at line 243 of file mcnpmit.cpp.
244 {
245
246 if( ( s.find( "Box" ) < 100 ) || ( s.find( "xyz" ) < 100 ) )
247 coord_system = CARTESIAN;
248 else if( s.find( "Cyl" ) < 100 )
249 coord_system = CYLINDRICAL;
250 else if( s.find( "Sph" ) < 100 )
251 coord_system = SPHERICAL;
252 else
253 return MCNP_FAILURE;
254
255 return MCNP_SUCCESS;
256 }
References CARTESIAN, coord_system, CYLINDRICAL, MCNP_FAILURE, MCNP_SUCCESS, and SPHERICAL.
Referenced by read_mcnpfile().
MCNPError McnpData::read_mcnpfile | ( | bool | skip_mesh | ) |
Definition at line 84 of file mcnpmit.cpp.
85 {
86
87 MCNPError result;
88 moab::ErrorCode MBresult;
89 moab::CartVect tvect;
90
91 std::vector< double > xvec[3];
92
93 // Open the file
94 std::ifstream mcnpfile;
95 mcnpfile.open( MCNP_filename.c_str() );
96 if( !mcnpfile )
97 {
98 std::cout << "Unable to open MCNP data file." << std::endl;
99 return MCNP_FAILURE;
100 }
101 std::cout << std::endl;
102 std::cout << "Reading MCNP input file..." << std::endl;
103
104 // Prepare for file reading ...
105 char line[10000];
106 int mode = 0; // Set the file reading mode to read proper data
107 int nv[3];
108
109 // Read in the file ...
110 while( !mcnpfile.eof() )
111 {
112
113 mcnpfile.getline( line, 10000 );
114 // std::cout << line << std::endl;
115
116 switch( mode )
117 {
118 case 0: // First line is a title
119 mode++;
120 break;
121 case 1: // Coordinate system
122 mode++;
123 result = read_coord_system( line );
124 if( result == MCNP_FAILURE ) return MCNP_FAILURE;
125 break;
126 case 2: // Rotation matrix
127 mode++;
128 for( int i = 0; i < 4; i++ )
129 {
130 mcnpfile.getline( line, 10000 );
131 result = read_rotation_matrix( line, i );
132 if( result == MCNP_FAILURE ) return MCNP_FAILURE;
133 }
134 if( skip_mesh ) return MCNP_SUCCESS;
135 break;
136 case 3: // Read in vertices and build elements
137 mode++;
138
139 for( int i = 0; i < 3; i++ )
140 {
141 // How many points in the x[i]-direction
142 nv[i] = how_many_numbers( line );
143 if( nv[i] <= 0 ) return MCNP_FAILURE;
144
145 // Get space and read in these points
146 result = read_numbers( line, nv[i], xvec[i] );
147 if( result == MCNP_FAILURE ) return MCNP_FAILURE;
148
149 // Update to the next line
150 mcnpfile.getline( line, 10000 );
151 }
152
153 // Make the elements and vertices
154 result = make_elements( xvec, nv );
155 if( result == MCNP_FAILURE ) return MCNP_FAILURE;
156 break;
157 case 4: // Read in tally data, make, and tag elements
158 mode++;
159 moab::EntityHandle elemhandle;
160
161 moab::EntityHandle vstart, vijk;
162 moab::EntityHandle connect[8];
163 // double d[3];
164
165 // vstart = MCNP_vertices.front();
166 vstart = *( vert_handles.begin() );
167
168 for( int i = 0; i < nv[0] - 1; i++ )
169 {
170 for( int j = 0; j < nv[1] - 1; j++ )
171 {
172 for( int k = 0; k < nv[2] - 1; k++ )
173 {
174 vijk = vstart + ( i + j * nv[0] + k * nv[0] * nv[1] );
175
176 // std::cout << vijk << std::endl;
177
178 connect[0] = vijk;
179 connect[1] = vijk + 1;
180 connect[2] = vijk + 1 + nv[0];
181 connect[3] = vijk + nv[0];
182 connect[4] = vijk + nv[0] * nv[1];
183 connect[5] = vijk + 1 + nv[0] * nv[1];
184 connect[6] = vijk + 1 + nv[0] + nv[0] * nv[1];
185 connect[7] = vijk + nv[0] + nv[0] * nv[1];
186
187 MBresult = MBI->create_element( moab::MBHEX, connect, 8, elemhandle );
188 if( MBresult != moab::MB_SUCCESS ) return MCNP_FAILURE;
189 elem_handles.insert( elemhandle );
190
191 mcnpfile.getline( line, 10000 );
192 result = extract_tally_data( line, elemhandle );
193 if( result == MCNP_FAILURE ) return MCNP_FAILURE;
194 }
195 }
196 }
197
198 /*
199 for (MBRange::iterator rit=vert_handles.begin(); rit !=
200 vert_handles.end(); ++rit) { std::cout << *rit << std::endl;
201 }
202
203
204 for (int i=0; i < nv[0]-1; i++) {
205 for (int j=0; j < nv[1]-1; j++) {
206 for (int k=0; k < nv[2]-1; k++) {
207 vijk = vstart + (i + j*nv[0] + k*nv[0]*nv[1]);
208 connect[0] = vijk;
209 connect[1] = vijk + 1;
210 connect[2] = vijk + 1 + nv[0];
211 connect[3] = vijk + nv[0];
212 connect[4] = vijk + nv[0]*nv[1];
213 connect[5] = vijk + 1 + nv[0]*nv[1];
214 connect[6] = vijk + 1 + nv[0] + nv[0]*nv[1];
215 connect[7] = vijk + nv[0] + nv[0]*nv[1];
216
217 MBresult = MBI->create_element(MBHEX, connect, 8,
218 elemhandle); if (MBresult != MB_SUCCESS) return MCNP_FAILURE;
219 elem_handles.insert(elemhandle);
220
221 mcnpfile.getline(line, 10000);
222 result = extract_tally_data(line, elemhandle);
223 if (result == MCNP_FAILURE) return MCNP_FAILURE;
224
225 }
226 }
227 }
228 */
229 break;
230 case 5: // Ckeck for weirdness at end of file
231 if( !mcnpfile.eof() ) return MCNP_FAILURE;
232 break;
233 }
234 }
235
236 std::cout << "SUCCESS! Read in " << elem_handles.size() << " elements!" << std::endl << std::endl;
237 // MCNP_vertices.clear();
238 vert_handles.clear();
239 MCNP_elems.clear();
240 return MCNP_SUCCESS;
241 }
References moab::Range::begin(), moab::Range::clear(), elem_handles, ErrorCode, extract_tally_data(), how_many_numbers(), moab::Range::insert(), make_elements(), MB_SUCCESS, MBHEX, MBI, MCNP_elems, MCNP_FAILURE, MCNP_filename, MCNP_SUCCESS, read_coord_system(), read_numbers(), read_rotation_matrix(), moab::Range::size(), and vert_handles.
MCNPError McnpData::read_rotation_matrix | ( | std::string | s, |
int | i | ||
) |
Definition at line 258 of file mcnpmit.cpp.
259 {
260
261 int fpos = 0;
262 MCNPError result;
263
264 for( int j = 0; j < 4; j++ )
265 {
266 result = next_number( s, rotation_matrix[4 * i + j], fpos );
267 if( result == MCNP_FAILURE ) return MCNP_FAILURE;
268 }
269
270 return MCNP_SUCCESS;
271 }
References MCNP_FAILURE, MCNP_SUCCESS, next_number(), and rotation_matrix.
Referenced by read_mcnpfile().
MCNPError McnpData::set_coord_system | ( | int | k | ) |
Definition at line 48 of file mcnpmit.cpp.
49 {
50 coord_system = k;
51 return MCNP_SUCCESS;
52 }
References coord_system, and MCNP_SUCCESS.
MCNPError McnpData::set_filename | ( | std::string | fname | ) |
Definition at line 73 of file mcnpmit.cpp.
74 {
75 MCNP_filename = fname;
76 return MCNP_SUCCESS;
77 }
References MCNP_filename, and MCNP_SUCCESS.
MCNPError McnpData::set_rotation_matrix | ( | double | r[16] | ) |
Definition at line 59 of file mcnpmit.cpp.
60 {
61 for( int i = 0; i < 16; i++ )
62 {
63 rotation_matrix[i] = r[i];
64 }
65 return MCNP_SUCCESS;
66 }
References MCNP_SUCCESS, and rotation_matrix.
MCNPError McnpData::transform_point | ( | double * | p, |
double * | r, | ||
int | csys, | ||
double * | rmat | ||
) |
Definition at line 425 of file mcnpmit.cpp.
426 {
427
428 double q[3];
429
430 // Apply the rotation matrix
431 for( unsigned int i = 0; i < 3; i++ )
432 {
433 q[i] = p[0] * rmat[4 * i] + p[1] * rmat[4 * i + 1] + p[2] * rmat[4 * i + 2] + rmat[4 * i + 3];
434 }
435
436 // Transform coordinate system
437 switch( csys )
438 {
439 case CARTESIAN:
440 r[0] = q[0];
441 r[1] = q[1];
442 r[2] = q[2]; // x, y, z
443 break;
444 case CYLINDRICAL:
445 r[0] = sqrt( q[0] * q[0] + q[1] * q[1] ); // r
446 r[1] = q[2]; // z
447 r[2] = c2pi * ( atan2( q[1], q[0] ) ); // theta (in rotations)
448 break;
449 case SPHERICAL:
450 return MCNP_FAILURE;
451 // break;
452 default:
453 return MCNP_FAILURE;
454 // break;
455 }
456
457 return MCNP_SUCCESS;
458 }
References c2pi, CARTESIAN, CYLINDRICAL, MCNP_FAILURE, MCNP_SUCCESS, and SPHERICAL.
moab::Tag McnpData::box_max_tag |
Definition at line 45 of file mcnpmit.hpp.
moab::Tag McnpData::box_min_tag |
Definition at line 45 of file mcnpmit.hpp.
int McnpData::coord_system |
Definition at line 35 of file mcnpmit.hpp.
Referenced by get_coord_system(), McnpData(), read_coord_system(), and set_coord_system().
moab::Range McnpData::elem_handles |
Definition at line 42 of file mcnpmit.hpp.
Referenced by read_mcnpfile().
std::vector< moab::EntityHandle > McnpData::MCNP_elems |
Definition at line 40 of file mcnpmit.hpp.
Referenced by read_mcnpfile().
std::string McnpData::MCNP_filename |
Definition at line 50 of file mcnpmit.hpp.
Referenced by get_filename(), read_mcnpfile(), and set_filename().
std::vector< moab::EntityHandle > McnpData::MCNP_vertices |
Definition at line 39 of file mcnpmit.hpp.
Referenced by ~McnpData().
moab::Tag McnpData::relerr_tag |
Definition at line 47 of file mcnpmit.hpp.
Referenced by extract_tally_data(), and initialize_tags().
double McnpData::rotation_matrix[16] |
Definition at line 36 of file mcnpmit.hpp.
Referenced by get_rotation_matrix(), McnpData(), read_rotation_matrix(), and set_rotation_matrix().
moab::Tag McnpData::tally_tag |
Definition at line 46 of file mcnpmit.hpp.
Referenced by extract_tally_data(), and initialize_tags().
moab::Range McnpData::vert_handles |
Definition at line 41 of file mcnpmit.hpp.
Referenced by make_elements(), and read_mcnpfile().