All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
tractionDisplacementCorrectionFvPatchVectorField.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2021 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
29 
30 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
31 
32 namespace Foam
33 {
34 
35 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
36 
39 (
40  const fvPatch& p,
41  const DimensionedField<vector, volMesh>& iF
42 )
43 :
44  fixedGradientFvPatchVectorField(p, iF),
45  traction_(p.size(), Zero),
46  pressure_(p.size(), 0.0)
47 {
48  fvPatchVectorField::operator=(patchInternalField());
49  gradient() = Zero;
50 }
51 
52 
53 tractionDisplacementCorrectionFvPatchVectorField::
54 tractionDisplacementCorrectionFvPatchVectorField
55 (
56  const fvPatch& p,
57  const DimensionedField<vector, volMesh>& iF,
58  const dictionary& dict
59 )
60 :
61  fixedGradientFvPatchVectorField(p, iF),
62  traction_("traction", dict, p.size()),
63  pressure_("pressure", dict, p.size())
64 {
65  fvPatchVectorField::operator=(patchInternalField());
66  gradient() = Zero;
67 }
68 
69 
70 tractionDisplacementCorrectionFvPatchVectorField::
71 tractionDisplacementCorrectionFvPatchVectorField
72 (
73  const tractionDisplacementCorrectionFvPatchVectorField& tdpvf,
74  const fvPatch& p,
75  const DimensionedField<vector, volMesh>& iF,
76  const fvPatchFieldMapper& mapper
77 )
78 :
79  fixedGradientFvPatchVectorField(tdpvf, p, iF, mapper),
80  traction_(mapper(tdpvf.traction_)),
81  pressure_(mapper(tdpvf.pressure_))
82 {}
83 
84 
85 tractionDisplacementCorrectionFvPatchVectorField::
86 tractionDisplacementCorrectionFvPatchVectorField
87 (
88  const tractionDisplacementCorrectionFvPatchVectorField& tdpvf,
89  const DimensionedField<vector, volMesh>& iF
90 )
91 :
92  fixedGradientFvPatchVectorField(tdpvf, iF),
93  traction_(tdpvf.traction_),
94  pressure_(tdpvf.pressure_)
95 {}
96 
97 
98 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
99 
100 void tractionDisplacementCorrectionFvPatchVectorField::autoMap
101 (
102  const fvPatchFieldMapper& m
103 )
104 {
105  fixedGradientFvPatchVectorField::autoMap(m);
106  m(traction_, traction_);
107  m(pressure_, pressure_);
108 }
109 
110 
111 // Reverse-map the given fvPatchField onto this fvPatchField
112 void tractionDisplacementCorrectionFvPatchVectorField::rmap
113 (
114  const fvPatchVectorField& ptf,
115  const labelList& addr
116 )
117 {
118  fixedGradientFvPatchVectorField::rmap(ptf, addr);
119 
120  const tractionDisplacementCorrectionFvPatchVectorField& dmptf =
121  refCast<const tractionDisplacementCorrectionFvPatchVectorField>(ptf);
122 
123  traction_.rmap(dmptf.traction_, addr);
124  pressure_.rmap(dmptf.pressure_, addr);
125 }
126 
127 
128 // Update the coefficients associated with the patch field
129 void tractionDisplacementCorrectionFvPatchVectorField::updateCoeffs()
130 {
131  if (updated())
132  {
133  return;
134  }
135 
136  const label patchi = patch().index();
137 
138  const solidDisplacementThermo& thermo =
139  db().lookupObject<solidDisplacementThermo>
140  (
142  );
143 
144  const scalarField& E = thermo.E(patchi);
145  const scalarField& nu = thermo.nu(patchi);
146 
147  const scalarField mu(E/(2.0*(1.0 + nu)));
148  const scalarField lambda
149  (
150  thermo.planeStress()
151  ? nu*E/((1 + nu)*(1 - nu))
152  : nu*E/((1 + nu)*(1 - 2*nu))
153  );
154 
155  const vectorField n(patch().nf());
156 
157  const fvPatchField<symmTensor>& sigmaD =
158  patch().lookupPatchField<volSymmTensorField, symmTensor>("sigmaD");
159 
160  const fvPatchField<tensor>& sigmaExp =
161  patch().lookupPatchField<volTensorField, tensor>("sigmaExp");
162 
163  gradient() =
164  (
165  (traction_ + pressure_*n) - (n & (sigmaD + sigmaExp))
166  )/(2.0*mu + lambda);
167 
168  fixedGradientFvPatchVectorField::updateCoeffs();
169 }
170 
171 
173 {
175  writeEntry(os, "traction", traction_);
176  writeEntry(os, "pressure", pressure_);
177  writeEntry(os, "value", *this);
178 }
179 
180 
181 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
182 
184 (
186  tractionDisplacementCorrectionFvPatchVectorField
187 );
188 
189 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
190 
191 } // End namespace Foam
192 
193 // ************************************************************************* //
#define makePatchTypeField(PatchTypeField, typePatchTypeField)
Definition: fvPatchField.H:626
GeometricField< symmTensor, fvPatchField, volMesh > volSymmTensorField
Definition: volFieldsFwd.H:61
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
fvPatchField< vector > fvPatchVectorField
dimensionedScalar lambda(laminarTransport.lookup("lambda"))
GeometricField< tensor, fvPatchField, volMesh > volTensorField
Definition: volFieldsFwd.H:62
volVectorField vectorField(fieldObject, mesh)
Macros for easy insertion into run-time selection tables.
SymmTensor< scalar > symmTensor
SymmTensor of scalars.
Definition: symmTensor.H:48
void write(std::ostream &os, const bool binary, List< floatScalar > &fField)
Write floats ascii or binary.
List< label > labelList
A List of labels.
Definition: labelList.H:56
static const zero Zero
Definition: zero.H:97
volScalarField scalarField(fieldObject, mesh)
const dimensionedScalar mu
Atomic mass unit.
tractionDisplacementCorrectionFvPatchVectorField(const fvPatch &, const DimensionedField< vector, volMesh > &)
Construct from patch and internal field.
void writeEntry(Ostream &os, const HashTable< T, Key, Hash > &ht)
Definition: HashTableIO.C:96
label n
Tensor< scalar > tensor
Tensor of scalars.
Definition: tensor.H:51
const word dictName("noiseDict")
Namespace for OpenFOAM.