www.mooseframework.org
PolycrystalVoronoiVoidIC.C
Go to the documentation of this file.
1 //* This file is part of the MOOSE framework
2 //* https://www.mooseframework.org
3 //*
4 //* All rights reserved, see COPYRIGHT for full restrictions
5 //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
6 //*
7 //* Licensed under LGPL 2.1, please see LICENSE for details
8 //* https://www.gnu.org/licenses/lgpl-2.1.html
9 
11 
12 // MOOSE includes
13 #include "MooseMesh.h"
14 #include "MooseVariable.h"
15 #include "DelimitedFileReader.h"
16 #include "GrainTrackerInterface.h"
17 #include "PolycrystalVoronoi.h"
18 
19 InputParameters
21 {
22  InputParameters params = ::validParams<MultiSmoothCircleIC>();
23 
24  params.addRequiredParam<unsigned int>("op_num", "Number of order parameters");
25 
26  params.addParam<bool>(
27  "columnar_3D", false, "3D microstructure will be columnar in the z-direction?");
28 
29  return params;
30 }
31 
33 
34 template <>
35 InputParameters
37 {
38  InputParameters params = PolycrystalVoronoiVoidIC::actionParameters();
39  MooseEnum structure_options("grains voids");
40  params.addRequiredParam<MooseEnum>("structure_type",
41  structure_options,
42  "Which structure type is being initialized, grains or voids");
43  params.addParam<unsigned int>("op_index",
44  0,
45  "The index for the current "
46  "order parameter, not needed if "
47  "structure_type = voids");
48  params.addRequiredParam<UserObjectName>(
49  "polycrystal_ic_uo", "UserObject for obtaining the polycrystal grain structure.");
50  params.addParam<FileName>(
51  "file_name",
52  "",
53  "File containing grain centroids, if file_name is provided, the centroids "
54  "from the file will be used.");
55  return params;
56 }
57 
58 PolycrystalVoronoiVoidIC::PolycrystalVoronoiVoidIC(const InputParameters & parameters)
59  : MultiSmoothCircleIC(parameters),
60  _structure_type(getParam<MooseEnum>("structure_type")),
61  _op_num(getParam<unsigned int>("op_num")),
62  _op_index(getParam<unsigned int>("op_index")),
63  _columnar_3D(getParam<bool>("columnar_3D")),
64  _poly_ic_uo(getUserObject<PolycrystalVoronoi>("polycrystal_ic_uo")),
65  _file_name(getParam<FileName>("file_name"))
66 {
67  if (_invalue < _outvalue)
68  mooseError("PolycrystalVoronoiVoidIC requires that the voids be "
69  "represented with invalue > outvalue");
70  if (_numbub == 0)
71  mooseError("PolycrystalVoronoiVoidIC requires numbub > 0. If you want no voids to "
72  "be "
73  "represented, use invalue = outvalue. In general, you should use "
74  "PolycrystalVoronoi to represent Voronoi grain structures without "
75  "voids.");
76 }
77 
78 void
80 {
81  if (_op_num <= _op_index)
82  mooseError("op_index is too large in CircleGrainVoidIC");
83 
84  // Obtain total number and centerpoints of the grains
87 
88  // Call initial setup from MultiSmoothCircleIC to create _centers and _radii
89  // for voids
91 }
92 
93 void
95 {
96  _centers.resize(_numbub);
97 
98  // This Code will place void center points on grain boundaries
99  for (unsigned int vp = 0; vp < _numbub; ++vp)
100  {
101  bool try_again;
102  unsigned int num_tries = 0;
103 
104  do
105  {
106  try_again = false;
107  num_tries++;
108 
109  if (num_tries > _max_num_tries)
110  mooseError("Too many tries of assigning void centers in "
111  "PolycrystalVoronoiVoidIC");
112 
113  Point rand_point;
114 
115  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
116  rand_point(i) = _bottom_left(i) + _range(i) * _random.rand(_tid);
117 
118  // Allow the vectors to be sorted based on their distance from the
119  // rand_point
120  std::vector<PolycrystalVoronoiVoidIC::DistancePoint> diff(_grain_num);
121 
122  for (unsigned int gr = 0; gr < _grain_num; ++gr)
123  {
124  diff[gr].d = _mesh.minPeriodicDistance(_var.number(), rand_point, _centerpoints[gr]);
125  diff[gr].gr = gr;
126  }
127 
128  std::sort(diff.begin(), diff.end(), _customLess);
129 
130  Point closest_point = _centerpoints[diff[0].gr];
131  Point next_closest_point = _centerpoints[diff[1].gr];
132 
133  // Find Slope of Line in the plane orthogonal to the diff_centerpoint
134  // vector
135  Point diff_centerpoints =
136  _mesh.minPeriodicVector(_var.number(), closest_point, next_closest_point);
137  Point diff_rand_center = _mesh.minPeriodicVector(_var.number(), closest_point, rand_point);
138  Point normal_vector = diff_centerpoints.cross(diff_rand_center);
139  Point slope = normal_vector.cross(diff_centerpoints);
140 
141  // Midpoint position vector between two center points
142  Point midpoint = closest_point + (0.5 * diff_centerpoints);
143 
144  // Solve for the scalar multiplier solution on the line
145  Real lambda = 0;
146  Point mid_rand_vector = _mesh.minPeriodicVector(_var.number(), midpoint, rand_point);
147 
148  Real slope_dot = slope * slope;
149  mooseAssert(slope_dot > 0, "The dot product of slope with itself is zero");
150  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
151  lambda += (mid_rand_vector(i) * slope(i)) / slope_dot;
152 
153  // Assigning points to vector
154  _centers[vp] = slope * lambda + midpoint;
155 
156  // Checking to see if points are in the domain ONLY WORKS FOR PERIODIC
157  for (unsigned int i = 0; i < LIBMESH_DIM; i++)
158  if ((_centers[vp](i) > _top_right(i)) || (_centers[vp](i) < _bottom_left(i)))
159  try_again = true;
160 
161  for (unsigned int i = 0; i < vp; ++i)
162  {
163  Real dist = _mesh.minPeriodicDistance(_var.number(), _centers[vp], _centers[i]);
164 
165  if (dist < _bubspac)
166  try_again = true;
167  }
168 
169  // Two algorithms are available for screening bubbles falling in grain
170  // interior. They produce
171  // nearly identical results.
172  // Here only one is listed. The other one is available upon request.
173 
174  // Use circle center for checking whether voids are at GBs
175  if (try_again == false)
176  {
177  Real min_rij_1, min_rij_2, rij, rij_diff_tol;
178 
179  min_rij_1 = _range.norm();
180  min_rij_2 = _range.norm();
181 
182  rij_diff_tol = 0.1 * _radius;
183 
184  for (unsigned int gr = 0; gr < _grain_num; ++gr)
185  {
186  rij = _mesh.minPeriodicDistance(_var.number(), _centers[vp], _centerpoints[gr]);
187 
188  if (rij < min_rij_1)
189  {
190  min_rij_2 = min_rij_1;
191  min_rij_1 = rij;
192  }
193  else if (rij < min_rij_2)
194  min_rij_2 = rij;
195  }
196 
197  if (std::abs(min_rij_1 - min_rij_2) > rij_diff_tol)
198  try_again = true;
199  }
200 
201  } while (try_again == true);
202  }
203 }
204 
205 Real
207 {
208  Real value = 0.0;
209 
210  // Determine value for voids
211  Real void_value = MultiSmoothCircleIC::value(p);
212 
213  // Determine value for grains
214  Real grain_value = _poly_ic_uo.getVariableValue(_op_index, p);
215 
216  switch (_structure_type)
217  {
218  case 0: // assigning values for grains (order parameters)
219  if (grain_value == 0) // Not in this grain
220  value = grain_value;
221  else // in this grain, but might be in a void
222  if (void_value == _outvalue) // Not in a void
223  value = grain_value;
224  else if (void_value > _outvalue && void_value < _invalue) // On void interface
225  value = grain_value * (_invalue - void_value) / (_invalue - _outvalue);
226  else if (void_value == _invalue) // In a void, so op = 0
227  value = 0.0;
228  break;
229 
230  case 1: // assigning values for voids (concentration)
231  value = void_value;
232  break;
233  }
234 
235  return value;
236 }
237 
240 {
242  RealGradient void_gradient = MultiSmoothCircleIC::gradient(p);
243 
244  // Order parameter assignment assumes zero gradient (sharp interface)
245  switch (_structure_type)
246  {
247  case 1: // assigning gradient for voids
248  gradient = void_gradient;
249  break;
250  }
251 
252  return gradient;
253 }
SmoothCircleBaseIC::_random
MooseRandom _random
Definition: SmoothCircleBaseIC.h:63
SmoothCircleBaseIC::_centers
std::vector< Point > _centers
Definition: SmoothCircleBaseIC.h:54
PolycrystalVoronoi::getNumGrains
virtual unsigned int getNumGrains() const override
Must be overridden by the deriving class to provide the number of grains in the polycrystal structure...
Definition: PolycrystalVoronoi.h:30
SmoothCircleBaseIC::value
virtual Real value(const Point &p)
Definition: SmoothCircleBaseIC.C:71
PolycrystalVoronoi::getVariableValue
virtual Real getVariableValue(unsigned int op_index, const Point &p) const override
Returns the variable value for a given op_index and mesh point.
Definition: PolycrystalVoronoi.C:94
PolycrystalVoronoiVoidIC::gradient
virtual RealGradient gradient(const Point &p) override
Definition: PolycrystalVoronoiVoidIC.C:239
libMesh::RealGradient
VectorValue< Real > RealGradient
Definition: GrainForceAndTorqueInterface.h:17
PolycrystalVoronoiVoidIC::actionParameters
static InputParameters actionParameters()
Definition: PolycrystalVoronoiVoidIC.C:20
MultiSmoothCircleIC::_max_num_tries
const unsigned int _max_num_tries
Definition: MultiSmoothCircleIC.h:41
MultiSmoothCircleIC::_bottom_left
Point _bottom_left
Definition: MultiSmoothCircleIC.h:47
PolycrystalVoronoiVoidIC
PolycrystalVoronoiVoidIC initializes either grain or void values for a voronoi tesselation with voids...
Definition: PolycrystalVoronoiVoidIC.h:28
PolycrystalVoronoiVoidIC::value
virtual Real value(const Point &p) override
Definition: PolycrystalVoronoiVoidIC.C:206
PolycrystalVoronoiVoidIC::_op_num
const unsigned int _op_num
Definition: PolycrystalVoronoiVoidIC.h:40
PolycrystalVoronoiVoidIC::initialSetup
virtual void initialSetup() override
Definition: PolycrystalVoronoiVoidIC.C:79
PolycrystalVoronoiVoidIC::computeCircleCenters
virtual void computeCircleCenters() override
Definition: PolycrystalVoronoiVoidIC.C:94
PolycrystalVoronoiVoidIC::_centerpoints
std::vector< Point > _centerpoints
Definition: PolycrystalVoronoiVoidIC.h:55
SmoothCircleBaseIC::_invalue
Real _invalue
Definition: SmoothCircleBaseIC.h:46
PolycrystalVoronoiVoidIC::_poly_ic_uo
const PolycrystalVoronoi & _poly_ic_uo
Definition: PolycrystalVoronoiVoidIC.h:45
PolycrystalVoronoiVoidIC.h
registerMooseObject
registerMooseObject("PhaseFieldApp", PolycrystalVoronoiVoidIC)
MultiSmoothCircleIC::_top_right
Point _top_right
Definition: MultiSmoothCircleIC.h:48
PolycrystalVoronoiVoidIC::_op_index
const unsigned int _op_index
Definition: PolycrystalVoronoiVoidIC.h:41
SmoothCircleBaseIC::_mesh
MooseMesh & _mesh
Definition: SmoothCircleBaseIC.h:44
SmoothCircleBaseIC::_outvalue
Real _outvalue
Definition: SmoothCircleBaseIC.h:47
GrainTrackerInterface.h
MultiSmoothCircleIC
MultismoothCircleIC creates multiple SmoothCircles (number = numbub) that are randomly positioned aro...
Definition: MultiSmoothCircleIC.h:27
PolycrystalVoronoiVoidIC::_structure_type
const MooseEnum _structure_type
Definition: PolycrystalVoronoiVoidIC.h:38
PolycrystalVoronoi
Definition: PolycrystalVoronoi.h:20
MultiSmoothCircleIC::_numbub
const unsigned int _numbub
Definition: MultiSmoothCircleIC.h:38
validParams< MultiSmoothCircleIC >
InputParameters validParams< MultiSmoothCircleIC >()
Definition: MultiSmoothCircleIC.C:20
PolycrystalVoronoiVoidIC::_customLess
struct PolycrystalVoronoiVoidIC::DistancePointComparator _customLess
PolycrystalVoronoi.h
SmoothCircleBaseIC::gradient
virtual RealGradient gradient(const Point &p)
Definition: SmoothCircleBaseIC.C:87
PolycrystalVoronoi::getGrainCenters
virtual std::vector< Point > getGrainCenters() const
Definition: PolycrystalVoronoi.h:31
validParams< PolycrystalVoronoiVoidIC >
InputParameters validParams< PolycrystalVoronoiVoidIC >()
Definition: PolycrystalVoronoiVoidIC.C:36
MultiSmoothCircleIC::_radius
const Real _radius
Definition: MultiSmoothCircleIC.h:43
MultiSmoothCircleIC::_range
Point _range
Definition: MultiSmoothCircleIC.h:49
PolycrystalVoronoiVoidIC::PolycrystalVoronoiVoidIC
PolycrystalVoronoiVoidIC(const InputParameters &parameters)
Definition: PolycrystalVoronoiVoidIC.C:58
PolycrystalVoronoiVoidIC::_grain_num
unsigned int _grain_num
Definition: PolycrystalVoronoiVoidIC.h:54
MultiSmoothCircleIC::_bubspac
const Real _bubspac
Definition: MultiSmoothCircleIC.h:39
MultiSmoothCircleIC::initialSetup
virtual void initialSetup() override
Definition: MultiSmoothCircleIC.C:53