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 "GeneralizedPlaneStrainOffDiagNOSPD.h"
11 : #include "MooseVariableScalar.h"
12 : #include "PeridynamicsMesh.h"
13 :
14 : registerMooseObject("PeridynamicsApp", GeneralizedPlaneStrainOffDiagNOSPD);
15 :
16 : InputParameters
17 83 : GeneralizedPlaneStrainOffDiagNOSPD::validParams()
18 : {
19 83 : InputParameters params = MechanicsBaseNOSPD::validParams();
20 83 : params.addClassDescription(
21 : "Class for calculating the off-diagonal Jacobian of the coupling between displacements (or "
22 : "temperature) with scalar out-of-plane strain for the generalized plane strain using the "
23 : "H1NOSPD formulation");
24 :
25 166 : params.addCoupledVar("scalar_out_of_plane_strain",
26 : "Scalar variable for strain in the out-of-plane direction");
27 :
28 83 : return params;
29 0 : }
30 :
31 62 : GeneralizedPlaneStrainOffDiagNOSPD::GeneralizedPlaneStrainOffDiagNOSPD(
32 62 : const InputParameters & parameters)
33 : : MechanicsBaseNOSPD(parameters),
34 62 : _scalar_out_of_plane_strain_var_num(coupledScalar("scalar_out_of_plane_strain"))
35 : {
36 : // Consistency check
37 62 : if (_disp_var.size() != 2)
38 0 : mooseError("GeneralizedPlaneStrain only works for two dimensional models!");
39 62 : }
40 :
41 : void
42 5880 : GeneralizedPlaneStrainOffDiagNOSPD::computeOffDiagJacobianScalar(unsigned int jvar_num)
43 : {
44 5880 : if (jvar_num == _scalar_out_of_plane_strain_var_num)
45 : {
46 5880 : prepare();
47 :
48 5880 : if (_var.number() == _disp_var[0]->number())
49 2646 : if (_use_full_jacobian)
50 1176 : computeDispFullOffDiagJacobianScalar(0, jvar_num);
51 : else
52 1470 : computeDispPartialOffDiagJacobianScalar(0, jvar_num);
53 3234 : else if (_var.number() == _disp_var[1]->number())
54 2646 : if (_use_full_jacobian)
55 1176 : computeDispFullOffDiagJacobianScalar(1, jvar_num);
56 : else
57 1470 : computeDispPartialOffDiagJacobianScalar(1, jvar_num);
58 588 : else if (_temp_coupled ? _var.number() == _temp_var->number() : 0)
59 588 : computeTempOffDiagJacobianScalar(jvar_num);
60 : }
61 5880 : }
62 :
63 : void
64 2940 : GeneralizedPlaneStrainOffDiagNOSPD::computeDispPartialOffDiagJacobianScalar(unsigned int component,
65 : unsigned int jvar_num)
66 : {
67 : // off-diagonal jacobian entries on the column and row corresponding to
68 : // scalar_out_of_plane_strain for coupling with displacements
69 :
70 2940 : prepareMatrixTag(_assembly, _var.number(), jvar_num, _ken);
71 2940 : prepareMatrixTag(_assembly, jvar_num, _var.number(), _kne);
72 2940 : MooseVariableScalar & jvar = _sys.getScalarVariable(_tid, jvar_num);
73 :
74 : // fill in the column corresponding to the scalar variable
75 2940 : std::vector<RankTwoTensor> dSdE33(_nnodes);
76 8820 : for (unsigned int nd = 0; nd < _nnodes; ++nd)
77 23520 : for (unsigned int i = 0; i < 3; ++i)
78 70560 : for (unsigned int j = 0; j < 3; ++j)
79 52920 : dSdE33[nd](i, j) = _Jacobian_mult[nd](i, j, 2, 2);
80 :
81 8820 : for (unsigned int i = 0; i < _nnodes; ++i)
82 11760 : for (unsigned int j = 0; j < jvar.order(); ++j)
83 8820 : _ken(i, j) += (i == j ? -1 : 1) *
84 11760 : (_multi[0] * (dSdE33[0] * _shape2[0].inverse()).row(component) +
85 17640 : _multi[1] * (dSdE33[1] * _shape2[1].inverse()).row(component)) *
86 5880 : _origin_vec * _bond_status;
87 :
88 2940 : _kne(0, 0) +=
89 2940 : computeDSDU(component, 0)(2, 2) * _node_vol[0] * _dg_vol_frac[0] * _bond_status; // node i
90 2940 : _kne(0, 1) +=
91 2940 : computeDSDU(component, 1)(2, 2) * _node_vol[1] * _dg_vol_frac[1] * _bond_status; // node j
92 2940 : accumulateTaggedLocalMatrix(_assembly, _var.number(), jvar_num, _ken);
93 2940 : accumulateTaggedLocalMatrix(_assembly, jvar_num, _var.number(), _kne);
94 2940 : }
95 :
96 : void
97 2352 : GeneralizedPlaneStrainOffDiagNOSPD::computeDispFullOffDiagJacobianScalar(unsigned int component,
98 : unsigned int jvar_num)
99 : {
100 : // LOCAL contribution
101 :
102 : // off-diagonal jacobian entries on the column and row corresponding to
103 : // scalar_out_of_plane_strain for coupling with displacements
104 2352 : prepareMatrixTag(_assembly, _var.number(), jvar_num, _ken);
105 2352 : prepareMatrixTag(_assembly, jvar_num, _var.number(), _kne);
106 2352 : MooseVariableScalar & jvar = _sys.getScalarVariable(_tid, jvar_num);
107 :
108 : // fill in the column corresponding to the scalar variable
109 2352 : std::vector<RankTwoTensor> dSdE33(_nnodes);
110 7056 : for (unsigned int nd = 0; nd < _nnodes; ++nd)
111 18816 : for (unsigned int i = 0; i < 3; ++i)
112 56448 : for (unsigned int j = 0; j < 3; ++j)
113 42336 : dSdE33[nd](i, j) = _Jacobian_mult[nd](i, j, 2, 2);
114 :
115 7056 : for (unsigned int i = 0; i < _nnodes; ++i)
116 9408 : for (unsigned int j = 0; j < jvar.order(); ++j)
117 7056 : _ken(i, j) += (i == j ? -1 : 1) *
118 9408 : (_multi[0] * (dSdE33[0] * _shape2[0].inverse()).row(component) +
119 14112 : _multi[1] * (dSdE33[1] * _shape2[1].inverse()).row(component)) *
120 4704 : _origin_vec * _bond_status;
121 :
122 : // fill in the row corresponding to the scalar variable
123 2352 : _kne(0, 0) +=
124 2352 : computeDSDU(component, 0)(2, 2) * _node_vol[0] * _dg_vol_frac[0] * _bond_status; // node i
125 2352 : _kne(0, 1) +=
126 2352 : computeDSDU(component, 1)(2, 2) * _node_vol[1] * _dg_vol_frac[1] * _bond_status; // node j
127 :
128 2352 : accumulateTaggedLocalMatrix(_assembly, _var.number(), jvar_num, _ken);
129 2352 : accumulateTaggedLocalMatrix(_assembly, jvar_num, _var.number(), _kne);
130 :
131 : // NONLOCAL contribution
132 :
133 : // fill in the row corresponding to the scalar variable
134 7056 : for (unsigned int nd = 0; nd < _nnodes; ++nd)
135 : {
136 : // calculation of jacobian contribution to current_node's neighbors
137 : // NOT including the contribution to nodes i and j, which is considered as local off-diagonal
138 4704 : std::vector<dof_id_type> ivardofs(_nnodes);
139 4704 : ivardofs[0] = _current_elem->node_ptr(nd)->dof_number(_sys.number(), _var.number(), 0);
140 4704 : std::vector<dof_id_type> neighbors = _pdmesh.getNeighbors(_current_elem->node_id(nd));
141 4704 : std::vector<dof_id_type> bonds = _pdmesh.getBonds(_current_elem->node_id(nd));
142 :
143 : dof_id_type nb_index =
144 4704 : std::find(neighbors.begin(), neighbors.end(), _current_elem->node_id(1 - nd)) -
145 4704 : neighbors.begin();
146 : std::vector<dof_id_type> dg_neighbors =
147 4704 : _pdmesh.getBondDeformationGradientNeighbors(_current_elem->node_id(nd), nb_index);
148 :
149 : Real vol_nb;
150 : RealGradient origin_vec_nb;
151 4704 : RankTwoTensor dFdUk, dPdUk;
152 :
153 63552 : for (unsigned int nb = 0; nb < dg_neighbors.size(); ++nb)
154 58848 : if (_bond_status_var->getElementalValue(_pdmesh.elemPtr(bonds[dg_neighbors[nb]])) > 0.5)
155 : {
156 56656 : Node * dgneighbor_nb = _pdmesh.nodePtr(neighbors[dg_neighbors[nb]]);
157 56656 : ivardofs[1] = dgneighbor_nb->dof_number(_sys.number(), _var.number(), 0);
158 56656 : vol_nb = _pdmesh.getNodeVolume(neighbors[dg_neighbors[nb]]);
159 :
160 : // obtain bond nb's origin vector
161 56656 : origin_vec_nb = _pdmesh.getNodeCoord(dgneighbor_nb->id()) -
162 113312 : _pdmesh.getNodeCoord(_current_elem->node_id(nd));
163 :
164 : dFdUk.zero();
165 169968 : for (unsigned int i = 0; i < _dim; ++i)
166 113312 : dFdUk(component, i) =
167 113312 : _horizon_radius[nd] / origin_vec_nb.norm() * origin_vec_nb(i) * vol_nb;
168 :
169 56656 : dFdUk *= _shape2[nd].inverse();
170 113312 : dPdUk = _Jacobian_mult[nd] * 0.5 * (dFdUk.transpose() + dFdUk);
171 :
172 56656 : _local_ke.resize(_ken.n(), _ken.m());
173 : _local_ke.zero();
174 56656 : _local_ke(0, 1) = dPdUk(2, 2) * _dg_vol_frac[nd] * _node_vol[nd] * _bond_status;
175 :
176 56656 : addJacobian(_assembly, _local_ke, jvar.dofIndices(), ivardofs, _var.scalingFactor());
177 : }
178 4704 : }
179 2352 : }
180 :
181 : void
182 588 : GeneralizedPlaneStrainOffDiagNOSPD::computeTempOffDiagJacobianScalar(unsigned int jvar_num)
183 : {
184 : // off-diagonal jacobian entries on the row corresponding to scalar_out_of_plane_strain for
185 : // coupling with temperature
186 588 : prepareMatrixTag(_assembly, jvar_num, _var.number(), _kne);
187 :
188 : // one-way coupling between the scalar_out_of_plane_strain and temperature. fill in the row
189 : // corresponding to the scalar_out_of_plane_strain
190 588 : std::vector<RankTwoTensor> dSdT(_nnodes);
191 1764 : for (unsigned int nd = 0; nd < _nnodes; ++nd)
192 2352 : for (unsigned int es = 0; es < _deigenstrain_dT.size(); ++es)
193 1176 : dSdT[nd] = -_Jacobian_mult[nd] * (*_deigenstrain_dT[es])[nd];
194 :
195 588 : _kne(0, 0) += dSdT[0](2, 2) * _dg_vol_frac[0] * _node_vol[0] * _bond_status; // node i
196 588 : _kne(0, 1) += dSdT[1](2, 2) * _dg_vol_frac[1] * _node_vol[1] * _bond_status; // node j
197 588 : accumulateTaggedLocalMatrix(_assembly, jvar_num, _var.number(), _kne);
198 588 : }
|