Line data Source code
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 :
10 : //* This file is part of the MOOSE framework
11 : //* https://www.mooseframework.org
12 : //*
13 : //* All rights reserved, see COPYRIGHT for full restrictions
14 : //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
15 : //*
16 : //* Licensed under LGPL 2.1, please see LICENSE for details
17 : //* https://www.gnu.org/licenses/lgpl-2.1.html
18 :
19 : #include "HomogenizationConstraint.h"
20 :
21 : #include "Function.h"
22 :
23 : registerMooseObject("TensorMechanicsApp", HomogenizationConstraint);
24 :
25 : InputParameters
26 148 : HomogenizationConstraint::validParams()
27 : {
28 148 : InputParameters params = ElementUserObject::validParams();
29 296 : params.addRequiredParam<MultiMooseEnum>(
30 : "constraint_types",
31 : Homogenization::constraintType,
32 : "Type of each constraint: strain, stress, or none. The types are specified in the "
33 : "column-major order, and there must be 9 entries in total.");
34 296 : params.addRequiredParam<std::vector<FunctionName>>(
35 : "targets", "Functions giving the targets to hit for constraint types that are not none.");
36 296 : params.addParam<bool>("large_kinematics", false, "Using large displacements?");
37 296 : params.addParam<std::string>("base_name", "Material property base name");
38 148 : return params;
39 0 : }
40 :
41 74 : HomogenizationConstraint::HomogenizationConstraint(const InputParameters & parameters)
42 : : ElementUserObject(parameters),
43 74 : _large_kinematics(getParam<bool>("large_kinematics")),
44 148 : _base_name(isParamValid("base_name") ? getParam<std::string>("base_name") + "_" : ""),
45 148 : _F(getMaterialPropertyByName<RankTwoTensor>(_base_name + "deformation_gradient")),
46 148 : _pk1_stress(getMaterialPropertyByName<RankTwoTensor>(_base_name + "pk1_stress")),
47 296 : _pk1_jacobian(getMaterialPropertyByName<RankFourTensor>(_base_name + "pk1_jacobian"))
48 : {
49 : // Constraint types
50 148 : auto types = getParam<MultiMooseEnum>("constraint_types");
51 74 : if (types.size() != Moose::dim * Moose::dim)
52 0 : mooseError("Number of constraint types must equal dim * dim. ", types.size(), " are provided.");
53 :
54 : // Targets to hit
55 74 : const std::vector<FunctionName> & fnames = getParam<std::vector<FunctionName>>("targets");
56 :
57 : // Prepare the constraint map
58 : unsigned int fcount = 0;
59 296 : for (const auto j : make_range(Moose::dim))
60 888 : for (const auto i : make_range(Moose::dim))
61 : {
62 666 : const auto idx = i + Moose::dim * j;
63 666 : const auto ctype = static_cast<Homogenization::ConstraintType>(types.get(idx));
64 666 : if (ctype != Homogenization::ConstraintType::None)
65 : {
66 381 : const Function * const f = &getFunctionByName(fnames[fcount++]);
67 381 : _cmap[{i, j}] = {ctype, f};
68 : }
69 : }
70 74 : }
71 :
72 : void
73 12356 : HomogenizationConstraint::initialize()
74 : {
75 : _residual.zero();
76 12356 : _jacobian.zero();
77 12356 : }
78 :
79 : void
80 801654 : HomogenizationConstraint::execute()
81 : {
82 7044766 : for (_qp = 0; _qp < _qrule->n_points(); _qp++)
83 : {
84 6243112 : Real dV = _JxW[_qp] * _coord[_qp];
85 6243112 : _residual += computeResidual() * dV;
86 12486224 : _jacobian += computeJacobian() * dV;
87 : }
88 801654 : }
89 :
90 : void
91 0 : HomogenizationConstraint::threadJoin(const UserObject & y)
92 : {
93 : const HomogenizationConstraint & other = static_cast<const HomogenizationConstraint &>(y);
94 0 : _residual += other._residual;
95 0 : _jacobian += other._jacobian;
96 0 : }
97 :
98 : void
99 12356 : HomogenizationConstraint::finalize()
100 : {
101 12356 : std::vector<Real> residual(&_residual(0, 0), &_residual(0, 0) + RankTwoTensor::N2);
102 12356 : std::vector<Real> jacobian(&_jacobian(0, 0, 0, 0), &_jacobian(0, 0, 0, 0) + RankFourTensor::N4);
103 :
104 : gatherSum(residual);
105 : gatherSum(jacobian);
106 :
107 12356 : std::copy(residual.begin(), residual.end(), &_residual(0, 0));
108 12356 : std::copy(jacobian.begin(), jacobian.end(), &_jacobian(0, 0, 0, 0));
109 12356 : }
110 :
111 : RankTwoTensor
112 6243112 : HomogenizationConstraint::computeResidual()
113 : {
114 6243112 : RankTwoTensor res;
115 :
116 51307280 : for (auto && [indices, constraint] : _cmap)
117 : {
118 : auto && [i, j] = indices;
119 : auto && [ctype, ctarget] = constraint;
120 :
121 45064168 : if (_large_kinematics)
122 : {
123 24445608 : if (ctype == Homogenization::ConstraintType::Stress)
124 8073632 : res(i, j) = _pk1_stress[_qp](i, j) - ctarget->value(_t, _q_point[_qp]);
125 16371976 : else if (ctype == Homogenization::ConstraintType::Strain)
126 28650120 : res(i, j) = _F[_qp](i, j) - (Real(i == j) + ctarget->value(_t, _q_point[_qp]));
127 : else
128 0 : mooseError("Unknown constraint type in the integral!");
129 : }
130 : else
131 : {
132 20618560 : if (ctype == Homogenization::ConstraintType::Stress)
133 10358336 : res(i, j) = _pk1_stress[_qp](i, j) - ctarget->value(_t, _q_point[_qp]);
134 10260224 : else if (ctype == Homogenization::ConstraintType::Strain)
135 10260224 : res(i, j) = 0.5 * (_F[_qp](i, j) + _F[_qp](j, i)) -
136 15364256 : (Real(i == j) + ctarget->value(_t, _q_point[_qp]));
137 : else
138 0 : mooseError("Unknown constraint type in the integral!");
139 : }
140 : }
141 :
142 6243112 : return res;
143 : }
144 :
145 : RankFourTensor
146 6243112 : HomogenizationConstraint::computeJacobian()
147 : {
148 6243112 : RankFourTensor res;
149 :
150 51307280 : for (auto && [indices1, constraint1] : _cmap)
151 : {
152 : auto && [i, j] = indices1;
153 : auto && ctype = constraint1.first;
154 386693840 : for (auto && indices2 : _cmap)
155 : {
156 : auto && [a, b] = indices2.first;
157 341629672 : if (ctype == Homogenization::ConstraintType::Stress)
158 134065024 : res(i, j, a, b) = _pk1_jacobian[_qp](i, j, a, b);
159 207564648 : else if (ctype == Homogenization::ConstraintType::Strain)
160 : {
161 207564648 : if (_large_kinematics)
162 276156232 : res(i, j, a, b) = Real(i == a && j == b);
163 : else
164 168485216 : res(i, j, a, b) = 0.5 * Real(i == a && j == b) + 0.5 * Real(i == b && j == a);
165 : }
166 : else
167 0 : mooseError("Unknown constraint type in Jacobian calculator!");
168 : }
169 : }
170 :
171 6243112 : return res;
172 : }
|