Commit b6fbff52 authored by mattijs's avatar mattijs
Browse files

Merge branch 'master' of /home/dm4/OpenFOAM/OpenFOAM-dev

parents 0752591f 2d773bc0
......@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
......@@ -85,7 +85,8 @@ int main(int argc, char *argv[])
for (int corr=1; corr<=1; corr++)
{
volScalarField rAU("Dp", 1.0/UEqn.A());
volScalarField rAU(1.0/UEqn.A());
surfaceScalarField Dp("Dp", fvc::interpolate(rAU));
volVectorField HbyA("HbyA", U);
HbyA = rAU*UEqn.H();
......@@ -93,12 +94,12 @@ int main(int argc, char *argv[])
(
"phiHbyA",
(fvc::interpolate(HbyA) & mesh.Sf())
+ fvc::ddtPhiCorr(rAU, U, phi)
+ Dp*fvc::ddtCorr(U, phi)
);
fvScalarMatrix pEqn
(
fvm::laplacian(rAU, p) == fvc::div(phiHbyA)
fvm::laplacian(Dp, p) == fvc::div(phiHbyA)
);
pEqn.solve();
......
rho = thermo.rho();
volScalarField rAU(1.0/UEqn.A());
volVectorField HbyA("HbyA", U);
HbyA = invA & UEqn.H();
......@@ -11,9 +12,9 @@ if (pimple.transonic())
"phid",
fvc::interpolate(psi)
*(
(fvc::interpolate(HbyA) & mesh.Sf())
+ fvc::ddtPhiCorr(rAU, rho, U, phi)
)
(fvc::interpolate(rho*HbyA) & mesh.Sf())
+ fvc::interpolate(rho*rAU)*fvc::ddtCorr(rho, U, phi)
)/fvc::interpolate(rho)
);
while (pimple.correctNonOrthogonal())
......@@ -38,10 +39,9 @@ else
surfaceScalarField phiHbyA
(
"phiHbyA",
fvc::interpolate(rho)*
(
(fvc::interpolate(HbyA) & mesh.Sf())
+ fvc::ddtPhiCorr(rAU, rho, U, phi)
(fvc::interpolate(rho*HbyA) & mesh.Sf())
+ fvc::interpolate(rho*rAU)*fvc::ddtCorr(rho, U, phi)
)
);
......
rho = thermo.rho();
volScalarField rAU(1.0/UEqn.A());
surfaceScalarField Dp("Dp", fvc::interpolate(rho*rAU));
volVectorField HbyA("HbyA", U);
HbyA = rAU*UEqn.H();
......@@ -11,9 +13,9 @@ if (pimple.transonic())
"phid",
fvc::interpolate(psi)
*(
(fvc::interpolate(HbyA) & mesh.Sf())
+ fvc::ddtPhiCorr(rAU, rho, U, phi)
)
(fvc::interpolate(rho*HbyA) & mesh.Sf())
+ Dp*fvc::ddtCorr(rho, U, phi)
)/fvc::interpolate(rho)
);
fvOptions.makeRelative(fvc::interpolate(psi), phid);
......@@ -24,7 +26,7 @@ if (pimple.transonic())
(
fvm::ddt(psi, p)
+ fvm::div(phid, p)
- fvm::laplacian(rho*rAU, p)
- fvm::laplacian(Dp, p)
==
fvOptions(psi, p, rho.name())
);
......@@ -44,10 +46,9 @@ else
surfaceScalarField phiHbyA
(
"phiHbyA",
fvc::interpolate(rho)
*(
(fvc::interpolate(HbyA) & mesh.Sf())
+ fvc::ddtPhiCorr(rAU, rho, U, phi)
(
(fvc::interpolate(rho*HbyA) & mesh.Sf())
+ Dp*fvc::ddtCorr(rho, U, phi)
)
);
......@@ -59,7 +60,7 @@ else
(
fvm::ddt(psi, p)
+ fvc::div(phiHbyA)
- fvm::laplacian(rho*rAU, p)
- fvm::laplacian(Dp, p)
==
fvOptions(psi, p, rho.name())
);
......
rho = thermo.rho();
volScalarField rAU(1.0/UEqn.A());
surfaceScalarField Dp("Dp", fvc::interpolate(rho*rAU));
volVectorField HbyA("HbyA", U);
HbyA = rAU*UEqn.H();
......@@ -10,7 +12,13 @@ if (pimple.transonic())
(
"phid",
fvc::interpolate(psi)
*((fvc::interpolate(HbyA) & mesh.Sf()) - fvc::meshPhi(rho, U))
*(
(
(fvc::interpolate(rho*HbyA) & mesh.Sf())
//***HGW + Dp*fvc::ddtCorr(rho, U, phi)
)/fvc::interpolate(rho)
- fvc::meshPhi(rho, U)
)
);
fvOptions.makeRelative(fvc::interpolate(psi), phid);
......@@ -41,8 +49,11 @@ else
surfaceScalarField phiHbyA
(
"phiHbyA",
fvc::interpolate(rho)
*((fvc::interpolate(HbyA) & mesh.Sf()) - fvc::meshPhi(rho, U))
(
(fvc::interpolate(rho*HbyA) & mesh.Sf())
//***HGW + Dp*fvc::ddtCorr(rho, U, phi)
)
- fvc::interpolate(rho)*fvc::meshPhi(rho, U)
);
fvOptions.makeRelative(fvc::interpolate(rho), phiHbyA);
......
rho = thermo.rho();
volScalarField rAU(1.0/UEqn.A());
surfaceScalarField rhorAUf("Dp", fvc::interpolate(rho*rAU));
surfaceScalarField Dp("Dp", fvc::interpolate(rho*rAU));
volVectorField HbyA("HbyA", U);
HbyA = rAU*UEqn.H();
phi.boundaryField() =
fvc::interpolate(rho.boundaryField()*U.boundaryField())
& mesh.Sf().boundaryField();
surfaceScalarField phig(-rhorAUf*ghf*fvc::snGrad(rho)*mesh.magSf());
surfaceScalarField phig(-Dp*ghf*fvc::snGrad(rho)*mesh.magSf());
surfaceScalarField phiHbyA
(
"phiHbyA",
fvc::interpolate(rho)
*(
(fvc::interpolate(HbyA) & mesh.Sf())
+ fvc::ddtPhiCorr(rAU, rho, U, phi)
(
(fvc::interpolate(rho*HbyA) & mesh.Sf())
+ Dp*fvc::ddtCorr(rho, U, phi)
)
+ phig
);
......@@ -30,7 +30,7 @@ while (pimple.correctNonOrthogonal())
fvc::ddt(psi, rho)*gh
+ fvc::div(phiHbyA)
+ fvm::ddt(psi, p_rgh)
- fvm::laplacian(rhorAUf, p_rgh)
- fvm::laplacian(Dp, p_rgh)
==
parcels.Srho()
+ surfaceFilm.Srho()
......@@ -44,7 +44,7 @@ while (pimple.correctNonOrthogonal())
if (pimple.finalNonOrthogonalIter())
{
phi = phiHbyA + p_rghEqn.flux();
U = HbyA + rAU*fvc::reconstruct((p_rghEqn.flux() + phig)/rhorAUf);
U = HbyA + rAU*fvc::reconstruct((p_rghEqn.flux() + phig)/Dp);
U.correctBoundaryConditions();
fvOptions.correct(U);
}
......
rho = thermo.rho();
volScalarField rAU(1.0/UEqn.A());
surfaceScalarField Dp("Dp", fvc::interpolate(rho*rAU));
volVectorField HbyA("HbyA", U);
HbyA = rAU*UEqn.H();
......@@ -11,9 +13,9 @@ if (pimple.transonic())
"phid",
fvc::interpolate(psi)
*(
(fvc::interpolate(HbyA) & mesh.Sf())
+ fvc::ddtPhiCorr(rAU, rho, U, phi)
)
(fvc::interpolate(rho*HbyA) & mesh.Sf())
+ Dp*fvc::ddtCorr(rho, U, phi)
)/fvc::interpolate(rho)
);
fvOptions.makeRelative(fvc::interpolate(psi), phid);
......@@ -44,10 +46,9 @@ else
surfaceScalarField phiHbyA
(
"phiHbyA",
fvc::interpolate(rho)
*(
(fvc::interpolate(HbyA) & mesh.Sf())
+ fvc::ddtPhiCorr(rAU, rho, U, phi)
(
(fvc::interpolate(rho*HbyA) & mesh.Sf())
+ Dp*fvc::ddtCorr(rho, U, phi)
)
);
......
......@@ -6,20 +6,19 @@
thermo.rho() -= psi*p;
volScalarField rAU(1.0/UEqn.A());
surfaceScalarField rhorAUf("Dp", fvc::interpolate(rho*rAU));
surfaceScalarField Dp("Dp", fvc::interpolate(rho*rAU));
volVectorField HbyA("HbyA", U);
HbyA = rAU*UEqn.H();
surfaceScalarField phig(-rhorAUf*ghf*fvc::snGrad(rho)*mesh.magSf());
surfaceScalarField phig(-Dp*ghf*fvc::snGrad(rho)*mesh.magSf());
surfaceScalarField phiHbyA
(
"phiHbyA",
fvc::interpolate(rho)
*(
(fvc::interpolate(U) & mesh.Sf())
+ fvc::ddtPhiCorr(rAU, rho, U, phi)
(
(fvc::interpolate(rho*HbyA) & mesh.Sf())
+ Dp*fvc::ddtCorr(rho, U, phi)
)
+ phig
);
......@@ -39,7 +38,7 @@
fvScalarMatrix p_rghEqn
(
p_rghDDtEqn
- fvm::laplacian(rhorAUf, p_rgh)
- fvm::laplacian(Dp, p_rgh)
);
fvOptions.constrain(p_rghEqn);
......@@ -56,7 +55,7 @@
// Correct the momentum source with the pressure gradient flux
// calculated from the relaxed pressure
U = HbyA + rAU*fvc::reconstruct((phig + p_rghEqn.flux())/rhorAUf);
U = HbyA + rAU*fvc::reconstruct((phig + p_rghEqn.flux())/Dp);
U.correctBoundaryConditions();
fvOptions.correct(U);
K = 0.5*magSqr(U);
......
......@@ -6,6 +6,8 @@
thermo.rho() -= psi*p;
volScalarField rAU(1.0/UEqn.A());
surfaceScalarField Dp("Dp", fvc::interpolate(rho*rAU));
volVectorField HbyA("HbyA", U);
HbyA = rAU*UEqn.H();
......@@ -14,8 +16,10 @@
surfaceScalarField phiHbyA
(
"phiHbyA",
(fvc::interpolate(HbyA) & mesh.Sf())
+ fvc::ddtPhiCorr(rAU, rho, U, phi)
(
(fvc::interpolate(rho*HbyA) & mesh.Sf())
+ Dp*fvc::ddtCorr(rho, U, phi)
)/fvc::interpolate(rho)
);
fvOptions.makeRelative(phiHbyA);
......@@ -35,7 +39,7 @@
fvScalarMatrix pEqn
(
pDDtEqn
- fvm::laplacian(rho*rAU, p)
- fvm::laplacian(Dp, p)
==
fvOptions(psi, p, rho.name())
);
......@@ -55,10 +59,9 @@
surfaceScalarField phiHbyA
(
"phiHbyA",
fvc::interpolate(rho)
*(
(fvc::interpolate(HbyA) & mesh.Sf())
+ fvc::ddtPhiCorr(rAU, rho, U, phi)
(
(fvc::interpolate(rho*HbyA) & mesh.Sf())
+ Dp*fvc::ddtCorr(rho, U, phi)
)
);
......@@ -77,7 +80,7 @@
fvScalarMatrix pEqn
(
pDDtEqn
- fvm::laplacian(rho*rAU, p)
- fvm::laplacian(Dp, p)
);
fvOptions.constrain(pEqn);
......
......@@ -4,6 +4,8 @@ rho = min(rho, rhoMax);
rho.relax();
volScalarField rAU(1.0/UEqn().A());
surfaceScalarField Dp("Dp", fvc::interpolate(rho*rAU));
volVectorField HbyA("HbyA", U);
HbyA = rAU*UEqn().H();
......@@ -19,15 +21,13 @@ if (pimple.transonic())
"phid",
fvc::interpolate(psi)
*(
(fvc::interpolate(HbyA) & mesh.Sf())
+ fvc::ddtPhiCorr(rAU, rho, U, phi)
)
(fvc::interpolate(rho*HbyA) & mesh.Sf())
+ Dp*fvc::ddtCorr(rho, U, phi)
)/fvc::interpolate(rho)
);
fvOptions.makeRelative(fvc::interpolate(psi), phid);
volScalarField Dp("Dp", rho*rAU);
while (pimple.correctNonOrthogonal())
{
fvScalarMatrix pEqn
......@@ -54,10 +54,9 @@ else
surfaceScalarField phiHbyA
(
"phiHbyA",
fvc::interpolate(rho)
*(
(fvc::interpolate(HbyA) & mesh.Sf())
+ fvc::ddtPhiCorr(rAU, rho, U, phi)
(
(fvc::interpolate(rho*HbyA) & mesh.Sf())
+ Dp*fvc::ddtCorr(rho, U, phi)
)
);
......
......@@ -4,6 +4,8 @@ rho = min(rho, rhoMax);
rho.relax();
volScalarField rAU(1.0/UEqn().A());
surfaceScalarField Dp("Dp", fvc::interpolate(rho*rAU));
volVectorField HbyA("HbyA", U);
HbyA = rAU*UEqn().H();
......@@ -19,15 +21,13 @@ if (pimple.transonic())
"phid",
fvc::interpolate(psi)
*(
(fvc::interpolate(HbyA) & mesh.Sf())
+ fvc::ddtPhiCorr(rAU, rho, U, phiAbs)
)
(fvc::interpolate(rho*HbyA) & mesh.Sf())
+ Dp*fvc::ddtCorr(rho, U, phiAbs)
)/fvc::interpolate(rho)
);
fvOptions.makeRelative(fvc::interpolate(psi), phid);
volScalarField Dp("Dp", rho*rAU);
while (pimple.correctNonOrthogonal())
{
fvScalarMatrix pEqn
......@@ -54,18 +54,12 @@ else
surfaceScalarField phiHbyA
(
"phiHbyA",
fvc::interpolate(rho)
*(
(fvc::interpolate(HbyA) & mesh.Sf())
- fvc::meshPhi(rho, U)
+ fvc::ddtPhiCorr(rAU, rho, U, phiAbs)
)
(fvc::interpolate(rho*HbyA) & mesh.Sf())
+ Dp*fvc::ddtCorr(rho, U, phiAbs)
);
fvOptions.makeRelative(fvc::interpolate(rho), phiHbyA);
volScalarField Dp("Dp", rho*rAU);
while (pimple.correctNonOrthogonal())
{
// Pressure corrector
......
......@@ -59,7 +59,7 @@ int main(int argc, char *argv[])
#include "CourantNo.H"
#include "setInitialDeltaT.H"
// Create old-time absolute flux for ddtPhiCorr
// Create old-time absolute flux for ddtCorr
surfaceScalarField phiAbs("phiAbs", phi);
......@@ -75,7 +75,7 @@ int main(int argc, char *argv[])
// Make the fluxes absolute before mesh-motion
fvc::makeAbsolute(phi, rho, U);
// Update absolute flux for ddtPhiCorr
// Update absolute flux for ddtCorr
phiAbs = phi;
#include "setDeltaT.H"
......
......@@ -20,9 +20,9 @@ if (pimple.transonic())
"phid",
fvc::interpolate(psi)
*(
(fvc::interpolate(HbyA) & mesh.Sf())
+ fvc::ddtPhiCorr(rAU, rho, U, phi)
)
(fvc::interpolate(rho*HbyA) & mesh.Sf())
+ fvc::interpolate(rho*rAU)*fvc::ddtCorr(rho, U, phi)
)/fvc::interpolate(rho)
);
fvOptions.makeRelative(fvc::interpolate(psi), phid);
......@@ -64,10 +64,9 @@ else
surfaceScalarField phiHbyA
(
"phiHbyA",
fvc::interpolate(rho)
*(
(fvc::interpolate(HbyA) & mesh.Sf())
+ fvc::ddtPhiCorr(rAU, rho, U, phi)
(
(fvc::interpolate(rho*HbyA) & mesh.Sf())
+ fvc::interpolate(rho*rAU)*fvc::ddtCorr(rho, U, phi)
)
);
......
......@@ -12,7 +12,9 @@
surfaceScalarField phid
(
"phid",
fvc::interpolate(psi)*(fvc::interpolate(HbyA) & mesh.Sf())
fvc::interpolate(psi)
*(fvc::interpolate(rho*HbyA) & mesh.Sf())
/fvc::interpolate(rho)
);
fvOptions.makeRelative(fvc::interpolate(psi), phid);
......@@ -47,7 +49,7 @@
surfaceScalarField phiHbyA
(
"phiHbyA",
fvc::interpolate(rho)*(fvc::interpolate(HbyA) & mesh.Sf())
fvc::interpolate(rho*HbyA) & mesh.Sf()
);
fvOptions.makeRelative(fvc::interpolate(rho), phiHbyA);
......
......@@ -13,7 +13,9 @@ if (simple.transonic())
surfaceScalarField phid
(
"phid",
fvc::interpolate(psi)*(fvc::interpolate(HbyA) & mesh.Sf())
fvc::interpolate(psi)
*(fvc::interpolate(rho*HbyA) & mesh.Sf())
/fvc::interpolate(rho)
);
surfaceScalarField phic
......@@ -57,7 +59,7 @@ else
surfaceScalarField phiHbyA
(
"phiHbyA",
fvc::interpolate(rho)*(fvc::interpolate(HbyA) & mesh.Sf())
fvc::interpolate(rho*HbyA) & mesh.Sf()
);
closedVolume = adjustPhi(phiHbyA, U, p);
......
rho = thermo.rho();
volScalarField rAU(1.0/UEqn.A());
surfaceScalarField Dp("Dp", fvc::interpolate(rho*rAU));
volVectorField HbyA("HbyA", U);
HbyA = rAU*UEqn.H();
surfaceScalarField phid
(
"phid",
fvc::interpolate(psi)
*(
(fvc::interpolate(HbyA) & mesh.Sf())
+ fvc::ddtPhiCorr(rAU, rho, U, phi)
)
fvc::interpolate(psi)*
(
(mesh.Sf() & fvc::interpolate(rho*HbyA))
+ Dp*fvc::ddtCorr(rho, U, phi)
)/fvc::interpolate(rho)
);
volScalarField Dp("Dp", rho*rAU);
// Non-orthogonal pressure corrector loop
while (pimple.correctNonOrthogonal())
{
......
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2013 OpenFOAM Foundation
\\/ 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
createUf
Description
Creates and initialises the velocity velocity field Uf.
\*---------------------------------------------------------------------------*/
#ifndef createUf_H
#define createUf_H
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Info<< "Reading/calculating face velocity rhoUf\n" << endl;
surfaceVectorField rhoUf
(
IOobject
(
"rhoUf",
runTime.timeName(),
mesh,
IOobject::READ_IF_PRESENT,
IOobject::AUTO_WRITE
),