www.mooseframework.org
ComputeGrainForceAndTorque.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 
25 ComputeGrainForceAndTorque::ComputeGrainForceAndTorque(const InputParameters & parameters)
26  : DerivativeMaterialInterface<ShapeElementUserObject>(parameters),
28  _c_name(getVar("c", 0)->name()),
29  _c_var(coupled("c")),
30  _dF_name(getParam<MaterialPropertyName>("force_density")),
31  _dF(getMaterialPropertyByName<std::vector<RealGradient>>(_dF_name)),
32  _dFdc(
33  getMaterialPropertyByName<std::vector<RealGradient>>(propertyNameFirst(_dF_name, _c_name))),
34  _op_num(coupledComponents("etas")),
35  _grain_tracker(getUserObject<GrainTrackerInterface>("grain_data")),
36  _vals_var(_op_num),
37  _vals_name(_op_num),
38  _dFdgradeta(_op_num)
39 {
40  for (unsigned int i = 0; i < _op_num; ++i)
41  {
42  _vals_var[i] = coupled("etas", i);
43  _vals_name[i] = getVar("etas", i)->name();
44  _dFdgradeta[i] =
45  &getMaterialPropertyByName<std::vector<Real>>(propertyNameFirst(_dF_name, _vals_name[i]));
46  }
47 }
48 
49 void
51 {
53  _ncomp = 6 * _grain_num;
54 
55  _force_values.resize(_grain_num);
56  _torque_values.resize(_grain_num);
57  _force_torque_store.assign(_ncomp, 0.0);
58 
59  if (_fe_problem.currentlyComputingJacobian())
60  {
61  _total_dofs = _subproblem.es().n_dofs();
64 
65  for (unsigned int i = 0; i < _op_num; ++i)
66  _force_torque_eta_jacobian_store[i].assign(_ncomp * _total_dofs, 0.0);
67  }
68 }
69 
70 void
72 {
73  const auto & op_to_grains = _grain_tracker.getVarToFeatureVector(_current_elem->id());
74 
75  for (unsigned int i = 0; i < _grain_num; ++i)
76  for (unsigned int j = 0; j < _op_num; ++j)
77  if (i == op_to_grains[j])
78  {
79  const auto centroid = _grain_tracker.getGrainCentroid(i);
80  for (_qp = 0; _qp < _qrule->n_points(); ++_qp)
81  if (_dF[_qp][j](0) != 0.0 || _dF[_qp][j](1) != 0.0 || _dF[_qp][j](2) != 0.0)
82  {
83  const RealGradient compute_torque =
84  _JxW[_qp] * _coord[_qp] * (_current_elem->centroid() - centroid).cross(_dF[_qp][j]);
85  _force_torque_store[6 * i + 0] += _JxW[_qp] * _coord[_qp] * _dF[_qp][j](0);
86  _force_torque_store[6 * i + 1] += _JxW[_qp] * _coord[_qp] * _dF[_qp][j](1);
87  _force_torque_store[6 * i + 2] += _JxW[_qp] * _coord[_qp] * _dF[_qp][j](2);
88  _force_torque_store[6 * i + 3] += compute_torque(0);
89  _force_torque_store[6 * i + 4] += compute_torque(1);
90  _force_torque_store[6 * i + 5] += compute_torque(2);
91  }
92  }
93 }
94 
95 void
97 {
98  const auto & op_to_grains = _grain_tracker.getVarToFeatureVector(_current_elem->id());
99 
100  if (jvar == _c_var)
101  for (unsigned int i = 0; i < _grain_num; ++i)
102  for (unsigned int j = 0; j < _op_num; ++j)
103  if (i == op_to_grains[j])
104  {
105  const auto centroid = _grain_tracker.getGrainCentroid(i);
106  for (_qp = 0; _qp < _qrule->n_points(); ++_qp)
107  if (_dFdc[_qp][j](0) != 0.0 || _dFdc[_qp][j](1) != 0.0 || _dFdc[_qp][j](2) != 0.0)
108  {
109  const Real factor = _JxW[_qp] * _coord[_qp] * _phi[_j][_qp];
110  const RealGradient compute_torque_jacobian_c =
111  factor * (_current_elem->centroid() - centroid).cross(_dFdc[_qp][j]);
112  _force_torque_c_jacobian_store[(6 * i + 0) * _total_dofs + _j_global] +=
113  factor * _dFdc[_qp][j](0);
114  _force_torque_c_jacobian_store[(6 * i + 1) * _total_dofs + _j_global] +=
115  factor * _dFdc[_qp][j](1);
116  _force_torque_c_jacobian_store[(6 * i + 2) * _total_dofs + _j_global] +=
117  factor * _dFdc[_qp][j](2);
118  _force_torque_c_jacobian_store[(6 * i + 3) * _total_dofs + _j_global] +=
119  compute_torque_jacobian_c(0);
120  _force_torque_c_jacobian_store[(6 * i + 4) * _total_dofs + _j_global] +=
121  compute_torque_jacobian_c(1);
122  _force_torque_c_jacobian_store[(6 * i + 5) * _total_dofs + _j_global] +=
123  compute_torque_jacobian_c(2);
124  }
125  }
126 
127  for (unsigned int i = 0; i < _op_num; ++i)
128  if (jvar == _vals_var[i])
129  for (unsigned int j = 0; j < _grain_num; ++j)
130  for (unsigned int k = 0; k < _op_num; ++k)
131  if (j == op_to_grains[k])
132  {
133  const auto centroid = _grain_tracker.getGrainCentroid(j);
134  for (_qp = 0; _qp < _qrule->n_points(); ++_qp)
135  if ((*_dFdgradeta[i])[_qp][j] != 0.0)
136  {
137  const Real factor = _JxW[_qp] * _coord[_qp] * (*_dFdgradeta[i])[_qp][k];
138  const RealGradient compute_torque_jacobian_eta =
139  factor * (_current_elem->centroid() - centroid).cross(_grad_phi[_j][_qp]);
140  _force_torque_eta_jacobian_store[i][(6 * j + 0) * _total_dofs + _j_global] +=
141  factor * _grad_phi[_j][_qp](0);
142  _force_torque_eta_jacobian_store[i][(6 * j + 1) * _total_dofs + _j_global] +=
143  factor * _grad_phi[_j][_qp](1);
144  _force_torque_eta_jacobian_store[i][(6 * j + 2) * _total_dofs + _j_global] +=
145  factor * _grad_phi[_j][_qp](2);
146  _force_torque_eta_jacobian_store[i][(6 * j + 3) * _total_dofs + _j_global] +=
147  compute_torque_jacobian_eta(0);
148  _force_torque_eta_jacobian_store[i][(6 * j + 4) * _total_dofs + _j_global] +=
149  compute_torque_jacobian_eta(1);
150  _force_torque_eta_jacobian_store[i][(6 * j + 5) * _total_dofs + _j_global] +=
151  compute_torque_jacobian_eta(2);
152  }
153  }
154 }
155 
156 void
158 {
159  gatherSum(_force_torque_store);
160  for (unsigned int i = 0; i < _grain_num; ++i)
161  {
162  _force_values[i](0) = _force_torque_store[6 * i + 0];
163  _force_values[i](1) = _force_torque_store[6 * i + 1];
164  _force_values[i](2) = _force_torque_store[6 * i + 2];
165  _torque_values[i](0) = _force_torque_store[6 * i + 3];
166  _torque_values[i](1) = _force_torque_store[6 * i + 4];
167  _torque_values[i](2) = _force_torque_store[6 * i + 5];
168  }
169 
170  if (_fe_problem.currentlyComputingJacobian())
171  {
173  for (unsigned int i = 0; i < _op_num; ++i)
174  gatherSum(_force_torque_eta_jacobian_store[i]);
175  }
176 }
177 
178 void
180 {
181  const ComputeGrainForceAndTorque & pps = static_cast<const ComputeGrainForceAndTorque &>(y);
182  for (unsigned int i = 0; i < _ncomp; ++i)
184  if (_fe_problem.currentlyComputingJacobian())
185  {
186  for (unsigned int i = 0; i < _ncomp * _total_dofs; ++i)
188  for (unsigned int i = 0; i < _op_num; ++i)
189  for (unsigned int j = 0; j < _ncomp * _total_dofs; ++j)
191  }
192 }
193 
194 const std::vector<RealGradient> &
196 {
197  return _force_values;
198 }
199 
200 const std::vector<RealGradient> &
202 {
203  return _torque_values;
204 }
205 
206 const std::vector<Real> &
208 {
210 }
211 const std::vector<std::vector<Real>> &
213 {
215 }
const GrainTrackerInterface & _grain_tracker
provide UserObject for calculating grain volumes and centers
ComputeGrainForceAndTorque(const InputParameters &parameters)
virtual Point getGrainCentroid(unsigned int grain_id) const =0
Returns the centroid for the given grain number.
std::vector< Real > _force_torque_c_jacobian_store
vector storing jacobian of grain force and torque values
This class defines the interface for the GrainTracking objects.
This class provides interface for extracting the forces and torques computed in other UserObjects...
std::vector< RealGradient > _force_values
providing grain forces, torques and their jacobians w. r. t c
virtual std::size_t getTotalFeatureCount() const =0
Returns a number large enough to contain the largest ID for all grains in use.
virtual void executeJacobian(unsigned int jvar)
std::vector< std::vector< Real > > _force_torque_eta_jacobian_store
const MaterialProperty< std::vector< RealGradient > > & _dF
virtual void threadJoin(const UserObject &y)
std::vector< RealGradient > _torque_values
std::vector< const MaterialProperty< std::vector< Real > > * > _dFdgradeta
InputParameters validParams< ComputeGrainForceAndTorque >()
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.
This class is here to get the force and torque acting on a grain.
virtual const std::vector< std::vector< Real > > & getForceEtaJacobians() const
virtual const std::vector< RealGradient > & getTorqueValues() const
virtual const std::vector< RealGradient > & getForceValues() const
std::vector< unsigned int > _vals_var
MaterialPropertyName _dF_name
material property that provides force density
unsigned int _op_num
no. of order parameters
const MaterialProperty< std::vector< RealGradient > > & _dFdc
material property that provides jacobian of force density with respect to c
std::vector< VariableName > _vals_name
virtual const std::vector< Real > & getForceCJacobians() const
std::vector< Real > _force_torque_store
vector storing grain force and torque values