Commit 6ff341b5 authored by sergio's avatar sergio
Browse files

Adding alphaEqn.H with interpolation method.

Adding special alphaCourantNo for overlaping
Adding bounded term to UEq.H for overInterDyMFoam
Changing to NO_WRITE for the cellMask field
Changing twoSimpleRotors tutorial to open domain
parent ad8dfd4a
......@@ -3,6 +3,7 @@
fvVectorMatrix UEqn
(
fvm::ddt(rho, U) + fvm::div(rhoPhi, U)
- fvm::Sp(fvc::ddt(rho) + fvc::div(rhoPhi), U)
+ MRF.DDt(rho, U)
+ turbulence->divDevRhoReff(rho, U)
==
......
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2011-2014 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
alphaCourantNo
Description
Calculates and outputs the mean and maximum Courant Numbers.
\*---------------------------------------------------------------------------*/
scalar maxAlphaCo
(
readScalar(runTime.controlDict().lookup("maxAlphaCo"))
);
scalar alphaCoNum = 0.0;
scalar meanAlphaCoNum = 0.0;
if (mesh.nInternalFaces())
{
surfaceScalarField phiMask(localMin<scalar>(mesh).interpolate(cellMask));
scalarField sumPhi
(
mixture.nearInterface()().internalField()
*fvc::surfaceSum(mag(phiMask*phi))().internalField()
);
alphaCoNum = 0.5*gMax(sumPhi/mesh.V().field())*runTime.deltaTValue();
meanAlphaCoNum =
0.5*(gSum(sumPhi)/gSum(mesh.V().field()))*runTime.deltaTValue();
}
Info<< "Interface Courant Number mean: " << meanAlphaCoNum
<< " max: " << alphaCoNum << endl;
// ************************************************************************* //
{
word alphaScheme("div(phi,alpha)");
word alpharScheme("div(phirb,alpha)");
// Standard face-flux compression coefficient
surfaceScalarField phic(mixture.cAlpha()*mag(phi/mesh.magSf()));
// Add the optional isotropic compression contribution
if (icAlpha > 0)
{
phic *= (1.0 - icAlpha);
phic += (mixture.cAlpha()*icAlpha)*fvc::interpolate(mag(U));
}
surfaceScalarField::Boundary& phicBf =
phic.boundaryFieldRef();
// Do not compress interface at non-coupled boundary faces
// (inlets, outlets etc.)
forAll(phic.boundaryField(), patchi)
{
fvsPatchScalarField& phicp = phicBf[patchi];
if (!phicp.coupled())
{
phicp == 0;
}
}
tmp<surfaceScalarField> tphiAlpha;
if (MULESCorr)
{
mesh.interpolate(alpha1);
fvScalarMatrix alpha1Eqn
(
fv::EulerDdtScheme<scalar>(mesh).fvmDdt(alpha1)
+ fv::gaussConvectionScheme<scalar>
(
mesh,
phi,
upwind<scalar>(mesh, phi)
).fvmDiv(phi, alpha1)
);
alpha1Eqn.solve();
tmp<surfaceScalarField> tphiAlphaUD(alpha1Eqn.flux());
alphaPhi = tphiAlphaUD();
if (alphaApplyPrevCorr && tphiAlphaCorr0.valid())
{
MULES::correct(alpha1, alphaPhi, tphiAlphaCorr0.ref(), 1, 0);
alphaPhi += tphiAlphaCorr0();
}
// Cache the upwind-flux
tphiAlphaCorr0 = tphiAlphaUD;
alpha2 = 1.0 - alpha1;
mixture.correct();
}
for (int aCorr=0; aCorr<nAlphaCorr; aCorr++)
{
surfaceScalarField phir(phic*mixture.nHatf());
tmp<surfaceScalarField> tphiAlphaUn
(
fvc::flux
(
phi,
alpha1,
alphaScheme
)
+ fvc::flux
(
-fvc::flux(-phir, alpha2, alpharScheme),
alpha1,
alpharScheme
)
);
if (MULESCorr)
{
tmp<surfaceScalarField> tphiAlphaCorr(tphiAlphaUn() - alphaPhi);
volScalarField alpha10("alpha10", alpha1);
MULES::correct(alpha1, tphiAlphaUn(), tphiAlphaCorr.ref(), 1, 0);
mesh.interpolate(alpha1);
// Under-relax the correction for all but the 1st corrector
if (aCorr == 0)
{
alphaPhi += tphiAlphaCorr();
}
else
{
alpha1 = 0.5*alpha1 + 0.5*alpha10;
alphaPhi += 0.5*tphiAlphaCorr();
}
}
else
{
alphaPhi = tphiAlphaUn;
MULES::explicitSolve(alpha1, phi, alphaPhi, 1, 0);
mesh.interpolate(alpha1);
}
alpha2 = 1.0 - alpha1;
mixture.correct();
}
rhoPhi = alphaPhi*(rho1 - rho2) + phi*rho2;
if (alphaApplyPrevCorr && MULESCorr)
{
tphiAlphaCorr0 = alphaPhi - tphiAlphaCorr0;
}
Info<< "Phase-1 volume fraction = "
<< alpha1.weightedAverage(mesh.Vsc()).value()
<< " Min(alpha1) = " << min(alpha1).value()
<< " Max(alpha1) = " << max(alpha1).value()
<< endl;
}
......@@ -41,7 +41,7 @@ volScalarField cellMask
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::AUTO_WRITE
IOobject::NO_WRITE
),
mesh,
dimensionedScalar("cellMask", dimless, 1.0),
......
......@@ -16,7 +16,7 @@ FoamFile
dimensions [0 1 -1 0 0 0 0];
internalField uniform (0 0 0);
internalField uniform (0.1 0 0);
boundaryField
{
......@@ -36,16 +36,16 @@ boundaryField
value uniform (0 0 0);
}
// left1
// {
// type pressureInletOutletVelocity;
// value uniform (0 0 0);
// }
// left1
// {
// type fixedValue;
// value $internalField;
// }
outlet
{
type pressureInletOutletVelocity;
value uniform (0 0 0);
}
inlet
{
type fixedValue;
value $internalField;
}
//
// right1
// {
......
......@@ -28,26 +28,15 @@ boundaryField
type zeroGradient;
}
// left1
// {
// type uniformTotalPressure;
// pressure table
// (
// (0 10)
// (1 40)
// );
// p0 40; // only used for restarts
// U U;
// phi phi;
// rho none;
// psi none;
// gamma 1;
// value uniform 40;
// }
// left1
// {
// type zeroGradient;
// }
outlet
{
type fixedValue;
value uniform 0;
}
inlet
{
type zeroGradient;
}
//
// right1
// {
......
......@@ -94,7 +94,23 @@ boundary
(
(3 7 6 2)
(1 5 4 0)
);
}
inlet
{
type patch;
faces
(
(0 4 7 3)
);
}
outlet
{
type patch;
faces
(
(2 6 5 1)
);
}
......
......@@ -24,7 +24,7 @@ startTime 0;
stopAt endTime;
endTime 0.06;
endTime 0.6;
deltaT 0.00025;
......
......@@ -20,7 +20,7 @@ libs ("liboverset.so");
application overInterDyMFoam ;
startFrom startTime;
startFrom latestTime;//startTime;
startTime 0.0;
......@@ -32,7 +32,7 @@ deltaT 0.001;
writeControl adjustableRunTime;
writeInterval 0.05;
writeInterval 0.1;
purgeWrite 0;
......
......@@ -80,8 +80,8 @@ solvers
PIMPLE
{
momentumPredictor no;
nOuterCorrectors 3;
nCorrectors 1;
nOuterCorrectors 2;
nCorrectors 2;
nNonOrthogonalCorrectors 0;
ddtCorr yes;
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment