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 : #include "MechanicsFiniteStrainBaseNOSPD.h"
11 :
12 : InputParameters
13 17 : MechanicsFiniteStrainBaseNOSPD::validParams()
14 : {
15 17 : InputParameters params = MechanicsBaseNOSPD::validParams();
16 17 : params.addClassDescription("Base class for kernels of the stabilized non-ordinary "
17 : "state-based peridynamic correspondence models");
18 :
19 34 : params.addParam<std::vector<MaterialPropertyName>>(
20 : "eigenstrain_names",
21 : {},
22 : "List of eigenstrains to be coupled in non-ordinary state-based mechanics kernels");
23 :
24 17 : return params;
25 0 : }
26 :
27 12 : MechanicsFiniteStrainBaseNOSPD::MechanicsFiniteStrainBaseNOSPD(const InputParameters & parameters)
28 : : MechanicsBaseNOSPD(parameters),
29 12 : _dgrad_old(getMaterialPropertyOld<RankTwoTensor>("deformation_gradient")),
30 24 : _E_inc(getMaterialProperty<RankTwoTensor>("strain_increment")),
31 36 : _R_inc(getMaterialProperty<RankTwoTensor>("rotation_increment"))
32 : {
33 12 : }
34 :
35 : RankTwoTensor
36 36720 : MechanicsFiniteStrainBaseNOSPD::computeDSDU(unsigned int component, unsigned int nd)
37 : {
38 : // compute the derivative of stress w.r.t the solution components for finite strain
39 36720 : RankTwoTensor dSdU;
40 :
41 : // fetch the derivative of stress w.r.t the Fhat
42 36720 : RankFourTensor DSDFhat = computeDSDFhat(nd);
43 :
44 : // third calculate derivative of Fhat w.r.t solution components
45 36720 : RankTwoTensor Tp3;
46 36720 : if (component == 0)
47 18360 : Tp3 = _dgrad_old[nd].inverse() * _ddgraddu[nd];
48 18360 : else if (component == 1)
49 18360 : Tp3 = _dgrad_old[nd].inverse() * _ddgraddv[nd];
50 0 : else if (component == 2)
51 0 : Tp3 = _dgrad_old[nd].inverse() * _ddgraddw[nd];
52 :
53 : // assemble the fetched and calculated quantities to form the derivative of Cauchy stress w.r.t
54 : // solution components
55 36720 : dSdU = DSDFhat * Tp3;
56 :
57 36720 : return dSdU;
58 : }
59 :
60 : RankFourTensor
61 36720 : MechanicsFiniteStrainBaseNOSPD::computeDSDFhat(unsigned int nd)
62 : {
63 : // compute the derivative of stress w.r.t the Fhat for finite strain
64 36720 : RankTwoTensor I(RankTwoTensor::initIdentity);
65 36720 : RankFourTensor dSdFhat;
66 36720 : dSdFhat.zero();
67 :
68 : // first calculate the derivative of incremental Cauchy stress w.r.t the inverse of Fhat
69 : // Reference: M. M. Rashid (1993), Incremental Kinematics for finite element applications, IJNME
70 36720 : RankTwoTensor S_inc = _Jacobian_mult[nd] * _E_inc[nd];
71 36720 : RankFourTensor Tp1;
72 36720 : Tp1.zero();
73 146880 : for (unsigned int i = 0; i < 3; ++i)
74 440640 : for (unsigned int j = 0; j < 3; ++j)
75 1321920 : for (unsigned int k = 0; k < 3; ++k)
76 3965760 : for (unsigned int l = 0; l < 3; ++l)
77 11897280 : for (unsigned int m = 0; m < 3; ++m)
78 35691840 : for (unsigned int n = 0; n < 3; ++n)
79 107075520 : for (unsigned int r = 0; r < 3; ++r)
80 80306640 : Tp1(i, j, k, l) +=
81 80306640 : S_inc(m, n) *
82 80306640 : (_R_inc[nd](j, n) * (0.5 * I(k, m) * I(i, l) - I(m, l) * _R_inc[nd](i, k) +
83 80306640 : 0.5 * _R_inc[nd](i, k) * _R_inc[nd](m, l)) +
84 80306640 : _R_inc[nd](i, m) * (0.5 * I(k, n) * I(j, l) - I(n, l) * _R_inc[nd](j, k) +
85 80306640 : 0.5 * _R_inc[nd](j, k) * _R_inc[nd](n, l))) -
86 80306640 : _R_inc[nd](l, m) * _R_inc[nd](i, n) * _R_inc[nd](j, r) *
87 80306640 : _Jacobian_mult[nd](n, r, m, k);
88 :
89 : // second calculate derivative of inverse of Fhat w.r.t Fhat
90 : // d(inv(Fhat)_kl)/dFhat_mn = - inv(Fhat)_km * inv(Fhat)_nl
91 : // the bases are gk, gl, gm, gn, indictates the inverse rather than the inverse transpose
92 :
93 36720 : RankFourTensor Tp2;
94 36720 : Tp2.zero();
95 36720 : RankTwoTensor invFhat = (_dgrad[nd] * _dgrad_old[nd].inverse()).inverse();
96 146880 : for (unsigned int k = 0; k < 3; ++k)
97 440640 : for (unsigned int l = 0; l < 3; ++l)
98 1321920 : for (unsigned int m = 0; m < 3; ++m)
99 3965760 : for (unsigned int n = 0; n < 3; ++n)
100 2974320 : Tp2(k, l, m, n) += -invFhat(k, m) * invFhat(n, l);
101 :
102 : // assemble two calculated quantities to form the derivative of Cauchy stress w.r.t
103 : // Fhat
104 36720 : dSdFhat = Tp1 * Tp2;
105 :
106 36720 : return dSdFhat;
107 : }
108 :
109 : Real
110 36720 : MechanicsFiniteStrainBaseNOSPD::computeDJDU(unsigned int component, unsigned int nd)
111 : {
112 : // for finite formulation, compute the derivative of determinant of deformation gradient w.r.t the
113 : // solution components
114 : // dJ / du = dJ / dF_ij * dF_ij / du = J * inv(F)_ji * dF_ij / du
115 :
116 : Real dJdU = 0.0;
117 36720 : RankTwoTensor invF = _dgrad[nd].inverse();
118 36720 : Real detF = _dgrad[nd].det();
119 146880 : for (unsigned int i = 0; i < 3; ++i)
120 440640 : for (unsigned int j = 0; j < 3; ++j)
121 : {
122 330480 : if (component == 0)
123 165240 : dJdU += detF * invF(j, i) * _ddgraddu[nd](i, j);
124 165240 : else if (component == 1)
125 165240 : dJdU += detF * invF(j, i) * _ddgraddv[nd](i, j);
126 0 : else if (component == 2)
127 0 : dJdU += detF * invF(j, i) * _ddgraddw[nd](i, j);
128 : }
129 :
130 36720 : return dJdU;
131 : }
132 :
133 : RankTwoTensor
134 36720 : MechanicsFiniteStrainBaseNOSPD::computeDinvFTDU(unsigned int component, unsigned int nd)
135 : {
136 : // for finite formulation, compute the derivative of transpose of inverse of deformation gradient
137 : // w.r.t the solution components
138 : // d(inv(F)_ji)/du = d(inv(F)_ji)/dF_kl * dF_kl/du = - inv(F)_jk * inv(F)_li * dF_kl/du
139 : // the bases are gi, gj, gk, gl, indictates the inverse transpose rather than the inverse
140 :
141 36720 : RankTwoTensor dinvFTdU;
142 : dinvFTdU.zero();
143 36720 : RankTwoTensor invF = _dgrad[nd].inverse();
144 36720 : if (component == 0)
145 : {
146 18360 : dinvFTdU(0, 1) =
147 18360 : _ddgraddu[nd](0, 2) * _dgrad[nd](2, 1) - _ddgraddu[nd](0, 1) * _dgrad[nd](2, 2);
148 18360 : dinvFTdU(0, 2) =
149 18360 : _ddgraddu[nd](0, 1) * _dgrad[nd](1, 2) - _ddgraddu[nd](0, 2) * _dgrad[nd](1, 1);
150 18360 : dinvFTdU(1, 1) =
151 18360 : _ddgraddu[nd](0, 0) * _dgrad[nd](2, 2) - _ddgraddu[nd](0, 2) * _dgrad[nd](2, 0);
152 18360 : dinvFTdU(1, 2) =
153 18360 : _ddgraddu[nd](0, 2) * _dgrad[nd](1, 0) - _ddgraddu[nd](0, 0) * _dgrad[nd](1, 2);
154 18360 : dinvFTdU(2, 1) =
155 18360 : _ddgraddu[nd](0, 1) * _dgrad[nd](2, 0) - _ddgraddu[nd](0, 0) * _dgrad[nd](2, 1);
156 18360 : dinvFTdU(2, 2) =
157 18360 : _ddgraddu[nd](0, 0) * _dgrad[nd](1, 1) - _ddgraddu[nd](0, 1) * _dgrad[nd](1, 0);
158 : }
159 18360 : else if (component == 1)
160 : {
161 18360 : dinvFTdU(0, 0) =
162 18360 : _ddgraddv[nd](1, 1) * _dgrad[nd](2, 2) - _ddgraddv[nd](1, 2) * _dgrad[nd](2, 1);
163 18360 : dinvFTdU(0, 2) =
164 18360 : _ddgraddv[nd](1, 2) * _dgrad[nd](0, 1) - _ddgraddv[nd](0, 2) * _dgrad[nd](1, 1);
165 18360 : dinvFTdU(1, 0) =
166 18360 : _ddgraddv[nd](1, 2) * _dgrad[nd](2, 0) - _ddgraddv[nd](1, 0) * _dgrad[nd](2, 2);
167 18360 : dinvFTdU(1, 2) =
168 18360 : _ddgraddv[nd](1, 0) * _dgrad[nd](0, 2) - _ddgraddv[nd](1, 2) * _dgrad[nd](0, 0);
169 18360 : dinvFTdU(2, 0) =
170 18360 : _ddgraddv[nd](1, 0) * _dgrad[nd](2, 1) - _ddgraddv[nd](1, 1) * _dgrad[nd](2, 0);
171 18360 : dinvFTdU(2, 2) =
172 18360 : _ddgraddv[nd](1, 1) * _dgrad[nd](0, 0) - _ddgraddv[nd](1, 0) * _dgrad[nd](0, 1);
173 : }
174 0 : else if (component == 2)
175 : {
176 0 : dinvFTdU(0, 0) =
177 0 : _ddgraddw[nd](2, 2) * _dgrad[nd](1, 1) - _ddgraddw[nd](2, 1) * _dgrad[nd](1, 2);
178 0 : dinvFTdU(0, 1) =
179 0 : _ddgraddw[nd](2, 1) * _dgrad[nd](0, 2) - _ddgraddw[nd](2, 2) * _dgrad[nd](0, 1);
180 0 : dinvFTdU(1, 0) =
181 0 : _ddgraddw[nd](2, 0) * _dgrad[nd](1, 2) - _ddgraddw[nd](2, 2) * _dgrad[nd](1, 0);
182 0 : dinvFTdU(1, 1) =
183 0 : _ddgraddw[nd](2, 2) * _dgrad[nd](0, 0) - _ddgraddw[nd](2, 0) * _dgrad[nd](0, 2);
184 0 : dinvFTdU(2, 0) =
185 0 : _ddgraddw[nd](2, 1) * _dgrad[nd](1, 0) - _ddgraddw[nd](2, 0) * _dgrad[nd](1, 1);
186 0 : dinvFTdU(2, 1) =
187 0 : _ddgraddw[nd](2, 0) * _dgrad[nd](0, 1) - _ddgraddw[nd](2, 1) * _dgrad[nd](0, 0);
188 : }
189 :
190 36720 : dinvFTdU /= _dgrad[nd].det();
191 146880 : for (unsigned int i = 0; i < 3; ++i)
192 440640 : for (unsigned int j = 0; j < 3; ++j)
193 1321920 : for (unsigned int k = 0; k < 3; ++k)
194 3965760 : for (unsigned int l = 0; l < 3; ++l)
195 : {
196 2974320 : if (component == 0)
197 1487160 : dinvFTdU(i, j) -= invF(i, j) * invF(l, k) * _ddgraddu[nd](k, l);
198 1487160 : else if (component == 1)
199 1487160 : dinvFTdU(i, j) -= invF(i, j) * invF(l, k) * _ddgraddv[nd](k, l);
200 0 : else if (component == 2)
201 0 : dinvFTdU(i, j) -= invF(i, j) * invF(l, k) * _ddgraddw[nd](k, l);
202 : }
203 :
204 36720 : return dinvFTdU;
205 : }
|
