www.mooseframework.org
ComputeReducedOrderEigenstrain.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 /****************************************************************/
7 
9 
10 #include "MooseMesh.h"
11 #include "libmesh/quadrature.h"
12 
13 template <>
14 InputParameters
16 {
17  InputParameters params = validParams<ComputeEigenstrainBase>();
18  params.addRequiredParam<std::vector<MaterialPropertyName>>(
19  "input_eigenstrain_names", "List of eigenstrains to be applied in this strain calculation");
20  return params;
21 }
22 
24  : ComputeEigenstrainBase(parameters),
25  _input_eigenstrain_names(
26  getParam<std::vector<MaterialPropertyName>>("input_eigenstrain_names")),
27  _eigenstrains(_input_eigenstrain_names.size()),
28  _eigenstrains_old(_input_eigenstrain_names.size()),
29  _subproblem(*parameters.get<SubProblem *>("_subproblem")),
30  _ncols(1 + _subproblem.mesh().dimension()),
31  _second_order(_subproblem.mesh().hasSecondOrderElements()),
32  _eigsum(),
33  _A(),
34  _b(6),
35  _AT(),
36  _ATb(_ncols),
37  _x(6, DenseVector<Real>(_ncols)),
38  _vals(6),
39  _adjusted_eigenstrain()
40 {
41  for (unsigned int i = 0; i < _eigenstrains.size(); ++i)
42  {
44  _eigenstrains[i] = &getMaterialProperty<RankTwoTensor>(_input_eigenstrain_names[i]);
46  _eigenstrains_old[i] = &getMaterialPropertyOld<RankTwoTensor>(_input_eigenstrain_names[i]);
47  }
48 }
49 
50 void
52 {
53  _eigenstrain[_qp].zero();
54 }
55 
56 void
58 {
60 
62 
63  Material::computeProperties();
64 }
65 
66 void
68 {
69  if (_second_order)
70  {
71  for (unsigned i = 0; i < 6; ++i)
72  {
73  _vals[i] = _x[i](0);
74  for (unsigned j = 1; j < _ncols; ++j)
75  _vals[i] += _x[i](j) * _q_point[_qp](j - 1);
76  }
77  _adjusted_eigenstrain.fillFromInputVector(_vals);
78  }
80 }
81 
82 void
84 {
85  _eigsum.resize(_qrule->n_points());
86  for (_qp = 0; _qp < _qrule->n_points(); ++_qp)
87  {
88  _eigsum[_qp].zero();
89  for (auto es : _eigenstrains)
90  _eigsum[_qp] += (*es)[_qp];
91  }
92 }
93 
94 void
96 {
97  // The eigenstrains can either be constant in an element or linear in x, y, z
98  // If constant, do volume averaging.
99  if (!_second_order || _qrule->n_points() == 1)
100  {
101  // Volume average
102  _adjusted_eigenstrain.zero();
103  Real vol = 0.0;
104  for (unsigned qp = 0; qp < _qrule->n_points(); ++qp)
105  {
106  _adjusted_eigenstrain += _eigsum[qp] * _JxW[qp] * _coord[qp];
107  vol += _JxW[qp] * _coord[qp];
108  }
109  _adjusted_eigenstrain /= vol;
110  }
111  else
112  {
113  // 1 x y z
114 
115  // Six sets (one for each unique component of eigen) of qp data
116  _A.resize(_qrule->n_points(), _ncols);
117  for (auto && b : _b)
118  b.resize(_qrule->n_points());
119 
120  for (unsigned qp = 0; qp < _qrule->n_points(); ++qp)
121  {
122  _A(qp, 0) = 1.0;
123  for (unsigned j = 1; j < _ncols; ++j)
124  _A(qp, j) = _q_point[qp](j - 1);
125 
126  _b[0](qp) = _eigsum[qp](0, 0);
127  _b[1](qp) = _eigsum[qp](1, 1);
128  _b[2](qp) = _eigsum[qp](2, 2);
129  _b[3](qp) = _eigsum[qp](1, 2);
130  _b[4](qp) = _eigsum[qp](0, 2);
131  _b[5](qp) = _eigsum[qp](0, 1);
132  }
133 
134  _A.get_transpose(_AT);
135  _A.left_multiply(_AT);
136  for (unsigned i = 0; i < 6; ++i)
137  {
138  _AT.vector_mult(_ATb, _b[i]);
139 
140  _A.cholesky_solve(_ATb, _x[i]);
141  }
142  }
143 }
InputParameters validParams< ComputeEigenstrainBase >()
std::vector< RankTwoTensor > _eigsum
The sum of all eigenstrains at each integration point.
std::vector< const MaterialProperty< RankTwoTensor > * > _eigenstrains_old
void prepareEigenstrain()
Compute either the volume average or linear eigenstrain field in an element.
RankTwoTensor _adjusted_eigenstrain
Filled with _vals and subracted from strain.
const unsigned _ncols
Number of columns in A matrix (1 plus mesh dimension)
std::vector< DenseVector< Real > > _b
The b array holding the unique eigenstrain components for each integration point. ...
DenseMatrix< Real > _A
The (num q points x ncols) array for the least squares. Holds 1, xcoor, ycoor, zcoor.
std::string _base_name
Base name prepended to material property name.
DenseVector< Real > _ATb
Transpose of A times b.
MaterialProperty< RankTwoTensor > & _eigenstrain
Stores the current total eigenstrain.
ComputeReducedOrderEigenstrain(const InputParameters &parameters)
const bool _second_order
Whether the mesh is made of second order elements.
std::vector< const MaterialProperty< RankTwoTensor > * > _eigenstrains
void sumEigenstrain()
Add contributions from every eigenstrain at each integration point.
std::vector< DenseVector< Real > > _x
The solution vector for each unique component of the adjusted eigenstrain.
InputParameters validParams< ComputeReducedOrderEigenstrain >()
DenseMatrix< Real > _AT
Transpose of A.
bool _incremental_form
Whether the eigenstrain model should compute the total or incremental eigenstrain.
std::vector< Real > _vals
Vector to hold the adjusted strain as computed with _x.
ComputeEigenstrainBase is the base class for eigenstrain tensors.
std::vector< MaterialPropertyName > _input_eigenstrain_names
void computeQpEigenstrain()
Compute the eigenstrain and store in _eigenstrain.