applications/solvers/stressAnalysis/solidDisplacementFoam/solidDisplacementFoam.C

   1 /*---------------------------------------------------------------------------*\
   2   =========                 |
   3   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   4    \\    /   O peration     |
   5     \\  /    A nd           | Copyright (C) 1991-2010 OpenCFD Ltd.
   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 "Switch.H"
  40
  41 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
  42
  43 int main(int argc, char *argv[])
  44 {
  45     #include "setRootCase.H"
  46
  47     #include "createTime.H"
  48     #include "createMesh.H"
  49     #include "readMechanicalProperties.H"
  50     #include "readThermalProperties.H"
  51     #include "readSolidDisplacementFoamControls.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
  62         #include "readSolidDisplacementFoamControls.H"
  63
  64         int iCorr = 0;
  65         scalar initialResidual = 0;
  66
  67         do
  68         {
  69             if (thermalStress)
  70             {
  71                 volScalarField& T = Tptr();
  72                 solve
  73                 (
  74                     fvm::ddt(T) == fvm::laplacian(DT, T)
  75                 );
  76             }
  77
  78             {
  79                 fvVectorMatrix DEqn
  80                 (
  81                     fvm::d2dt2(D)
  82                  ==
  83                     fvm::laplacian(2*mu + lambda, D, "laplacian(DD,D)")
  84                   + divSigmaExp
  85                 );
  86
  87                 if (thermalStress)
  88                 {
  89                     const volScalarField& T = Tptr();
  90                     DEqn += fvc::grad(threeKalpha*T);
  91                 }
  92
  93                 //DEqn.setComponentReference(1, 0, vector::X, 0);
  94                 //DEqn.setComponentReference(1, 0, vector::Z, 0);
  95
  96                 initialResidual = DEqn.solve().initialResidual();
  97
  98                 if (!compactNormalStress)
  99                 {
 100                     divSigmaExp = fvc::div(DEqn.flux());
 101                 }
 102             }
 103
 104             {
 105                 volTensorField gradD = fvc::grad(D);
 106                 sigmaD = mu*twoSymm(gradD) + (lambda*I)*tr(gradD);
 107
 108                 if (compactNormalStress)
 109                 {
 110                     divSigmaExp = fvc::div
 111                     (
 112                         sigmaD - (2*mu + lambda)*gradD,
 113                         "div(sigmaD)"
 114                     );
 115                 }
 116                 else
 117                 {
 118                     divSigmaExp += fvc::div(sigmaD);
 119                 }
 120             }
 121
 122         } while (initialResidual > convergenceTolerance && ++iCorr < nCorr);
 123
 124         #include "calculateStress.H"
 125
 126         Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
 127             << "  ClockTime = " << runTime.elapsedClockTime() << " s"
 128             << nl << endl;
 129     }
 130
 131     Info<< "End\n" << endl;
 132
 133     return 0;
 134 }
 135
 136
 137 // ************************************************************************* //