www.mooseframework.org
ComputeExternalGrainForceAndTorque.C
Go to the documentation of this file.
1 /****************************************************************/
2 /* MOOSE - Multiphysics Object Oriented Simulation Environment */
3 /* */
4 /* All contents are licensed under LGPL V2.1 */
5 /* See LICENSE for full restrictions */
6 /****************************************************************/
9 
10 #include "libmesh/quadrature.h"
11 
12 template <>
13 InputParameters
15 {
16  InputParameters params = validParams<ShapeElementUserObject>();
17  params.addClassDescription("Userobject for calculating force and torque acting on a grain");
18  params.addParam<MaterialPropertyName>("force_density", "force_density", "Force density material");
19  params.addParam<UserObjectName>("grain_data", "center of mass of grains");
20  params.addCoupledVar("c", "Concentration field");
21  params.addCoupledVar("etas", "Array of coupled order parameters");
22  return params;
23 }
24 
26  const InputParameters & parameters)
27  : DerivativeMaterialInterface<ShapeElementUserObject>(parameters),
29  _c_name(getVar("c", 0)->name()),
30  _c_var(coupled("c")),
31  _dF_name(getParam<MaterialPropertyName>("force_density")),
32  _dF(getMaterialPropertyByName<std::vector<RealGradient>>(_dF_name)),
33  _dFdc(
34  getMaterialPropertyByName<std::vector<RealGradient>>(propertyNameFirst(_dF_name, _c_name))),
35  _op_num(coupledComponents("etas")),
36  _grain_tracker(getUserObject<GrainTrackerInterface>("grain_data")),
37  _vals_var(_op_num),
38  _vals_name(_op_num),
39  _dFdeta(_op_num)
40 {
41  for (unsigned int i = 0; i < _op_num; ++i)
42  {
43  _vals_var[i] = coupled("etas", i);
44  _vals_name[i] = getVar("etas", i)->name();
45  _dFdeta[i] = &getMaterialPropertyByName<std::vector<RealGradient>>(
46  propertyNameFirst(_dF_name, _vals_name[i]));
47  }
48 }
49 
50 void
52 {
54  _ncomp = 6 * _grain_num;
55 
56  _force_values.resize(_grain_num);
57  _torque_values.resize(_grain_num);
58  _force_torque_store.assign(_ncomp, 0.0);
59 
60  if (_fe_problem.currentlyComputingJacobian())
61  {
62  _total_dofs = _subproblem.es().n_dofs();
65 
66  for (unsigned int i = 0; i < _op_num; ++i)
67  _force_torque_eta_jacobian_store[i].assign(_ncomp * _total_dofs, 0.0);
68  }
69 }
70 
71 void
73 {
74  const auto & op_to_grains = _grain_tracker.getVarToFeatureVector(_current_elem->id());
75 
76  for (unsigned int i = 0; i < _grain_num; ++i)
77  for (unsigned int j = 0; j < _op_num; ++j)
78  if (i == op_to_grains[j])
79  {
80  const auto centroid = _grain_tracker.getGrainCentroid(i);
81  for (_qp = 0; _qp < _qrule->n_points(); ++_qp)
82  if (_dF[_qp][j](0) != 0.0 || _dF[_qp][j](1) != 0.0 || _dF[_qp][j](2) != 0.0)
83  {
84  const RealGradient compute_torque =
85  _JxW[_qp] * _coord[_qp] * (_current_elem->centroid() - centroid).cross(_dF[_qp][j]);
86  _force_torque_store[6 * i + 0] += _JxW[_qp] * _coord[_qp] * _dF[_qp][j](0);
87  _force_torque_store[6 * i + 1] += _JxW[_qp] * _coord[_qp] * _dF[_qp][j](1);
88  _force_torque_store[6 * i + 2] += _JxW[_qp] * _coord[_qp] * _dF[_qp][j](2);
89  _force_torque_store[6 * i + 3] += compute_torque(0);
90  _force_torque_store[6 * i + 4] += compute_torque(1);
91  _force_torque_store[6 * i + 5] += compute_torque(2);
92  }
93  }
94 }
95 
96 void
98 {
99  const auto & op_to_grains = _grain_tracker.getVarToFeatureVector(_current_elem->id());
100 
101  if (jvar == _c_var)
102  for (unsigned int i = 0; i < _grain_num; ++i)
103  for (unsigned int j = 0; j < _op_num; ++j)
104  if (i == op_to_grains[j])
105  {
106  const auto centroid = _grain_tracker.getGrainCentroid(i);
107  for (_qp = 0; _qp < _qrule->n_points(); ++_qp)
108  if (_dFdc[_qp][j](0) != 0.0 || _dFdc[_qp][j](1) != 0.0 || _dFdc[_qp][j](2) != 0.0)
109  {
110  const Real factor = _JxW[_qp] * _coord[_qp] * _phi[_j][_qp];
111  const RealGradient compute_torque_jacobian_c =
112  factor * (_current_elem->centroid() - centroid).cross(_dFdc[_qp][j]);
113  _force_torque_c_jacobian_store[(6 * i + 0) * _total_dofs + _j_global] +=
114  factor * _dFdc[_qp][j](0);
115  _force_torque_c_jacobian_store[(6 * i + 1) * _total_dofs + _j_global] +=
116  factor * _dFdc[_qp][j](1);
117  _force_torque_c_jacobian_store[(6 * i + 2) * _total_dofs + _j_global] +=
118  factor * _dFdc[_qp][j](2);
119  _force_torque_c_jacobian_store[(6 * i + 3) * _total_dofs + _j_global] +=
120  compute_torque_jacobian_c(0);
121  _force_torque_c_jacobian_store[(6 * i + 4) * _total_dofs + _j_global] +=
122  compute_torque_jacobian_c(1);
123  _force_torque_c_jacobian_store[(6 * i + 5) * _total_dofs + _j_global] +=
124  compute_torque_jacobian_c(2);
125  }
126  }
127 
128  for (unsigned int i = 0; i < _op_num; ++i)
129  if (jvar == _vals_var[i])
130  for (unsigned int j = 0; j < _grain_num; ++j)
131  for (unsigned int k = 0; k < _op_num; ++k)
132  if (j == op_to_grains[k])
133  {
134  const auto centroid = _grain_tracker.getGrainCentroid(j);
135  for (_qp = 0; _qp < _qrule->n_points(); ++_qp)
136  if ((*_dFdeta[i])[_qp][j](0) != 0.0 || (*_dFdeta[i])[_qp][j](1) != 0.0 ||
137  (*_dFdeta[i])[_qp][j](2) != 0.0)
138  {
139  const Real factor = _JxW[_qp] * _coord[_qp] * _phi[_j][_qp];
140  const RealGradient compute_torque_jacobian_eta =
141  factor * (_current_elem->centroid() - centroid).cross((*_dFdeta[i])[_qp][k]);
142  _force_torque_eta_jacobian_store[i][(6 * j + 0) * _total_dofs + _j_global] +=
143  factor * (*_dFdeta[i])[_qp][k](0);
144  _force_torque_eta_jacobian_store[i][(6 * j + 1) * _total_dofs + _j_global] +=
145  factor * (*_dFdeta[i])[_qp][k](1);
146  _force_torque_eta_jacobian_store[i][(6 * j + 2) * _total_dofs + _j_global] +=
147  factor * (*_dFdeta[i])[_qp][k](2);
148  _force_torque_eta_jacobian_store[i][(6 * j + 3) * _total_dofs + _j_global] +=
149  compute_torque_jacobian_eta(0);
150  _force_torque_eta_jacobian_store[i][(6 * j + 4) * _total_dofs + _j_global] +=
151  compute_torque_jacobian_eta(1);
152  _force_torque_eta_jacobian_store[i][(6 * j + 5) * _total_dofs + _j_global] +=
153  compute_torque_jacobian_eta(2);
154  }
155  }
156 }
157 
158 void
160 {
161  gatherSum(_force_torque_store);
162  for (unsigned int i = 0; i < _grain_num; ++i)
163  {
164  _force_values[i](0) = _force_torque_store[6 * i + 0];
165  _force_values[i](1) = _force_torque_store[6 * i + 1];
166  _force_values[i](2) = _force_torque_store[6 * i + 2];
167  _torque_values[i](0) = _force_torque_store[6 * i + 3];
168  _torque_values[i](1) = _force_torque_store[6 * i + 4];
169  _torque_values[i](2) = _force_torque_store[6 * i + 5];
170  }
171 
172  if (_fe_problem.currentlyComputingJacobian())
173  {
175  for (unsigned int i = 0; i < _op_num; ++i)
176  gatherSum(_force_torque_eta_jacobian_store[i]);
177  }
178 }
179 
180 void
182 {
184  static_cast<const ComputeExternalGrainForceAndTorque &>(y);
185  for (unsigned int i = 0; i < _ncomp; ++i)
187  if (_fe_problem.currentlyComputingJacobian())
188  {
189  for (unsigned int i = 0; i < _ncomp * _total_dofs; ++i)
191  for (unsigned int i = 0; i < _op_num; ++i)
192  for (unsigned int j = 0; j < _ncomp * _total_dofs; ++j)
194  }
195 }
196 
197 const std::vector<RealGradient> &
199 {
200  return _force_values;
201 }
202 
203 const std::vector<RealGradient> &
205 {
206  return _torque_values;
207 }
208 
209 const std::vector<Real> &
211 {
213 }
214 const std::vector<std::vector<Real>> &
216 {
218 }
virtual Point getGrainCentroid(unsigned int grain_id) const =0
Returns the centroid for the given grain number.
const GrainTrackerInterface & _grain_tracker
provide UserObject for calculating grain volumes and centers
This class defines the interface for the GrainTracking objects.
const MaterialProperty< std::vector< RealGradient > > & _dFdc
material property that provides jacobian of force density with respect to c
std::vector< const MaterialProperty< std::vector< RealGradient > > * > _dFdeta
This class provides interface for extracting the forces and torques computed in other UserObjects...
virtual std::size_t getTotalFeatureCount() const =0
Returns a number large enough to contain the largest ID for all grains in use.
std::vector< RealGradient > _force_values
providing grain forces, torques and their jacobians w. r. t c
virtual const std::vector< std::vector< Real > > & getForceEtaJacobians() const
This class is here to get the force and torque acting on a grain.
ComputeExternalGrainForceAndTorque(const InputParameters &parameters)
std::vector< std::vector< Real > > _force_torque_eta_jacobian_store
MaterialPropertyName _dF_name
material property that provides force density
const unsigned int _op_num
no. of order parameters
std::vector< Real > _force_torque_c_jacobian_store
vector storing jacobian of grain force and torque values
virtual const std::vector< unsigned int > & getVarToFeatureVector(dof_id_type elem_id) const =0
Returns a list of active unique feature ids for a particular element.
virtual const std::vector< RealGradient > & getForceValues() const
virtual const std::vector< RealGradient > & getTorqueValues() const
const MaterialProperty< std::vector< RealGradient > > & _dF
InputParameters validParams< ComputeExternalGrainForceAndTorque >()
virtual const std::vector< Real > & getForceCJacobians() const
std::vector< Real > _force_torque_store
vector storing grain force and torque values