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

PolycrystalVoronoiVoidIC initializes either grain or void values for a voronoi tesselation with voids distributed along the grain boundaries. More...

#include <PolycrystalVoronoiVoidIC.h>

Inheritance diagram for PolycrystalVoronoiVoidIC:
[legend]

Classes

struct  DistancePoint
 Type for distance and point. More...
 
struct  DistancePointComparator
 Sorts the temp_centerpoints into order of magnitude. More...
 

Public Member Functions

 PolycrystalVoronoiVoidIC (const InputParameters &parameters)
 
virtual void initialSetup () override
 

Static Public Member Functions

static InputParameters actionParameters ()
 

Protected Types

enum  ProfileType { ProfileType::COS, ProfileType::TANH }
 

Protected Member Functions

virtual void computeCircleCenters () override
 
virtual Real value (const Point &p) override
 
virtual RealGradient gradient (const Point &p) override
 
virtual Real grainValueCalc (const Point &p)
 
virtual void computeGrainCenters ()
 
virtual void computeCircleRadii () override
 
virtual Real computeCircleValue (const Point &p, const Point &center, const Real &radius)
 
virtual RealGradient computeCircleGradient (const Point &p, const Point &center, const Real &radius)
 

Protected Attributes

const MooseEnum _structure_type
 
const unsigned int _op_num
 
unsigned int _grain_num
 
const unsigned int _op_index
 
const unsigned int _rand_seed
 
const bool _columnar_3D
 
const FileName _file_name
 
std::vector< Point > _centerpoints
 
std::vector< unsigned int > _assigned_op
 
struct PolycrystalVoronoiVoidIC::DistancePointComparator _customLess
 
const unsigned int _numbub
 
const Real _bubspac
 
const unsigned int _max_num_tries
 
const Real _radius
 
const Real _radius_variation
 
const MooseEnum _radius_variation_type
 
Point _bottom_left
 
Point _top_right
 
Point _range
 
MooseMesh & _mesh
 
Real _invalue
 
Real _outvalue
 
Real _int_width
 
bool _3D_spheres
 
bool _zero_gradient
 
unsigned int _num_dim
 
std::vector< Point > _centers
 
std::vector< Real > _radii
 
enum SmoothCircleBaseIC::ProfileType _profile
 
MooseRandom _random
 

Detailed Description

PolycrystalVoronoiVoidIC initializes either grain or void values for a voronoi tesselation with voids distributed along the grain boundaries.

Definition at line 27 of file PolycrystalVoronoiVoidIC.h.

Member Enumeration Documentation

◆ ProfileType

enum SmoothCircleBaseIC::ProfileType
strongprotectedinherited
Enumerator
COS 
TANH 

Definition at line 58 of file SmoothCircleBaseIC.h.

59  {
60  COS,
61  TANH
62  } _profile;
enum SmoothCircleBaseIC::ProfileType _profile

Constructor & Destructor Documentation

◆ PolycrystalVoronoiVoidIC()

PolycrystalVoronoiVoidIC::PolycrystalVoronoiVoidIC ( const InputParameters &  parameters)

Definition at line 59 of file PolycrystalVoronoiVoidIC.C.

60  : MultiSmoothCircleIC(parameters),
61  _structure_type(getParam<MooseEnum>("structure_type")),
62  _op_num(getParam<unsigned int>("op_num")),
63  _grain_num(getParam<unsigned int>("grain_num")),
64  _op_index(getParam<unsigned int>("op_index")),
65  _rand_seed(getParam<unsigned int>("rand_seed")),
66  _columnar_3D(getParam<bool>("columnar_3D")),
67  _file_name(getParam<FileName>("file_name"))
68 {
69  if (_invalue < _outvalue)
70  mooseError("PolycrystalVoronoiVoidIC requires that the voids be "
71  "represented with invalue > outvalue");
72  if (_numbub == 0)
73  mooseError("PolycrystalVoronoiVoidIC requires numbub > 0. If you want no voids to "
74  "be "
75  "represented, use invalue = outvalue. In general, you should use "
76  "PolycrystalVoronoi to represent Voronoi grain structures without "
77  "voids.");
78 
79  if (_file_name == "" && _grain_num == 0)
80  mooseError("grain_num must be provided if the grain centriods are not read from a file");
81 
82  if (_file_name != "" && _grain_num > 0)
83  mooseWarning("grain_num is ignored and will be determined from the file.");
84 }
const unsigned int _numbub
MultiSmoothCircleIC(const InputParameters &parameters)

Member Function Documentation

◆ actionParameters()

InputParameters PolycrystalVoronoiVoidIC::actionParameters ( )
static

Definition at line 18 of file PolycrystalVoronoiVoidIC.C.

Referenced by validParams< PolycrystalVoronoiVoidIC >(), and validParams< PolycrystalVoronoiVoidICAction >().

19 {
20  InputParameters params = validParams<MultiSmoothCircleIC>();
21 
22  params.addRequiredParam<unsigned int>("op_num", "Number of order parameters");
23  params.addParam<unsigned int>(
24  "grain_num", 0, "Number of grains being represented by the order parameters");
25 
26  params.addParam<unsigned int>("rand_seed", 12444, "The random seed");
27 
28  params.addParam<bool>(
29  "columnar_3D", false, "3D microstructure will be columnar in the z-direction?");
30 
31  return params;
32 }
InputParameters validParams< MultiSmoothCircleIC >()

◆ computeCircleCenters()

void PolycrystalVoronoiVoidIC::computeCircleCenters ( )
overrideprotectedvirtual

Reimplemented from MultiSmoothCircleIC.

Definition at line 110 of file PolycrystalVoronoiVoidIC.C.

111 {
112  _centers.resize(_numbub);
113 
114  // This Code will place void center points on grain boundaries
115  for (unsigned int vp = 0; vp < _numbub; ++vp)
116  {
117  bool try_again;
118  unsigned int num_tries = 0;
119 
120  do
121  {
122  try_again = false;
123  num_tries++;
124 
125  if (num_tries > _max_num_tries)
126  mooseError("Too many tries of assigning void centers in "
127  "PolycrystalVoronoiVoidIC");
128 
129  Point rand_point;
130 
131  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
132  rand_point(i) = _bottom_left(i) + _range(i) * MooseRandom::rand();
133 
134  // Allow the vectors to be sorted based on their distance from the
135  // rand_point
136  std::vector<PolycrystalVoronoiVoidIC::DistancePoint> diff(_grain_num);
137 
138  for (unsigned int gr = 0; gr < _grain_num; ++gr)
139  {
140  diff[gr].d = _mesh.minPeriodicDistance(_var.number(), rand_point, _centerpoints[gr]);
141  diff[gr].gr = gr;
142  }
143 
144  std::sort(diff.begin(), diff.end(), _customLess);
145 
146  Point closest_point = _centerpoints[diff[0].gr];
147  Point next_closest_point = _centerpoints[diff[1].gr];
148 
149  // Find Slope of Line in the plane orthogonal to the diff_centerpoint
150  // vector
151  Point diff_centerpoints =
152  _mesh.minPeriodicVector(_var.number(), closest_point, next_closest_point);
153  Point diff_rand_center = _mesh.minPeriodicVector(_var.number(), closest_point, rand_point);
154  Point normal_vector = diff_centerpoints.cross(diff_rand_center);
155  Point slope = normal_vector.cross(diff_centerpoints);
156 
157  // Midpoint position vector between two center points
158  Point midpoint = closest_point + (0.5 * diff_centerpoints);
159 
160  // Solve for the scalar multiplier solution on the line
161  Real lambda = 0;
162  Point mid_rand_vector = _mesh.minPeriodicVector(_var.number(), midpoint, rand_point);
163 
164  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
165  lambda += (mid_rand_vector(i) * slope(i)) /
166  (slope(0) * slope(0) + slope(1) * slope(1) + slope(2) * slope(2));
167 
168  // Assigning points to vector
169  _centers[vp] = slope * lambda + midpoint;
170 
171  // Checking to see if points are in the domain ONLY WORKS FOR PERIODIC
172  for (unsigned int i = 0; i < LIBMESH_DIM; i++)
173  if ((_centers[vp](i) > _top_right(i)) || (_centers[vp](i) < _bottom_left(i)))
174  try_again = true;
175 
176  for (unsigned int i = 0; i < vp; ++i)
177  {
178  Real dist = _mesh.minPeriodicDistance(_var.number(), _centers[vp], _centers[i]);
179 
180  if (dist < _bubspac)
181  try_again = true;
182  }
183 
184  // Two algorithms are available for screening bubbles falling in grain
185  // interior. They produce
186  // nearly identical results.
187  // Here only one is listed. The other one is available upon request.
188 
189  // Use circle center for checking whether voids are at GBs
190  if (try_again == false)
191  {
192  Real min_rij_1, min_rij_2, rij, rij_diff_tol;
193 
194  min_rij_1 = _range.norm();
195  min_rij_2 = _range.norm();
196 
197  rij_diff_tol = 0.1 * _radius;
198 
199  for (unsigned int gr = 0; gr < _grain_num; ++gr)
200  {
201  rij = _mesh.minPeriodicDistance(_var.number(), _centers[vp], _centerpoints[gr]);
202 
203  if (rij < min_rij_1)
204  {
205  min_rij_2 = min_rij_1;
206  min_rij_1 = rij;
207  }
208  else if (rij < min_rij_2)
209  min_rij_2 = rij;
210  }
211 
212  if (std::abs(min_rij_1 - min_rij_2) > rij_diff_tol)
213  try_again = true;
214  }
215 
216  } while (try_again == true);
217  }
218 }
struct PolycrystalVoronoiVoidIC::DistancePointComparator _customLess
std::vector< Point > _centers
const unsigned int _max_num_tries
std::vector< Point > _centerpoints
const unsigned int _numbub

◆ computeCircleGradient()

RealGradient SmoothCircleBaseIC::computeCircleGradient ( const Point &  p,
const Point &  center,
const Real &  radius 
)
protectedvirtualinherited

Definition at line 150 of file SmoothCircleBaseIC.C.

Referenced by SmoothCircleBaseIC::gradient().

153 {
154  Point l_center = center;
155  Point l_p = p;
156  if (!_3D_spheres) // Create 3D cylinders instead of spheres
157  {
158  l_p(2) = 0.0;
159  l_center(2) = 0.0;
160  }
161  // Compute the distance between the current point and the center
162  Real dist = _mesh.minPeriodicDistance(_var.number(), l_p, l_center);
163 
164  // early return if we are probing the center of the circle
165  if (dist == 0.0)
166  return 0.0;
167 
168  Real DvalueDr = 0.0;
169  switch (_profile)
170  {
171  case ProfileType::COS:
172  if (dist < radius + _int_width / 2.0 && dist > radius - _int_width / 2.0)
173  {
174  const Real int_pos = (dist - radius + _int_width / 2.0) / _int_width;
175  const Real Dint_posDr = 1.0 / _int_width;
176  DvalueDr = Dint_posDr * (_invalue - _outvalue) *
177  (-std::sin(int_pos * libMesh::pi) * libMesh::pi) / 2.0;
178  }
179  break;
180 
181  case ProfileType::TANH:
182  DvalueDr = -(_invalue - _outvalue) * 0.5 / _int_width * libMesh::pi *
183  (1.0 - Utility::pow<2>(std::tanh(4.0 * (radius - dist) / _int_width)));
184  break;
185 
186  default:
187  mooseError("Internal error.");
188  }
189 
190  return _mesh.minPeriodicVector(_var.number(), center, p) * (DvalueDr / dist);
191 }
enum SmoothCircleBaseIC::ProfileType _profile

◆ computeCircleRadii()

void MultiSmoothCircleIC::computeCircleRadii ( )
overrideprotectedvirtualinherited

Implements SmoothCircleBaseIC.

Definition at line 72 of file MultiSmoothCircleIC.C.

73 {
74  _radii.resize(_numbub);
75 
76  for (unsigned int i = 0; i < _numbub; i++)
77  {
78  // Vary bubble radius
79  switch (_radius_variation_type)
80  {
81  case 0: // Uniform distribution
82  _radii[i] = _radius * (1.0 + (1.0 - 2.0 * _random.rand(_tid)) * _radius_variation);
83  break;
84  case 1: // Normal distribution
85  _radii[i] = _random.randNormal(_tid, _radius, _radius_variation);
86  break;
87  case 2: // No variation
88  _radii[i] = _radius;
89  }
90 
91  _radii[i] = std::max(_radii[i], 0.0);
92  }
93 }
std::vector< Real > _radii
const MooseEnum _radius_variation_type
const unsigned int _numbub

◆ computeCircleValue()

Real SmoothCircleBaseIC::computeCircleValue ( const Point &  p,
const Point &  center,
const Real &  radius 
)
protectedvirtualinherited

Reimplemented in RndSmoothCircleIC.

Definition at line 110 of file SmoothCircleBaseIC.C.

Referenced by SmoothCircleBaseIC::gradient(), and SmoothCircleBaseIC::value().

111 {
112  Point l_center = center;
113  Point l_p = p;
114  if (!_3D_spheres) // Create 3D cylinders instead of spheres
115  {
116  l_p(2) = 0.0;
117  l_center(2) = 0.0;
118  }
119  // Compute the distance between the current point and the center
120  Real dist = _mesh.minPeriodicDistance(_var.number(), l_p, l_center);
121 
122  switch (_profile)
123  {
124  case ProfileType::COS:
125  {
126  // Return value
127  Real value = _outvalue; // Outside circle
128 
129  if (dist <= radius - _int_width / 2.0) // Inside circle
130  value = _invalue;
131  else if (dist < radius + _int_width / 2.0) // Smooth interface
132  {
133  Real int_pos = (dist - radius + _int_width / 2.0) / _int_width;
134  value = _outvalue + (_invalue - _outvalue) * (1.0 + std::cos(int_pos * libMesh::pi)) / 2.0;
135  }
136  return value;
137  }
138 
139  case ProfileType::TANH:
140  return (_invalue - _outvalue) * 0.5 *
141  (std::tanh(libMesh::pi * (radius - dist) / _int_width) + 1.0) +
142  _outvalue;
143 
144  default:
145  mooseError("Internal error.");
146  }
147 }
virtual Real value(const Point &p)
enum SmoothCircleBaseIC::ProfileType _profile

◆ computeGrainCenters()

void PolycrystalVoronoiVoidIC::computeGrainCenters ( )
protectedvirtual

Definition at line 294 of file PolycrystalVoronoiVoidIC.C.

Referenced by initialSetup().

295 {
296  if (!_file_name.empty())
297  {
298  MooseUtils::DelimitedFileReader txt_reader(_file_name, &_communicator);
299 
300  txt_reader.read();
301  std::vector<std::string> col_names = txt_reader.getNames();
302  std::vector<std::vector<Real>> data = txt_reader.getData();
303 
304  _grain_num = data[0].size();
305  _centerpoints.resize(_grain_num);
306 
307  for (unsigned int i = 0; i < col_names.size(); ++i)
308  {
309  // Check vector lengths
310  if (data[i].size() != _grain_num)
311  paramError("Columns in ", _file_name, " do not have uniform lengths.");
312  }
313 
314  for (unsigned int grain = 0; grain < _grain_num; ++grain)
315  {
316  _centerpoints[grain](0) = data[0][grain];
317  _centerpoints[grain](1) = data[1][grain];
318  if (col_names.size() > 2)
319  _centerpoints[grain](2) = data[2][grain];
320  else
321  _centerpoints[grain](2) = 0.0;
322  }
323  }
324  else
325  {
326  _centerpoints.resize(_grain_num);
327  // Randomly generate the centers of the individual grains represented by the
328  // Voronoi tessellation
329  for (unsigned int grain = 0; grain < _grain_num; grain++)
330  {
331  for (unsigned int i = 0; i < LIBMESH_DIM; i++)
332  _centerpoints[grain](i) = _bottom_left(i) + _range(i) * MooseRandom::rand();
333 
334  if (_columnar_3D)
335  _centerpoints[grain](2) = _bottom_left(2) + _range(2) * 0.5;
336  }
337  }
338 
339  if (_op_num > _grain_num)
340  mooseError("ERROR in PolycrystalVoronoiVoidIC: Number of order parameters "
341  "(op_num) can't be "
342  "larger than the number of grains (grain_num)");
343 
344  _assigned_op.resize(_grain_num);
345 
346  // Assign grains to specific order parameters in a way that maximizes the
347  // distance
349 }
std::vector< Point > _centerpoints
std::vector< unsigned int > assignPointsToVariables(const std::vector< Point > &centerpoints, const Real op_num, const MooseMesh &mesh, const MooseVariable &var)
std::vector< unsigned int > _assigned_op

◆ gradient()

RealGradient PolycrystalVoronoiVoidIC::gradient ( const Point &  p)
overrideprotectedvirtual

Reimplemented from SmoothCircleBaseIC.

Definition at line 254 of file PolycrystalVoronoiVoidIC.C.

255 {
256  RealGradient gradient;
257  RealGradient void_gradient = MultiSmoothCircleIC::gradient(p);
258 
259  // Order parameter assignment assumes zero gradient (sharp interface)
260  switch (_structure_type)
261  {
262  case 1: // assigning gradient for voids
263  gradient = void_gradient;
264  break;
265  }
266 
267  return gradient;
268 }
virtual RealGradient gradient(const Point &p)
virtual RealGradient gradient(const Point &p) override

◆ grainValueCalc()

Real PolycrystalVoronoiVoidIC::grainValueCalc ( const Point &  p)
protectedvirtual

Definition at line 271 of file PolycrystalVoronoiVoidIC.C.

Referenced by value().

272 {
273  Real val = 0.0;
274 
275  unsigned int min_index =
277 
278  // If the current order parameter index (_op_index) is equal to the min_index,
279  // set the value to
280  // 1.0
281  if (_assigned_op[min_index] == _op_index)
282  val = 1.0;
283 
284  if (val > 1.0)
285  val = 1.0;
286 
287  if (val < 0.0)
288  val = 0.0;
289 
290  return val;
291 }
std::vector< Point > _centerpoints
unsigned int assignPointToGrain(const Point &p, const std::vector< Point > &centerpoints, const MooseMesh &mesh, const MooseVariable &var, const Real maxsize)
std::vector< unsigned int > _assigned_op

◆ initialSetup()

void PolycrystalVoronoiVoidIC::initialSetup ( )
overridevirtual

Reimplemented from MultiSmoothCircleIC.

Definition at line 87 of file PolycrystalVoronoiVoidIC.C.

88 {
89  if (_op_num <= _op_index)
90  mooseError("op_index is too large in CircleGrainVoidIC");
91 
92  MooseRandom::seed(getParam<unsigned int>("rand_seed"));
93  // Set up domain bounds with mesh tools
94  for (unsigned int i = 0; i < LIBMESH_DIM; i++)
95  {
96  _bottom_left(i) = _mesh.getMinInDimension(i);
97  _top_right(i) = _mesh.getMaxInDimension(i);
98  }
100 
101  // Create _centerpoints and _assigned_op vectors
103 
104  // Call initial setup from MultiSmoothCircleIC to create _centers and _radii
105  // for voids
107 }
virtual void initialSetup() override

◆ value()

Real PolycrystalVoronoiVoidIC::value ( const Point &  p)
overrideprotectedvirtual

Reimplemented from SmoothCircleBaseIC.

Definition at line 221 of file PolycrystalVoronoiVoidIC.C.

222 {
223  Real value = 0.0;
224 
225  // Determine value for voids
226  Real void_value = MultiSmoothCircleIC::value(p);
227 
228  // Determine value for grains
229  Real grain_value = grainValueCalc(p);
230 
231  switch (_structure_type)
232  {
233  case 0: // assigning values for grains (order parameters)
234  if (grain_value == 0) // Not in this grain
235  value = grain_value;
236  else // in this grain, but might be in a void
237  if (void_value == _outvalue) // Not in a void
238  value = grain_value;
239  else if (void_value > _outvalue && void_value < _invalue) // On void interface
240  value = 1.0 - (void_value - _outvalue) / (_invalue - _outvalue);
241  else if (void_value == _invalue) // In a void, so op = 0
242  value = 0.0;
243  break;
244 
245  case 1: // assigning values for voids (concentration)
246  value = void_value;
247  break;
248  }
249 
250  return value;
251 }
virtual Real value(const Point &p)
virtual Real value(const Point &p) override
virtual Real grainValueCalc(const Point &p)

Member Data Documentation

◆ _3D_spheres

bool SmoothCircleBaseIC::_3D_spheres
protectedinherited

◆ _assigned_op

std::vector<unsigned int> PolycrystalVoronoiVoidIC::_assigned_op
protected

Definition at line 58 of file PolycrystalVoronoiVoidIC.h.

Referenced by computeGrainCenters(), and grainValueCalc().

◆ _bottom_left

Point MultiSmoothCircleIC::_bottom_left
protectedinherited

◆ _bubspac

const Real MultiSmoothCircleIC::_bubspac
protectedinherited

◆ _centerpoints

std::vector<Point> PolycrystalVoronoiVoidIC::_centerpoints
protected

◆ _centers

std::vector<Point> SmoothCircleBaseIC::_centers
protectedinherited

◆ _columnar_3D

const bool PolycrystalVoronoiVoidIC::_columnar_3D
protected

Definition at line 45 of file PolycrystalVoronoiVoidIC.h.

Referenced by computeGrainCenters().

◆ _customLess

struct PolycrystalVoronoiVoidIC::DistancePointComparator PolycrystalVoronoiVoidIC::_customLess
protected

Referenced by computeCircleCenters().

◆ _file_name

const FileName PolycrystalVoronoiVoidIC::_file_name
protected

Definition at line 47 of file PolycrystalVoronoiVoidIC.h.

Referenced by computeGrainCenters(), and PolycrystalVoronoiVoidIC().

◆ _grain_num

unsigned int PolycrystalVoronoiVoidIC::_grain_num
protected

◆ _int_width

Real SmoothCircleBaseIC::_int_width
protectedinherited

◆ _invalue

Real SmoothCircleBaseIC::_invalue
protectedinherited

◆ _max_num_tries

const unsigned int MultiSmoothCircleIC::_max_num_tries
protectedinherited

◆ _mesh

MooseMesh& SmoothCircleBaseIC::_mesh
protectedinherited

◆ _num_dim

unsigned int SmoothCircleBaseIC::_num_dim
protectedinherited

Definition at line 53 of file SmoothCircleBaseIC.h.

◆ _numbub

const unsigned int MultiSmoothCircleIC::_numbub
protectedinherited

◆ _op_index

const unsigned int PolycrystalVoronoiVoidIC::_op_index
protected

Definition at line 41 of file PolycrystalVoronoiVoidIC.h.

Referenced by grainValueCalc(), and initialSetup().

◆ _op_num

const unsigned int PolycrystalVoronoiVoidIC::_op_num
protected

Definition at line 39 of file PolycrystalVoronoiVoidIC.h.

Referenced by computeGrainCenters(), and initialSetup().

◆ _outvalue

Real SmoothCircleBaseIC::_outvalue
protectedinherited

◆ _profile

enum SmoothCircleBaseIC::ProfileType SmoothCircleBaseIC::_profile
protectedinherited

◆ _radii

std::vector<Real> SmoothCircleBaseIC::_radii
protectedinherited

◆ _radius

const Real MultiSmoothCircleIC::_radius
protectedinherited

◆ _radius_variation

const Real MultiSmoothCircleIC::_radius_variation
protectedinherited

◆ _radius_variation_type

const MooseEnum MultiSmoothCircleIC::_radius_variation_type
protectedinherited

◆ _rand_seed

const unsigned int PolycrystalVoronoiVoidIC::_rand_seed
protected

Definition at line 43 of file PolycrystalVoronoiVoidIC.h.

◆ _random

MooseRandom SmoothCircleBaseIC::_random
protectedinherited

◆ _range

Point MultiSmoothCircleIC::_range
protectedinherited

◆ _structure_type

const MooseEnum PolycrystalVoronoiVoidIC::_structure_type
protected

Definition at line 37 of file PolycrystalVoronoiVoidIC.h.

Referenced by gradient(), and value().

◆ _top_right

Point MultiSmoothCircleIC::_top_right
protectedinherited

◆ _zero_gradient

bool SmoothCircleBaseIC::_zero_gradient
protectedinherited

Definition at line 51 of file SmoothCircleBaseIC.h.

Referenced by SmoothCircleBaseIC::gradient().


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