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

sonicDyMFoam.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     sonicDyMFoam
00026 
00027 Description
00028     Transient solver for trans-sonic/supersonic, laminar or turbulent flow
00029     of a compressible gas with mesh motion.
00030 
00031 Usage
00032     - sonicDyMFoam [OPTION]
00033 
00034     @param -case <dir> \n
00035     Specify the case directory
00036 
00037     @param -parallel \n
00038     Run the case in parallel
00039 
00040     @param -help \n
00041     Display short usage message
00042 
00043     @param -doc \n
00044     Display Doxygen documentation page
00045 
00046     @param -srcDoc \n
00047     Display source code
00048 
00049 \*---------------------------------------------------------------------------*/
00050 
00051 #include <finiteVolume/fvCFD.H>
00052 #include <basicThermophysicalModels/basicPsiThermo.H>
00053 #include <compressibleTurbulenceModel/turbulenceModel.H>
00054 #include <dynamicMesh/motionSolver.H>
00055 
00056 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00057 
00058 int main(int argc, char *argv[])
00059 {
00060     #include <OpenFOAM/setRootCase.H>
00061     #include <OpenFOAM/createTime.H>
00062     #include <OpenFOAM/createMesh.H>
00063     #include "createFields.H"
00064     #include <finiteVolume/initContinuityErrs.H>
00065 
00066     autoPtr<Foam::motionSolver> motionPtr = motionSolver::New(mesh);
00067 
00068     // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00069 
00070     Info<< "\nStarting time loop\n" << endl;
00071 
00072     while (runTime.loop())
00073     {
00074         Info<< "Time = " << runTime.timeName() << nl << endl;
00075 
00076         #include <finiteVolume/readPISOControls.H>
00077         #include <finiteVolume/compressibleCourantNo.H>
00078 
00079         mesh.movePoints(motionPtr->newPoints());
00080 
00081         #include <finiteVolume/rhoEqn.H>
00082 
00083         fvVectorMatrix UEqn
00084         (
00085             fvm::ddt(rho, U)
00086           + fvm::div(phi, U)
00087           + turbulence->divDevRhoReff(U)
00088         );
00089 
00090         solve(UEqn == -fvc::grad(p));
00091 
00092         #include "../sonicFoam/eEqn.H"
00093 
00094 
00095         // --- PISO loop
00096 
00097         for (int corr=0; corr<nCorr; corr++)
00098         {
00099             U = UEqn.H()/UEqn.A();
00100 
00101             surfaceScalarField phid
00102             (
00103                 "phid",
00104                 fvc::interpolate(psi)
00105                *(
00106                     (fvc::interpolate(U) & mesh.Sf()) - fvc::meshPhi(rho, U)
00107                 )
00108             );
00109 
00110             for (int nonOrth=0; nonOrth<=nNonOrthCorr; nonOrth++)
00111             {
00112                 fvScalarMatrix pEqn
00113                 (
00114                     fvm::ddt(psi, p)
00115                   + fvm::div(phid, p)
00116                   - fvm::laplacian(rho/UEqn.A(), p)
00117                 );
00118 
00119                 pEqn.solve();
00120 
00121                 phi = pEqn.flux();
00122             }
00123 
00124             #include "compressibleContinuityErrs.H"
00125 
00126             U -= fvc::grad(p)/UEqn.A();
00127             U.correctBoundaryConditions();
00128         }
00129 
00130         DpDt = fvc::DDt
00131         (
00132             surfaceScalarField
00133             (
00134                 "phiU",
00135                 phi/fvc::interpolate(rho) + fvc::meshPhi(rho, U)
00136             ),
00137             p
00138         );
00139 
00140         turbulence->correct();
00141 
00142         rho = psi*p;
00143 
00144         runTime.write();
00145 
00146         Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
00147             << "  ClockTime = " << runTime.elapsedClockTime() << " s"
00148             << nl << endl;
00149     }
00150 
00151     Info<< "End\n" << endl;
00152 
00153     return 0;
00154 }
00155 
00156 
00157 // ************************ vim: set sw=4 sts=4 et: ************************ //
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines