www.mooseframework.org
Classes | Public Member Functions | Static Public Member Functions | 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 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 ()
 
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
 
const unsigned int _grain_num
 
const unsigned int _op_index
 
const unsigned int _rand_seed
 
const bool _columnar_3D
 
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
 
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 24 of file PolycrystalVoronoiVoidIC.h.

Constructor & Destructor Documentation

PolycrystalVoronoiVoidIC::PolycrystalVoronoiVoidIC ( const InputParameters &  parameters)

Definition at line 47 of file PolycrystalVoronoiVoidIC.C.

48  : MultiSmoothCircleIC(parameters),
49  _structure_type(getParam<MooseEnum>("structure_type")),
50  _op_num(getParam<unsigned int>("op_num")),
51  _grain_num(getParam<unsigned int>("grain_num")),
52  _op_index(getParam<unsigned int>("op_index")),
53  _rand_seed(getParam<unsigned int>("rand_seed")),
54  _columnar_3D(getParam<bool>("columnar_3D"))
55 {
56  if (_invalue < _outvalue)
57  mooseError("PolycrystalVoronoiVoidIC requires that the voids be "
58  "represented with invalue > outvalue");
59  if (_numbub == 0)
60  mooseError("PolycrystalVoronoiVoidIC requires numbub > 0. If you want no voids to "
61  "be "
62  "represented, use invalue = outvalue. In general, you should use "
63  "PolycrystalReducedIC to represent Voronoi grain structures without "
64  "voids.");
65 }
const unsigned int _numbub
MultiSmoothCircleIC(const InputParameters &parameters)

Member Function Documentation

InputParameters PolycrystalVoronoiVoidIC::actionParameters ( )
static

Definition at line 14 of file PolycrystalVoronoiVoidIC.C.

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

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

Reimplemented from MultiSmoothCircleIC.

Definition at line 91 of file PolycrystalVoronoiVoidIC.C.

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

Definition at line 127 of file SmoothCircleBaseIC.C.

Referenced by SmoothCircleBaseIC::gradient().

130 {
131  Point l_center = center;
132  Point l_p = p;
133  if (!_3D_spheres) // Create 3D cylinders instead of spheres
134  {
135  l_p(2) = 0.0;
136  l_center(2) = 0.0;
137  }
138  // Compute the distance between the current point and the center
139  Real dist = _mesh.minPeriodicDistance(_var.number(), l_p, l_center);
140 
141  Real DvalueDr = 0.0;
142 
143  if (dist < radius + _int_width / 2.0 && dist > radius - _int_width / 2.0)
144  {
145  Real int_pos = (dist - radius + _int_width / 2.0) / _int_width;
146  Real Dint_posDr = 1.0 / _int_width;
147  DvalueDr = Dint_posDr * (_invalue - _outvalue) *
148  (-std::sin(int_pos * libMesh::pi) * libMesh::pi) / 2.0;
149  }
150 
151  // Set gradient over the smooth interface
152  if (dist != 0.0)
153  return _mesh.minPeriodicVector(_var.number(), center, p) * (DvalueDr / dist);
154  else
155  return 0.0;
156 }
void MultiSmoothCircleIC::computeCircleRadii ( )
protectedvirtualinherited

Implements SmoothCircleBaseIC.

Definition at line 68 of file MultiSmoothCircleIC.C.

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

Reimplemented in RndSmoothCircleIC.

Definition at line 100 of file SmoothCircleBaseIC.C.

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

101 {
102  Point l_center = center;
103  Point l_p = p;
104  if (!_3D_spheres) // Create 3D cylinders instead of spheres
105  {
106  l_p(2) = 0.0;
107  l_center(2) = 0.0;
108  }
109  // Compute the distance between the current point and the center
110  Real dist = _mesh.minPeriodicDistance(_var.number(), l_p, l_center);
111 
112  // Return value
113  Real value = _outvalue; // Outside circle
114 
115  if (dist <= radius - _int_width / 2.0) // Inside circle
116  value = _invalue;
117  else if (dist < radius + _int_width / 2.0) // Smooth interface
118  {
119  Real int_pos = (dist - radius + _int_width / 2.0) / _int_width;
120  value = _outvalue + (_invalue - _outvalue) * (1.0 + std::cos(int_pos * libMesh::pi)) / 2.0;
121  }
122 
123  return value;
124 }
virtual Real value(const Point &p)
void PolycrystalVoronoiVoidIC::computeGrainCenters ( )
protectedvirtual

Definition at line 275 of file PolycrystalVoronoiVoidIC.C.

Referenced by initialSetup().

276 {
277  if (_op_num > _grain_num)
278  mooseError("ERROR in PolycrystalVoronoiVoidIC: Number of order parameters "
279  "(op_num) can't be "
280  "larger than the number of grains (grain_num)");
281 
282  // Initialize vectors
283  _centerpoints.resize(_grain_num);
284  _assigned_op.resize(_grain_num);
285 
286  // Randomly generate the centers of the individual grains represented by the
287  // Voronoi tessellation
288  for (unsigned int grain = 0; grain < _grain_num; grain++)
289  {
290  for (unsigned int i = 0; i < LIBMESH_DIM; i++)
291  _centerpoints[grain](i) = _bottom_left(i) + _range(i) * MooseRandom::rand();
292 
293  if (_columnar_3D)
294  _centerpoints[grain](2) = _bottom_left(2) + _range(2) * 0.5;
295  }
296 
297  // Assign grains to specific order parameters in a way that maximizes the
298  // distance
300 }
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
RealGradient PolycrystalVoronoiVoidIC::gradient ( const Point &  p)
overrideprotectedvirtual

Reimplemented from SmoothCircleBaseIC.

Definition at line 235 of file PolycrystalVoronoiVoidIC.C.

236 {
237  RealGradient gradient;
238  RealGradient void_gradient = MultiSmoothCircleIC::gradient(p);
239 
240  // Order parameter assignment assumes zero gradient (sharp interface)
241  switch (_structure_type)
242  {
243  case 1: // assigning gradient for voids
244  gradient = void_gradient;
245  break;
246  }
247 
248  return gradient;
249 }
virtual RealGradient gradient(const Point &p)
virtual RealGradient gradient(const Point &p) override
Real PolycrystalVoronoiVoidIC::grainValueCalc ( const Point &  p)
protectedvirtual

Definition at line 252 of file PolycrystalVoronoiVoidIC.C.

Referenced by value().

253 {
254  Real val = 0.0;
255 
256  unsigned int min_index =
258 
259  // If the current order parameter index (_op_index) is equal to the min_index,
260  // set the value to
261  // 1.0
262  if (_assigned_op[min_index] == _op_index)
263  val = 1.0;
264 
265  if (val > 1.0)
266  val = 1.0;
267 
268  if (val < 0.0)
269  val = 0.0;
270 
271  return val;
272 }
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
void PolycrystalVoronoiVoidIC::initialSetup ( )
overridevirtual

Reimplemented from MultiSmoothCircleIC.

Definition at line 68 of file PolycrystalVoronoiVoidIC.C.

69 {
70  if (_op_num <= _op_index)
71  mooseError("op_index is too large in CircleGrainVoidIC");
72 
73  MooseRandom::seed(getParam<unsigned int>("rand_seed"));
74  // Set up domain bounds with mesh tools
75  for (unsigned int i = 0; i < LIBMESH_DIM; i++)
76  {
77  _bottom_left(i) = _mesh.getMinInDimension(i);
78  _top_right(i) = _mesh.getMaxInDimension(i);
79  }
81 
82  // Create _centerpoints and _assigned_op vectors
84 
85  // Call initial setup from MultiSmoothCircleIC to create _centers and _radii
86  // for voids
88 }
virtual void initialSetup()
Real PolycrystalVoronoiVoidIC::value ( const Point &  p)
overrideprotectedvirtual

Reimplemented from SmoothCircleBaseIC.

Definition at line 202 of file PolycrystalVoronoiVoidIC.C.

203 {
204  Real value = 0.0;
205 
206  // Determine value for voids
207  Real void_value = MultiSmoothCircleIC::value(p);
208 
209  // Determine value for grains
210  Real grain_value = grainValueCalc(p);
211 
212  switch (_structure_type)
213  {
214  case 0: // assigning values for grains (order parameters)
215  if (grain_value == 0) // Not in this grain
216  value = grain_value;
217  else // in this grain, but might be in a void
218  if (void_value == _outvalue) // Not in a void
219  value = grain_value;
220  else if (void_value > _outvalue && void_value < _invalue) // On void interface
221  value = 1.0 - (void_value - _outvalue) / (_invalue - _outvalue);
222  else if (void_value == _invalue) // In a void, so op = 0
223  value = 0.0;
224  break;
225 
226  case 1: // assigning values for voids (concentration)
227  value = void_value;
228  break;
229  }
230 
231  return value;
232 }
virtual Real value(const Point &p)
virtual Real value(const Point &p) override
virtual Real grainValueCalc(const Point &p)

Member Data Documentation

bool SmoothCircleBaseIC::_3D_spheres
protectedinherited
std::vector<unsigned int> PolycrystalVoronoiVoidIC::_assigned_op
protected

Definition at line 53 of file PolycrystalVoronoiVoidIC.h.

Referenced by computeGrainCenters(), and grainValueCalc().

Point MultiSmoothCircleIC::_bottom_left
protectedinherited
const Real MultiSmoothCircleIC::_bubspac
protectedinherited
std::vector<Point> PolycrystalVoronoiVoidIC::_centerpoints
protected
std::vector<Point> SmoothCircleBaseIC::_centers
protectedinherited
const bool PolycrystalVoronoiVoidIC::_columnar_3D
protected

Definition at line 42 of file PolycrystalVoronoiVoidIC.h.

Referenced by computeGrainCenters().

struct PolycrystalVoronoiVoidIC::DistancePointComparator PolycrystalVoronoiVoidIC::_customLess
protected
const unsigned int PolycrystalVoronoiVoidIC::_grain_num
protected

Definition at line 37 of file PolycrystalVoronoiVoidIC.h.

Referenced by computeCircleCenters(), and computeGrainCenters().

Real SmoothCircleBaseIC::_int_width
protectedinherited
Real SmoothCircleBaseIC::_invalue
protectedinherited
const unsigned int MultiSmoothCircleIC::_max_num_tries
protectedinherited
MooseMesh& SmoothCircleBaseIC::_mesh
protectedinherited
unsigned int SmoothCircleBaseIC::_num_dim
protectedinherited

Definition at line 50 of file SmoothCircleBaseIC.h.

const unsigned int MultiSmoothCircleIC::_numbub
protectedinherited
const unsigned int PolycrystalVoronoiVoidIC::_op_index
protected

Definition at line 38 of file PolycrystalVoronoiVoidIC.h.

Referenced by grainValueCalc(), and initialSetup().

const unsigned int PolycrystalVoronoiVoidIC::_op_num
protected

Definition at line 36 of file PolycrystalVoronoiVoidIC.h.

Referenced by computeGrainCenters(), and initialSetup().

Real SmoothCircleBaseIC::_outvalue
protectedinherited
std::vector<Real> SmoothCircleBaseIC::_radii
protectedinherited
const Real MultiSmoothCircleIC::_radius
protectedinherited
const Real MultiSmoothCircleIC::_radius_variation
protectedinherited
const MooseEnum MultiSmoothCircleIC::_radius_variation_type
protectedinherited
const unsigned int PolycrystalVoronoiVoidIC::_rand_seed
protected

Definition at line 40 of file PolycrystalVoronoiVoidIC.h.

MooseRandom SmoothCircleBaseIC::_random
protectedinherited
Point MultiSmoothCircleIC::_range
protectedinherited
const MooseEnum PolycrystalVoronoiVoidIC::_structure_type
protected

Definition at line 34 of file PolycrystalVoronoiVoidIC.h.

Referenced by gradient(), and value().

Point MultiSmoothCircleIC::_top_right
protectedinherited
bool SmoothCircleBaseIC::_zero_gradient
protectedinherited

Definition at line 48 of file SmoothCircleBaseIC.h.

Referenced by SmoothCircleBaseIC::gradient().


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