1 /**
2  * MOAB, a Mesh-Oriented datABase, is a software component for creating,
3  * storing and accessing finite element mesh data.
4  *
5  * Copyright 2004 Sandia Corporation.  Under the terms of Contract
6  * DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government
7  * retains certain rights in this software.
8  *
9  * This library is free software; you can redistribute it and/or
10  * modify it under the terms of the GNU Lesser General Public
11  * License as published by the Free Software Foundation; either
12  * version 2.1 of the License, or (at your option) any later version.
13  *
14  */
15 
16 //----------------------------------------------------------------------
17 // Filename : ReadMCNP5.hpp
18 // Purpose  : Read a meshtal file created by MCNP5 into MOAB
19 // Creator  : Brandon Smith
20 // Date     : 07/2009
21 //----------------------------------------------------------------------
22 
23 /**
24  * Data structure of MCNP5 data created by this reader:
25  *
26  * each file_meshset contains
27  *   DATA_AND_TIME_TAG
28  *   TITLE_TAG
29  *   NPS_TAG
30  *   each tally_meshset contains
31  *     TALLY_NUMBER_TAG
32  *     TALLY_COMMENT_TAG
33  *     TALLY_PARTICLE_TAG
34  *     TALLY_COORD_SYS_TAG
35  *     each mesh element contains
36  *       TALLY_TAG
37  *       ERROR_TAG
38  */
39 
40 #include "moab/Interface.hpp"
41 #include "moab/ReaderIface.hpp"
42 #include <iostream>
43 #include <fstream>
44 #include <sstream>
45 #include <vector>
46 
47 namespace moab {
48 
49 class ReadUtilIface;
50 
51 class ReadMCNP5 : public ReaderIface
52 {
53 
54 public:
55   // factory method
56   static ReaderIface* factory( Interface* );
57 
58   ErrorCode load_file( const char* file_name,
59                        const EntityHandle* file_set,
60                        const FileOptions& opts,
61                        const SubsetList* subset_list = 0,
62                        const Tag* file_id_tag = 0 );
63 
64   ErrorCode read_tag_values( const char* file_name,
65                              const char* tag_name,
66                              const FileOptions& opts,
67                              std::vector<int>& tag_values_out,
68                              const SubsetList* subset_list = 0 );
69 
70   // constructor
71   ReadMCNP5(Interface* impl = NULL);
72 
73   // destructor
74   virtual ~ReadMCNP5();
75 
76 protected:
77 
78 private:
79   // constants
80   static const double PI;
81   static const double C2PI;
82   static const double CPI;
83 
84   enum coordinate_system { NO_SYSTEM,
85                            CARTESIAN,
86                            CYLINDRICAL,
87                            SPHERICAL };
88   enum particle { NEUTRON,
89                   PHOTON,
90                   ELECTRON };
91 
92   // read mesh interface
93   ReadUtilIface* readMeshIface;
94 
95   // MOAB Interface
96   Interface* MBI;
97 
98   const Tag* fileIDTag;
99   int nodeId, elemId;
100 
101   // reads the meshtal file
102   ErrorCode load_one_file( const char           *fname,
103                              const EntityHandle *input_meshset,
104                              const FileOptions    &options,
105                              const bool           average );
106 
107   ErrorCode create_tags( Tag &date_and_time_tag,
108                            Tag &title_tag,
109                            Tag &nps_tag,
110                            Tag &tally_number_tag,
111                            Tag &tally_comment_tag,
112                            Tag &tally_particle_tag,
113                            Tag &tally_coord_sys_tag,
114                            Tag &tally_tag,
115                            Tag &error_tag );
116 
117   ErrorCode read_file_header( std::fstream      &file,
118                                 bool              debug,
119                                 char              date_and_time[100],
120                                 char              title[100],
121                                 unsigned long int &nps );
122 
123   ErrorCode set_header_tags( EntityHandle             output_meshset,
124                                         char              date_and_time[100],
125                                         char              title[100],
126                                         unsigned long int nps,
127                                         Tag             data_and_time_tag,
128                                         Tag             title_tag,
129                                         Tag             nps_tag );
130 
131   ErrorCode read_tally_header( std::fstream &file,
132                                  bool         debug,
133                                  unsigned int &tally_number,
134                                  char         tally_comment[100],
135                                  particle     &tally_particle );
136 
137   ErrorCode get_tally_particle( std::string a,
138                                   bool        debug,
139                                   particle    &tally_particle );
140 
141   ErrorCode read_mesh_planes( std::fstream         &file,
142                                 bool                 debug,
143                                 std::vector<double>  planes[3],
144                                 coordinate_system    &coord_sys);
145 
146   ErrorCode get_mesh_plane( std::istringstream  &ss,
147                               bool                debug,
148                               std::vector<double> &plane);
149 
150   ErrorCode read_element_values_and_errors( std::fstream        &file,
151                                               bool                debug,
152                                               std::vector<double> planes[3],
153                                               unsigned int        n_chopped_x0_planes,
154                                               unsigned int        n_chopped_x2_planes,
155                                               particle            tally_particle,
156                                               double              values[],
157                                               double              errors[] );
158 
159   ErrorCode set_tally_tags( EntityHandle    tally_meshset,
160                               unsigned int      tally_number,
161                               char              tally_comment[100],
162                               particle          tally_particle,
163                               coordinate_system tally_coord_sys,
164                               Tag             tally_number_tag,
165                               Tag             tally_comment_tag,
166                               Tag             tally_particle_tag,
167                               Tag             tally_coord_sys_tag );
168 
169   ErrorCode create_vertices( std::vector<double> planes[3],
170                                bool                debug,
171                                EntityHandle      &start_vert,
172                                coordinate_system   coord_sys,
173                                EntityHandle      tally_meshset );
174 
175   ErrorCode create_elements( bool                debug,
176                                std::vector<double> planes[3],
177                                unsigned int        n_chopped_x0_planes,
178                                unsigned int        n_chopped_x2_planes,
179                                EntityHandle      start_vert,
180                                double              values[],
181                                double              errors[],
182                                Tag               tally_tag,
183                                Tag               error_tag,
184                                EntityHandle      tally_meshset,
185                                coordinate_system   tally_coord_sys );
186 
187   ErrorCode average_with_existing_tally( bool              debug,
188                                            unsigned long int &new_nps,
189                                            unsigned long int nps,
190                                            unsigned int      tally_number,
191                                            Tag             tally_number_tag,
192                                            Tag             nps_tag,
193                                            Tag             tally_tag,
194                                            Tag             error_tag,
195                                            double            values[],
196                                            double            errors[],
197                                            unsigned int      n_elements );
198 
199   ErrorCode transform_point_to_cartesian(double *in,
200                                            double *out,
201                                            coordinate_system coord_sys);
202 
203   ErrorCode average_tally_values(const unsigned long int nps0,
204                                    const unsigned long int nps1,
205                                    double                  *values0,
206                                    const double            *values1,
207                                    double                  *errors0,
208                                    const double            *errors1,
209                                    const unsigned long int n_values);
210 };
211 
212 } // namespace moab
213