www.mooseframework.org
Public Member Functions | Static Public Member Functions | Protected Member Functions | Protected Attributes | List of all members
EBSDReader Class Reference

A GeneralUserObject that reads an EBSD file and stores the centroid data in a data structure which indexes on element centroids. More...

#include <EBSDReader.h>

Inheritance diagram for EBSDReader:
[legend]

Public Member Functions

 EBSDReader (const InputParameters &params)
 
virtual ~EBSDReader ()
 
virtual void readFile ()
 
virtual void initialize ()
 
virtual void execute ()
 
virtual void finalize ()
 
const EBSDPointDatagetData (const Point &p) const
 Get the requested type of data at the point p. More...
 
const EBSDAvgData & getAvgData (unsigned int i) const
 Get the requested type of average data for (global) grain number i. More...
 
const EBSDAvgData & getAvgData (unsigned int phase, unsigned int local_id) const
 Get the requested type of average data for a given phase and (local) grain. More...
 
virtual const EulerAnglesgetEulerAngles (unsigned int) const
 EulerAngleProvider interface implementation to fetch a triplet of Euler angles. More...
 
virtual unsigned int getGrainNum () const
 Return the total number of grains. More...
 
virtual unsigned int getPhaseNum () const
 Return the total number of phases. More...
 
unsigned int getGrainNum (unsigned int phase) const
 Return the number of grains in a given phase. More...
 
unsigned int getFeatureID (unsigned int phase, unsigned int local_id) const
 Return the EBSD feature id for a given phase and phase (local) grain number. More...
 
unsigned int getFeatureID (unsigned int global_id) const
 Return the EBSD feature id for a given (global) grain number. More...
 
unsigned int getGlobalID (unsigned int phase, unsigned int local_id) const
 Return the (global) grain id for a given phase and (local) grain number. More...
 
unsigned int getGlobalID (unsigned int feature_id) const
 Return the (global) grain id for a given phase and (local) grain number. More...
 
MooseSharedPointer< EBSDPointDataFunctorgetPointDataAccessFunctor (const MooseEnum &field_name) const
 Factory function to return a point functor specified by name. More...
 
MooseSharedPointer< EBSDAvgDataFunctorgetAvgDataAccessFunctor (const MooseEnum &field_name) const
 Factory function to return a average functor specified by name. More...
 
const std::map< dof_id_type, std::vector< Real > > & getNodeToGrainWeightMap () const
 Returns a map consisting of the node index followd by a vector of all grain weights for that node. More...
 
const std::map< dof_id_type, std::vector< Real > > & getNodeToPhaseWeightMap () const
 Returns a map consisting of the node index followd by a vector of all phase weights for that node. More...
 
void meshChanged ()
 Maps need to be updated when the mesh changes. More...
 

Static Public Member Functions

static MooseEnum getPointDataFieldType ()
 
static MooseEnum getAvgDataFieldType ()
 

Protected Member Functions

unsigned indexFromPoint (const Point &p) const
 Computes a global index in the _data array given an input centroid point. More...
 
unsigned indexFromIndex (unsigned int var) const
 Transfer the index into the _avg_data array from given index. More...
 
void buildNodeWeightMaps ()
 Build grain and phase weight maps. More...
 

Protected Attributes

unsigned int _custom_columns
 number of additional custom data columns More...
 
std::vector< EBSDPointData_data
 Logically three-dimensional data indexed by geometric points in a 1D vector. More...
 
std::vector< EBSDAvgData > _avg_data
 Averages by (global) grain ID. More...
 
std::vector< EulerAngles_avg_angles
 Euler Angles by (global) grain ID. More...
 
std::map< unsigned int, unsigned int > _global_id_map
 map from feature_id to global_id More...
 
std::vector< std::vector< unsigned int > > _global_id
 global ID for given phases and grains More...
 
std::map< dof_id_type, std::vector< Real > > _node_to_grain_weight_map
 Map of grain weights per node. More...
 
std::map< dof_id_type, std::vector< Real > > _node_to_phase_weight_map
 Map of phase weights per node. More...
 
const int & _time_step
 current timestep. Maps are only rebuild on mesh change during time step zero More...
 
unsigned int _mesh_dimension
 Dimension of the problem domain. More...
 
unsigned _nx
 The number of values in the x, y and z directions. More...
 
unsigned _ny
 
unsigned _nz
 
Real _dx
 The spacing of the values in x, y and z directions. More...
 
Real _dy
 
Real _dz
 
Real _minx
 Grid origin. More...
 
Real _miny
 
Real _minz
 
Real _maxx
 Maximum grid extent. More...
 
Real _maxy
 
Real _maxz
 
MooseMesh & _mesh
 MooseMesh Variables. More...
 
NonlinearSystemBase & _nl
 
unsigned int _grain_num
 Variables needed to determine reduced order parameter values. More...
 
Point _bottom_left
 
Point _top_right
 
Point _range
 

Detailed Description

A GeneralUserObject that reads an EBSD file and stores the centroid data in a data structure which indexes on element centroids.

Grains are indexed through multiple schemes:

Phases are referred to using the numbers in the EBSD data file. In case the phase number in the data file starts at 1 the phase 0 will simply contain no grains.

Definition at line 34 of file EBSDReader.h.

Constructor & Destructor Documentation

EBSDReader::EBSDReader ( const InputParameters &  params)

Definition at line 28 of file EBSDReader.C.

29  : EulerAngleProvider(params),
30  _mesh(_fe_problem.mesh()),
31  _nl(_fe_problem.getNonlinearSystemBase()),
32  _grain_num(0),
33  _custom_columns(getParam<unsigned int>("custom_columns")),
34  _time_step(_fe_problem.timeStep()),
35  _mesh_dimension(_mesh.dimension()),
36  _nx(0),
37  _ny(0),
38  _nz(0),
39  _dx(0.),
40  _dy(0.),
41  _dz(0.)
42 {
43  readFile();
44 }
unsigned int _custom_columns
number of additional custom data columns
Definition: EBSDReader.h:136
unsigned _nz
Definition: EBSDReader.h:166
unsigned int _grain_num
Variables needed to determine reduced order parameter values.
Definition: EBSDReader.h:129
unsigned _nx
The number of values in the x, y and z directions.
Definition: EBSDReader.h:166
MooseMesh & _mesh
MooseMesh Variables.
Definition: EBSDReader.h:124
const int & _time_step
current timestep. Maps are only rebuild on mesh change during time step zero
Definition: EBSDReader.h:160
unsigned _ny
Definition: EBSDReader.h:166
Real _dx
The spacing of the values in x, y and z directions.
Definition: EBSDReader.h:169
EulerAngleProvider(const InputParameters &parameters)
unsigned int _mesh_dimension
Dimension of the problem domain.
Definition: EBSDReader.h:163
NonlinearSystemBase & _nl
Definition: EBSDReader.h:125
virtual void readFile()
Definition: EBSDReader.C:47
EBSDReader::~EBSDReader ( )
virtual

Definition at line 221 of file EBSDReader.C.

221 {}

Member Function Documentation

void EBSDReader::buildNodeWeightMaps ( )
protected

Build grain and phase weight maps.

Definition at line 335 of file EBSDReader.C.

Referenced by meshChanged(), and readFile().

336 {
337  // Import nodeToElemMap from MooseMesh for current node
338  // This map consists of the node index followed by a vector of element indices that are associated
339  // with that node
340  const std::map<dof_id_type, std::vector<dof_id_type>> & node_to_elem_map =
341  _mesh.nodeToActiveSemilocalElemMap();
342  libMesh::MeshBase & mesh = _mesh.getMesh();
343 
344  // Loop through each node in mesh and calculate eta values for each grain associated with the node
345  MeshBase::const_node_iterator ni = mesh.active_nodes_begin();
346  const MeshBase::const_node_iterator nend = mesh.active_nodes_end();
347  for (; ni != nend; ++ni)
348  {
349  // Get node_id
350  const dof_id_type node_id = (*ni)->id();
351 
352  // Initialize map entries for current node
353  _node_to_grain_weight_map[node_id].assign(getGrainNum(), 0.0);
354  _node_to_phase_weight_map[node_id].assign(getPhaseNum(), 0.0);
355 
356  // Loop through element indices associated with the current node and record weighted eta value
357  // in new map
358  const auto & node_to_elem_pair = node_to_elem_map.find(node_id);
359  if (node_to_elem_pair != node_to_elem_map.end())
360  {
361  unsigned int n_elems =
362  node_to_elem_pair->second
363  .size(); // n_elems can range from 1 to 4 for 2D and 1 to 8 for 3D problems
364 
365  for (unsigned int ne = 0; ne < n_elems; ++ne)
366  {
367  // Current element index
368  unsigned int elem_id = (node_to_elem_pair->second)[ne];
369 
370  // Retrieve EBSD grain number for the current element index
371  const Elem * elem = mesh.elem(elem_id);
372  const EBSDReader::EBSDPointData & d = getData(elem->centroid());
373 
374  // get the (global) grain ID for the EBSD feature ID
375  const unsigned int global_id = getGlobalID(d._feature_id);
376 
377  // Calculate eta value and add to map
378  _node_to_grain_weight_map[node_id][global_id] += 1.0 / n_elems;
379  _node_to_phase_weight_map[node_id][d._phase] += 1.0 / n_elems;
380  }
381  }
382  }
383 }
std::map< dof_id_type, std::vector< Real > > _node_to_phase_weight_map
Map of phase weights per node.
Definition: EBSDReader.h:157
unsigned int getGlobalID(unsigned int phase, unsigned int local_id) const
Return the (global) grain id for a given phase and (local) grain number.
Definition: EBSDReader.h:93
MooseMesh & _mesh
MooseMesh Variables.
Definition: EBSDReader.h:124
unsigned int _feature_id
EBSD feature id, (gklobal) grain number, symmetry, and phase data.
Per element EBSD data point.
const EBSDPointData & getData(const Point &p) const
Get the requested type of data at the point p.
Definition: EBSDReader.C:224
virtual unsigned int getPhaseNum() const
Return the total number of phases.
Definition: EBSDReader.h:74
virtual unsigned int getGrainNum() const
Return the total number of grains.
Definition: EBSDReader.C:248
std::map< dof_id_type, std::vector< Real > > _node_to_grain_weight_map
Map of grain weights per node.
Definition: EBSDReader.h:154
virtual void EBSDReader::execute ( )
inlinevirtual

Definition at line 43 of file EBSDReader.h.

43 {}
virtual void EBSDReader::finalize ( )
inlinevirtual

Definition at line 44 of file EBSDReader.h.

44 {}
const EBSDReader::EBSDAvgData & EBSDReader::getAvgData ( unsigned int  i) const

Get the requested type of average data for (global) grain number i.

Definition at line 230 of file EBSDReader.C.

Referenced by finalize(), PolycrystalEBSD::getGrainsBasedOnPoint(), and EBSDReaderAvgDataAux::precalculateValue().

231 {
232  return _avg_data[indexFromIndex(var)];
233 }
std::vector< EBSDAvgData > _avg_data
Averages by (global) grain ID.
Definition: EBSDReader.h:142
unsigned indexFromIndex(unsigned int var) const
Transfer the index into the _avg_data array from given index.
Definition: EBSDReader.C:289
const EBSDReader::EBSDAvgData & EBSDReader::getAvgData ( unsigned int  phase,
unsigned int  local_id 
) const

Get the requested type of average data for a given phase and (local) grain.

Definition at line 242 of file EBSDReader.C.

243 {
244  return _avg_data[indexFromIndex(_global_id[phase][local_id])];
245 }
std::vector< EBSDAvgData > _avg_data
Averages by (global) grain ID.
Definition: EBSDReader.h:142
std::vector< std::vector< unsigned int > > _global_id
global ID for given phases and grains
Definition: EBSDReader.h:151
unsigned indexFromIndex(unsigned int var) const
Transfer the index into the _avg_data array from given index.
Definition: EBSDReader.C:289
MooseSharedPointer< EBSDAccessFunctors::EBSDAvgDataFunctor > EBSDReader::getAvgDataAccessFunctor ( const MooseEnum &  field_name) const

Factory function to return a average functor specified by name.

Definition at line 432 of file EBSDReader.C.

Referenced by getGlobalID().

433 {
434  EBSDAvgDataFunctor * ret_val = NULL;
435 
436  switch (field_name)
437  {
438  case 0: // phi1
439  ret_val = new EBSDAvgDataPhi1();
440  break;
441  case 1: // phi
442  ret_val = new EBSDAvgDataPhi();
443  break;
444  case 2: // phi2
445  ret_val = new EBSDAvgDataPhi2();
446  break;
447  case 3: // phase
448  ret_val = new EBSDAvgDataPhase();
449  break;
450  case 4: // symmetry
451  ret_val = new EBSDAvgDataSymmetry();
452  break;
453  case 5: // local_id
454  ret_val = new EBSDAvgDataLocalID();
455  break;
456  case 6: // feature_id
457  ret_val = new EBSDAvgDataFeatureID();
458  break;
459  default:
460  {
461  // check for custom columns
462  for (unsigned int i = 0; i < _custom_columns; ++i)
463  if (field_name == "CUSTOM" + Moose::stringify(i))
464  {
465  ret_val = new EBSDAvgDataCustom(i);
466  break;
467  }
468  }
469  }
470 
471  // If ret_val was not set by any of the above cases, throw an error.
472  if (!ret_val)
473  mooseError("Error: Please input supported EBSD_param");
474 
475  // If we made it here, wrap up the the ret_val in a
476  // MooseSharedPointer and ship it off.
477  return MooseSharedPointer<EBSDAvgDataFunctor>(ret_val);
478 }
unsigned int _custom_columns
number of additional custom data columns
Definition: EBSDReader.h:136
MooseEnum EBSDAccessFunctors::getAvgDataFieldType ( )
staticinherited

Definition at line 10 of file EBSDAccessFunctors.C.

Referenced by validParams< EBSDReaderAvgDataAux >().

11 {
12  return MooseEnum("phi1 phi phi2 phase symmetry local_id feature_id", "", true);
13 }
const EBSDReader::EBSDPointData & EBSDReader::getData ( const Point &  p) const

Get the requested type of data at the point p.

Definition at line 224 of file EBSDReader.C.

Referenced by buildNodeWeightMaps(), finalize(), PolycrystalEBSD::getGrainsBasedOnPoint(), and EBSDReaderPointDataAux::precalculateValue().

225 {
226  return _data[indexFromPoint(p)];
227 }
unsigned indexFromPoint(const Point &p) const
Computes a global index in the _data array given an input centroid point.
Definition: EBSDReader.C:260
std::vector< EBSDPointData > _data
Logically three-dimensional data indexed by geometric points in a 1D vector.
Definition: EBSDReader.h:139
const EulerAngles & EBSDReader::getEulerAngles ( unsigned int  var) const
virtual

EulerAngleProvider interface implementation to fetch a triplet of Euler angles.

Implements EulerAngleProvider.

Definition at line 236 of file EBSDReader.C.

Referenced by finalize().

237 {
238  return _avg_angles[indexFromIndex(var)];
239 }
std::vector< EulerAngles > _avg_angles
Euler Angles by (global) grain ID.
Definition: EBSDReader.h:145
unsigned indexFromIndex(unsigned int var) const
Transfer the index into the _avg_data array from given index.
Definition: EBSDReader.C:289
unsigned int EBSDReader::getFeatureID ( unsigned int  phase,
unsigned int  local_id 
) const
inline

Return the EBSD feature id for a given phase and phase (local) grain number.

Definition at line 82 of file EBSDReader.h.

83  {
84  return _avg_data[_global_id[phase][local_id]]._feature_id;
85  }
std::vector< EBSDAvgData > _avg_data
Averages by (global) grain ID.
Definition: EBSDReader.h:142
std::vector< std::vector< unsigned int > > _global_id
global ID for given phases and grains
Definition: EBSDReader.h:151
unsigned int EBSDReader::getFeatureID ( unsigned int  global_id) const
inline

Return the EBSD feature id for a given (global) grain number.

Definition at line 87 of file EBSDReader.h.

88  {
89  return _avg_data[global_id]._feature_id;
90  }
std::vector< EBSDAvgData > _avg_data
Averages by (global) grain ID.
Definition: EBSDReader.h:142
unsigned int EBSDReader::getGlobalID ( unsigned int  phase,
unsigned int  local_id 
) const
inline

Return the (global) grain id for a given phase and (local) grain number.

Definition at line 93 of file EBSDReader.h.

Referenced by buildNodeWeightMaps(), PolycrystalEBSD::getGrainsBasedOnPoint(), and PolycrystalEBSD::getNodalVariableValue().

94  {
95  return _global_id[phase][local_id];
96  }
std::vector< std::vector< unsigned int > > _global_id
global ID for given phases and grains
Definition: EBSDReader.h:151
unsigned int EBSDReader::getGlobalID ( unsigned int  feature_id) const

Return the (global) grain id for a given phase and (local) grain number.

Definition at line 318 of file EBSDReader.C.

319 {
320  auto it = _global_id_map.find(feature_id);
321  if (it == _global_id_map.end())
322  mooseError("Invalid Feature ID");
323  return it->second;
324 }
std::map< unsigned int, unsigned int > _global_id_map
map from feature_id to global_id
Definition: EBSDReader.h:148
unsigned int EBSDReader::getGrainNum ( ) const
virtual

Return the total number of grains.

Implements EulerAngleProvider.

Definition at line 248 of file EBSDReader.C.

Referenced by buildNodeWeightMaps(), finalize(), PolycrystalEBSD::getNumGrains(), and getPhaseNum().

249 {
250  return _grain_num;
251 }
unsigned int _grain_num
Variables needed to determine reduced order parameter values.
Definition: EBSDReader.h:129
unsigned int EBSDReader::getGrainNum ( unsigned int  phase) const

Return the number of grains in a given phase.

Definition at line 254 of file EBSDReader.C.

255 {
256  return _global_id[phase].size();
257 }
std::vector< std::vector< unsigned int > > _global_id
global ID for given phases and grains
Definition: EBSDReader.h:151
const std::map< dof_id_type, std::vector< Real > > & EBSDReader::getNodeToGrainWeightMap ( ) const

Returns a map consisting of the node index followd by a vector of all grain weights for that node.

Needed by ReconVarIC

Definition at line 306 of file EBSDReader.C.

Referenced by getGlobalID().

307 {
309 }
std::map< dof_id_type, std::vector< Real > > _node_to_grain_weight_map
Map of grain weights per node.
Definition: EBSDReader.h:154
const std::map< dof_id_type, std::vector< Real > > & EBSDReader::getNodeToPhaseWeightMap ( ) const

Returns a map consisting of the node index followd by a vector of all phase weights for that node.

Needed by ReconPhaseVarIC

Definition at line 312 of file EBSDReader.C.

Referenced by getGlobalID().

313 {
315 }
std::map< dof_id_type, std::vector< Real > > _node_to_phase_weight_map
Map of phase weights per node.
Definition: EBSDReader.h:157
virtual unsigned int EBSDReader::getPhaseNum ( ) const
inlinevirtual

Return the total number of phases.

Definition at line 74 of file EBSDReader.h.

Referenced by buildNodeWeightMaps().

74 { return _global_id.size(); }
std::vector< std::vector< unsigned int > > _global_id
global ID for given phases and grains
Definition: EBSDReader.h:151
MooseSharedPointer< EBSDAccessFunctors::EBSDPointDataFunctor > EBSDReader::getPointDataAccessFunctor ( const MooseEnum &  field_name) const

Factory function to return a point functor specified by name.

Definition at line 386 of file EBSDReader.C.

Referenced by getGlobalID().

387 {
388  EBSDPointDataFunctor * ret_val = NULL;
389 
390  switch (field_name)
391  {
392  case 0: // phi1
393  ret_val = new EBSDPointDataPhi1();
394  break;
395  case 1: // phi
396  ret_val = new EBSDPointDataPhi();
397  break;
398  case 2: // phi2
399  ret_val = new EBSDPointDataPhi2();
400  break;
401  case 3: // grain
402  ret_val = new EBSDPointDataFeatureID();
403  break;
404  case 4: // phase
405  ret_val = new EBSDPointDataPhase();
406  break;
407  case 5: // symmetry
408  ret_val = new EBSDPointDataSymmetry();
409  break;
410  default:
411  {
412  // check for custom columns
413  for (unsigned int i = 0; i < _custom_columns; ++i)
414  if (field_name == "CUSTOM" + Moose::stringify(i))
415  {
416  ret_val = new EBSDPointDataCustom(i);
417  break;
418  }
419  }
420  }
421 
422  // If ret_val was not set by any of the above cases, throw an error.
423  if (!ret_val)
424  mooseError("Error: Please input supported EBSD_param");
425 
426  // If we made it here, wrap up the the ret_val in a
427  // MooseSharedPointer and ship it off.
428  return MooseSharedPointer<EBSDPointDataFunctor>(ret_val);
429 }
unsigned int _custom_columns
number of additional custom data columns
Definition: EBSDReader.h:136
MooseEnum EBSDAccessFunctors::getPointDataFieldType ( )
staticinherited

Definition at line 4 of file EBSDAccessFunctors.C.

Referenced by validParams< EBSDReaderPointDataAux >().

5 {
6  return MooseEnum("phi1 phi phi2 feature_id phase symmetry", "", true);
7 }
unsigned int EBSDReader::indexFromIndex ( unsigned int  var) const
protected

Transfer the index into the _avg_data array from given index.

Definition at line 289 of file EBSDReader.C.

Referenced by getAvgData(), and getEulerAngles().

290 {
291 
292  // Transfer the index into the _avg_data array.
293  unsigned avg_index = var;
294 
295  // Don't access out of range!
296  if (avg_index >= _avg_data.size())
297  mooseError("Error! Index out of range in EBSDReader::indexFromIndex(), index: ",
298  avg_index,
299  " size: ",
300  _avg_data.size());
301 
302  return avg_index;
303 }
std::vector< EBSDAvgData > _avg_data
Averages by (global) grain ID.
Definition: EBSDReader.h:142
unsigned int EBSDReader::indexFromPoint ( const Point &  p) const
protected

Computes a global index in the _data array given an input centroid point.

Definition at line 260 of file EBSDReader.C.

Referenced by getData(), and readFile().

261 {
262  // Don't assume an ordering on the input data, use the (x, y,
263  // z) values of this centroid to determine the index.
264  unsigned int x_index, y_index, z_index, global_index;
265 
266  x_index = (unsigned int)((p(0) - _minx) / _dx);
267  y_index = (unsigned int)((p(1) - _miny) / _dy);
268 
269  if (_mesh_dimension == 3)
270  {
271  z_index = (unsigned int)((p(2) - _minz) / _dz);
272  global_index = z_index * _ny;
273  }
274  else
275  global_index = 0;
276 
277  // Compute the global index into the _data array. This stores points
278  // in a [z][y][x] ordering.
279  global_index = (global_index + y_index) * _nx + x_index;
280 
281  // Don't access out of range!
282  mooseAssert(global_index < _data.size(),
283  "global_index " << global_index << " points out of _data range: " << _data.size());
284 
285  return global_index;
286 }
Real _minx
Grid origin.
Definition: EBSDReader.h:172
Real _miny
Definition: EBSDReader.h:172
unsigned _nx
The number of values in the x, y and z directions.
Definition: EBSDReader.h:166
Real _minz
Definition: EBSDReader.h:172
unsigned _ny
Definition: EBSDReader.h:166
Real _dx
The spacing of the values in x, y and z directions.
Definition: EBSDReader.h:169
unsigned int _mesh_dimension
Dimension of the problem domain.
Definition: EBSDReader.h:163
std::vector< EBSDPointData > _data
Logically three-dimensional data indexed by geometric points in a 1D vector.
Definition: EBSDReader.h:139
virtual void EBSDReader::initialize ( )
inlinevirtual

Definition at line 42 of file EBSDReader.h.

42 {}
void EBSDReader::meshChanged ( )

Maps need to be updated when the mesh changes.

Definition at line 327 of file EBSDReader.C.

Referenced by getGlobalID().

328 {
329  // maps are only rebuild for use in initial conditions, which happens in time step zero
330  if (_time_step == 0)
332 }
const int & _time_step
current timestep. Maps are only rebuild on mesh change during time step zero
Definition: EBSDReader.h:160
void buildNodeWeightMaps()
Build grain and phase weight maps.
Definition: EBSDReader.C:335
void EBSDReader::readFile ( )
virtual

Definition at line 47 of file EBSDReader.C.

Referenced by EBSDReader().

48 {
49  // No need to re-read data upon recovery
50  if (_app.isRecovering())
51  return;
52 
53  // Fetch and check mesh
54  EBSDMesh * mesh = dynamic_cast<EBSDMesh *>(&_mesh);
55  if (mesh == NULL)
56  mooseError("Please use an EBSDMesh in your simulation.");
57 
58  std::ifstream stream_in(mesh->getEBSDFilename().c_str());
59  if (!stream_in)
60  mooseError("Can't open EBSD file: ", mesh->getEBSDFilename());
61 
62  const EBSDMesh::EBSDMeshGeometry & g = mesh->getEBSDGeometry();
63 
64  // Copy file header data from the EBSDMesh
65  _dx = g.d[0];
66  _nx = g.n[0];
67  _minx = g.min[0];
68  _maxx = _minx + _dx * _nx;
69 
70  _dy = g.d[1];
71  _ny = g.n[1];
72  _miny = g.min[1];
73  _maxy = _miny + _dy * _ny;
74 
75  _dz = g.d[2];
76  _nz = g.n[2];
77  _minz = g.min[2];
78  _maxz = _minz + _dz * _nz;
79 
80  // Resize the _data array
81  unsigned total_size = g.dim < 3 ? _nx * _ny : _nx * _ny * _nz;
82  _data.resize(total_size);
83 
84  std::string line;
85  while (std::getline(stream_in, line))
86  {
87  if (line.find("#") != 0)
88  {
89  // Temporary variables to read in on each line
90  EBSDPointData d;
91  Real x, y, z;
92 
93  std::istringstream iss(line);
94  iss >> d._phi1 >> d._Phi >> d._phi2 >> x >> y >> z >> d._feature_id >> d._phase >>
95  d._symmetry;
96 
97  // Transform angles to degrees
98  d._phi1 *= 180.0 / libMesh::pi;
99  d._Phi *= 180.0 / libMesh::pi;
100  d._phi2 *= 180.0 / libMesh::pi;
101 
102  // Custom columns
103  d._custom.resize(_custom_columns);
104  for (unsigned int i = 0; i < _custom_columns; ++i)
105  if (!(iss >> d._custom[i]))
106  mooseError("Unable to read in EBSD custom data column #", i);
107 
108  if (x < _minx || y < _miny || x > _maxx || y > _maxy ||
109  (g.dim == 3 && (z < _minz || z > _maxz)))
110  mooseError("EBSD Data ouside of the domain declared in the header ([",
111  _minx,
112  ':',
113  _maxx,
114  "], [",
115  _miny,
116  ':',
117  _maxy,
118  "], [",
119  _minz,
120  ':',
121  _maxz,
122  "]) dim=",
123  g.dim,
124  "\n",
125  line);
126 
127  d._p = Point(x, y, z);
128 
129  // determine number of grains in the dataset
130  if (_global_id_map.find(d._feature_id) == _global_id_map.end())
131  _global_id_map[d._feature_id] = _grain_num++;
132 
133  unsigned int global_index = indexFromPoint(Point(x, y, z));
134  _data[global_index] = d;
135  }
136  }
137  stream_in.close();
138 
139  // Resize the variables
140  _avg_data.resize(_grain_num);
141  _avg_angles.resize(_grain_num);
142 
143  // clear the averages
144  for (unsigned int i = 0; i < _grain_num; ++i)
145  {
146  EBSDAvgData & a = _avg_data[i];
147  a._symmetry = a._phase = a._n = 0;
148  a._p = 0.0;
149  a._custom.assign(_custom_columns, 0.0);
150 
151  EulerAngles & b = _avg_angles[i];
152  b.phi1 = b.Phi = b.phi2 = 0.0;
153  }
154 
155  // Iterate through data points to get average variable values for each grain
156  for (auto & j : _data)
157  {
158  EBSDAvgData & a = _avg_data[_global_id_map[j._feature_id]];
159  EulerAngles & b = _avg_angles[_global_id_map[j._feature_id]];
160 
161  // use Eigen::Quaternion<Real> here?
162  b.phi1 += j._phi1;
163  b.Phi += j._Phi;
164  b.phi2 += j._phi2;
165 
166  if (a._n == 0)
167  a._phase = j._phase;
168  else if (a._phase != j._phase)
169  mooseError("An EBSD feature needs to have a uniform phase.");
170 
171  if (a._n == 0)
172  a._symmetry = j._symmetry;
173  else if (a._symmetry != j._symmetry)
174  mooseError("An EBSD feature needs to have a uniform symmetry parameter.");
175 
176  for (unsigned int i = 0; i < _custom_columns; ++i)
177  a._custom[i] += j._custom[i];
178 
179  // store the feature (or grain) ID
180  a._feature_id = j._feature_id;
181 
182  a._p += j._p;
183  a._n++;
184  }
185 
186  for (unsigned int i = 0; i < _grain_num; ++i)
187  {
188  EBSDAvgData & a = _avg_data[i];
189  EulerAngles & b = _avg_angles[i];
190 
191  if (a._n == 0)
192  continue;
193 
194  // TODO: need better way to average angles
195  b.phi1 /= Real(a._n);
196  b.Phi /= Real(a._n);
197  b.phi2 /= Real(a._n);
198 
199  // link the EulerAngles into the EBSDAvgData for access via the functors
200  a._angles = &b;
201 
202  if (a._phase >= _global_id.size())
203  _global_id.resize(a._phase + 1);
204 
205  // The averaged per grain data locally contains the phase id, local id, and
206  // original feature id. It is stored contiguously indexed by global id.
207  a._local_id = _global_id[a._phase].size();
208  _global_id[a._phase].push_back(i);
209 
210  a._p *= 1.0 / Real(a._n);
211 
212  for (unsigned int i = 0; i < _custom_columns; ++i)
213  a._custom[i] /= Real(a._n);
214  }
215 
216  // Build maps to indicate the weights with which grain and phase data
217  // from the surrounding elements contributes to a node fo IC purposes
219 }
unsigned int _custom_columns
number of additional custom data columns
Definition: EBSDReader.h:136
Real _minx
Grid origin.
Definition: EBSDReader.h:172
Real _maxx
Maximum grid extent.
Definition: EBSDReader.h:175
unsigned indexFromPoint(const Point &p) const
Computes a global index in the _data array given an input centroid point.
Definition: EBSDReader.C:260
unsigned _nz
Definition: EBSDReader.h:166
Real _miny
Definition: EBSDReader.h:172
unsigned int _grain_num
Variables needed to determine reduced order parameter values.
Definition: EBSDReader.h:129
unsigned _nx
The number of values in the x, y and z directions.
Definition: EBSDReader.h:166
MooseMesh & _mesh
MooseMesh Variables.
Definition: EBSDReader.h:124
Real _minz
Definition: EBSDReader.h:172
std::map< unsigned int, unsigned int > _global_id_map
map from feature_id to global_id
Definition: EBSDReader.h:148
unsigned _ny
Definition: EBSDReader.h:166
Real _dx
The spacing of the values in x, y and z directions.
Definition: EBSDReader.h:169
std::vector< EBSDAvgData > _avg_data
Averages by (global) grain ID.
Definition: EBSDReader.h:142
Real _maxz
Definition: EBSDReader.h:175
Euler angle triplet.
Definition: EulerAngles.h:20
Real _maxy
Definition: EBSDReader.h:175
std::vector< EulerAngles > _avg_angles
Euler Angles by (global) grain ID.
Definition: EBSDReader.h:145
void buildNodeWeightMaps()
Build grain and phase weight maps.
Definition: EBSDReader.C:335
const EBSDMeshGeometry & getEBSDGeometry() const
Definition: EBSDMesh.h:43
Mesh generated from parameters.
Definition: EBSDMesh.h:22
const std::string & getEBSDFilename() const
Definition: EBSDMesh.h:44
std::vector< std::vector< unsigned int > > _global_id
global ID for given phases and grains
Definition: EBSDReader.h:151
std::vector< EBSDPointData > _data
Logically three-dimensional data indexed by geometric points in a 1D vector.
Definition: EBSDReader.h:139

Member Data Documentation

std::vector<EulerAngles> EBSDReader::_avg_angles
protected

Euler Angles by (global) grain ID.

Definition at line 145 of file EBSDReader.h.

Referenced by getEulerAngles(), and readFile().

std::vector<EBSDAvgData> EBSDReader::_avg_data
protected

Averages by (global) grain ID.

Definition at line 142 of file EBSDReader.h.

Referenced by getAvgData(), getFeatureID(), indexFromIndex(), and readFile().

Point EBSDReader::_bottom_left
protected

Definition at line 130 of file EBSDReader.h.

unsigned int EBSDReader::_custom_columns
protected

number of additional custom data columns

Definition at line 136 of file EBSDReader.h.

Referenced by getAvgDataAccessFunctor(), getPointDataAccessFunctor(), and readFile().

std::vector<EBSDPointData> EBSDReader::_data
protected

Logically three-dimensional data indexed by geometric points in a 1D vector.

Definition at line 139 of file EBSDReader.h.

Referenced by getData(), indexFromPoint(), and readFile().

Real EBSDReader::_dx
protected

The spacing of the values in x, y and z directions.

Definition at line 169 of file EBSDReader.h.

Referenced by indexFromPoint(), and readFile().

Real EBSDReader::_dy
protected

Definition at line 169 of file EBSDReader.h.

Referenced by indexFromPoint(), and readFile().

Real EBSDReader::_dz
protected

Definition at line 169 of file EBSDReader.h.

Referenced by indexFromPoint(), and readFile().

std::vector<std::vector<unsigned int> > EBSDReader::_global_id
protected

global ID for given phases and grains

Definition at line 151 of file EBSDReader.h.

Referenced by getAvgData(), getFeatureID(), getGlobalID(), getGrainNum(), getPhaseNum(), and readFile().

std::map<unsigned int, unsigned int> EBSDReader::_global_id_map
protected

map from feature_id to global_id

Definition at line 148 of file EBSDReader.h.

Referenced by getGlobalID(), and readFile().

unsigned int EBSDReader::_grain_num
protected

Variables needed to determine reduced order parameter values.

Definition at line 129 of file EBSDReader.h.

Referenced by getGrainNum(), and readFile().

Real EBSDReader::_maxx
protected

Maximum grid extent.

Definition at line 175 of file EBSDReader.h.

Referenced by readFile().

Real EBSDReader::_maxy
protected

Definition at line 175 of file EBSDReader.h.

Referenced by readFile().

Real EBSDReader::_maxz
protected

Definition at line 175 of file EBSDReader.h.

Referenced by readFile().

MooseMesh& EBSDReader::_mesh
protected

MooseMesh Variables.

Definition at line 124 of file EBSDReader.h.

Referenced by buildNodeWeightMaps(), and readFile().

unsigned int EBSDReader::_mesh_dimension
protected

Dimension of the problem domain.

Definition at line 163 of file EBSDReader.h.

Referenced by indexFromPoint().

Real EBSDReader::_minx
protected

Grid origin.

Definition at line 172 of file EBSDReader.h.

Referenced by indexFromPoint(), and readFile().

Real EBSDReader::_miny
protected

Definition at line 172 of file EBSDReader.h.

Referenced by indexFromPoint(), and readFile().

Real EBSDReader::_minz
protected

Definition at line 172 of file EBSDReader.h.

Referenced by indexFromPoint(), and readFile().

NonlinearSystemBase& EBSDReader::_nl
protected

Definition at line 125 of file EBSDReader.h.

std::map<dof_id_type, std::vector<Real> > EBSDReader::_node_to_grain_weight_map
protected

Map of grain weights per node.

Definition at line 154 of file EBSDReader.h.

Referenced by buildNodeWeightMaps(), and getNodeToGrainWeightMap().

std::map<dof_id_type, std::vector<Real> > EBSDReader::_node_to_phase_weight_map
protected

Map of phase weights per node.

Definition at line 157 of file EBSDReader.h.

Referenced by buildNodeWeightMaps(), and getNodeToPhaseWeightMap().

unsigned EBSDReader::_nx
protected

The number of values in the x, y and z directions.

Definition at line 166 of file EBSDReader.h.

Referenced by indexFromPoint(), and readFile().

unsigned EBSDReader::_ny
protected

Definition at line 166 of file EBSDReader.h.

Referenced by indexFromPoint(), and readFile().

unsigned EBSDReader::_nz
protected

Definition at line 166 of file EBSDReader.h.

Referenced by readFile().

Point EBSDReader::_range
protected

Definition at line 132 of file EBSDReader.h.

const int& EBSDReader::_time_step
protected

current timestep. Maps are only rebuild on mesh change during time step zero

Definition at line 160 of file EBSDReader.h.

Referenced by meshChanged().

Point EBSDReader::_top_right
protected

Definition at line 131 of file EBSDReader.h.


The documentation for this class was generated from the following files: