Merge branch 'filmPyr'

This commit is contained in:
andy
2011-04-12 14:01:03 +01:00
328 changed files with 47532 additions and 35716 deletions

View File

@ -4,13 +4,23 @@
+ fvm::div(phi, U)
+ turbulence->divDevRhoReff(U)
==
rho.dimensionedInternalField()*g
+ parcels.SU(U)
parcels.SU(U)
);
UEqn.relax();
if (momentumPredictor)
{
solve(UEqn == -fvc::grad(p));
solve
(
UEqn
==
fvc::reconstruct
(
(
- ghf*fvc::snGrad(rho)
- fvc::snGrad(p_rgh)
)*mesh.magSf()
)
);
}

View File

@ -5,7 +5,7 @@ tmp<fv::convectionScheme<scalar> > mvConvection
mesh,
fields,
phi,
mesh.divScheme("div(phi,Yi_h)")
mesh.divScheme("div(phi,Yi_hs)")
)
);
@ -23,7 +23,7 @@ tmp<fv::convectionScheme<scalar> > mvConvection
(
fvm::ddt(rho, Yi)
+ mvConvection->fvmDiv(phi, Yi)
- fvm::laplacian(turbulence->muEff(), Yi)
- fvm::laplacian(turbulence->alphaEff(), Yi)
==
parcels.SYi(i, Yi)
+ surfaceFilm.Srho(i)

View File

@ -1,3 +1,4 @@
if (chemistry.chemistry())
{
Info << "Solving chemistry" << endl;

View File

@ -15,11 +15,6 @@
const word inertSpecie(thermo.lookup("inertSpecie"));
volScalarField& p = thermo.p();
volScalarField& hs = thermo.hs();
const volScalarField& T = thermo.T();
const volScalarField& psi = thermo.psi();
Info<< "Creating field rho\n" << endl;
volScalarField rho
(
@ -34,6 +29,11 @@
thermo.rho()
);
volScalarField& p = thermo.p();
volScalarField& hs = thermo.hs();
const volScalarField& T = thermo.T();
const volScalarField& psi = thermo.psi();
Info<< "\nReading field U\n" << endl;
volVectorField U
(
@ -84,6 +84,28 @@
fvc::DDt(surfaceScalarField("phiU", phi/fvc::interpolate(rho)), p)
);
Info<< "Calculating field g.h\n" << endl;
volScalarField gh("gh", g & mesh.C());
surfaceScalarField ghf("gh", g & mesh.Cf());
volScalarField p_rgh
(
IOobject
(
"p_rgh",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);
// Force p_rgh to be consistent with p
p_rgh = p - rho*gh;
multivariateSurfaceInterpolationScheme<scalar>::fieldTable fields;
forAll(Y, i)

View File

@ -19,4 +19,6 @@
thermo.correct();
radiation->correct();
Info<< "min/max(T) = " << min(T).value() << ", " << max(T).value() << endl;
}

View File

@ -1,74 +1,58 @@
rho = thermo.rho();
volScalarField rAU(1.0/UEqn.A());
surfaceScalarField rhorAUf(rAU.name(), fvc::interpolate(rho*rAU));
U = rAU*UEqn.H();
if (transonic)
surfaceScalarField phiU
(
fvc::interpolate(rho)
*(
(fvc::interpolate(U) & mesh.Sf())
+ fvc::ddtPhiCorr(rAU, rho, U, phi)
)
);
phi = phiU - rhorAUf*ghf*fvc::snGrad(rho)*mesh.magSf();
for (int nonOrth=0; nonOrth<=nNonOrthCorr; nonOrth++)
{
surfaceScalarField phid
fvScalarMatrix p_rghEqn
(
"phid",
fvc::interpolate(psi)
*(
(fvc::interpolate(U) & mesh.Sf())
+ fvc::ddtPhiCorr(rAU, rho, U, phi)
fvc::ddt(psi, rho)*gh
+ fvc::div(phi)
+ fvm::ddt(psi, p_rgh)
- fvm::laplacian(rhorAUf, p_rgh)
==
parcels.Srho()
+ surfaceFilm.Srho()
);
p_rghEqn.solve
(
mesh.solver
(
p_rgh.select
(
pimpleCorr.finalIter()
&& corr == nCorr-1
&& nonOrth == nNonOrthCorr
)
)
);
for (int nonOrth=0; nonOrth<=nNonOrthCorr; nonOrth++)
if (nonOrth == nNonOrthCorr)
{
fvScalarMatrix pEqn
(
fvm::ddt(psi, p)
+ fvm::div(phid, p)
- fvm::laplacian(rho*rAU, p)
==
parcels.Srho()
+ surfaceFilm.Srho()
);
pEqn.solve();
if (nonOrth == nNonOrthCorr)
{
phi == pEqn.flux();
}
phi += p_rghEqn.flux();
}
}
else
{
phi =
fvc::interpolate(rho)
*(
(fvc::interpolate(U) & mesh.Sf())
+ fvc::ddtPhiCorr(rAU, rho, U, phi)
);
for (int nonOrth=0; nonOrth<=nNonOrthCorr; nonOrth++)
{
fvScalarMatrix pEqn
(
fvm::ddt(psi, p)
+ fvc::div(phi)
- fvm::laplacian(rho*rAU, p)
==
parcels.Srho()
+ surfaceFilm.Srho()
);
pEqn.solve();
if (nonOrth == nNonOrthCorr)
{
phi += pEqn.flux();
}
}
}
p = p_rgh + rho*gh;
#include "rhoEqn.H"
#include "compressibleContinuityErrs.H"
U -= rAU*fvc::grad(p);
U += rAU*fvc::reconstruct((phi - phiU)/rhorAUf);
U.correctBoundaryConditions();
DpDt = fvc::DDt(surfaceScalarField("phiU", phi/fvc::interpolate(rho)), p);

View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2009-2010 OpenCFD Ltd.
\\ / A nd | Copyright (C) 2009-2011 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -39,6 +39,7 @@ Description
#include "chemistrySolver.H"
#include "radiationModel.H"
#include "SLGThermo.H"
#include "pimpleLoop.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@ -66,8 +67,9 @@ int main(int argc, char *argv[])
while (runTime.run())
{
#include "readTimeControls.H"
#include "readPISOControls.H"
#include "readPIMPLEControls.H"
#include "compressibleCourantNo.H"
#include "setMultiRegionDeltaT.H"
#include "setDeltaT.H"
runTime++;
@ -84,20 +86,22 @@ int main(int argc, char *argv[])
#include "rhoEqn.H"
// --- PIMPLE loop
for (int ocorr=1; ocorr<=nOuterCorr; ocorr++)
for
(
pimpleLoop pimpleCorr(mesh, nOuterCorr);
pimpleCorr.loop();
pimpleCorr++
)
{
#include "UEqn.H"
#include "YEqn.H"
#include "hsEqn.H"
// --- PISO loop
for (int corr=1; corr<=nCorr; corr++)
{
#include "hsEqn.H"
#include "pEqn.H"
}
Info<< "T gas min/max = " << min(T).value() << ", "
<< max(T).value() << endl;
}
turbulence->correct();

View File

@ -0,0 +1,57 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2010-2011 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
OpenFOAM is free software: you can redistribute it and/or modify it
under the terms of the GNU General Public License as published by
the Free Software Foundation, either version 3 of the License, or
(at your option) any later version.
OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
for more details.
You should have received a copy of the GNU General Public License
along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
Global
setMultiRegionDeltaT
Description
Reset the timestep to maintain a constant maximum Courant numbers.
Reduction of time-step is immediate, but increase is damped to avoid
unstable oscillations.
\*---------------------------------------------------------------------------*/
if (adjustTimeStep)
{
if (CoNum == -GREAT)
{
CoNum = SMALL;
}
const scalar TFactorFluid = maxCo/(CoNum + SMALL);
const scalar TFactorFilm = maxCo/(surfaceFilm.CourantNumber() + SMALL);
const scalar dt0 = runTime.deltaTValue();
runTime.setDeltaT
(
min
(
dt0*min(min(TFactorFluid, TFactorFilm), 1.2),
maxDeltaT
)
);
}
// ************************************************************************* //