solidDisplacementFoam.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-2018 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 Application
25  solidDisplacementFoam
26 
27 Description
28  Transient segregated finite-volume solver of linear-elastic,
29  small-strain deformation of a solid body, with optional thermal
30  diffusion and thermal stresses.
31 
32  Simple linear elasticity structural analysis code.
33  Solves for the displacement vector field D, also generating the
34  stress tensor field sigma.
35 
36 \*---------------------------------------------------------------------------*/
37 
38 #include "fvCFD.H"
39 #include "fvOptions.H"
40 #include "Switch.H"
41 
42 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
43 
44 int main(int argc, char *argv[])
45 {
46  #include "postProcess.H"
47 
48  #include "setRootCaseLists.H"
49  #include "createTime.H"
50  #include "createMesh.H"
51  #include "createControls.H"
52  #include "createFields.H"
53 
54  // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
55 
56  Info<< "\nCalculating displacement field\n" << endl;
57 
58  while (runTime.loop())
59  {
60  Info<< "Iteration: " << runTime.value() << nl << endl;
61 
63 
64  int iCorr = 0;
65  scalar initialResidual = 0;
66 
67  do
68  {
69  if (thermalStress)
70  {
71  volScalarField& T = Tptr();
72  fvScalarMatrix TEqn
73  (
74  fvm::ddt(T) == fvm::laplacian(DT, T) + fvOptions(T)
75  );
76 
77  fvOptions.constrain(TEqn);
78 
79  TEqn.solve();
80 
81  fvOptions.correct(T);
82  }
83 
84  {
85  fvVectorMatrix DEqn
86  (
87  fvm::d2dt2(D)
88  ==
89  fvm::laplacian(2*mu + lambda, D, "laplacian(DD,D)")
90  + divSigmaExp
91  + fvOptions.d2dt2(D)
92  );
93 
94  if (thermalStress)
95  {
96  const volScalarField& T = Tptr();
97  DEqn += fvc::grad(threeKalpha*T);
98  }
99 
100  fvOptions.constrain(DEqn);
101 
102  initialResidual = DEqn.solve().max().initialResidual();
103 
104  if (!compactNormalStress)
105  {
106  divSigmaExp = fvc::div(DEqn.flux());
107  }
108  }
109 
110  {
111  volTensorField gradD(fvc::grad(D));
112  sigmaD = mu*twoSymm(gradD) + (lambda*I)*tr(gradD);
113 
115  {
116  divSigmaExp = fvc::div
117  (
118  sigmaD - (2*mu + lambda)*gradD,
119  "div(sigmaD)"
120  );
121  }
122  else
123  {
124  divSigmaExp += fvc::div(sigmaD);
125  }
126  }
127 
128  } while (initialResidual > convergenceTolerance && ++iCorr < nCorr);
129 
130  #include "calculateStress.H"
131 
132  Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
133  << " ClockTime = " << runTime.elapsedClockTime() << " s"
134  << nl << endl;
135  }
136 
137  Info<< "End\n" << endl;
138 
139  return 0;
140 }
141 
142 
143 // ************************************************************************* //
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcGrad.C:52
Info<< "Reading field D\"<< endl;volVectorField D(IOobject("D", runTime.timeName(), mesh, IOobject::MUST_READ, IOobject::AUTO_WRITE), mesh);autoPtr< volScalarField > Tptr(nullptr)
fvMatrix< scalar > fvScalarMatrix
Definition: fvMatricesFwd.H:42
tmp< GeometricField< Type, fvPatchField, volMesh > > d2dt2(const GeometricField< Type, fvPatchField, volMesh > &vf)
Definition: fvcD2dt2.C:45
fv::options & fvOptions
GeometricField< tensor, fvPatchField, volMesh > volTensorField
Definition: volFieldsFwd.H:59
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcDiv.C:47
Switch compactNormalStress(stressControl.lookup("compactNormalStress"))
engineTime & runTime
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:256
tmp< GeometricField< Type, fvPatchField, volMesh > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvcLaplacian.C:45
dimensionedSymmTensor twoSymm(const dimensionedSymmTensor &dt)
tmp< GeometricField< Type, fvPatchField, volMesh > > ddt(const dimensioned< Type > dt, const fvMesh &mesh)
Definition: fvcDdt.C:45
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
static const Identity< scalar > I
Definition: Identity.H:93
dimensionedScalar tr(const dimensionedSphericalTensor &dt)
static const char nl
Definition: Ostream.H:265
const dimensionedScalar mu
Atomic mass unit.
int nCorr
Definition: createControls.H:3
fvMatrix< vector > fvVectorMatrix
Definition: fvMatricesFwd.H:45
dimensionedScalar lambda(laminarTransport.lookup("lambda"))
messageStream Info
Execute application functionObjects to post-process existing results.