OpenFOAM: API Guide
v2112
The open source CFD toolbox
solidDisplacementFoam.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) 2011-2016 OpenFOAM Foundation
9
-------------------------------------------------------------------------------
10
License
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
Application
27
solidDisplacementFoam
28
29
Group
30
grpStressAnalysisSolvers
31
32
Description
33
Transient segregated finite-volume solver of linear-elastic,
34
small-strain deformation of a solid body, with optional thermal
35
diffusion and thermal stresses.
36
37
Simple linear elasticity structural analysis code.
38
Solves for the displacement vector field D, also generating the
39
stress tensor field sigma.
40
41
\*---------------------------------------------------------------------------*/
42
43
#include "
fvCFD.H
"
44
45
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
46
47
int
main(
int
argc,
char
*argv[])
48
{
49
argList::addNote
50
(
51
"Transient segregated finite-volume solver of linear-elastic,"
52
" small-strain deformation of a solid body, with optional thermal"
53
" diffusion and thermal stresses"
54
);
55
56
#include "
postProcess.H
"
57
58
#include "
addCheckCaseOptions.H
"
59
#include "
setRootCaseLists.H
"
60
#include "
createTime.H
"
61
#include "
createMesh.H
"
62
#include "createControls.H"
63
#include "createFields.H"
64
65
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
66
67
Info
<<
"\nCalculating displacement field\n"
<<
endl
;
68
69
while
(
runTime
.loop())
70
{
71
Info
<<
"Iteration: "
<<
runTime
.value() <<
nl
<<
endl
;
72
73
#include "
readSolidDisplacementFoamControls.H
"
74
75
int
iCorr = 0;
76
scalar initialResidual = 0;
77
78
do
79
{
80
if
(thermalStress)
81
{
82
volScalarField
&
T
=
Tptr
();
83
solve
84
(
85
fvm::ddt(
T
) == fvm::laplacian(DT,
T
)
86
);
87
}
88
89
{
90
fvVectorMatrix
DEqn
91
(
92
fvm::d2dt2(
D
)
93
==
94
fvm::laplacian(2*
mu
+
lambda
,
D
,
"laplacian(DD,D)"
)
95
+ divSigmaExp
96
);
97
98
if
(thermalStress)
99
{
100
const
volScalarField
&
T
=
Tptr
();
101
DEqn += fvc::grad(threeKalpha*
T
);
102
}
103
104
//DEqn.setComponentReference(1, 0, vector::X, 0);
105
//DEqn.setComponentReference(1, 0, vector::Z, 0);
106
107
initialResidual = DEqn.solve().max().initialResidual();
108
109
if
(!
compactNormalStress
)
110
{
111
divSigmaExp = fvc::div(DEqn.flux());
112
}
113
}
114
115
{
116
volTensorField
gradD(fvc::grad(
D
));
117
sigmaD =
mu
*
twoSymm
(gradD) + (
lambda
*I)*
tr
(gradD);
118
119
if
(
compactNormalStress
)
120
{
121
divSigmaExp = fvc::div
122
(
123
sigmaD - (2*
mu
+
lambda
)*gradD,
124
"div(sigmaD)"
125
);
126
}
127
else
128
{
129
divSigmaExp += fvc::div(sigmaD);
130
}
131
}
132
133
}
while
(initialResidual >
convergenceTolerance
&& ++iCorr <
nCorr
);
134
135
#include "calculateStress.H"
136
137
runTime
.printExecutionTime(Info);
138
}
139
140
Info
<<
"End\n"
<<
endl
;
141
142
return
0;
143
}
144
145
146
// ************************************************************************* //
addCheckCaseOptions.H
Required Classes.
T
const volScalarField & T
Definition:
createFieldRefs.H:2
mu
const volScalarField & mu
Definition:
createFieldRefs.H:4
runTime
engineTime & runTime
Definition:
createEngineTime.H:13
createMesh.H
Required Variables.
createTime.H
nCorr
const int nCorr
Definition:
readFluidMultiRegionPIMPLEControls.H:3
fvCFD.H
Foam::tr
dimensionedScalar tr(const dimensionedSphericalTensor &dt)
Definition:
dimensionedSphericalTensor.C:51
Foam::volScalarField
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition:
volFieldsFwd.H:82
Foam::twoSymm
dimensionedSymmTensor twoSymm(const dimensionedSymmTensor &dt)
Definition:
dimensionedSymmTensor.C:95
Foam::Info
messageStream Info
Information stream (stdout output on master, null elsewhere)
Foam::endl
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition:
Ostream.H:372
Foam::volTensorField
GeometricField< tensor, fvPatchField, volMesh > volTensorField
Definition:
volFieldsFwd.H:87
Foam::fvVectorMatrix
fvMatrix< vector > fvVectorMatrix
Definition:
fvMatricesFwd.H:46
Foam::nl
constexpr char nl
The newline '\n' character (0x0a)
Definition:
Ostream.H:53
postProcess.H
Execute application functionObjects to post-process existing results.
readSolidDisplacementFoamControls.H
compactNormalStress
compactNormalStress
Definition:
readSolidDisplacementFoamControls.H:3
convergenceTolerance
convergenceTolerance
Definition:
readSolidDisplacementFoamControls.H:2
setRootCaseLists.H
D
const dimensionedScalar & D
Definition:
solveBulkSurfactant.H:4
solve
CEqn solve()
lambda
dimensionedScalar lambda("lambda", dimTime/sqr(dimLength), laminarTransport)
Tptr
Info<< "Reading field D\n"<< endl;volVectorField D(IOobject("D", runTime.timeName(), mesh, IOobject::MUST_READ, IOobject::AUTO_WRITE), mesh);autoPtr< volScalarField > Tptr
Definition:
createFields.H:19
applications
solvers
stressAnalysis
solidDisplacementFoam
solidDisplacementFoam.C
Generated by
1.9.5
OPENFOAM® is a registered
trademark
of OpenCFD Ltd.