FreeFOAM The Cross-Platform CFD Toolkit
Hosted by SourceForge:
Get FreeFOAM at SourceForge.net.
            Fast, secure and Free Open Source software downloads

pimpleDyMFoam.C

Go to the documentation of this file.
00001 /*---------------------------------------------------------------------------*\
00002   =========                 |
00003   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
00004    \\    /   O peration     |
00005     \\  /    A nd           | Copyright (C) 1991-2010 OpenCFD Ltd.
00006      \\/     M anipulation  |
00007 -------------------------------------------------------------------------------
00008 License
00009     This file is part of OpenFOAM.
00010 
00011     OpenFOAM is free software: you can redistribute it and/or modify it
00012     under the terms of the GNU General Public License as published by
00013     the Free Software Foundation, either version 3 of the License, or
00014     (at your option) any later version.
00015 
00016     OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
00017     ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
00018     FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
00019     for more details.
00020 
00021     You should have received a copy of the GNU General Public License
00022     along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
00023 
00024 Application
00025     pimpleDyMFoam.C
00026 
00027 Description
00028     Transient solver for incompressible, flow of Newtonian fluids
00029     on a moving mesh using the PIMPLE (merged PISO-SIMPLE) algorithm.
00030 
00031     Turbulence modelling is generic, i.e. laminar, RAS or LES may be selected.
00032 
00033 Usage
00034     - PimpleDyMFoam [OPTION]
00035 
00036     @param -case <dir> \n
00037     Specify the case directory
00038 
00039     @param -parallel \n
00040     Run the case in parallel
00041 
00042     @param -help \n
00043     Display short usage message
00044 
00045     @param -doc \n
00046     Display Doxygen documentation page
00047 
00048     @param -srcDoc \n
00049     Display source code
00050 
00051 \*---------------------------------------------------------------------------*/
00052 
00053 #include <finiteVolume/fvCFD.H>
00054 #include <incompressibleTransportModels/singlePhaseTransportModel.H>
00055 #include <incompressibleTurbulenceModel/turbulenceModel.H>
00056 #include <dynamicFvMesh/dynamicFvMesh.H>
00057 
00058 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00059 
00060 int main(int argc, char *argv[])
00061 {
00062     #include <OpenFOAM/setRootCase.H>
00063 
00064     #include <OpenFOAM/createTime.H>
00065     #include <dynamicFvMesh/createDynamicFvMesh.H>
00066     #include <finiteVolume/readPIMPLEControls.H>
00067     #include <finiteVolume/initContinuityErrs.H>
00068     #include "createFields.H"
00069     #include <finiteVolume/readTimeControls.H>
00070 
00071     // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00072 
00073     Info<< "\nStarting time loop\n" << endl;
00074 
00075     while (runTime.run())
00076     {
00077         #include "readControls.H"
00078         #include <finiteVolume/CourantNo.H>
00079 
00080         // Make the fluxes absolute
00081         fvc::makeAbsolute(phi, U);
00082 
00083         #include <finiteVolume/setDeltaT.H>
00084 
00085         runTime++;
00086 
00087         Info<< "Time = " << runTime.timeName() << nl << endl;
00088 
00089         mesh.update();
00090 
00091         if (mesh.changing() && correctPhi)
00092         {
00093             #include "correctPhi.H"
00094         }
00095 
00096         // Make the fluxes relative to the mesh motion
00097         fvc::makeRelative(phi, U);
00098 
00099         if (mesh.changing() && checkMeshCourantNo)
00100         {
00101             #include <dynamicFvMesh/meshCourantNo.H>
00102         }
00103 
00104         // --- PIMPLE loop
00105         for (int ocorr=0; ocorr<nOuterCorr; ocorr++)
00106         {
00107             if (nOuterCorr != 1)
00108             {
00109                 p.storePrevIter();
00110             }
00111 
00112             #include "UEqn.H"
00113 
00114             // --- PISO loop
00115             for (int corr=0; corr<nCorr; corr++)
00116             {
00117                 rAU = 1.0/UEqn.A();
00118 
00119                 U = rAU*UEqn.H();
00120                 phi = (fvc::interpolate(U) & mesh.Sf());
00121 
00122                 if (p.needReference())
00123                 {
00124                     fvc::makeRelative(phi, U);
00125                     adjustPhi(phi, U, p);
00126                     fvc::makeAbsolute(phi, U);
00127                 }
00128 
00129                 for (int nonOrth=0; nonOrth<=nNonOrthCorr; nonOrth++)
00130                 {
00131                     fvScalarMatrix pEqn
00132                     (
00133                         fvm::laplacian(rAU, p) == fvc::div(phi)
00134                     );
00135 
00136                     pEqn.setReference(pRefCell, pRefValue);
00137 
00138                     if
00139                     (
00140                         ocorr == nOuterCorr-1
00141                      && corr == nCorr-1
00142                      && nonOrth == nNonOrthCorr)
00143                     {
00144                         pEqn.solve(mesh.solver(p.name() + "Final"));
00145                     }
00146                     else
00147                     {
00148                         pEqn.solve(mesh.solver(p.name()));
00149                     }
00150 
00151                     if (nonOrth == nNonOrthCorr)
00152                     {
00153                         phi -= pEqn.flux();
00154                     }
00155                 }
00156 
00157                 #include <finiteVolume/continuityErrs.H>
00158 
00159                 // Explicitly relax pressure for momentum corrector
00160                 if (ocorr != nOuterCorr-1)
00161                 {
00162                     p.relax();
00163                 }
00164 
00165                 // Make the fluxes relative to the mesh motion
00166                 fvc::makeRelative(phi, U);
00167 
00168                 U -= rAU*fvc::grad(p);
00169                 U.correctBoundaryConditions();
00170             }
00171 
00172             turbulence->correct();
00173         }
00174 
00175         runTime.write();
00176 
00177         Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
00178             << "  ClockTime = " << runTime.elapsedClockTime() << " s"
00179             << nl << endl;
00180     }
00181 
00182     Info<< "End\n" << endl;
00183 
00184     return 0;
00185 }
00186 
00187 
00188 // ************************ vim: set sw=4 sts=4 et: ************************ //
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines