Fixed URL for libccmio-2.6.1 (bug report #5 by Thomas Oliveira)
[foam-extend-3.2.git] / applications / solvers / solidMechanics / deprecatedSolvers / newContactStressFoam / newContactStressFoam.C
blob893eb86419699dda836c8b19d772ac5ae83d0719
1 /*---------------------------------------------------------------------------*\
2   =========                 |
3   \\      /  F ield         | foam-extend: Open Source CFD
4    \\    /   O peration     | Version:     3.2
5     \\  /    A nd           | Web:         http://www.foam-extend.org
6      \\/     M anipulation  | For copyright notice see file Copyright
7 -------------------------------------------------------------------------------
8 License
9     This file is part of foam-extend.
11     foam-extend is free software: you can redistribute it and/or modify it
12     under the terms of the GNU General Public License as published by the
13     Free Software Foundation, either version 3 of the License, or (at your
14     option) any later version.
16     foam-extend is distributed in the hope that it will be useful, but
17     WITHOUT ANY WARRANTY; without even the implied warranty of
18     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
19     General Public License for more details.
21     You should have received a copy of the GNU General Public License
22     along with foam-extend.  If not, see <http://www.gnu.org/licenses/>.
24 Application
25     stressedFoam
27 Description
28     Transient/steady-state solver for solid bodies in contact.
30     Solves for the displacement vector field U, also generating the
31     stress tensor field sigma.
33 \*---------------------------------------------------------------------------*/
35 #include "fvCFD.H"
36 #include "rheologyModel.H"
37 #include "contactProblem.H"
38 #include "componentReferenceList.H"
40 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
42 int main(int argc, char *argv[])
45 #   include "setRootCase.H"
47 #   include "createTime.H"
48 #   include "createMesh.H"
49 #   include "createFields.H"
51 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
53     Info<< "\nCalculating displacement field\n" << endl;
55     volScalarField rho = rheology.rho();
57     // Force n-sqaured projection
58 //     polyPatch::setNSquaredProjection(true);
60     while (runTime.loop())
61     {
62         Info<< "Iteration: " << runTime.timeName() << nl << endl;
64 #       include "readStressedFoamControls.H"
66         volScalarField mu = rheology.mu();
67         volScalarField lambda = rheology.lambda();
69         int iCorr=0;
70         scalar initialResidual=0;
72         contact.correct();
74         do
75         {
76             fvVectorMatrix UEqn
77             (
78                 fvm::d2dt2(rho, U)
79               ==
80                 fvm::laplacian(2*mu + lambda, U, "laplacian(DU,U)")
82               + fvc::div
83                 (
84                     mu*gradU.T() + lambda*(I*tr(gradU)) - (mu + lambda)*gradU,
85                     "div(sigma)"
86                 )
87             );
89 #           include "setComponentReference.H"
91             initialResidual = UEqn.solve().initialResidual();
93             gradU = fvc::grad(U);
95 #           include "calculateSigma.H"
97             rheology.correct();
99             rho = rheology.rho();
100             mu = rheology.mu();
101             lambda = rheology.lambda();
103         } while (initialResidual > convergenceTolerance && ++iCorr < nCorr);
105 #       include "calculateStress.H"
106 #       include "calculateContactArea.H"
108         Info<< "ExecutionTime = "
109             << runTime.elapsedCpuTime()
110             << " s\n\n" << endl;
111     }
113     Info<< "End\n" << endl;
115     return(0);
119 // ************************************************************************* //