|
moab
|
#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 15 of file mcnpmit.hpp.
Definition at line 20 of file mcnpmit.cpp.
{
// Default value for coordinate system
coord_system = 0;
// Default rotation matrix is identity matrix
for (int i = 0; i < 4; i++) {
for (int j = 0; j < 4; j++) {
if (i == j)
rotation_matrix[4*i + j] = 1;
else
rotation_matrix[4*i + j] = 0;
}
}
}
Definition at line 38 of file mcnpmit.cpp.
{
// Vertices and elements
MCNP_vertices.clear();
}
| MCNPError McnpData::extract_tally_data | ( | std::string | s, |
| moab::EntityHandle | handle | ||
| ) |
Definition at line 321 of file mcnpmit.cpp.
{
int fpos = 0;
double d = 0;
MCNPError result;
moab::ErrorCode MBresult;
// Discard first three lines
for (int i = 0; i < 3; i++) {
result = next_number(s, d, fpos);
if (result == MCNP_FAILURE) return MCNP_FAILURE;
}
// Need to read in tally entry and tag ...
result = next_number(s, d, fpos);
if (result == MCNP_FAILURE) return MCNP_FAILURE;
MBresult = MBI -> tag_set_data(tally_tag, &handle, 1, &d);
if (MBresult != moab::MB_SUCCESS) return MCNP_FAILURE;
// Need to read in relative error entry and tag ...
result = next_number(s, d, fpos);
if (result == MCNP_FAILURE) return MCNP_FAILURE;
MBresult = MBI -> tag_set_data(relerr_tag, &handle, 1, &d);
if (MBresult != moab::MB_SUCCESS) return MCNP_FAILURE;
return MCNP_SUCCESS;
}
| int McnpData::get_coord_system | ( | ) |
Definition at line 51 of file mcnpmit.cpp.
{
return coord_system;
}
| std::string McnpData::get_filename | ( | ) |
Definition at line 72 of file mcnpmit.cpp.
{
return MCNP_filename;
}
| double * McnpData::get_rotation_matrix | ( | ) |
Definition at line 62 of file mcnpmit.cpp.
{
return rotation_matrix;
}
Definition at line 312 of file mcnpmit.cpp.
{
MBI->tag_get_handle(TALLY_TAG, 1, moab::MB_TYPE_DOUBLE, tally_tag, moab::MB_TAG_DENSE|moab::MB_TAG_CREAT);
MBI->tag_get_handle(ERROR_TAG, 1, moab::MB_TYPE_DOUBLE, relerr_tag, moab::MB_TAG_DENSE|moab::MB_TAG_CREAT);
return MCNP_SUCCESS;
}
| MCNPError McnpData::make_adjacencies | ( | int * | ) |
| MCNPError McnpData::make_elements | ( | std::vector< double > | x[3], |
| int * | n | ||
| ) |
Definition at line 262 of file mcnpmit.cpp.
{
// double v[3];
// MBEntityHandle dumhandle;
// MBEntityHandle vstart, vijk;
unsigned int num_verts = n[0]*n[1]*n[2];
double *coords;
coords = new double [ 3 * num_verts ];
/*
// Enter the vertices ...
for (int k=0; k < n[2]; k++) {
v[2] = x[2].at(k);
for (int j=0; j < n[1]; j++) {
v[1] = x[1].at(j);
for (int i=0; i < n[0]; i++) {
v[0] = x[0].at(i);
MBresult = MBI->create_vertex(v, dumhandle);
if (MBresult != MB_SUCCESS) return MCNP_FAILURE;
MCNP_vertices.push_back(dumhandle);
}
}
}
*/
// Enter the vertices ...
for (int k=0; k < n[2]; k++) {
for (int j=0; j < n[1]; j++) {
for (int i=0; i < n[0]; i++) {
unsigned int ijk = 3*(k*n[0]*n[1] + j*n[0] + i);
coords[ ijk ] = x[0][i];
coords[ ijk+1 ] = x[1][j];
coords[ ijk+2 ] = x[2][k];
// std::cout << coords[ijk] << " " << coords[ijk+1] << " "
// << coords[ijk+2] << std::endl;
// MCNP_vertices.push_back(dumhandle);
}
}
}
MBI->create_vertices(coords, num_verts, vert_handles);
delete coords;
return MCNP_SUCCESS;
}
| MCNPError McnpData::read_coord_system | ( | std::string | s | ) |
Definition at line 234 of file mcnpmit.cpp.
{
if ((s.find("Box") < 100) || (s.find("xyz") < 100))
coord_system = CARTESIAN;
else if (s.find("Cyl") < 100)
coord_system = CYLINDRICAL;
else if (s.find("Sph") < 100)
coord_system = SPHERICAL;
else
return MCNP_FAILURE;
return MCNP_SUCCESS;
}
| MCNPError McnpData::read_mcnpfile | ( | bool | skip_mesh | ) |
Definition at line 78 of file mcnpmit.cpp.
{
MCNPError result;
moab::ErrorCode MBresult;
moab::CartVect tvect;
std::vector<double> xvec[3];
// Open the file
std::ifstream mcnpfile;
mcnpfile.open( MCNP_filename.c_str() );
if (!mcnpfile) {
std::cout << "Unable to open MCNP data file." << std::endl;
return MCNP_FAILURE;
}
std::cout << std::endl;
std::cout << "Reading MCNP input file..." << std::endl;
// Prepare for file reading ...
char line[10000];
int mode = 0; // Set the file reading mode to read proper data
int nv[3];
// Read in the file ...
while (! mcnpfile.eof() ) {
mcnpfile.getline(line, 10000);
// std::cout << line << std::endl;
switch(mode) {
case 0: // First line is a title
mode++;
break;
case 1: // Coordinate system
mode++;
result = read_coord_system(line);
if (result == MCNP_FAILURE)
return MCNP_FAILURE;
break;
case 2: // Rotation matrix
mode++;
for (int i = 0; i < 4; i++) {
mcnpfile.getline(line, 10000);
result = read_rotation_matrix(line, i);
if (result == MCNP_FAILURE)
return MCNP_FAILURE;
}
if (skip_mesh) return MCNP_SUCCESS;
break;
case 3: // Read in vertices and build elements
mode++;
for (int i = 0; i < 3; i++) {
// How many points in the x[i]-direction
nv[i] = how_many_numbers(line);
if (nv[i] <= 0) return MCNP_FAILURE;
// Get space and read in these points
result = read_numbers(line , nv[i], xvec[i]);
if (result == MCNP_FAILURE) return MCNP_FAILURE;
// Update to the next line
mcnpfile.getline(line, 10000);
}
// Make the elements and vertices
result = make_elements(xvec, nv);
if (result == MCNP_FAILURE) return MCNP_FAILURE;
break;
case 4: // Read in tally data, make, and tag elements
mode++;
moab::EntityHandle elemhandle;
moab::EntityHandle vstart, vijk;
moab::EntityHandle connect[8];
// double d[3];
// vstart = MCNP_vertices.front();
vstart = *(vert_handles.begin());
for (int i=0; i < nv[0]-1; i++) {
for (int j=0; j < nv[1]-1; j++) {
for (int k=0; k < nv[2]-1; k++) {
vijk = vstart + (i + j*nv[0] + k*nv[0]*nv[1]);
//std::cout << vijk << std::endl;
connect[0] = vijk;
connect[1] = vijk + 1;
connect[2] = vijk + 1 + nv[0];
connect[3] = vijk + nv[0];
connect[4] = vijk + nv[0]*nv[1];
connect[5] = vijk + 1 + nv[0]*nv[1];
connect[6] = vijk + 1 + nv[0] + nv[0]*nv[1];
connect[7] = vijk + nv[0] + nv[0]*nv[1];
MBresult = MBI->create_element(moab::MBHEX, connect, 8, elemhandle);
if (MBresult != moab::MB_SUCCESS) return MCNP_FAILURE;
elem_handles.insert(elemhandle);
mcnpfile.getline(line, 10000);
result = extract_tally_data(line, elemhandle);
if (result == MCNP_FAILURE) return MCNP_FAILURE;
}
}
}
/*
for (MBRange::iterator rit=vert_handles.begin(); rit != vert_handles.end(); rit++) {
std::cout << *rit << std::endl;
}
for (int i=0; i < nv[0]-1; i++) {
for (int j=0; j < nv[1]-1; j++) {
for (int k=0; k < nv[2]-1; k++) {
vijk = vstart + (i + j*nv[0] + k*nv[0]*nv[1]);
connect[0] = vijk;
connect[1] = vijk + 1;
connect[2] = vijk + 1 + nv[0];
connect[3] = vijk + nv[0];
connect[4] = vijk + nv[0]*nv[1];
connect[5] = vijk + 1 + nv[0]*nv[1];
connect[6] = vijk + 1 + nv[0] + nv[0]*nv[1];
connect[7] = vijk + nv[0] + nv[0]*nv[1];
MBresult = MBI->create_element(MBHEX, connect, 8, elemhandle);
if (MBresult != MB_SUCCESS) return MCNP_FAILURE;
elem_handles.insert(elemhandle);
mcnpfile.getline(line, 10000);
result = extract_tally_data(line, elemhandle);
if (result == MCNP_FAILURE) return MCNP_FAILURE;
}
}
}
*/
break;
case 5: // Ckeck for weirdness at end of file
if (! mcnpfile.eof() ) return MCNP_FAILURE;
break;
}
}
std::cout << "SUCCESS! Read in " << elem_handles.size()
<< " elements!" << std::endl << std::endl;
// MCNP_vertices.clear();
vert_handles.clear();
MCNP_elems.clear();
return MCNP_SUCCESS;
}
| MCNPError McnpData::read_rotation_matrix | ( | std::string | s, |
| int | i | ||
| ) |
Definition at line 248 of file mcnpmit.cpp.
{
int fpos = 0;
MCNPError result;
for (int j = 0; j < 4; j++) {
result = next_number(s, rotation_matrix[4*i+j], fpos);
if (result == MCNP_FAILURE)
return MCNP_FAILURE;
}
return MCNP_SUCCESS;
}
| MCNPError McnpData::set_coord_system | ( | int | k | ) |
Definition at line 47 of file mcnpmit.cpp.
{
coord_system = k;
return MCNP_SUCCESS;
}
| MCNPError McnpData::set_filename | ( | std::string | fname | ) |
Definition at line 68 of file mcnpmit.cpp.
{
MCNP_filename = fname;
return MCNP_SUCCESS;
}
| MCNPError McnpData::set_rotation_matrix | ( | double | r[16] | ) |
Definition at line 56 of file mcnpmit.cpp.
{
for (int i = 0; i < 16; i++) {
rotation_matrix[i] = r[i];
}
return MCNP_SUCCESS;
}
| MCNPError McnpData::transform_point | ( | double * | p, |
| double * | r, | ||
| int | csys, | ||
| double * | rmat | ||
| ) |
Definition at line 406 of file mcnpmit.cpp.
{
double q[3];
// Apply the rotation matrix
for (unsigned int i=0; i < 3; i++) {
q[i] = p[0] * rmat[4*i ] + p[1] * rmat[4*i+1]
+ p[2] * rmat[4*i+2] + rmat[4*i+3];
}
// Transform coordinate system
switch( csys ) {
case CARTESIAN :
r[0] = q[0]; r[1] = q[1]; r[2] = q[2]; // x, y, z
break;
case CYLINDRICAL :
r[0] = sqrt( q[0]*q[0] + q[1]*q[1] ); // r
r[1] = q[2]; // z
r[2] = c2pi * ( atan2( q[1], q[0] ) ); // theta (in rotations)
break;
case SPHERICAL :
return MCNP_FAILURE;
break;
default :
return MCNP_FAILURE;
break;
}
return MCNP_SUCCESS;
}
Definition at line 33 of file mcnpmit.hpp.
Definition at line 33 of file mcnpmit.hpp.
Definition at line 23 of file mcnpmit.hpp.
Definition at line 30 of file mcnpmit.hpp.
| std::vector<moab::EntityHandle> McnpData::MCNP_elems |
Definition at line 28 of file mcnpmit.hpp.
| std::string McnpData::MCNP_filename |
Definition at line 38 of file mcnpmit.hpp.
| std::vector<moab::EntityHandle> McnpData::MCNP_vertices |
Definition at line 27 of file mcnpmit.hpp.
Definition at line 35 of file mcnpmit.hpp.
| double McnpData::rotation_matrix[16] |
Definition at line 24 of file mcnpmit.hpp.
Definition at line 34 of file mcnpmit.hpp.
Definition at line 29 of file mcnpmit.hpp.
1.7.6.1