KirchhoffShell.C
Go to the documentation of this file.
1/*---------------------------------------------------------------------------*\
2 ========= |
3 \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4 \\ / O peration |
5 \\ / A nd | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2019-2021 OpenCFD Ltd.
9-------------------------------------------------------------------------------
10License
11 This file is part of OpenFOAM.
12
13 OpenFOAM is free software: you can redistribute it and/or modify it
14 under the terms of the GNU General Public License as published by
15 the Free Software Foundation, either version 3 of the License, or
16 (at your option) any later version.
17
18 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21 for more details.
22
23 You should have received a copy of the GNU General Public License
24 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25
26\*---------------------------------------------------------------------------*/
27
28#include "KirchhoffShell.H"
30#include "fvPatchFields.H"
32#include "subCycle.H"
33
34// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
35
36namespace Foam
37{
38namespace regionModels
39{
40
41// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
42
44
46
47// * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
48
50{
51 this->solution().readEntry("nNonOrthCorr", nNonOrthCorr_);
52 return true;
53}
54
55// * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
56
58{
59 if (debug)
60 {
62 }
63
64 const Time& time = primaryMesh().time();
65
66 areaScalarField solidMass(rho()*h_);
67 areaScalarField solidD(D()/solidMass);
68
69 // Save old times
72
73 if (nSubCycles_ > 1)
74 {
75 // Restore the oldTime in sub-cycling
76 w_.oldTime() = w0_;
77 w_.oldTime().oldTime() = w00_;
80 }
81
82 for
83 (
84 subCycleTime wSubCycle
85 (
86 const_cast<Time&>(time),
88 );
89 !(++wSubCycle).end();
90 )
91 {
92
95
97 (
99 + f1_*fam::ddt(w_)
100 - f0_*sqrt(solidD)*fac::ddt(laplaceW_)
101 + solidD*(laplace2W_ + f2_*fac::ddt(laplace2W_))
102 ==
103 ps_/solidMass
104 + faOptions()(solidMass, w_, dimLength/sqr(dimTime))
105 );
106
107 faOptions().constrain(wEqn);
108
109 wEqn.solve();
110
111 if (wSubCycle.index() >= wSubCycle.nSubCycles())
112 {
113 // Cache oldTimes inside the sub-cycling
114 w0_ = w_.oldTime();
115 w00_ = w_.oldTime().oldTime();
118
119 // Update shell acceleration
120 a_ = fac::d2dt2(w_);
121 }
122 }
123
124 Info<< "ws_vibrationShell: "
125 << "min = " << min(w_).value() << ", "
126 << "max = " << max(w_).value() << endl;
127
128 // Restore old time in main time
129 w_.oldTime() = w0;
130 w_.oldTime().oldTime() = w00;
131
133}
134
135
136// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
137
139(
140 const word& modelType,
141 const fvPatch& patch,
142 const dictionary& dict
143)
144:
145 vibrationShellModel(modelType, patch, dict),
146 f0_("f0", dimless, dict),
147 f1_("f1", inv(dimTime), dict),
148 f2_("f2", dimTime, dict),
149 nNonOrthCorr_(1),
150 nSubCycles_(1),
151 ps_
152 (
154 (
155 "ps_" + regionName_,
156 primaryMesh().time().timeName(),
157 primaryMesh(),
158 IOobject::READ_IF_PRESENT,
159 IOobject::NO_WRITE
160 ),
161 regionMesh(),
163 ),
164 h_
165 (
167 (
168 "h_" + regionName_,
169 primaryMesh().time().timeName(),
170 primaryMesh(),
171 IOobject::MUST_READ,
172 IOobject::AUTO_WRITE
173 ),
174 regionMesh()
175 ),
176 laplaceW_
177 (
179 (
180 "laplaceW_" + regionName_,
181 primaryMesh().time().timeName(),
182 primaryMesh(),
183 IOobject::NO_READ,
184 IOobject::NO_WRITE
185 ),
186 regionMesh(),
188 ),
189 laplace2W_
190 (
192 (
193 "laplace2W_" + regionName_,
194 primaryMesh().time().timeName(),
195 primaryMesh(),
196 IOobject::NO_READ,
197 IOobject::NO_WRITE
198 ),
199 regionMesh(),
201 ),
202 w0_
203 (
205 (
206 "w0_" + regionName_,
207 primaryMesh().time().timeName(),
208 primaryMesh(),
209 IOobject::NO_READ,
210 IOobject::NO_WRITE
211 ),
212 regionMesh(),
214 ),
215 w00_
216 (
218 (
219 "w00_" + regionName_,
220 primaryMesh().time().timeName(),
221 primaryMesh(),
222 IOobject::NO_READ,
223 IOobject::NO_WRITE
224 ),
225 regionMesh(),
227 ),
228 laplaceW0_
229 (
231 (
232 "laplaceW0_" + regionName_,
233 primaryMesh().time().timeName(),
234 primaryMesh(),
235 IOobject::NO_READ,
236 IOobject::NO_WRITE
237 ),
238 regionMesh(),
240 ),
241 laplace2W0_
242 (
244 (
245 "laplace2W0_" + regionName_,
246 primaryMesh().time().timeName(),
247 primaryMesh(),
248 IOobject::NO_READ,
249 IOobject::NO_WRITE
250 ),
251 regionMesh(),
253 )
254{
255 init(dict);
256}
257
258// * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
259
261{}
262
263
265{
266 nNonOrthCorr_ = solution().get<label>("nNonOrthCorr");
267 nSubCycles_ = solution().get<label>("nSubCycles");
268
269 for (int nonOrth=0; nonOrth<=nNonOrthCorr_; nonOrth++)
270 {
272 }
273}
274
275
277{
278 const dimensionedScalar E("E", dimForce/dimArea , solid().E());
279 const dimensionedScalar nu("nu", dimless, solid().nu());
280
281 return tmp<areaScalarField>(E*pow3(h_)/(12*(1 - sqr(nu))));
282}
283
284
286{
288 (
290 (
292 (
293 "rhos",
295 primaryMesh(),
298 false
299 ),
300 regionMesh(),
302 zeroGradientFaPatchScalarField::typeName
303 )
304 );
305}
306
308{}
309
310// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
311
312} // End namespace regionModels
313} // End namespace Foam
314
315// ************************************************************************* //
Macros for easy insertion into run-time selection tables.
#define addToRunTimeSelectionTable(baseType, thisType, argNames)
Add to construction table with typeName as the key.
#define w0
Definition: blockCreate.C:33
const GeometricField< Type, PatchField, GeoMesh > & oldTime() const
Return old time field.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:170
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:80
const iterator & end()
End of list for forward iterators.
Definition: UILList.H:414
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:126
T get(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
bool readEntry(const word &keyword, T &val, enum keyType::option matchOpt=keyType::REGEX, bool mandatory=true) const
A special matrix type and solver, designed for finite area solutions of scalar equations....
Definition: faMatrix.H:76
SolverPerformance< Type > solve(const dictionary &)
Solve returning the solution statistics.
Definition: faMatrixSolve.C:55
void correct(GeometricField< Type, faPatchField, areaMesh > &field)
Apply correction to field.
void constrain(faMatrix< Type > &eqn)
Apply constraints to equation.
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:290
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:71
areaScalarField h_
Thickness [m].
areaScalarField laplace2W0_
Cache laplace2.oldTime() in sub-cycling.
const tmp< areaScalarField > rho() const
Return density [Kg/m3].
areaScalarField ps_
External surface source [Pa].
areaScalarField laplaceW_
Laplace of the displacement.
void solveDisplacement()
Solve energy equation.
areaScalarField w0_
Cache w.oldTime() in sub-cycling.
areaScalarField w00_
Cache w.oldTime.oldTime() in sub-cycling.
virtual void preEvolveRegion()
Pre-evolve thermal baffle.
virtual void info()
Provide some feedback.
areaScalarField laplaceW0_
Cache laplaceW.oldTime() in sub-cycling.
areaScalarField laplace2W_
Laplace of the Laplace for the displacement.
label nNonOrthCorr_
Number of non orthogonal correctors.
const tmp< areaScalarField > D() const
Return stiffness.
virtual void evolveRegion()
Evolve the thermal baffle.
const Time & time() const
Return the reference to the time database.
const dictionary & solution() const
Return the solution dictionary.
const fvMesh & primaryMesh() const
Return the reference to the primary mesh database.
const faMesh & regionMesh() const
Return the region mesh database.
Foam::fa::options & faOptions() noexcept
Return faOptions.
areaScalarField w_
Shell displacement.
areaScalarField a_
Shell acceleration.
const solidProperties & solid() const noexcept
Return solid properties.
A class for managing sub-cycling times.
Definition: subCycleTime.H:51
A class for managing temporary objects.
Definition: tmp.H:65
A class for handling words, derived from Foam::string.
Definition: word.H:68
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
Definition: className.H:121
word timeName
Definition: getTimeIndex.H:3
#define InfoInFunction
Report an information message using Foam::Info.
tmp< GeometricField< Type, faPatchField, areaMesh > > laplacian(const GeometricField< Type, faPatchField, areaMesh > &vf, const word &name)
Definition: facLaplacian.C:47
tmp< GeometricField< Type, faPatchField, areaMesh > > d2dt2(const dimensioned< Type > dt, const faMesh &mesh)
Definition: facD2dt2.C:46
tmp< GeometricField< Type, faPatchField, areaMesh > > ddt(const dimensioned< Type > dt, const faMesh &mesh)
Definition: facDdt.C:47
tmp< faMatrix< Type > > d2dt2(const GeometricField< Type, faPatchField, areaMesh > &vf)
Definition: famD2dt2.C:47
tmp< faMatrix< Type > > ddt(const GeometricField< Type, faPatchField, areaMesh > &vf)
Definition: famDdt.C:48
Namespace for OpenFOAM.
const dimensionSet dimPressure
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:47
const dimensionSet dimless
Dimensionless.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:52
dimensionedScalar pow3(const dimensionedScalar &ds)
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:53
const dimensionSet dimArea(sqr(dimLength))
Definition: dimensionSets.H:59
messageStream Info
Information stream (stdout output on master, null elsewhere)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
dimensionedScalar sqrt(const dimensionedScalar &ds)
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:33
const dimensionSet dimForce
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
dimensionedSphericalTensor inv(const dimensionedSphericalTensor &dt)
const dimensionSet dimDensity
volScalarField & nu
dictionary dict