Line data Source code
1 : //* This file is part of the MOOSE framework
2 : //* https://mooseframework.inl.gov
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://mooseframework.inl.gov
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("SolidMechanicsApp", HomogenizationConstraint);
24 :
25 : InputParameters
26 296 : HomogenizationConstraint::validParams()
27 : {
28 296 : InputParameters params = ElementUserObject::validParams();
29 296 : params.addClassDescription("Calculate element contribution to the homogenization constraint "
30 : "depending on the homogenization constraint type.");
31 592 : params.addRequiredParam<MultiMooseEnum>(
32 : "constraint_types",
33 : Homogenization::constraintType,
34 : "Type of each constraint: strain, stress, or none. The types are specified in the "
35 : "column-major order, and there must be 9 entries in total.");
36 592 : params.addRequiredParam<std::vector<FunctionName>>(
37 : "targets", "Functions giving the targets to hit for constraint types that are not none.");
38 592 : params.addParam<bool>("large_kinematics", false, "Using large displacements?");
39 592 : params.addParam<std::string>("base_name", "Material property base name");
40 296 : return params;
41 0 : }
42 :
43 148 : HomogenizationConstraint::HomogenizationConstraint(const InputParameters & parameters)
44 : : ElementUserObject(parameters),
45 148 : _large_kinematics(getParam<bool>("large_kinematics")),
46 296 : _base_name(isParamValid("base_name") ? getParam<std::string>("base_name") + "_" : ""),
47 296 : _F(getMaterialPropertyByName<RankTwoTensor>(_base_name + "deformation_gradient")),
48 296 : _pk1_stress(getMaterialPropertyByName<RankTwoTensor>(_base_name + "pk1_stress")),
49 592 : _pk1_jacobian(getMaterialPropertyByName<RankFourTensor>(_base_name + "pk1_jacobian"))
50 : {
51 : // Constraint types
52 296 : auto types = getParam<MultiMooseEnum>("constraint_types");
53 148 : if (types.size() != Moose::dim * Moose::dim)
54 0 : mooseError("Number of constraint types must equal dim * dim. ", types.size(), " are provided.");
55 :
56 : // Targets to hit
57 148 : const std::vector<FunctionName> & fnames = getParam<std::vector<FunctionName>>("targets");
58 :
59 : // Prepare the constraint map
60 : unsigned int fcount = 0;
61 592 : for (const auto j : make_range(Moose::dim))
62 1776 : for (const auto i : make_range(Moose::dim))
63 : {
64 1332 : const auto idx = i + Moose::dim * j;
65 1332 : const auto ctype = static_cast<Homogenization::ConstraintType>(types.get(idx));
66 1332 : if (ctype != Homogenization::ConstraintType::None)
67 : {
68 762 : const Function * const f = &getFunctionByName(fnames[fcount++]);
69 762 : _cmap[{i, j}] = {ctype, f};
70 : }
71 : }
72 148 : }
73 :
74 : void
75 24126 : HomogenizationConstraint::initialize()
76 : {
77 : _residual.zero();
78 24126 : _jacobian.zero();
79 24126 : }
80 :
81 : void
82 1574252 : HomogenizationConstraint::execute()
83 : {
84 13889628 : for (_qp = 0; _qp < _qrule->n_points(); _qp++)
85 : {
86 12315376 : Real dV = _JxW[_qp] * _coord[_qp];
87 12315376 : _residual += computeResidual() * dV;
88 24630752 : _jacobian += computeJacobian() * dV;
89 : }
90 1574252 : }
91 :
92 : void
93 0 : HomogenizationConstraint::threadJoin(const UserObject & y)
94 : {
95 : const HomogenizationConstraint & other = static_cast<const HomogenizationConstraint &>(y);
96 0 : _residual += other._residual;
97 0 : _jacobian += other._jacobian;
98 0 : }
99 :
100 : void
101 24126 : HomogenizationConstraint::finalize()
102 : {
103 24126 : std::vector<Real> residual(&_residual(0, 0), &_residual(0, 0) + RankTwoTensor::N2);
104 24126 : std::vector<Real> jacobian(&_jacobian(0, 0, 0, 0), &_jacobian(0, 0, 0, 0) + RankFourTensor::N4);
105 :
106 : gatherSum(residual);
107 : gatherSum(jacobian);
108 :
109 24126 : std::copy(residual.begin(), residual.end(), &_residual(0, 0));
110 24126 : std::copy(jacobian.begin(), jacobian.end(), &_jacobian(0, 0, 0, 0));
111 24126 : }
112 :
113 : RankTwoTensor
114 12315376 : HomogenizationConstraint::computeResidual()
115 : {
116 12315376 : RankTwoTensor res;
117 :
118 101410784 : for (auto && [indices, constraint] : _cmap)
119 : {
120 : auto && [i, j] = indices;
121 : auto && [ctype, ctarget] = constraint;
122 :
123 89095408 : if (_large_kinematics)
124 : {
125 48306096 : if (ctype == Homogenization::ConstraintType::Stress)
126 15962016 : res(i, j) = _pk1_stress[_qp](i, j) - ctarget->value(_t, _q_point[_qp]);
127 32344080 : else if (ctype == Homogenization::ConstraintType::Strain)
128 56598544 : res(i, j) = _F[_qp](i, j) - (Real(i == j) + ctarget->value(_t, _q_point[_qp]));
129 : else
130 0 : mooseError("Unknown constraint type in the integral!");
131 : }
132 : else
133 : {
134 40789312 : if (ctype == Homogenization::ConstraintType::Stress)
135 20469408 : res(i, j) = _pk1_stress[_qp](i, j) - ctarget->value(_t, _q_point[_qp]);
136 20319904 : else if (ctype == Homogenization::ConstraintType::Strain)
137 20319904 : res(i, j) = 0.5 * (_F[_qp](i, j) + _F[_qp](j, i)) -
138 30440096 : (Real(i == j) + ctarget->value(_t, _q_point[_qp]));
139 : else
140 0 : mooseError("Unknown constraint type in the integral!");
141 : }
142 : }
143 :
144 12315376 : return res;
145 : }
146 :
147 : RankFourTensor
148 12315376 : HomogenizationConstraint::computeJacobian()
149 : {
150 12315376 : RankFourTensor res;
151 :
152 101410784 : for (auto && [indices1, constraint1] : _cmap)
153 : {
154 : auto && [i, j] = indices1;
155 : auto && ctype = constraint1.first;
156 765062624 : for (auto && indices2 : _cmap)
157 : {
158 : auto && [a, b] = indices2.first;
159 675967216 : if (ctype == Homogenization::ConstraintType::Stress)
160 265277632 : res(i, j, a, b) = _pk1_jacobian[_qp](i, j, a, b);
161 410689584 : else if (ctype == Homogenization::ConstraintType::Strain)
162 : {
163 410689584 : if (_large_kinematics)
164 545991440 : res(i, j, a, b) = Real(i == a && j == b);
165 : else
166 334045856 : res(i, j, a, b) = 0.5 * Real(i == a && j == b) + 0.5 * Real(i == b && j == a);
167 : }
168 : else
169 0 : mooseError("Unknown constraint type in Jacobian calculator!");
170 : }
171 : }
172 :
173 12315376 : return res;
174 : }
|
