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 37 of file EBSDReader.h.

Constructor & Destructor Documentation

◆ EBSDReader()

EBSDReader::EBSDReader ( const InputParameters &  params)

Definition at line 32 of file EBSDReader.C.

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

◆ ~EBSDReader()

EBSDReader::~EBSDReader ( )
virtual

Definition at line 225 of file EBSDReader.C.

225 {}

Member Function Documentation

◆ buildNodeWeightMaps()

void EBSDReader::buildNodeWeightMaps ( )
protected

Build grain and phase weight maps.

Definition at line 343 of file EBSDReader.C.

Referenced by meshChanged(), and readFile().

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

◆ execute()

virtual void EBSDReader::execute ( )
inlinevirtual

Definition at line 46 of file EBSDReader.h.

46 {}

◆ finalize()

virtual void EBSDReader::finalize ( )
inlinevirtual

Definition at line 47 of file EBSDReader.h.

47 {}

◆ getAvgData() [1/2]

const EBSDReader::EBSDAvgData & EBSDReader::getAvgData ( unsigned int  i) const

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

Definition at line 234 of file EBSDReader.C.

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

235 {
236  return _avg_data[indexFromIndex(var)];
237 }
std::vector< EBSDAvgData > _avg_data
Averages 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:297

◆ getAvgData() [2/2]

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 246 of file EBSDReader.C.

247 {
248  return _avg_data[indexFromIndex(_global_id[phase][local_id])];
249 }
std::vector< EBSDAvgData > _avg_data
Averages 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:297
std::vector< std::vector< unsigned int > > _global_id
global ID for given phases and grains
Definition: EBSDReader.h:154

◆ getAvgDataAccessFunctor()

MooseSharedPointer< EBSDAccessFunctors::EBSDAvgDataFunctor > EBSDReader::getAvgDataAccessFunctor ( const MooseEnum &  field_name) const

Factory function to return a average functor specified by name.

Definition at line 438 of file EBSDReader.C.

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

◆ getAvgDataFieldType()

MooseEnum EBSDAccessFunctors::getAvgDataFieldType ( )
staticinherited

Definition at line 18 of file EBSDAccessFunctors.C.

Referenced by validParams< EBSDReaderAvgDataAux >().

19 {
20  return MooseEnum("phi1 phi phi2 phase symmetry local_id feature_id", "", true);
21 }

◆ getData()

const EBSDReader::EBSDPointData & EBSDReader::getData ( const Point &  p) const

Get the requested type of data at the point p.

Definition at line 228 of file EBSDReader.C.

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

229 {
230  return _data[indexFromPoint(p)];
231 }
unsigned indexFromPoint(const Point &p) const
Computes a global index in the _data array given an input centroid point.
Definition: EBSDReader.C:264
std::vector< EBSDPointData > _data
Logically three-dimensional data indexed by geometric points in a 1D vector.
Definition: EBSDReader.h:142

◆ getEulerAngles()

const EulerAngles & EBSDReader::getEulerAngles ( unsigned int  var) const
virtual

EulerAngleProvider interface implementation to fetch a triplet of Euler angles.

Implements EulerAngleProvider.

Definition at line 240 of file EBSDReader.C.

241 {
242  return _avg_angles[indexFromIndex(var)];
243 }
unsigned indexFromIndex(unsigned int var) const
Transfer the index into the _avg_data array from given index.
Definition: EBSDReader.C:297
std::vector< EulerAngles > _avg_angles
Euler Angles by (global) grain ID.
Definition: EBSDReader.h:148

◆ getFeatureID() [1/2]

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 85 of file EBSDReader.h.

86  {
87  return _avg_data[_global_id[phase][local_id]]._feature_id;
88  }
std::vector< EBSDAvgData > _avg_data
Averages by (global) grain ID.
Definition: EBSDReader.h:145
std::vector< std::vector< unsigned int > > _global_id
global ID for given phases and grains
Definition: EBSDReader.h:154

◆ getFeatureID() [2/2]

unsigned int EBSDReader::getFeatureID ( unsigned int  global_id) const
inline

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

Definition at line 90 of file EBSDReader.h.

91  {
92  return _avg_data[global_id]._feature_id;
93  }
std::vector< EBSDAvgData > _avg_data
Averages by (global) grain ID.
Definition: EBSDReader.h:145

◆ getGlobalID() [1/2]

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 96 of file EBSDReader.h.

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

97  {
98  return _global_id[phase][local_id];
99  }
std::vector< std::vector< unsigned int > > _global_id
global ID for given phases and grains
Definition: EBSDReader.h:154

◆ getGlobalID() [2/2]

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 326 of file EBSDReader.C.

327 {
328  auto it = _global_id_map.find(feature_id);
329  if (it == _global_id_map.end())
330  mooseError("Invalid Feature ID");
331  return it->second;
332 }
std::map< unsigned int, unsigned int > _global_id_map
map from feature_id to global_id
Definition: EBSDReader.h:151

◆ getGrainNum() [1/2]

unsigned int EBSDReader::getGrainNum ( ) const
virtual

Return the total number of grains.

Implements EulerAngleProvider.

Definition at line 252 of file EBSDReader.C.

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

253 {
254  return _grain_num;
255 }
unsigned int _grain_num
Variables needed to determine reduced order parameter values.
Definition: EBSDReader.h:132

◆ getGrainNum() [2/2]

unsigned int EBSDReader::getGrainNum ( unsigned int  phase) const

Return the number of grains in a given phase.

Definition at line 258 of file EBSDReader.C.

259 {
260  return _global_id[phase].size();
261 }
std::vector< std::vector< unsigned int > > _global_id
global ID for given phases and grains
Definition: EBSDReader.h:154

◆ getNodeToGrainWeightMap()

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 314 of file EBSDReader.C.

315 {
317 }
std::map< dof_id_type, std::vector< Real > > _node_to_grain_weight_map
Map of grain weights per node.
Definition: EBSDReader.h:157

◆ getNodeToPhaseWeightMap()

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 320 of file EBSDReader.C.

321 {
323 }
std::map< dof_id_type, std::vector< Real > > _node_to_phase_weight_map
Map of phase weights per node.
Definition: EBSDReader.h:160

◆ getPhaseNum()

virtual unsigned int EBSDReader::getPhaseNum ( ) const
inlinevirtual

Return the total number of phases.

Definition at line 77 of file EBSDReader.h.

Referenced by buildNodeWeightMaps().

77 { return _global_id.size(); }
std::vector< std::vector< unsigned int > > _global_id
global ID for given phases and grains
Definition: EBSDReader.h:154

◆ getPointDataAccessFunctor()

MooseSharedPointer< EBSDAccessFunctors::EBSDPointDataFunctor > EBSDReader::getPointDataAccessFunctor ( const MooseEnum &  field_name) const

Factory function to return a point functor specified by name.

Definition at line 392 of file EBSDReader.C.

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

◆ getPointDataFieldType()

MooseEnum EBSDAccessFunctors::getPointDataFieldType ( )
staticinherited

Definition at line 12 of file EBSDAccessFunctors.C.

Referenced by validParams< EBSDReaderPointDataAux >().

13 {
14  return MooseEnum("phi1 phi phi2 feature_id phase symmetry", "", true);
15 }

◆ indexFromIndex()

unsigned int EBSDReader::indexFromIndex ( unsigned int  var) const
protected

Transfer the index into the _avg_data array from given index.

Definition at line 297 of file EBSDReader.C.

Referenced by getAvgData(), and getEulerAngles().

298 {
299 
300  // Transfer the index into the _avg_data array.
301  unsigned avg_index = var;
302 
303  // Don't access out of range!
304  if (avg_index >= _avg_data.size())
305  mooseError("Error! Index out of range in EBSDReader::indexFromIndex(), index: ",
306  avg_index,
307  " size: ",
308  _avg_data.size());
309 
310  return avg_index;
311 }
std::vector< EBSDAvgData > _avg_data
Averages by (global) grain ID.
Definition: EBSDReader.h:145

◆ indexFromPoint()

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 264 of file EBSDReader.C.

Referenced by getData(), and readFile().

265 {
266  // Don't assume an ordering on the input data, use the (x, y,
267  // z) values of this centroid to determine the index.
268  unsigned int x_index, y_index, z_index, global_index;
269 
270  x_index = (unsigned int)((p(0) - _minx) / _dx);
271  y_index = (unsigned int)((p(1) - _miny) / _dy);
272  if (p(0) <= _minx || p(0) >= _maxx || p(1) <= _miny || p(1) >= _maxy)
273  mooseError("Data points must be on the interior of the mesh elements. In EBSDReader ", name());
274  if (_mesh_dimension == 3)
275  {
276  z_index = (unsigned int)((p(2) - _minz) / _dz);
277  global_index = z_index * _ny;
278  if (p(2) <= _minz || p(2) >= _maxz)
279  mooseError("Data points must be on the interior of the mesh elements. In EBSDReader ",
280  name());
281  }
282  else
283  global_index = 0;
284 
285  // Compute the global index into the _data array. This stores points
286  // in a [z][y][x] ordering.
287  global_index = (global_index + y_index) * _nx + x_index;
288 
289  // Don't access out of range!
290  mooseAssert(global_index < _data.size(),
291  "global_index " << global_index << " points out of _data range: " << _data.size());
292 
293  return global_index;
294 }
Real _minx
Grid origin.
Definition: EBSDReader.h:175
Real _maxx
Maximum grid extent.
Definition: EBSDReader.h:178
Real _miny
Definition: EBSDReader.h:175
unsigned _nx
The number of values in the x, y and z directions.
Definition: EBSDReader.h:169
Real _minz
Definition: EBSDReader.h:175
const std::string name
Definition: Setup.h:22
unsigned _ny
Definition: EBSDReader.h:169
Real _dx
The spacing of the values in x, y and z directions.
Definition: EBSDReader.h:172
unsigned int _mesh_dimension
Dimension of the problem domain.
Definition: EBSDReader.h:166
Real _maxz
Definition: EBSDReader.h:178
Real _maxy
Definition: EBSDReader.h:178
std::vector< EBSDPointData > _data
Logically three-dimensional data indexed by geometric points in a 1D vector.
Definition: EBSDReader.h:142

◆ initialize()

virtual void EBSDReader::initialize ( )
inlinevirtual

Definition at line 45 of file EBSDReader.h.

45 {}

◆ meshChanged()

void EBSDReader::meshChanged ( )

Maps need to be updated when the mesh changes.

Definition at line 335 of file EBSDReader.C.

336 {
337  // maps are only rebuild for use in initial conditions, which happens in time step zero
338  if (_time_step == 0)
340 }
const int & _time_step
current timestep. Maps are only rebuild on mesh change during time step zero
Definition: EBSDReader.h:163
void buildNodeWeightMaps()
Build grain and phase weight maps.
Definition: EBSDReader.C:343

◆ readFile()

void EBSDReader::readFile ( )
virtual

Definition at line 51 of file EBSDReader.C.

Referenced by EBSDReader().

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

Member Data Documentation

◆ _avg_angles

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

Euler Angles by (global) grain ID.

Definition at line 148 of file EBSDReader.h.

Referenced by getEulerAngles(), and readFile().

◆ _avg_data

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

Averages by (global) grain ID.

Definition at line 145 of file EBSDReader.h.

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

◆ _bottom_left

Point EBSDReader::_bottom_left
protected

Definition at line 133 of file EBSDReader.h.

◆ _custom_columns

unsigned int EBSDReader::_custom_columns
protected

number of additional custom data columns

Definition at line 139 of file EBSDReader.h.

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

◆ _data

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

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

Definition at line 142 of file EBSDReader.h.

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

◆ _dx

Real EBSDReader::_dx
protected

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

Definition at line 172 of file EBSDReader.h.

Referenced by indexFromPoint(), and readFile().

◆ _dy

Real EBSDReader::_dy
protected

Definition at line 172 of file EBSDReader.h.

Referenced by indexFromPoint(), and readFile().

◆ _dz

Real EBSDReader::_dz
protected

Definition at line 172 of file EBSDReader.h.

Referenced by indexFromPoint(), and readFile().

◆ _global_id

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

global ID for given phases and grains

Definition at line 154 of file EBSDReader.h.

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

◆ _global_id_map

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

map from feature_id to global_id

Definition at line 151 of file EBSDReader.h.

Referenced by getGlobalID(), and readFile().

◆ _grain_num

unsigned int EBSDReader::_grain_num
protected

Variables needed to determine reduced order parameter values.

Definition at line 132 of file EBSDReader.h.

Referenced by getGrainNum(), and readFile().

◆ _maxx

Real EBSDReader::_maxx
protected

Maximum grid extent.

Definition at line 178 of file EBSDReader.h.

Referenced by indexFromPoint(), and readFile().

◆ _maxy

Real EBSDReader::_maxy
protected

Definition at line 178 of file EBSDReader.h.

Referenced by indexFromPoint(), and readFile().

◆ _maxz

Real EBSDReader::_maxz
protected

Definition at line 178 of file EBSDReader.h.

Referenced by indexFromPoint(), and readFile().

◆ _mesh

MooseMesh& EBSDReader::_mesh
protected

MooseMesh Variables.

Definition at line 127 of file EBSDReader.h.

Referenced by buildNodeWeightMaps(), and readFile().

◆ _mesh_dimension

unsigned int EBSDReader::_mesh_dimension
protected

Dimension of the problem domain.

Definition at line 166 of file EBSDReader.h.

Referenced by indexFromPoint().

◆ _minx

Real EBSDReader::_minx
protected

Grid origin.

Definition at line 175 of file EBSDReader.h.

Referenced by indexFromPoint(), and readFile().

◆ _miny

Real EBSDReader::_miny
protected

Definition at line 175 of file EBSDReader.h.

Referenced by indexFromPoint(), and readFile().

◆ _minz

Real EBSDReader::_minz
protected

Definition at line 175 of file EBSDReader.h.

Referenced by indexFromPoint(), and readFile().

◆ _nl

NonlinearSystemBase& EBSDReader::_nl
protected

Definition at line 128 of file EBSDReader.h.

◆ _node_to_grain_weight_map

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

Map of grain weights per node.

Definition at line 157 of file EBSDReader.h.

Referenced by buildNodeWeightMaps(), and getNodeToGrainWeightMap().

◆ _node_to_phase_weight_map

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

Map of phase weights per node.

Definition at line 160 of file EBSDReader.h.

Referenced by buildNodeWeightMaps(), and getNodeToPhaseWeightMap().

◆ _nx

unsigned EBSDReader::_nx
protected

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

Definition at line 169 of file EBSDReader.h.

Referenced by indexFromPoint(), and readFile().

◆ _ny

unsigned EBSDReader::_ny
protected

Definition at line 169 of file EBSDReader.h.

Referenced by indexFromPoint(), and readFile().

◆ _nz

unsigned EBSDReader::_nz
protected

Definition at line 169 of file EBSDReader.h.

Referenced by readFile().

◆ _range

Point EBSDReader::_range
protected

Definition at line 135 of file EBSDReader.h.

◆ _time_step

const int& EBSDReader::_time_step
protected

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

Definition at line 163 of file EBSDReader.h.

Referenced by meshChanged().

◆ _top_right

Point EBSDReader::_top_right
protected

Definition at line 134 of file EBSDReader.h.


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