diff --git a/applications/solvers/multiphase/MPPICInterFoam/Allwclean b/applications/solvers/multiphase/MPPICInterFoam/Allwclean new file mode 100755 index 0000000000..d6607c3db0 --- /dev/null +++ b/applications/solvers/multiphase/MPPICInterFoam/Allwclean @@ -0,0 +1,8 @@ +#!/bin/sh +cd ${0%/*} || exit 1 # run from this directory +set -x + +wclean libso CompressibleTwoPhaseMixtureTurbulenceModels +wclean + +# ----------------------------------------------------------------- end-of-file diff --git a/applications/solvers/multiphase/MPPICInterFoam/Allwmake b/applications/solvers/multiphase/MPPICInterFoam/Allwmake new file mode 100755 index 0000000000..110d4b093a --- /dev/null +++ b/applications/solvers/multiphase/MPPICInterFoam/Allwmake @@ -0,0 +1,8 @@ +#!/bin/sh +cd ${0%/*} || exit 1 # run from this directory +set -x + +wmake libso CompressibleTwoPhaseMixtureTurbulenceModels +wmake + +# ----------------------------------------------------------------- end-of-file diff --git a/applications/solvers/multiphase/MPPICInterFoam/CompressibleTwoPhaseMixtureTurbulenceModels/CompressibleTwoPhaseMixtureTurbulenceModels.C b/applications/solvers/multiphase/MPPICInterFoam/CompressibleTwoPhaseMixtureTurbulenceModels/CompressibleTwoPhaseMixtureTurbulenceModels.C new file mode 100644 index 0000000000..6124bf2522 --- /dev/null +++ b/applications/solvers/multiphase/MPPICInterFoam/CompressibleTwoPhaseMixtureTurbulenceModels/CompressibleTwoPhaseMixtureTurbulenceModels.C @@ -0,0 +1,64 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2016 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 . + +\*---------------------------------------------------------------------------*/ + +#include "PhaseCompressibleTurbulenceModel.H" +#include "immiscibleIncompressibleTwoPhaseMixture.H" +#include "addToRunTimeSelectionTable.H" +#include "makeTurbulenceModel.H" + +#include "laminar.H" +#include "turbulentTransportModel.H" +#include "LESModel.H" + +makeBaseTurbulenceModel +( + volScalarField, + volScalarField, + compressibleTurbulenceModel, + PhaseCompressibleTurbulenceModel, + immiscibleIncompressibleTwoPhaseMixture +); + +#define makeRASModel(Type) \ + makeTemplatedTurbulenceModel \ + (immiscibleIncompressibleTwoPhaseMixturePhaseCompressibleTurbulenceModel, RAS, Type) + +#define makeLESModel(Type) \ + makeTemplatedTurbulenceModel \ + (immiscibleIncompressibleTwoPhaseMixturePhaseCompressibleTurbulenceModel, LES, Type) + +#include "kEpsilon.H" +makeRASModel(kEpsilon); + +#include "Smagorinsky.H" +makeLESModel(Smagorinsky); + +#include "kEqn.H" +makeLESModel(kEqn); + +#include "kOmega.H" +makeRASModel(kOmega); + +// ************************************************************************* // diff --git a/applications/solvers/multiphase/MPPICInterFoam/CompressibleTwoPhaseMixtureTurbulenceModels/Make/files b/applications/solvers/multiphase/MPPICInterFoam/CompressibleTwoPhaseMixtureTurbulenceModels/Make/files new file mode 100644 index 0000000000..de1a922fbc --- /dev/null +++ b/applications/solvers/multiphase/MPPICInterFoam/CompressibleTwoPhaseMixtureTurbulenceModels/Make/files @@ -0,0 +1,3 @@ +CompressibleTwoPhaseMixtureTurbulenceModels.C + +LIB = $(FOAM_LIBBIN)/libCompressibleTwoPhaseMixtureTurbulenceModels diff --git a/applications/solvers/multiphase/MPPICInterFoam/CompressibleTwoPhaseMixtureTurbulenceModels/Make/options b/applications/solvers/multiphase/MPPICInterFoam/CompressibleTwoPhaseMixtureTurbulenceModels/Make/options new file mode 100644 index 0000000000..03cad33b12 --- /dev/null +++ b/applications/solvers/multiphase/MPPICInterFoam/CompressibleTwoPhaseMixtureTurbulenceModels/Make/options @@ -0,0 +1,15 @@ +EXE_INC = \ + -I$(LIB_SRC)/finiteVolume/lnInclude \ + -I$(LIB_SRC)/meshTools/lnInclude \ + -I$(LIB_SRC)/transportModels/compressible/lnInclude \ + -I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \ + -I$(LIB_SRC)/transportModels \ + -I$(LIB_SRC)/transportModels/immiscibleIncompressibleTwoPhaseMixture/lnInclude \ + -I$(LIB_SRC)/transportModels/twoPhaseMixture/lnInclude \ + -I$(LIB_SRC)/transportModels/incompressible/lnInclude \ + -I$(LIB_SRC)/transportModels/interfaceProperties/lnInclude \ + -I$(LIB_SRC)/TurbulenceModels/turbulenceModels/lnInclude \ + -I$(LIB_SRC)/TurbulenceModels/incompressible/lnInclude \ + -I$(LIB_SRC)/TurbulenceModels/compressible/lnInclude \ + -I$(LIB_SRC)/TurbulenceModels/phaseCompressible/lnInclude \ + -I$(LIB_SRC)/TurbulenceModels/phaseIncompressible/lnInclude diff --git a/applications/solvers/multiphase/MPPICInterFoam/MPPICInterFoam.C b/applications/solvers/multiphase/MPPICInterFoam/MPPICInterFoam.C new file mode 100644 index 0000000000..35797278d7 --- /dev/null +++ b/applications/solvers/multiphase/MPPICInterFoam/MPPICInterFoam.C @@ -0,0 +1,167 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2016 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 . + +Application + MPPICInterFoam + +Description + Solver for 2 incompressible, isothermal immiscible fluids using a VOF + (volume of fluid) phase-fraction based interface capturing approach. + The momentum and other fluid properties are of the "mixture" and a single + momentum equation is solved. + + It includes MRF and MPPIC clouds + + Turbulence modelling is generic, i.e. laminar, RAS or LES may be selected. + +\*---------------------------------------------------------------------------*/ + +#include "fvCFD.H" +#include "CMULES.H" +#include "subCycle.H" + +#include "immiscibleIncompressibleTwoPhaseMixture.H" +#include "PhaseCompressibleTurbulenceModel.H" +#include "pimpleControl.H" +#include "fvOptions.H" +#include "CorrectPhi.H" +#include "fixedFluxPressureFvPatchScalarField.H" + +#include "basicKinematicMPPICCloud.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +int main(int argc, char *argv[]) +{ + #include "setRootCase.H" + #include "createTime.H" + #include "createMesh.H" + + pimpleControl pimple(mesh); + + #include "createTimeControls.H" + #include "initContinuityErrs.H" + #include "createFields.H" + #include "createMRF.H" + #include "createFvOptions.H" + #include "correctPhi.H" + + turbulence->validate(); + + + #include "readTimeControls.H" + #include "CourantNo.H" + #include "setInitialDeltaT.H" + + // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + + Info<< "\nStarting time loop\n" << endl; + + while (runTime.run()) + { + #include "readTimeControls.H" + #include "CourantNo.H" + #include "alphaCourantNo.H" + #include "setDeltaT.H" + + runTime++; + + Info<< "Time = " << runTime.timeName() << nl << endl; + + Info<< "Evolving " << kinematicCloud.name() << endl; + + kinematicCloud.evolve(); + + // Update continuous phase volume fraction field + alphac = max(1.0 - kinematicCloud.theta(), alphacMin); + alphac.correctBoundaryConditions(); + + Info<< "Continous phase-1 volume fraction = " + << alphac.weightedAverage(mesh.Vsc()).value() + << " Min(alphac) = " << min(alphac).value() + << " Max(alphac) = " << max(alphac).value() + << endl; + + alphacf = fvc::interpolate(alphac); + alphaRhoPhic = alphacf*rhoPhi; + alphaPhic = alphacf*phi; + alphacRho = alphac*rho; + + fvVectorMatrix cloudSU(kinematicCloud.SU(U)); + volVectorField cloudVolSUSu + ( + IOobject + ( + "cloudVolSUSu", + runTime.timeName(), + mesh + ), + mesh, + dimensionedVector + ( + "0", + cloudSU.dimensions()/dimVolume, + vector::zero + ), + zeroGradientFvPatchVectorField::typeName + ); + + cloudVolSUSu.internalField() = -cloudSU.source()/mesh.V(); + cloudVolSUSu.correctBoundaryConditions(); + + cloudSU.source() = vector::zero; + + // --- Pressure-velocity PIMPLE corrector loop + while (pimple.loop()) + { + #include "alphaControls.H" + #include "alphaEqnSubCycle.H" + + #include "UEqn.H" + + // --- Pressure corrector loop + while (pimple.correct()) + { + #include "pEqn.H" + } + + if (pimple.turbCorr()) + { + turbulence->correct(); + } + } + + runTime.write(); + + Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s" + << " ClockTime = " << runTime.elapsedClockTime() << " s" + << nl << endl; + } + + Info<< "End\n" << endl; + + return 0; +} + + +// ************************************************************************* // diff --git a/applications/solvers/multiphase/MPPICInterFoam/Make/files b/applications/solvers/multiphase/MPPICInterFoam/Make/files new file mode 100644 index 0000000000..cbedde119e --- /dev/null +++ b/applications/solvers/multiphase/MPPICInterFoam/Make/files @@ -0,0 +1,3 @@ +MPPICInterFoam.C + +EXE = $(FOAM_APPBIN)/MPPICInterFoam diff --git a/applications/solvers/multiphase/MPPICInterFoam/Make/options b/applications/solvers/multiphase/MPPICInterFoam/Make/options new file mode 100644 index 0000000000..8ecd0b4a1e --- /dev/null +++ b/applications/solvers/multiphase/MPPICInterFoam/Make/options @@ -0,0 +1,42 @@ +interFoamPath = $(FOAM_SOLVERS)/multiphase/interFoam + +EXE_INC = \ + -I./IncompressibleTwoPhaseMixtureTurbulenceModels/lnInclude \ + -I$(interFoamPath) \ + -I$(LIB_SRC)/finiteVolume/lnInclude \ + -I$(LIB_SRC)/fvOptions/lnInclude \ + -I$(LIB_SRC)/meshTools/lnInclude \ + -I$(LIB_SRC)/sampling/lnInclude \ + -I$(LIB_SRC)/lagrangian/basic/lnInclude \ + -I$(LIB_SRC)/lagrangian/intermediate/lnInclude \ + -I$(LIB_SRC)/thermophysicalModels/specie/lnInclude \ + -I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \ + -I$(LIB_SRC)/thermophysicalModels/reactionThermo/lnInclude \ + -I$(LIB_SRC)/thermophysicalModels/radiation/lnInclude \ + -I$(LIB_SRC)/transportModels \ + -I$(LIB_SRC)/transportModels/immiscibleIncompressibleTwoPhaseMixture/lnInclude \ + -I$(LIB_SRC)/transportModels/twoPhaseMixture/lnInclude \ + -I$(LIB_SRC)/transportModels/incompressible/lnInclude \ + -I$(LIB_SRC)/transportModels/interfaceProperties/lnInclude \ + -I$(LIB_SRC)/TurbulenceModels/turbulenceModels/lnInclude \ + -I$(LIB_SRC)/TurbulenceModels/compressible/lnInclude \ + -I$(LIB_SRC)/TurbulenceModels/phaseCompressible/lnInclude \ + -I$(LIB_SRC)/regionModels/regionModel/lnInclude \ + -I$(LIB_SRC)/regionModels/surfaceFilmModels/lnInclude + + +EXE_LIBS = \ + -lfiniteVolume \ + -lmeshTools \ + -llagrangian \ + -llagrangianIntermediate \ + -lthermophysicalFunctions \ + -lspecie \ + -lincompressibleTransportModels \ + -limmiscibleIncompressibleTwoPhaseMixture \ + -lturbulenceModels \ + -lsampling \ + -lregionModels \ + -lsurfaceFilmModels \ + -lfvOptions \ + -lCompressibleTwoPhaseMixtureTurbulenceModels diff --git a/applications/solvers/multiphase/MPPICInterFoam/UEqn.H b/applications/solvers/multiphase/MPPICInterFoam/UEqn.H new file mode 100644 index 0000000000..e402e3a5d5 --- /dev/null +++ b/applications/solvers/multiphase/MPPICInterFoam/UEqn.H @@ -0,0 +1,50 @@ +MRF.correctBoundaryVelocity(U); + +fvVectorMatrix UEqn +( + fvm::ddt(alphacRho, U) + + MRF.DDt(alphacRho, U) + - fvm::Sp(fvc::ddt(rho) + fvc::div(rhoPhi), U) + + fvm::div(rhoPhi, U) + + turbulence->divDevRhoReff(U) + == + fvOptions(rho, U) + + cloudSU +); + +UEqn.relax(); +fvOptions.constrain(UEqn); + +volScalarField rAUc(1.0/UEqn.A()); +surfaceScalarField rAUcf(fvc::interpolate(rAUc)); + + +surfaceScalarField phicForces +( + (fvc::interpolate(rAUc*cloudVolSUSu) & mesh.Sf()) +); + + +if (pimple.momentumPredictor()) +{ + solve + ( + UEqn + == + fvc::reconstruct + ( + phicForces/rAUcf + + + ( + fvc::interpolate + ( + mixture.sigmaK() + )*fvc::snGrad(alpha1) + - ghf*fvc::snGrad(rho) + - fvc::snGrad(p_rgh) + ) * mesh.magSf() + ) + ); + + fvOptions.correct(U); +} diff --git a/applications/solvers/multiphase/MPPICInterFoam/alphaEqn.H b/applications/solvers/multiphase/MPPICInterFoam/alphaEqn.H new file mode 100644 index 0000000000..a177326504 --- /dev/null +++ b/applications/solvers/multiphase/MPPICInterFoam/alphaEqn.H @@ -0,0 +1,163 @@ +{ + word alphaScheme("div(phi,alpha)"); + word alpharScheme("div(phirb,alpha)"); + + // Standard face-flux compression coefficient + surfaceScalarField phic + ( + mixture.cAlpha()*mag(alphaPhic/mesh.magSf()) + ); + + // Add the optional isotropic compression contribution + if (icAlpha > 0) + { + phic *= (1.0 - icAlpha); + phic += (mixture.cAlpha()*icAlpha)*fvc::interpolate(mag(U)); + } + + // Do not compress interface at non-coupled boundary faces + // (inlets, outlets etc.) + forAll(phic.boundaryField(), patchi) + { + fvsPatchScalarField& phicp = phic.boundaryField()[patchi]; + + if (!phicp.coupled()) + { + phicp == 0; + } + } + + tmp tphiAlpha; + + if (MULESCorr) + { + fvScalarMatrix alpha1Eqn + ( + fv::EulerDdtScheme(mesh).fvmDdt(alphac, alpha1) + + fv::gaussConvectionScheme + ( + mesh, + alphaPhic, + upwind(mesh, alphaPhic) + ).fvmDiv(alphaPhic, alpha1) + - fvm::Sp(fvc::ddt(alphac) + fvc::div(alphaPhic), alpha1) + ); + + alpha1Eqn.solve(); + + Info<< "Phase-1 volume fraction = " + << alpha1.weightedAverage(mesh.Vsc()).value() + << " Min(alpha1) = " << min(alpha1).value() + << " Max(alpha1) = " << max(alpha1).value() + << endl; + + tmp tphiAlphaUD(alpha1Eqn.flux()); + alphaPhi = tphiAlphaUD(); + + if (alphaApplyPrevCorr && tphiAlphaCorr0.valid()) + { + Info<< "Applying the previous iteration compression flux" << endl; + + MULES::correct + ( + alphac, + alpha1, + alphaPhi, + tphiAlphaCorr0.ref(), + zeroField(), zeroField(), + 1, 0 + ); + + alphaPhi += tphiAlphaCorr0(); + } + + // Cache the upwind-flux + tphiAlphaCorr0 = tphiAlphaUD; + + alpha2 = 1.0 - alpha1; + + mixture.correct(); + } + + for (int aCorr=0; aCorr tphiAlphaUn + ( + fvc::flux + ( + alphaPhic, + alpha1, + alphaScheme + ) + + fvc::flux + ( + -fvc::flux(-phir, alpha2, alpharScheme), + alpha1, + alpharScheme + ) + ); + + if (MULESCorr) + { + tmp tphiAlphaCorr(tphiAlphaUn() - alphaPhi); + volScalarField alpha10("alpha10", alpha1); + + //MULES::correct(alpha1, tphiAlphaUn(), tphiAlphaCorr(), 1, 0); + + MULES::correct + ( + alphac, + alpha1, + tphiAlphaUn(), + tphiAlphaCorr.ref(), + zeroField(), zeroField(), + 1, 0 + ); + + // 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 + ( + alphac, + alpha1, + alphaPhic, + alphaPhi, + zeroField(), zeroField(), + 1, 0 + ); + + } + + alpha2 = 1.0 - alpha1; + + mixture.correct(); + } + + rhoPhi = alphaPhi*(rho1 - rho2) + alphaPhic*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; +} diff --git a/applications/solvers/multiphase/MPPICInterFoam/alphaEqnSubCycle.H b/applications/solvers/multiphase/MPPICInterFoam/alphaEqnSubCycle.H new file mode 100644 index 0000000000..c0d3c8e43a --- /dev/null +++ b/applications/solvers/multiphase/MPPICInterFoam/alphaEqnSubCycle.H @@ -0,0 +1,34 @@ +if (nAlphaSubCycles > 1) +{ + dimensionedScalar totalDeltaT = runTime.deltaT(); + surfaceScalarField rhoPhiSum + ( + IOobject + ( + "rhoPhiSum", + runTime.timeName(), + mesh + ), + mesh, + dimensionedScalar("0", rhoPhi.dimensions(), 0) + ); + + for + ( + subCycle alphaSubCycle(alpha1, nAlphaSubCycles); + !(++alphaSubCycle).end(); + ) + { + #include "alphaEqn.H" + rhoPhiSum += (runTime.deltaT()/totalDeltaT)*rhoPhi; + } + + rhoPhi = rhoPhiSum; +} +else +{ + #include "alphaEqn.H" +} + +rho == alpha1*rho1 + alpha2*rho2; +mu = mixture.mu(); diff --git a/applications/solvers/multiphase/MPPICInterFoam/continuityErrs.H b/applications/solvers/multiphase/MPPICInterFoam/continuityErrs.H new file mode 100644 index 0000000000..cb806bed7b --- /dev/null +++ b/applications/solvers/multiphase/MPPICInterFoam/continuityErrs.H @@ -0,0 +1,48 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2015 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 . + +Global + continuityErrs + +Description + Calculates and prints the continuity errors. + +\*---------------------------------------------------------------------------*/ + +{ + volScalarField contErr(fvc::ddt(alphac) + fvc::div(alphacf*phi)); + + scalar sumLocalContErr = runTime.deltaTValue()* + mag(contErr)().weightedAverage(mesh.V()).value(); + + scalar globalContErr = runTime.deltaTValue()* + contErr.weightedAverage(mesh.V()).value(); + cumulativeContErr += globalContErr; + + Info<< "time step continuity errors : sum local = " << sumLocalContErr + << ", global = " << globalContErr + << ", cumulative = " << cumulativeContErr + << endl; +} + +// ************************************************************************* // diff --git a/applications/solvers/multiphase/MPPICInterFoam/correctPhi.H b/applications/solvers/multiphase/MPPICInterFoam/correctPhi.H new file mode 100644 index 0000000000..e54630529a --- /dev/null +++ b/applications/solvers/multiphase/MPPICInterFoam/correctPhi.H @@ -0,0 +1,11 @@ +CorrectPhi +( + U, + phi, + p_rgh, + dimensionedScalar("rAUf", dimTime/rho.dimensions(), 1), + geometricZeroField(), + pimple +); + +#include "continuityErrs.H" \ No newline at end of file diff --git a/applications/solvers/multiphase/MPPICInterFoam/createFields.H b/applications/solvers/multiphase/MPPICInterFoam/createFields.H new file mode 100644 index 0000000000..a7d456f532 --- /dev/null +++ b/applications/solvers/multiphase/MPPICInterFoam/createFields.H @@ -0,0 +1,219 @@ + Info<< "Reading field p_rgh\n" << endl; + volScalarField p_rgh + ( + IOobject + ( + "p_rgh", + runTime.timeName(), + mesh, + IOobject::MUST_READ, + IOobject::AUTO_WRITE + ), + mesh + ); + + Info<< "Reading field U\n" << endl; + volVectorField U + ( + IOobject + ( + "U", + runTime.timeName(), + mesh, + IOobject::MUST_READ, + IOobject::AUTO_WRITE + ), + mesh + ); + + #include "createPhi.H" + + Info<< "Reading transportProperties\n" << endl; + immiscibleIncompressibleTwoPhaseMixture mixture(U, phi); + + volScalarField& alpha1(mixture.alpha1()); + volScalarField& alpha2(mixture.alpha2()); + + const dimensionedScalar& rho1 = mixture.rho1(); + const dimensionedScalar& rho2 = mixture.rho2(); + + // Need to store rho for ddt(rho, U) + volScalarField rho + ( + IOobject + ( + "rho", + runTime.timeName(), + mesh, + IOobject::READ_IF_PRESENT, + IOobject::AUTO_WRITE + ), + alpha1*rho1 + alpha2*rho2 + ); + rho.oldTime(); + + // Need to store mu as incompressibleTwoPhaseMixture does not store it + volScalarField mu + ( + IOobject + ( + "mu", + runTime.timeName(), + mesh, + IOobject::READ_IF_PRESENT + ), + mixture.mu(), + calculatedFvPatchScalarField::typeName + ); + + + // Mass flux + surfaceScalarField rhoPhi + ( + IOobject + ( + "rhoPhi", + runTime.timeName(), + mesh, + IOobject::NO_READ, + IOobject::NO_WRITE + ), + fvc::interpolate(rho)*phi + ); + + + #include "readGravitationalAcceleration.H" + + #include "readhRef.H" + #include "gh.H" + + volScalarField p + ( + IOobject + ( + "p", + runTime.timeName(), + mesh, + IOobject::NO_READ, + IOobject::AUTO_WRITE + ), + p_rgh + rho*gh + ); + + label pRefCell = 0; + scalar pRefValue = 0.0; + setRefCell + ( + p, + p_rgh, + mesh.solutionDict().subDict("PIMPLE"), + pRefCell, + pRefValue + ); + + if (p_rgh.needReference()) + { + p += dimensionedScalar + ( + "p", + p.dimensions(), + pRefValue - getRefCellValue(p, pRefCell) + ); + p_rgh = p - rho*gh; + } + + mesh.setFluxRequired(p_rgh.name()); + mesh.setFluxRequired(alpha1.name()); + + // MULES flux from previous time-step + surfaceScalarField alphaPhi + ( + IOobject + ( + "alphaPhi", + runTime.timeName(), + mesh, + IOobject::READ_IF_PRESENT, + IOobject::AUTO_WRITE + ), + phi*fvc::interpolate(alpha1) + ); + + + tmp tphiAlphaCorr0; + + // alphac must be constructed before the cloud + // so that the drag-models can find it + volScalarField alphac + ( + IOobject + ( + "alphac", + runTime.timeName(), + mesh, + IOobject::READ_IF_PRESENT, + IOobject::AUTO_WRITE + ), + mesh, + dimensionedScalar("0", dimless, 0), + zeroGradientFvPatchScalarField::typeName + ); + alphac.oldTime(); + + volScalarField alphacRho(alphac*rho); + alphacRho.oldTime(); + + Info<< "Constructing kinematicCloud " << endl; + basicKinematicMPPICCloud kinematicCloud + ( + "kinematicCloud", + rho, + U, + mu, + g + ); + + // Particle fraction upper limit + scalar alphacMin + ( + 1.0 + - readScalar + ( + kinematicCloud.particleProperties().subDict("constantProperties") + .lookup("alphaMax") + ) + ); + + // Update alphac from the particle locations + alphac = max(1.0 - kinematicCloud.theta(), alphacMin); + alphac.correctBoundaryConditions(); + + surfaceScalarField alphacf("alphacf", fvc::interpolate(alphac)); + + // Phase mass flux + surfaceScalarField alphaRhoPhic("alphaRhoPhic", alphacf*rhoPhi); + + // Volumetric phase flux + surfaceScalarField alphaPhic("alphaPhic", alphacf*phi); + + autoPtr + < + PhaseCompressibleTurbulenceModel + < + immiscibleIncompressibleTwoPhaseMixture + > + >turbulence + ( + PhaseCompressibleTurbulenceModel + < + immiscibleIncompressibleTwoPhaseMixture + >::New + ( + alphac, + rho, + U, + alphaRhoPhic, + rhoPhi, + mixture + ) + ); diff --git a/applications/solvers/multiphase/MPPICInterFoam/pEqn.H b/applications/solvers/multiphase/MPPICInterFoam/pEqn.H new file mode 100644 index 0000000000..2e1330b3b8 --- /dev/null +++ b/applications/solvers/multiphase/MPPICInterFoam/pEqn.H @@ -0,0 +1,73 @@ +{ + volVectorField HbyA(constrainHbyA(rAUc*UEqn.H(), U, p_rgh)); + + surfaceScalarField phiHbyA + ( + "phiHbyA", + fvc::flux(HbyA) + + alphacf*fvc::interpolate(rho*rAUc)*fvc::ddtCorr(U, phi) + ); + + MRF.makeRelative(phiHbyA); + adjustPhi(phiHbyA, U, p_rgh); + + surfaceScalarField phig + ( + phicForces + + ( + fvc::interpolate(mixture.sigmaK())*fvc::snGrad(alpha1) + - ghf*fvc::snGrad(rho) + )*rAUcf*mesh.magSf() + ); + + phiHbyA += phig; + + // Update the pressure BCs to ensure flux consistency + constrainPressure(p_rgh, U, phiHbyA, rAUcf, MRF); + + while (pimple.correctNonOrthogonal()) + { + surfaceScalarField Dp("Dp", alphacf*rAUcf); + + fvScalarMatrix p_rghEqn + ( + fvm::laplacian(Dp, p_rgh) + == + fvc::ddt(alphac) + fvc::div(alphacf*phiHbyA) + ); + + p_rghEqn.setReference(pRefCell, getRefCellValue(p_rgh, pRefCell)); + + p_rghEqn.solve(mesh.solver(p_rgh.select(pimple.finalInnerIter()))); + + if (pimple.finalNonOrthogonalIter()) + { + phi = phiHbyA - p_rghEqn.flux()/alphacf; + + p_rgh.relax(); + + U = + HbyA + + rAUc*fvc::reconstruct((phig - p_rghEqn.flux()/alphacf)/rAUcf); + + U.correctBoundaryConditions(); + + fvOptions.correct(U); + } + } + + #include "continuityErrs.H" + + p == p_rgh + rho*gh; + + if (p_rgh.needReference()) + { + p += dimensionedScalar + ( + "p", + p.dimensions(), + pRefValue - getRefCellValue(p, pRefCell) + ); + p_rgh = p - rho*gh; + } +} diff --git a/applications/solvers/multiphase/interCondensingEvaporatingFoam/Allwclean b/applications/solvers/multiphase/interCondensingEvaporatingFoam/Allwclean new file mode 100755 index 0000000000..8717ef8953 --- /dev/null +++ b/applications/solvers/multiphase/interCondensingEvaporatingFoam/Allwclean @@ -0,0 +1,8 @@ +#!/bin/sh +cd ${0%/*} || exit 1 # Run from this directory +set -x + +wclean libso temperaturePhaseChangeTwoPhaseMixtures +wclean + +#------------------------------------------------------------------------------ diff --git a/applications/solvers/multiphase/interCondensingEvaporatingFoam/Allwmake b/applications/solvers/multiphase/interCondensingEvaporatingFoam/Allwmake new file mode 100755 index 0000000000..523c10c33c --- /dev/null +++ b/applications/solvers/multiphase/interCondensingEvaporatingFoam/Allwmake @@ -0,0 +1,12 @@ +#!/bin/sh +cd ${0%/*} || exit 1 # Run from this directory + +# Parse arguments for library compilation +targetType=libso +. $WM_PROJECT_DIR/wmake/scripts/AllwmakeParseArguments +set -x + +wmake $targetType temperaturePhaseChangeTwoPhaseMixtures +wmake + +#------------------------------------------------------------------------------ diff --git a/applications/solvers/multiphase/interCondensingEvaporatingFoam/Make/files b/applications/solvers/multiphase/interCondensingEvaporatingFoam/Make/files new file mode 100644 index 0000000000..7e0882ec64 --- /dev/null +++ b/applications/solvers/multiphase/interCondensingEvaporatingFoam/Make/files @@ -0,0 +1,3 @@ +interCondensatingEvaporatingFoam.C + +EXE = $(FOAM_APPBIN)/interCondensatingEvaporatingFoam diff --git a/applications/solvers/multiphase/interCondensingEvaporatingFoam/Make/options b/applications/solvers/multiphase/interCondensingEvaporatingFoam/Make/options new file mode 100644 index 0000000000..373084db4b --- /dev/null +++ b/applications/solvers/multiphase/interCondensingEvaporatingFoam/Make/options @@ -0,0 +1,30 @@ +interPhaseChangePath = $(FOAM_SOLVERS)/multiphase/interPhaseChangeFoam + +EXE_INC = \ + -ItemperaturePhaseChangeTwoPhaseMixtures/lnInclude \ + -I$(interPhaseChangePath) \ + -I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \ + -I$(LIB_SRC)/transportModels/twoPhaseMixture/lnInclude \ + -I$(LIB_SRC)/transportModels \ + -I$(LIB_SRC)/transportModels/incompressible/lnInclude \ + -I$(LIB_SRC)/transportModels/interfaceProperties/lnInclude \ + -I$(LIB_SRC)/TurbulenceModels/turbulenceModels/lnInclude \ + -I$(LIB_SRC)/TurbulenceModels/incompressible/lnInclude \ + -I$(LIB_SRC)/finiteVolume/lnInclude \ + -I$(LIB_SRC)/meshTools/lnInclude \ + -I$(LIB_SRC)/fvOptions/lnInclude\ + -I$(LIB_SRC)/sampling/lnInclude + +EXE_LIBS = \ + -lphaseTemperatureChangeTwoPhaseMixtures \ + -ltwoPhaseMixture \ + -linterfaceProperties \ + -ltwoPhaseProperties \ + -lincompressibleTransportModels \ + -lturbulenceModels \ + -lincompressibleTurbulenceModels \ + -lfiniteVolume \ + -lmeshTools \ + -lfvOptions \ + -lsampling \ + -lfluidThermophysicalModels diff --git a/applications/solvers/multiphase/interCondensingEvaporatingFoam/UEqn.H b/applications/solvers/multiphase/interCondensingEvaporatingFoam/UEqn.H new file mode 100644 index 0000000000..1e27597292 --- /dev/null +++ b/applications/solvers/multiphase/interCondensingEvaporatingFoam/UEqn.H @@ -0,0 +1,26 @@ + fvVectorMatrix UEqn + ( + fvm::ddt(rho, U) + + fvm::div(rhoPhi, U) + - fvm::Sp(fvc::ddt(rho) + fvc::div(rhoPhi), U) + + turbulence->divDevRhoReff(rho, U) + ); + + UEqn.relax(); + + if (pimple.momentumPredictor()) + { + solve + ( + UEqn + == + fvc::reconstruct + ( + ( + interface.surfaceTensionForce() + - ghf*fvc::snGrad(rho) + - fvc::snGrad(p_rgh) + ) * mesh.magSf() + ) + ); + } diff --git a/applications/solvers/multiphase/interCondensingEvaporatingFoam/continuityError.H b/applications/solvers/multiphase/interCondensingEvaporatingFoam/continuityError.H new file mode 100644 index 0000000000..9a33fb896d --- /dev/null +++ b/applications/solvers/multiphase/interCondensingEvaporatingFoam/continuityError.H @@ -0,0 +1,12 @@ +volScalarField contErr(fvc::ddt(rho) + fvc::div(rhoPhi)); + +scalar sumLocalContErr = runTime.deltaTValue()* + mag(contErr)().weightedAverage(mesh.V()).value(); + +scalar globalContErr = runTime.deltaTValue()* + contErr.weightedAverage(mesh.V()).value(); + +Info<< "time step continuity errors : sum local = " << sumLocalContErr + << ", global = " << globalContErr + << endl; + diff --git a/applications/solvers/multiphase/interCondensingEvaporatingFoam/createFields.H b/applications/solvers/multiphase/interCondensingEvaporatingFoam/createFields.H new file mode 100644 index 0000000000..20d3677bfc --- /dev/null +++ b/applications/solvers/multiphase/interCondensingEvaporatingFoam/createFields.H @@ -0,0 +1,158 @@ + Info<< "Reading field p_rgh\n" << endl; + volScalarField p_rgh + ( + IOobject + ( + "p_rgh", + runTime.timeName(), + mesh, + IOobject::MUST_READ, + IOobject::AUTO_WRITE + ), + mesh + ); + + Info<< "Reading field U\n" << endl; + volVectorField U + ( + IOobject + ( + "U", + runTime.timeName(), + mesh, + IOobject::MUST_READ, + IOobject::AUTO_WRITE + ), + mesh + ); + + #include "createPhi.H" + + // Create p before the thermo + volScalarField p + ( + IOobject + ( + "p", + runTime.timeName(), + mesh, + IOobject::NO_READ, + IOobject::AUTO_WRITE + ), + p_rgh + ); + + // Creating e based thermo + autoPtr thermo; + thermo.set(new twoPhaseMixtureEThermo(U, phi)); + + // Create mixture and + Info<< "Creating temperaturePhaseChangeTwoPhaseMixture\n" << endl; + autoPtr mixture = + temperaturePhaseChangeTwoPhaseMixture::New(thermo(), mesh); + + + volScalarField& T = thermo->T(); + volScalarField& e = thermo->he(); + + // Correct e from T and alpha + thermo->correct(); + + volScalarField& alpha1(thermo->alpha1()); + volScalarField& alpha2(thermo->alpha2()); + + const dimensionedScalar& rho1 = thermo->rho1(); + const dimensionedScalar& rho2 = thermo->rho2(); + + // Need to store rho for ddt(rho, U) + volScalarField rho + ( + IOobject + ( + "rho", + runTime.timeName(), + mesh, + IOobject::READ_IF_PRESENT, + IOobject::AUTO_WRITE + ), + alpha1*rho1 + alpha2*rho2, + alpha1.boundaryField().types() + ); + rho.oldTime(); + + + // Construct interface from alpha1 distribution + interfaceProperties interface + ( + alpha1, + U, + thermo->transportPropertiesDict() + ); + + // Construct incompressible turbulence model + autoPtr turbulence + ( + incompressible::turbulenceModel::New(U, phi, thermo()) + ); + + + Info<< "Calculating field g.h\n" << endl; + volScalarField gh("gh", g & mesh.C()); + surfaceScalarField ghf("ghf", g & mesh.Cf()); + + //Update p with rho + p = p_rgh + rho*gh; + + label pRefCell = 0; + scalar pRefValue = 0.0; + setRefCell + ( + p, + p_rgh, + pimple.dict(), + pRefCell, + pRefValue + ); + + if (p_rgh.needReference()) + { + p += dimensionedScalar + ( + "p", + p.dimensions(), + pRefValue - getRefCellValue(p, pRefCell) + ); + p_rgh = p - rho*gh; + } + + // Turbulent Prandtl number + dimensionedScalar Prt("Prt", dimless, thermo->transportPropertiesDict()); + + volScalarField kappaEff + ( + IOobject + ( + "kappaEff", + runTime.timeName(), + mesh, + IOobject::NO_READ, + IOobject::NO_WRITE + ), + thermo->kappa() + ); + + Info<< "Creating field kinetic energy K\n" << endl; + volScalarField K("K", 0.5*magSqr(U)); + + Info<< "Creating field pDivU\n" << endl; + volScalarField pDivU + ( + IOobject + ( + "pDivU", + runTime.timeName(), + mesh + ), + mesh, + dimensionedScalar("pDivU", p.dimensions()/dimTime, 0) + ); diff --git a/applications/solvers/multiphase/interCondensingEvaporatingFoam/eEqn.H b/applications/solvers/multiphase/interCondensingEvaporatingFoam/eEqn.H new file mode 100644 index 0000000000..291d574d8a --- /dev/null +++ b/applications/solvers/multiphase/interCondensingEvaporatingFoam/eEqn.H @@ -0,0 +1,31 @@ +{ + tmp tcp(thermo->Cp()); + const volScalarField& cp = tcp(); + + kappaEff = thermo->kappa() + rho*cp*turbulence->nut()/Prt; + + pDivU = dimensionedScalar("pDivU", p.dimensions()/dimTime, 0.0); + + if (thermo->pDivU()) + { + pDivU = (p*fvc::div(rhoPhi/fvc::interpolate(rho))); + } + + fvScalarMatrix eEqn + ( + fvm::ddt(rho, e) + + fvc::ddt(rho, K) + fvc::div(rhoPhi, K) + + fvm::div(rhoPhi, e) + - fvm::Sp(fvc::ddt(rho) + fvc::div(rhoPhi), e) + - fvm::laplacian(kappaEff/cp, e) + + pDivU + ); + + eEqn.relax(); + eEqn.solve(); + + thermo->correct(); + + Info<< "min/max(T) = " << min(T).value() << ", " + << max(T).value() <. + +Application + interCondensatingEvaporatingFoam + +Group + grpMultiphaseSolvers + +Description + Solver for 2 incompressible, non-isothermal immiscible fluids with + phase-change (evaporation-condensation) between a fluid and its vapour. + Uses a VOF (volume of fluid) phase-fraction based interface capturing + approach. + + The momentum, energy and other fluid properties are of the "mixture" and a + single momentum equation is solved. + + Turbulence modelling is generic, i.e. laminar, RAS or LES may be selected. + +\*---------------------------------------------------------------------------*/ + +#include "fvCFD.H" +#include "CMULES.H" +#include "subCycle.H" +#include "interfaceProperties.H" +#include "twoPhaseMixtureEThermo.H" +#include "temperaturePhaseChangeTwoPhaseMixture.H" +#include "turbulentTransportModel.H" +#include "turbulenceModel.H" +#include "pimpleControl.H" +#include "fvOptions.H" +#include "fixedFluxPressureFvPatchScalarField.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +int main(int argc, char *argv[]) +{ + #include "setRootCase.H" + #include "createTime.H" + #include "createMesh.H" + + pimpleControl pimple(mesh); + + #include "readGravitationalAcceleration.H" + #include "initContinuityErrs.H" + #include "createFields.H" + #include "createFvOptions.H" + #include "createTimeControls.H" + #include "CourantNo.H" + #include "setInitialDeltaT.H" + + // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + + Info<< "\nStarting time loop\n" << endl; + + while (runTime.run()) + { + #include "readTimeControls.H" + #include "CourantNo.H" + #include "setDeltaT.H" + + runTime++; + + Info<< "Time = " << runTime.timeName() << nl << endl; + + // --- Pressure-velocity PIMPLE corrector loop + while (pimple.loop()) + { + #include "alphaControls.H" + + surfaceScalarField rhoPhi + ( + IOobject + ( + "rhoPhi", + runTime.timeName(), + mesh + ), + mesh, + dimensionedScalar("0", dimMass/dimTime, 0) + ); + + #include "alphaEqnSubCycle.H" + mixture->correct(); + + #include "UEqn.H" + #include "eEqn.H" + + // --- Pressure corrector loop + while (pimple.correct()) + { + #include "pEqn.H" + } + + #include "continuityError.H" + + if (pimple.turbCorr()) + { + turbulence->correct(); + } + } + + runTime.write(); + + Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s" + << " ClockTime = " << runTime.elapsedClockTime() << " s" + << nl << endl; + } + + Info<< "End\n" << endl; + + return 0; +} + + +// ************************************************************************* // diff --git a/applications/solvers/multiphase/interCondensingEvaporatingFoam/pEqn.H b/applications/solvers/multiphase/interCondensingEvaporatingFoam/pEqn.H new file mode 100644 index 0000000000..c2b043c57b --- /dev/null +++ b/applications/solvers/multiphase/interCondensingEvaporatingFoam/pEqn.H @@ -0,0 +1,69 @@ +{ + volScalarField rAU("rAU", 1.0/UEqn.A()); + surfaceScalarField rAUf("rAUf", fvc::interpolate(rAU)); + + volVectorField HbyA(constrainHbyA(rAU*UEqn.H(), U, p_rgh)); + + surfaceScalarField phiHbyA + ( + "phiHbyA", + (fvc::interpolate(HbyA) & mesh.Sf()) + + fvc::interpolate(rho*rAU)*fvc::ddtCorr(U, phi) + ); + adjustPhi(phiHbyA, U, p_rgh); + + surfaceScalarField phig + ( + ( + interface.surfaceTensionForce() + - ghf*fvc::snGrad(rho) + )*rAUf*mesh.magSf() + ); + + phiHbyA += phig; + + // Update the pressure BCs to ensure flux consistency + constrainPressure(p_rgh, U, phiHbyA, rAUf); + + Pair > vDot = mixture->vDot(); + const volScalarField& vDotc = vDot[0](); + const volScalarField& vDotv = vDot[1](); + + while (pimple.correctNonOrthogonal()) + { + fvScalarMatrix p_rghEqn + ( + fvc::div(phiHbyA) + - fvm::laplacian(rAUf, p_rgh) + - (vDotc - vDotv) + ); + + p_rghEqn.setReference(pRefCell, pRefValue); + + p_rghEqn.solve(mesh.solver(p_rgh.select(pimple.finalInnerIter()))); + + if (pimple.finalNonOrthogonalIter()) + { + phi = phiHbyA + p_rghEqn.flux(); + + U = HbyA + rAU*fvc::reconstruct((phig + p_rghEqn.flux())/rAUf); + U.correctBoundaryConditions(); + fvOptions.correct(U); + K = 0.5*magSqr(U); + } + } + + + p == p_rgh + rho*gh; + + if (p_rgh.needReference()) + { + p += dimensionedScalar + ( + "p", + p.dimensions(), + pRefValue - getRefCellValue(p, pRefCell) + ); + p_rgh = p - rho*gh; + } +} diff --git a/applications/solvers/multiphase/interCondensingEvaporatingFoam/temperaturePhaseChangeTwoPhaseMixtures/Make/files b/applications/solvers/multiphase/interCondensingEvaporatingFoam/temperaturePhaseChangeTwoPhaseMixtures/Make/files new file mode 100644 index 0000000000..407c604655 --- /dev/null +++ b/applications/solvers/multiphase/interCondensingEvaporatingFoam/temperaturePhaseChangeTwoPhaseMixtures/Make/files @@ -0,0 +1,7 @@ +temperaturePhaseChangeTwoPhaseMixtures/newtemperaturePhaseChangeTwoPhaseMixture.C +temperaturePhaseChangeTwoPhaseMixtures/temperaturePhaseChangeTwoPhaseMixture.C +thermoIncompressibleTwoPhaseMixture/thermoIncompressibleTwoPhaseMixture.C +twoPhaseMixtureEThermo/twoPhaseMixtureEThermo.C +constant/constant.C + +LIB = $(FOAM_LIBBIN)/libphaseTemperatureChangeTwoPhaseMixtures diff --git a/applications/solvers/multiphase/interCondensingEvaporatingFoam/temperaturePhaseChangeTwoPhaseMixtures/Make/options b/applications/solvers/multiphase/interCondensingEvaporatingFoam/temperaturePhaseChangeTwoPhaseMixtures/Make/options new file mode 100644 index 0000000000..4a85a0ae9f --- /dev/null +++ b/applications/solvers/multiphase/interCondensingEvaporatingFoam/temperaturePhaseChangeTwoPhaseMixtures/Make/options @@ -0,0 +1,16 @@ + +EXE_INC = \ + -I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \ + -I$(LIB_SRC)/transportModels/twoPhaseMixture/lnInclude \ + -I$(LIB_SRC)/transportModels \ + -I$(LIB_SRC)/transportModels/interfaceProperties/lnInclude \ + -I$(LIB_SRC)/transportModels/incompressible/lnInclude \ + -I$(LIB_SRC)/finiteVolume/lnInclude + +LIB_LIBS = \ + -ltwoPhaseMixture \ + -linterfaceProperties \ + -ltwoPhaseProperties \ + -lincompressibleTransportModels \ + -lfiniteVolume \ + -lfluidThermophysicalModels diff --git a/applications/solvers/multiphase/interCondensingEvaporatingFoam/temperaturePhaseChangeTwoPhaseMixtures/constant/constant.C b/applications/solvers/multiphase/interCondensingEvaporatingFoam/temperaturePhaseChangeTwoPhaseMixtures/constant/constant.C new file mode 100644 index 0000000000..adc2598e86 --- /dev/null +++ b/applications/solvers/multiphase/interCondensingEvaporatingFoam/temperaturePhaseChangeTwoPhaseMixtures/constant/constant.C @@ -0,0 +1,174 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2016 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 . + +\*---------------------------------------------------------------------------*/ + +#include "constant.H" +#include "addToRunTimeSelectionTable.H" +#include "fvcGrad.H" +#include "twoPhaseMixtureEThermo.H" + +// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // + +namespace Foam +{ +namespace temperaturePhaseChangeTwoPhaseMixtures +{ + defineTypeNameAndDebug(constant, 0); + addToRunTimeSelectionTable + ( + temperaturePhaseChangeTwoPhaseMixture, + constant, + components + ); +} +} + +// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // + +Foam::temperaturePhaseChangeTwoPhaseMixtures::constant::constant +( + const thermoIncompressibleTwoPhaseMixture& mixture, + const fvMesh& mesh +) +: + temperaturePhaseChangeTwoPhaseMixture(mixture, mesh), + coeffC_(subDict(type() + "Coeffs").lookup("coeffC")), + coeffE_(subDict(type() + "Coeffs").lookup("coeffE")) +{ +} + + +// * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * // + +Foam::Pair > +Foam::temperaturePhaseChangeTwoPhaseMixtures::constant::mDotAlphal() const +{ + const volScalarField& T = mesh_.lookupObject("T"); + + const twoPhaseMixtureEThermo& thermo = + refCast + ( + mesh_.lookupObject(basicThermo::dictName) + ); + + const dimensionedScalar& TSat = thermo.TSat(); + + dimensionedScalar T0("0", dimTemperature, 0.0); + + return Pair > + ( + coeffC_*mixture_.rho2()*max(TSat - T, T0), + -coeffE_*mixture_.rho1()*max(T - TSat, T0) + ); +} + + +Foam::Pair > +Foam::temperaturePhaseChangeTwoPhaseMixtures::constant::mDot() const +{ + + volScalarField limitedAlpha1 + ( + min(max(mixture_.alpha1(), scalar(0)), scalar(1)) + ); + + volScalarField limitedAlpha2 + ( + min(max(mixture_.alpha2(), scalar(0)), scalar(1)) + ); + + const volScalarField& T = mesh_.lookupObject("T"); + + const twoPhaseMixtureEThermo& thermo = + refCast + ( + mesh_.lookupObject(basicThermo::dictName) + ); + + const dimensionedScalar& TSat = thermo.TSat(); + + dimensionedScalar T0("0", dimTemperature, 0.0); + + return Pair > + ( + coeffC_*mixture_.rho2()*limitedAlpha2*max(TSat - T, T0), + coeffE_*mixture_.rho1()*limitedAlpha1*max(T - TSat, T0) + ); +} + + +Foam::Pair > +Foam::temperaturePhaseChangeTwoPhaseMixtures::constant::mDotDeltaT() const +{ + volScalarField limitedAlpha1 + ( + min(max(mixture_.alpha1(), scalar(0)), scalar(1)) + ); + + volScalarField limitedAlpha2 + ( + min(max(mixture_.alpha2(), scalar(0)), scalar(1)) + ); + + const volScalarField& T = mesh_.lookupObject("T"); + + const twoPhaseMixtureEThermo& thermo = + refCast + ( + mesh_.lookupObject(basicThermo::dictName) + ); + + const dimensionedScalar& TSat = thermo.TSat(); + + + return Pair > + ( + coeffC_*mixture_.rho2()*limitedAlpha2*pos(TSat - T), + coeffE_*mixture_.rho1()*limitedAlpha1*pos(T - TSat) + ); +} + + +void Foam::temperaturePhaseChangeTwoPhaseMixtures::constant::correct() +{ +} + + +bool Foam::temperaturePhaseChangeTwoPhaseMixtures::constant::read() +{ + if (temperaturePhaseChangeTwoPhaseMixture::read()) + { + subDict(type() + "Coeffs").lookup("coeffC") >> coeffC_; + subDict(type() + "Coeffs").lookup("coeffE") >> coeffE_; + + return true; + } + else + { + return false; + } +} + + +// ************************************************************************* // diff --git a/applications/solvers/multiphase/interCondensingEvaporatingFoam/temperaturePhaseChangeTwoPhaseMixtures/constant/constant.H b/applications/solvers/multiphase/interCondensingEvaporatingFoam/temperaturePhaseChangeTwoPhaseMixtures/constant/constant.H new file mode 100644 index 0000000000..424a8a7c2a --- /dev/null +++ b/applications/solvers/multiphase/interCondensingEvaporatingFoam/temperaturePhaseChangeTwoPhaseMixtures/constant/constant.H @@ -0,0 +1,118 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2016 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 . + +Class + Foam::phaseChangeTwoPhaseMixtures::constant + +Description + constant condensation/saturation model. + + +SourceFiles + constant.C + +\*--------------------------------------------------------------------*/ + +#ifndef constant_H +#define constant_H + +#include "temperaturePhaseChangeTwoPhaseMixture.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ +namespace temperaturePhaseChangeTwoPhaseMixtures +{ + +/*--------------------------------------------------------------------*\ + Class constant +\*--------------------------------------------------------------------*/ + +class constant +: + public temperaturePhaseChangeTwoPhaseMixture +{ + // Private data + + //- Condensation coefficient [1/s/K] + dimensionedScalar coeffC_; + + //- Evaporation coefficient [1/s/K] + dimensionedScalar coeffE_; + + +public: + + //- Runtime type information + TypeName("constant"); + + + // Constructors + + //- Construct from components + constant + ( + const thermoIncompressibleTwoPhaseMixture& mixture, + const fvMesh& mesh + ); + + + //- Destructor + virtual ~constant() + {} + + + // Member Functions + + //- Return the mass condensation and vaporisation rates as a + // coefficient to multiply (1 - alphal) for the condensation rate + // and a coefficient to multiply alphal for the vaporisation rate + virtual Pair > mDotAlphal() const; + + //- Return the mass condensation and vaporisation rates as coefficients + virtual Pair > mDot() const; + + //- Return the mass condensation and vaporisation rates as a + // coefficient to multiply (Tsat - T) for the condensation rate + // and a coefficient to multiply (T - Tsat) for the vaporisation rate + virtual Pair > mDotDeltaT() const; + + //- Correct the constant phaseChange model + virtual void correct(); + + //- Read the transportProperties dictionary and update + virtual bool read(); +}; + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace temperaturePhaseChangeTwoPhaseMixtures +} // End namespace Foam + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* // diff --git a/applications/solvers/multiphase/interCondensingEvaporatingFoam/temperaturePhaseChangeTwoPhaseMixtures/temperaturePhaseChangeTwoPhaseMixtures/newtemperaturePhaseChangeTwoPhaseMixture.C b/applications/solvers/multiphase/interCondensingEvaporatingFoam/temperaturePhaseChangeTwoPhaseMixtures/temperaturePhaseChangeTwoPhaseMixtures/newtemperaturePhaseChangeTwoPhaseMixture.C new file mode 100644 index 0000000000..8dbb7c038c --- /dev/null +++ b/applications/solvers/multiphase/interCondensingEvaporatingFoam/temperaturePhaseChangeTwoPhaseMixtures/temperaturePhaseChangeTwoPhaseMixtures/newtemperaturePhaseChangeTwoPhaseMixture.C @@ -0,0 +1,81 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2016 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 . + +\*---------------------------------------------------------------------------*/ + +#include "temperaturePhaseChangeTwoPhaseMixture.H" +#include "basicThermo.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +Foam::autoPtr +Foam::temperaturePhaseChangeTwoPhaseMixture::New +( + const thermoIncompressibleTwoPhaseMixture& thermo, + const fvMesh& mesh +) +{ + IOdictionary phaseChangePropertiesDict + ( + IOobject + ( + "phaseChangeProperties", + mesh.time().constant(), + mesh, + IOobject::MUST_READ, + IOobject::NO_WRITE, + false + ) + ); + + word temperaturePhaseChangeTwoPhaseMixtureTypeName + ( + phaseChangePropertiesDict.lookup + ( + "phaseChangeTwoPhaseModel" + ) + ); + + Info<< "Selecting phaseChange model " + << temperaturePhaseChangeTwoPhaseMixtureTypeName << endl; + + componentsConstructorTable::iterator cstrIter = + componentsConstructorTablePtr_ + ->find(temperaturePhaseChangeTwoPhaseMixtureTypeName); + + if (cstrIter == componentsConstructorTablePtr_->end()) + { + FatalErrorInFunction + << "Unknown temperaturePhaseChangeTwoPhaseMixture type " + << temperaturePhaseChangeTwoPhaseMixtureTypeName << endl << endl + << "Valid temperaturePhaseChangeTwoPhaseMixtures are : " << endl + << componentsConstructorTablePtr_->sortedToc() + << exit(FatalError); + } + + return autoPtr + (cstrIter()(thermo, mesh)); +} + + +// ************************************************************************* // diff --git a/applications/solvers/multiphase/interCondensingEvaporatingFoam/temperaturePhaseChangeTwoPhaseMixtures/temperaturePhaseChangeTwoPhaseMixtures/temperaturePhaseChangeTwoPhaseMixture.C b/applications/solvers/multiphase/interCondensingEvaporatingFoam/temperaturePhaseChangeTwoPhaseMixtures/temperaturePhaseChangeTwoPhaseMixtures/temperaturePhaseChangeTwoPhaseMixture.C new file mode 100644 index 0000000000..9727f5978d --- /dev/null +++ b/applications/solvers/multiphase/interCondensingEvaporatingFoam/temperaturePhaseChangeTwoPhaseMixtures/temperaturePhaseChangeTwoPhaseMixtures/temperaturePhaseChangeTwoPhaseMixture.C @@ -0,0 +1,109 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2016 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 . + +\*---------------------------------------------------------------------------*/ + +#include "temperaturePhaseChangeTwoPhaseMixture.H" + +// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // + +namespace Foam +{ + defineTypeNameAndDebug(temperaturePhaseChangeTwoPhaseMixture, 0); + defineRunTimeSelectionTable + ( + temperaturePhaseChangeTwoPhaseMixture, + components + ); +} + +// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // + +Foam::temperaturePhaseChangeTwoPhaseMixture:: +temperaturePhaseChangeTwoPhaseMixture +( + const thermoIncompressibleTwoPhaseMixture& mixture, + const fvMesh& mesh +) +: + IOdictionary + ( + IOobject + ( + "phaseChangeProperties", + mesh.time().constant(), + mesh, + IOobject::MUST_READ_IF_MODIFIED, + IOobject::NO_WRITE + ) + ), + mixture_(mixture), + mesh_(mesh) +{} + + +// * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * // + +Foam::Pair > +Foam::temperaturePhaseChangeTwoPhaseMixture::vDotAlphal() const +{ + volScalarField alphalCoeff + ( + 1.0/mixture_.rho1() - mixture_.alpha1() + *(1.0/mixture_.rho1() - 1.0/mixture_.rho2()) + ); + + Pair > mDotAlphal = this->mDotAlphal(); + + return Pair > + ( + alphalCoeff*mDotAlphal[0], + alphalCoeff*mDotAlphal[1] + ); +} + + +Foam::Pair > +Foam::temperaturePhaseChangeTwoPhaseMixture::vDot() const +{ + dimensionedScalar pCoeff(1.0/mixture_.rho1() - 1.0/mixture_.rho2()); + Pair > mDot = this->mDot(); + + return Pair >(pCoeff*mDot[0], pCoeff*mDot[1]); +} + + +bool Foam::temperaturePhaseChangeTwoPhaseMixture::read() +{ + if (regIOobject::read()) + { + return true; + } + else + { + return false; + } +} + + +// ************************************************************************* // diff --git a/applications/solvers/multiphase/interCondensingEvaporatingFoam/temperaturePhaseChangeTwoPhaseMixtures/temperaturePhaseChangeTwoPhaseMixtures/temperaturePhaseChangeTwoPhaseMixture.H b/applications/solvers/multiphase/interCondensingEvaporatingFoam/temperaturePhaseChangeTwoPhaseMixtures/temperaturePhaseChangeTwoPhaseMixtures/temperaturePhaseChangeTwoPhaseMixture.H new file mode 100644 index 0000000000..9ecd5a544d --- /dev/null +++ b/applications/solvers/multiphase/interCondensingEvaporatingFoam/temperaturePhaseChangeTwoPhaseMixtures/temperaturePhaseChangeTwoPhaseMixtures/temperaturePhaseChangeTwoPhaseMixture.H @@ -0,0 +1,168 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2016 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 . + +Class + Foam::temperaturePhaseChangeTwoPhaseMixture + +Description + +SourceFiles + temperaturePhaseChangeTwoPhaseMixture.C + newtemperaturePhaseChangeTwoPhaseMixture.C + +\*---------------------------------------------------------------------------*/ + +#ifndef temperaturePhaseChangeTwoPhaseMixture_H +#define temperaturePhaseChangeTwoPhaseMixture_H + +#include "thermoIncompressibleTwoPhaseMixture.H" +#include "typeInfo.H" +#include "runTimeSelectionTables.H" +#include "volFields.H" +#include "dimensionedScalar.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + +/*---------------------------------------------------------------------------*\ + Class temperaturePhaseChangeTwoPhaseMixture Declaration +\*---------------------------------------------------------------------------*/ + +class temperaturePhaseChangeTwoPhaseMixture +: + public IOdictionary +{ + +protected: + + // Protected data + + //- Reference to the thermoIncompressibleTwoPhaseMixture + const thermoIncompressibleTwoPhaseMixture& mixture_; + + //- Reference to fvMesh + const fvMesh& mesh_; + + + // Private Member Functions + + //- Disallow copy construct + temperaturePhaseChangeTwoPhaseMixture + ( + const temperaturePhaseChangeTwoPhaseMixture& + ); + + //- Disallow default bitwise assignment + void operator=(const temperaturePhaseChangeTwoPhaseMixture&); + + +public: + + //- Runtime type information + TypeName("temperaturePhaseChangeTwoPhaseMixture"); + + + // Declare run-time constructor selection table + + declareRunTimeSelectionTable + ( + autoPtr, + temperaturePhaseChangeTwoPhaseMixture, + components, + ( + const thermoIncompressibleTwoPhaseMixture& mixture, + const fvMesh& mesh + ), + (mixture, mesh) + ); + + + // Selectors + + //- Return a reference to the selected phaseChange model + static autoPtr New + ( + const thermoIncompressibleTwoPhaseMixture& mixture, + const fvMesh& mesh + ); + + + // Constructors + + //- Construct from components + temperaturePhaseChangeTwoPhaseMixture + ( + const thermoIncompressibleTwoPhaseMixture& mixture, + const fvMesh& mesh + ); + + + //- Destructor + virtual ~temperaturePhaseChangeTwoPhaseMixture() + {} + + + // Member Functions + + + //- Return the mass condensation and vaporisation rates as a + // coefficient to multiply (1 - alphal) for the condensation rate + // and a coefficient to multiply alphal for the vaporisation rate + virtual Pair > mDotAlphal() const = 0; + + //- Return the mass condensation and vaporisation rates as coefficients + virtual Pair > mDot() const = 0; + + //- Return the mass condensation and vaporisation rates as a + // coefficient to multiply (Tsat - T) for the condensation rate + // and a coefficient to multiply (T - Tsat) for the vaporisation rate + virtual Pair > mDotDeltaT() const = 0; + + //- Return the volumetric condensation and vaporisation rates as a + // coefficient to multiply (1 - alphal) for the condensation rate + // and a coefficient to multiply alphal for the vaporisation rate + Pair > vDotAlphal() const; + + //- Return the volumetric condensation and vaporisation rates as + // coefficients + Pair > vDot() const; + + //- Correct the phaseChange model + virtual void correct() = 0; + + //- Read the transportProperties dictionary and update + virtual bool read(); +}; + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace Foam + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* // diff --git a/applications/solvers/multiphase/interCondensingEvaporatingFoam/temperaturePhaseChangeTwoPhaseMixtures/thermoIncompressibleTwoPhaseMixture/thermoIncompressibleTwoPhaseMixture.C b/applications/solvers/multiphase/interCondensingEvaporatingFoam/temperaturePhaseChangeTwoPhaseMixtures/thermoIncompressibleTwoPhaseMixture/thermoIncompressibleTwoPhaseMixture.C new file mode 100644 index 0000000000..6822e0bbc5 --- /dev/null +++ b/applications/solvers/multiphase/interCondensingEvaporatingFoam/temperaturePhaseChangeTwoPhaseMixtures/thermoIncompressibleTwoPhaseMixture/thermoIncompressibleTwoPhaseMixture.C @@ -0,0 +1,133 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2016 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 . + +\*---------------------------------------------------------------------------*/ + +#include "thermoIncompressibleTwoPhaseMixture.H" + + +// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // + +namespace Foam +{ + defineTypeNameAndDebug(thermoIncompressibleTwoPhaseMixture, 0); +} + + +// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // + +Foam::thermoIncompressibleTwoPhaseMixture::thermoIncompressibleTwoPhaseMixture +( + const volVectorField& U, + const surfaceScalarField& phi +) +: + incompressibleTwoPhaseMixture(U, phi), + + kappa1_ + ( + "kappa1", + dimEnergy/dimTime/dimLength/dimTemperature, + subDict(phase1Name_).lookup("kappa") + ), + kappa2_ + ( + "kappa2", + kappa1_.dimensions(), + subDict(phase2Name_).lookup("kappa") + ), + + Cp1_ + ( + "Cp1", + dimEnergy/dimTemperature/dimMass, + subDict(phase1Name_).lookup("Cp") + ), + Cp2_ + ( + "Cp2", + dimEnergy/dimTemperature/dimMass, + subDict(phase2Name_).lookup("Cp") + ), + + Cv1_ + ( + "Cv1", + dimEnergy/dimTemperature/dimMass, + subDict(phase1Name_).lookup("Cv") + ), + + Cv2_ + ( + "Cv2", + dimEnergy/dimTemperature/dimMass, + subDict(phase2Name_).lookup("Cv") + ), + + Hf1_ + ( + "Hf1", + dimEnergy/dimMass, + subDict(phase1Name_).lookup("hf") + ), + + Hf2_ + ( + "Hf2", + dimEnergy/dimMass, + subDict(phase2Name_).lookup("hf") + ) +{ + +} + + +// * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * // + + +bool Foam::thermoIncompressibleTwoPhaseMixture::read() +{ + if (incompressibleTwoPhaseMixture::read()) + { + subDict(phase1Name_).lookup("kappa") >> kappa1_; + subDict(phase2Name_).lookup("kappa") >> kappa2_; + + subDict(phase1Name_).lookup("Cp") >> Cp1_; + subDict(phase2Name_).lookup("Cp") >> Cp2_; + + subDict(phase1Name_).lookup("Cv") >> Cv1_; + subDict(phase2Name_).lookup("Cv") >> Cv2_; + + subDict(phase1Name_).lookup("hf") >> Hf1_; + subDict(phase2Name_).lookup("hf") >> Hf2_; + + return true; + } + else + { + return false; + } +} + + +// ************************************************************************* // diff --git a/applications/solvers/multiphase/interCondensingEvaporatingFoam/temperaturePhaseChangeTwoPhaseMixtures/thermoIncompressibleTwoPhaseMixture/thermoIncompressibleTwoPhaseMixture.H b/applications/solvers/multiphase/interCondensingEvaporatingFoam/temperaturePhaseChangeTwoPhaseMixtures/thermoIncompressibleTwoPhaseMixture/thermoIncompressibleTwoPhaseMixture.H new file mode 100644 index 0000000000..ea2bec1b5d --- /dev/null +++ b/applications/solvers/multiphase/interCondensingEvaporatingFoam/temperaturePhaseChangeTwoPhaseMixtures/thermoIncompressibleTwoPhaseMixture/thermoIncompressibleTwoPhaseMixture.H @@ -0,0 +1,156 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2016 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 . + +Class + Foam::thermoIncompressibleTwoPhaseMixture + +Description + A two-phase incompressible transportModel + +SourceFiles + thermoIncompressibleTwoPhaseMixture.C + +\*---------------------------------------------------------------------------*/ + +#ifndef thermoIncompressibleTwoPhaseMixture_H +#define thermoIncompressibleTwoPhaseMixture_H + +#include "incompressibleTwoPhaseMixture.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + +/*---------------------------------------------------------------------------*\ + Class thermoIncompressibleTwoPhaseMixture Declaration +\*---------------------------------------------------------------------------*/ + +class thermoIncompressibleTwoPhaseMixture +: + public incompressibleTwoPhaseMixture +{ +protected: + + // Protected data + + + //- Thermal variables + dimensionedScalar kappa1_; + dimensionedScalar kappa2_; + + dimensionedScalar Cp1_; + dimensionedScalar Cp2_; + + dimensionedScalar Cv1_; + dimensionedScalar Cv2_; + + dimensionedScalar Hf1_; + dimensionedScalar Hf2_; + + +public: + + TypeName("thermoIncompressibleTwoPhaseMixture"); + + + // Constructors + + //- Construct from U and phi + thermoIncompressibleTwoPhaseMixture + ( + const volVectorField& U, + const surfaceScalarField& phi + ); + + + //- Destructor + virtual ~thermoIncompressibleTwoPhaseMixture() + {} + + + // Access function + + //- Return const-access to phase2 kappa + const dimensionedScalar& kappa2() const + { + return kappa2_; + }; + + //- Return const-access to phase1 kappa + const dimensionedScalar& kappa1() const + { + return kappa1_; + }; + + //- Return const-access to phase1 Cp + const dimensionedScalar& Cp1() const + { + return Cp1_; + }; + + //- Return const-access to phase1 Cp + const dimensionedScalar& Cp2() const + { + return Cp2_; + }; + + //- Return const-access to phase1 Cv + const dimensionedScalar& Cv1() const + { + return Cv1_; + }; + + //- Return const-access to phase1 Cv + const dimensionedScalar& Cv2() const + { + return Cv2_; + }; + + //- Return latent heat for phase 1 + const dimensionedScalar& Hf1() const + { + return Hf1_; + }; + + //- Return latent heat for phase 2 + const dimensionedScalar& Hf2() const + { + return Hf2_; + }; + + //- Read base transportProperties dictionary + virtual bool read(); + +}; + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace Foam + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* // diff --git a/applications/solvers/multiphase/interCondensingEvaporatingFoam/temperaturePhaseChangeTwoPhaseMixtures/twoPhaseMixtureEThermo/twoPhaseMixtureEThermo.C b/applications/solvers/multiphase/interCondensingEvaporatingFoam/temperaturePhaseChangeTwoPhaseMixtures/twoPhaseMixtureEThermo/twoPhaseMixtureEThermo.C new file mode 100644 index 0000000000..19eb8ddca7 --- /dev/null +++ b/applications/solvers/multiphase/interCondensingEvaporatingFoam/temperaturePhaseChangeTwoPhaseMixtures/twoPhaseMixtureEThermo/twoPhaseMixtureEThermo.C @@ -0,0 +1,592 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2016 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 . + +\*---------------------------------------------------------------------------*/ + +#include "twoPhaseMixtureEThermo.H" + +#include "zeroGradientFvPatchFields.H" +#include "fixedEnergyFvPatchScalarField.H" +#include "gradientEnergyFvPatchScalarField.H" +#include "mixedEnergyFvPatchScalarField.H" +#include "twoPhaseMixtureEThermo.H" + +// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // + +namespace Foam +{ + defineTypeNameAndDebug(twoPhaseMixtureEThermo, 0); +} + +// * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * // + +void Foam::twoPhaseMixtureEThermo::eBoundaryCorrection(volScalarField& h) +{ + volScalarField::GeometricBoundaryField& hbf = h.boundaryField(); + + forAll(hbf, patchi) + { + if (isA(hbf[patchi])) + { + refCast(hbf[patchi]).gradient() + = hbf[patchi].fvPatchField::snGrad(); + } + else if (isA(hbf[patchi])) + { + refCast(hbf[patchi]).refGrad() + = hbf[patchi].fvPatchField::snGrad(); + } + } +} + +void Foam::twoPhaseMixtureEThermo::init() +{ + const volScalarField alpha1Rho1(alpha1()*rho1()); + const volScalarField alpha2Rho2(alpha2()*rho2()); + + e_ = + ( + (T_ - TSat_)*(alpha1Rho1*Cv1() + alpha2Rho2*Cv2()) + + (alpha1Rho1*Hf1() + alpha2Rho2*Hf2()) + ) + /(alpha1Rho1 + alpha2Rho2); + + e_.correctBoundaryConditions(); +} + +// * * * * * * * * * * * * * * * * Constructor * * * * * * * * * * * * * * * // + +Foam::twoPhaseMixtureEThermo::twoPhaseMixtureEThermo +( + const volVectorField& U, + const surfaceScalarField& phi +) +: + basicThermo(U.mesh(), word::null), + thermoIncompressibleTwoPhaseMixture(U, phi), + + e_ + ( + volScalarField + ( + IOobject + ( + "e", + U.mesh().time().timeName(), + U.mesh(), + IOobject::NO_READ, + IOobject::NO_WRITE + ), + U.mesh(), + dimensionedScalar("zero", dimEnergy/dimMass, 0.0), + heBoundaryTypes() + ) + ), + + TSat_ + ( + "TSat", + dimTemperature, + basicThermo::lookup("TSat") + ), + + pDivU_(basicThermo::lookupOrDefault("pDivU", true)) + +{ + // Initialise e + init(); +} + +// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * // + +Foam::twoPhaseMixtureEThermo::~twoPhaseMixtureEThermo() +{} + +// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // + +void Foam::twoPhaseMixtureEThermo::correct() +{ + incompressibleTwoPhaseMixture::correct(); + + const volScalarField alpha1Rho1(alpha1()*rho1()); + const volScalarField alpha2Rho2(alpha2()*rho2()); + + T_ = + ( + (e_*(alpha1Rho1 + alpha2Rho2)) + - (alpha1Rho1*Hf1() + alpha2Rho2*Hf2()) + ) + /(alpha1Rho1*Cv1() + alpha2Rho2*Cv2()) + + TSat_; + + T().correctBoundaryConditions(); +} + + +Foam::tmp Foam::twoPhaseMixtureEThermo::he +( + const volScalarField& p, + const volScalarField& T +) const +{ + const volScalarField alpha1Rho1(alpha1()*rho1()); + const volScalarField alpha2Rho2(alpha2()*rho2()); + + return + ( + (T - TSat_)*(alpha1Rho1*Cv1() + alpha2Rho2*Cv2()) + + (alpha1Rho1*Hf1() + alpha2Rho2*Hf2()) + ) + / (alpha1Rho1 + alpha2Rho2); +} + + +Foam::tmp Foam::twoPhaseMixtureEThermo::he +( + const scalarField& p, + const scalarField& T, + const labelList& cells +) const +{ + tmp the(new scalarField(T.size())); + scalarField& he = the.ref(); + + const volScalarField alpha1Rho1(alpha1()*rho1()); + const volScalarField alpha2Rho2(alpha2()*rho2()); + + forAll(T, i) + { + label celli = cells[i]; + he[i] = + ( + (T[i] - TSat_.value()) + *( + alpha1Rho1[celli]*Cv1().value() + + alpha2Rho2[celli]*Cv2().value() + ) + + ( + alpha1Rho1[celli]*Hf1().value() + + alpha2Rho2[celli]*Hf2().value() + ) + ) + / (alpha1Rho1[celli] + alpha2Rho2[celli]); + } + + return the; +} + + +Foam::tmp Foam::twoPhaseMixtureEThermo::he +( + const scalarField& p, + const scalarField& T, + const label patchI +) const +{ + const scalarField& alpha1p = alpha1().boundaryField()[patchI]; + const scalarField& alpha2p = alpha2().boundaryField()[patchI]; + + const scalarField& Tp = T_.boundaryField()[patchI]; + + return + ( + ( + (Tp - TSat_.value()) + *( + alpha1p*rho1().value()*Cv1().value() + + alpha2p*rho2().value()*Cv2().value() + ) + + ( + alpha1p*rho1().value()*Hf1().value() + + alpha2p*rho2().value()*Hf2().value() + ) + ) + / (alpha1p*rho1().value() + alpha2p*rho2().value()) + ); +} + + +Foam::tmp +Foam::twoPhaseMixtureEThermo::hc() const +{ + + const fvMesh& mesh = this->T_.mesh(); + + return tmp + ( + new volScalarField + ( + IOobject + ( + "hc", + mesh.time().timeName(), + mesh, + IOobject::NO_READ, + IOobject::AUTO_WRITE + ), + mesh, + dimensionedScalar("hc",Hf2() - Hf1()) + ) + ); +} + + +Foam::tmp Foam::twoPhaseMixtureEThermo::THE +( + const scalarField& h, + const scalarField& p, + const scalarField& T0, // starting temperature + const labelList& cells +) const +{ + NotImplemented; + return tmp(); +} + + +Foam::tmp Foam::twoPhaseMixtureEThermo::THE +( + const scalarField& h, + const scalarField& p, + const scalarField& T0, // starting temperature + const label patchi +) const +{ + NotImplemented; + return tmp(); +} + + +Foam::tmp +Foam::twoPhaseMixtureEThermo::Cp() const +{ + const volScalarField limitedAlpha1 + ( + min(max(alpha1_, scalar(0)), scalar(1)) + ); + + return tmp + ( + new volScalarField + ( + "cp", + limitedAlpha1*Cp1() + (scalar(1) - limitedAlpha1)*Cp2() + ) + ); +} + + +Foam::tmp Foam::twoPhaseMixtureEThermo::Cp +( + const scalarField& p, + const scalarField& T, + const label patchi +) const +{ + const volScalarField limitedAlpha1 + ( + min(max(alpha1_, scalar(0)), scalar(1)) + ); + + const scalarField& alpha1p = limitedAlpha1.boundaryField()[patchi]; + + return + ( + alpha1p*Cp1().value() + (scalar(1) - alpha1p)*Cp2().value() + ); +} + + +Foam::tmp +Foam::twoPhaseMixtureEThermo::rho() const +{ + const volScalarField limitedAlpha1 + ( + min(max(alpha1_, scalar(0)), scalar(1)) + ); + + return tmp + ( + new volScalarField + ( + "rho", + limitedAlpha1*rho1().value() + + (scalar(1) - limitedAlpha1)*rho2().value() + ) + ); +} + + +Foam::tmp +Foam::twoPhaseMixtureEThermo::rho +( + const label patchi +) const +{ + const volScalarField limitedAlpha1 + ( + min(max(alpha1_, scalar(0)), scalar(1)) + ); + + const scalarField& alpha1p = limitedAlpha1.boundaryField()[patchi]; + + return + ( + alpha1p*rho1().value() + (scalar(1.0) - alpha1p)*rho2().value() + ); +} + + +Foam::tmp +Foam::twoPhaseMixtureEThermo::Cv() const +{ + const volScalarField limitedAlpha1 + ( + min(max(alpha1_, scalar(0)), scalar(1)) + ); + + return tmp + ( + new volScalarField + ( + "cv", + limitedAlpha1*Cv1() + (scalar(1) - limitedAlpha1)*Cv2() + ) + ); +} + + +Foam::tmp +Foam::twoPhaseMixtureEThermo::Cv +( + const scalarField& p, + const scalarField& T, + const label patchi +) const +{ + const volScalarField limitedAlpha1 + ( + min(max(alpha1_, scalar(0)), scalar(1)) + ); + + const scalarField& alpha1p = limitedAlpha1.boundaryField()[patchi]; + + return + ( + alpha1p*Cv1().value() + (scalar(1) - alpha1p)*Cv2().value() + ); +} + + +Foam::tmp +Foam::twoPhaseMixtureEThermo::gamma() const +{ + return tmp + ( + (alpha1_*Cp1() + alpha2_*Cp2()) / (alpha1_*Cv1() + alpha2_*Cv2()) + ); +} + + +Foam::tmp Foam::twoPhaseMixtureEThermo::gamma +( + const scalarField& p, + const scalarField& T, + const label patchi +) const +{ + return + ( + gamma()().boundaryField()[patchi] + ); +} + + +Foam::tmp +Foam::twoPhaseMixtureEThermo::Cpv() const +{ + // This is an e thermo (Cpv = Cv) + return Cv(); +} + + +Foam::tmp Foam::twoPhaseMixtureEThermo::Cpv +( + const scalarField& p, + const scalarField& T, + const label patchI +) const +{ + // This is a e thermo (Cpv = Cv) + return Cv(p, T, patchI); +} + + +Foam::tmp +Foam::twoPhaseMixtureEThermo::CpByCpv() const +{ + NotImplemented; + return tmp(); +} + + +Foam::tmp +Foam::twoPhaseMixtureEThermo::CpByCpv +( + const scalarField& p, + const scalarField& T, + const label patchi +) const +{ + NotImplemented; + return tmp(); +} + + +Foam::tmp +Foam::twoPhaseMixtureEThermo::kappa() const +{ + const volScalarField limitedAlpha1 + ( + min(max(alpha1_, scalar(0)), scalar(1)) + ); + + return tmp + ( + new volScalarField + ( + "kappa", + limitedAlpha1*kappa1() + (scalar(1) - limitedAlpha1)*kappa2() + ) + ); +} + + +Foam::tmp +Foam::twoPhaseMixtureEThermo::kappa(const label patchi) const +{ + const volScalarField limitedAlpha1 + ( + min(max(alpha1_, scalar(0)), scalar(1)) + ); + + const scalarField& alpha1p = limitedAlpha1.boundaryField()[patchi]; + + return (alpha1p*kappa1().value() + (1 - alpha1p)*kappa2().value()); +} + + +Foam::tmp +Foam::twoPhaseMixtureEThermo::kappaEff +( + const volScalarField& kappat +) const +{ + tmp kappaEff(kappa() + kappat); + kappaEff.ref().rename("kappaEff"); + return kappaEff; +} + + +Foam::tmp +Foam::twoPhaseMixtureEThermo::kappaEff +( + const scalarField& kappat, + const label patchi +) const +{ + const volScalarField limitedAlpha1 + ( + min(max(alpha1_, scalar(0)), scalar(1)) + ); + + const scalarField& alpha1p = limitedAlpha1.boundaryField()[patchi]; + + return + (alpha1p*kappa1().value() + (1 - alpha1p)*kappa2().value()) + kappat; + +} + + +Foam::tmp +Foam::twoPhaseMixtureEThermo::alphaEff +( + const volScalarField& alphat +) const +{ + const volScalarField rho + ( + alpha1_*rho1() + (1.0 - alpha1_)*rho2() + ); + + + return (kappa()/Cp()/rho + alphat); +} + +Foam::tmp +Foam::twoPhaseMixtureEThermo::alphaEff +( + const scalarField& alphat, + const label patchi +) const +{ + const volScalarField limitedAlpha1 + ( + min(max(alpha1_, scalar(0)), scalar(1)) + ); + + const scalarField& alpha1p = limitedAlpha1.boundaryField()[patchi]; + + const scalarField rho + ( + alpha1p*rho1().value() + (1.0 - alpha1p)*rho2().value() + ); + + const scalarField kappa + ( + alpha1p*kappa1().value() + (1.0 - alpha1p)*kappa2().value() + ); + + const scalarField Cp + ( + alpha1p*Cp1().value() + (1.0 - alpha1p)*Cp2().value() + ); + + return kappa/Cp/rho + alphat; +} + +bool Foam::twoPhaseMixtureEThermo::read() +{ + if (basicThermo::read() && thermoIncompressibleTwoPhaseMixture::read()) + { + basicThermo::lookup("pDivU") >> pDivU_; + basicThermo::lookup("TSat") >> TSat_; + return true; + } + else + { + return false; + } +} + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // diff --git a/applications/solvers/multiphase/interCondensingEvaporatingFoam/temperaturePhaseChangeTwoPhaseMixtures/twoPhaseMixtureEThermo/twoPhaseMixtureEThermo.H b/applications/solvers/multiphase/interCondensingEvaporatingFoam/temperaturePhaseChangeTwoPhaseMixtures/twoPhaseMixtureEThermo/twoPhaseMixtureEThermo.H new file mode 100644 index 0000000000..8616713b19 --- /dev/null +++ b/applications/solvers/multiphase/interCondensingEvaporatingFoam/temperaturePhaseChangeTwoPhaseMixtures/twoPhaseMixtureEThermo/twoPhaseMixtureEThermo.H @@ -0,0 +1,313 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2016 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 . + +Class + Foam::twoPhaseMixtureEThermo + +Description + Two phases thermo Internal energy mixture Defined as: + e1 = Cv1(T - Tsat) + Hv1 + e2 = Cv2(T - Tsat) + Hv2 + e = (alpha1*rho1*e1 + alpha2*rho2*e2) / (alpha1*rho1 + alpha2*rho2) + +SourceFiles + twoPhaseMixtureEThermo.C + +\*---------------------------------------------------------------------------*/ + +#ifndef flashThermo_H +#define flashThermo_H + +#include "volFields.H" + +#include "basicThermo.H" +#include "thermoIncompressibleTwoPhaseMixture.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + +/*---------------------------------------------------------------------------*\ + Class twoPhaseMixtureEThermo Declaration +\*---------------------------------------------------------------------------*/ + +class twoPhaseMixtureEThermo +: + public basicThermo, + public thermoIncompressibleTwoPhaseMixture +{ + +protected: + + // Protected Data + + //- Internal energy field [J] + volScalarField e_; + + //- Saturation Temperature + dimensionedScalar TSat_; + + //- Should the p*div(U) term be included in the energy equation + Switch pDivU_; + + + // Protected Member Functions + + //- Correct boundary for e + void eBoundaryCorrection (volScalarField& h); + + //- Initialize + void init(); + + +public: + + TypeName("twoPhaseMixtureEThermo"); + + // Constructor + twoPhaseMixtureEThermo + ( + const volVectorField& U, + const surfaceScalarField& phi + ); + + + //- Destructor + virtual ~twoPhaseMixtureEThermo(); + + + // Member Functions + + + //- Return access to the inernal energy field [J/Kg] + virtual volScalarField& he() + { + return e_; + } + + //- Return access to the inernal energy field [J/Kg] + virtual const volScalarField& he() const + { + return e_; + } + + //- Enthalpy/Internal energy + // for given pressure and temperature [J/kg] + virtual tmp he + ( + const volScalarField& p, + const volScalarField& T + ) const; + + //- Enthalpy/Internal energy for cell-set [J/kg] + virtual tmp he + ( + const scalarField& p, + const scalarField& T, + const labelList& cells + ) const; + + //- Enthalpy/Internal energy for patch [J/kg] + virtual tmp he + ( + const scalarField& p, + const scalarField& T, + const label patchi + ) const; + + //- Chemical enthalpy [J/kg] + virtual tmp hc() const; + + //- Temperature from enthalpy/internal energy for cell-set + virtual tmp THE + ( + const scalarField& h, + const scalarField& p, + const scalarField& T0, // starting temperature + const labelList& cells + ) const; + + //- Temperature from enthalpy/internal energy for patch + virtual tmp THE + ( + const scalarField& h, + const scalarField& p, + const scalarField& T0, // starting temperature + const label patchi + ) const; + + + //- Return true if the equation of state is incompressible + // i.e. rho != f(p) + bool incompressible() const + { + return (true); + } + + //- Return true if the equation of state is isochoric + // i.e. rho = const + bool isochoric() const + { + return (false); + } + + //- Return rho of the mixture + virtual tmp rho() const; + + //- Return rho for patch + virtual tmp rho(const label patchi) const; + + //- Return Cp of the mixture + virtual tmp Cp() const; + + //- Heat capacity at constant pressure for patch [J/kg/K] + virtual tmp Cp + ( + const scalarField& p, + const scalarField& T, + const label patchi + ) const; + + //- Return Cv of the mixture + virtual tmp Cv() const; + + //- Heat capacity at constant volume for patch [J/kg/K] + virtual tmp Cv + ( + const scalarField& p, + const scalarField& T, + const label patchI + ) const; + + //- Gamma = Cp/Cv [] + virtual tmp gamma() const; + + //- Gamma = Cp/Cv for patch [] + virtual tmp gamma + ( + const scalarField& p, + const scalarField& T, + const label patchi + ) const; + + //- Heat capacity at constant pressure/volume [J/kg/K] + virtual tmp Cpv() const; + + //- Heat capacity at constant pressure/volume for patch [J/kg/K] + virtual tmp Cpv + ( + const scalarField& p, + const scalarField& T, + const label patchi + ) const; + + //- Heat capacity ratio [] + virtual tmp CpByCpv() const; + + //- Heat capacity ratio for patch [] + virtual tmp CpByCpv + ( + const scalarField& p, + const scalarField& T, + const label patchi + ) const; + + //- Thermal diffusivity for temperature of mixture [J/m/s/K] + virtual tmp kappa() const; + + //- Thermal diffusivity for temperature + // of mixture for patch [J/m/s/K] + virtual tmp kappa + ( + const label patchi + ) const; + + //- Effective thermal diffusivity for temperature + // of mixture [J/m/s/K] + virtual tmp kappaEff + ( + const volScalarField& + ) const; + + //- Effective thermal diffusivity for temperature + // of mixture for patch [J/m/s/K] + virtual tmp kappaEff + ( + const scalarField& alphat, + const label patchi + ) const; + + //- Effective thermal diffusivity of mixture [kg/m/s] + virtual tmp alphaEff + ( + const volScalarField& alphat + ) const; + + //- Effective thermal diffusivity of mixture for patch [kg/m/s] + virtual tmp alphaEff + ( + const scalarField& alphat, + const label patchi + ) const; + + + //- Correct the thermo fields + virtual void correct(); + + //- Read properties + virtual bool read(); + + + // Access to thermodynamic state variables + + //- Return const-access to the saturation temperature + const dimensionedScalar& TSat() const + { + return TSat_; + } + + //- Return transport properties dictionary + const incompressibleTwoPhaseMixture& transportPropertiesDict() + { + return *this; + } + + //- Should the dpdt term be included in the enthalpy equation + Switch pDivU() const + { + return pDivU_; + } + +}; + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace Foam + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* // diff --git a/applications/test/ensightFile/Make/files b/applications/test/ensightFile/Make/files new file mode 100644 index 0000000000..856556c316 --- /dev/null +++ b/applications/test/ensightFile/Make/files @@ -0,0 +1,3 @@ +Test-ensightFile.C + +EXE = $(FOAM_USER_APPBIN)/Test-ensightFile diff --git a/applications/test/ensightFile/Make/options b/applications/test/ensightFile/Make/options new file mode 100644 index 0000000000..a4f7656262 --- /dev/null +++ b/applications/test/ensightFile/Make/options @@ -0,0 +1,10 @@ +EXE_INC = \ + -I$(LIB_SRC)/conversion/lnInclude \ + -I$(LIB_SRC)/fileFormats/lnInclude \ + -I$(LIB_SRC)/meshTools/lnInclude + +EXE_LIBS = \ + -lconversion \ + -lfileFormats \ + -lmeshTools + diff --git a/applications/test/ensightFile/Test-ensightFile.C b/applications/test/ensightFile/Test-ensightFile.C new file mode 100644 index 0000000000..6921176bc2 --- /dev/null +++ b/applications/test/ensightFile/Test-ensightFile.C @@ -0,0 +1,77 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2016 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 . + +Application + Test-ensightFile + +Description + check cleanup of ensight file and variable names + +\*---------------------------------------------------------------------------*/ + +#include "argList.H" +#include "ensightFileName.H" +#include "ensightVarName.H" +#include "IOstreams.H" + +using namespace Foam; + +void printCleaning(const fileName& pathName) +{ + Info<< "input = " << pathName << nl; + Info<< "file = " << ensight::FileName(pathName) << nl; + Info<< "var = " << ensight::VarName(pathName) << nl; + Info<< nl; +} + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // +// Main program: + +int main(int argc, char *argv[]) +{ + argList::noBanner(); + argList::noParallel(); + argList::validArgs.insert("fileName .. fileNameN"); + + argList args(argc, argv, false, true); + + if (args.size() <= 1 && args.options().empty()) + { + args.printUsage(); + } + + fileName pathName; + + for (label argI=1; argI < args.size(); ++argI) + { + pathName = args[argI]; + printCleaning(pathName); + } + + Info<< "\nEnd\n" << endl; + return 0; +} + + +// ************************************************************************* // diff --git a/applications/utilities/miscellaneous/foamHelp/helpTypes/helpSolver/helpSolver.C b/applications/utilities/miscellaneous/foamHelp/helpTypes/helpSolver/helpSolver.C index 116be3f552..cac59967c3 100644 --- a/applications/utilities/miscellaneous/foamHelp/helpTypes/helpSolver/helpSolver.C +++ b/applications/utilities/miscellaneous/foamHelp/helpTypes/helpSolver/helpSolver.C @@ -3,7 +3,7 @@ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | \\ / A nd | Copyright (C) 2015 OpenFOAM Foundation - \\/ M anipulation | + \\/ M anipulation | Copyright (C) 2016 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -25,6 +25,7 @@ License #include "helpSolver.H" #include "addToRunTimeSelectionTable.H" +#include "fvMesh.H" // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // @@ -63,6 +64,11 @@ void Foam::helpTypes::helpSolver::init() helpType::init(); argList::validArgs.append("solver"); + argList::addBoolOption + ( + "read", + "read solver type from the system/controlDict" + ); } @@ -78,6 +84,11 @@ void Foam::helpTypes::helpSolver::execute { displayDoc(solver, ".*solvers/.*Foam/", true, "C"); } + else if (args.optionFound("read")) + { + mesh.time().controlDict().lookup("application") >> solver; + displayDoc(solver, ".*solvers/.*Foam/", true, "C"); + } else { displayDocOptions(".*solvers/.*Foam/", true, "C"); diff --git a/applications/utilities/postProcessing/dataConversion/foamToEnsightParts/checkHasMovingMesh.H b/applications/utilities/postProcessing/dataConversion/foamToEnsightParts/checkHasMovingMesh.H index 454cb3ebb6..d59f258d35 100644 --- a/applications/utilities/postProcessing/dataConversion/foamToEnsightParts/checkHasMovingMesh.H +++ b/applications/utilities/postProcessing/dataConversion/foamToEnsightParts/checkHasMovingMesh.H @@ -1,4 +1,5 @@ // check for "points" in all of the result directories +// - could restrict to the selected times bool hasMovingMesh = false; if (timeDirs.size() > 1) diff --git a/applications/utilities/postProcessing/dataConversion/foamToEnsightParts/ensightOutputCase.H b/applications/utilities/postProcessing/dataConversion/foamToEnsightParts/ensightOutputCase.H index 0d22159079..83bae0f842 100644 --- a/applications/utilities/postProcessing/dataConversion/foamToEnsightParts/ensightOutputCase.H +++ b/applications/utilities/postProcessing/dataConversion/foamToEnsightParts/ensightOutputCase.H @@ -8,7 +8,7 @@ if (timeDirs[0].value() < 0) } // the case file is always ASCII -Info << "write case: " << caseFileName.c_str() << endl; +Info<< "write case: " << caseFileName.c_str() << endl; OFstream caseFile(ensightDir/caseFileName, IOstream::ASCII); caseFile.setf(ios_base::left); @@ -20,7 +20,9 @@ caseFile << "FORMAT" << nl << setw(16) << "type:" << "ensight gold" << nl << nl; -if (hasMovingMesh) +// time-set for geometries +// TODO: split off into separate time-set, but need to verify ensight spec +if (geometryTimesUsed.size()) { caseFile << "GEOMETRY" << nl @@ -43,7 +45,7 @@ forAllConstIter(HashTable>, cloudFields, cloudIter) caseFile << setw(16) << "measured: 2" << fileName(dataMask/cloud::prefix/cloudName/"positions").c_str() - << nl; + << nl; } caseFile << nl << "VARIABLE" << nl; @@ -87,9 +89,8 @@ forAllConstIter(HashTable, volumeFields, fieldIter) } } + // TODO: allow similar/different time-steps for each cloud - - label cloudNo = 0; forAllConstIter(HashTable>, cloudFields, cloudIter) { @@ -135,7 +136,7 @@ forAllConstIter(HashTable>, cloudFields, cloudIter) // add time values caseFile << nl << "TIME" << nl; -// time set 1 - geometry and volume fields +// time set 1 - volume fields if (fieldTimesUsed.size()) { caseFile @@ -161,9 +162,9 @@ if (fieldTimesUsed.size()) count = 0; forAll(fieldTimesUsed, i) { + const label& index = fieldTimesUsed[i]; caseFile - << " " << setw(12) - << timeIndices[fieldTimesUsed[i]] + timeCorrection; + << " " << setw(12) << timeIndices[index] + timeCorrection; if (++count % 6 == 0) { @@ -173,6 +174,49 @@ if (fieldTimesUsed.size()) caseFile << nl << nl; } + +// time set 2 - geometry +// THIS NEEDS MORE CHECKING +#if 0 +if (geometryTimesUsed.size()) +{ + caseFile + << "time set: " << 2 << nl + << "number of steps: " << geometryTimesUsed.size() << nl + << "filename numbers:" << nl; + + label count = 0; + forAll(geometryTimesUsed, i) + { + caseFile + << " " << setw(12) << geometryTimesUsed[i]; + + if (++count % 6 == 0) + { + caseFile << nl; + } + } + + caseFile + << nl << "time values:" << nl; + + count = 0; + forAll(geometryTimesUsed, i) + { + const label& index = geometryTimesUsed[i]; + caseFile + << " " << setw(12) << timeIndices[index] + timeCorrection; + + if (++count % 6 == 0) + { + caseFile << nl; + } + } + caseFile << nl << nl; +} +#endif + +// time set - clouds // TODO: allow similar/different time-steps for each cloud cloudNo = 0; forAllConstIter(HashTable>, cloudTimesUsed, cloudIter) @@ -205,9 +249,9 @@ forAllConstIter(HashTable>, cloudTimesUsed, cloudIter) count = 0; forAll(timesUsed, i) { + const label& index = timesUsed[i]; caseFile - << " " << setw(12) - << timeIndices[timesUsed[i]] + timeCorrection; + << " " << setw(12) << timeIndices[index] + timeCorrection; if (++count % 6 == 0) { diff --git a/applications/utilities/postProcessing/dataConversion/foamToEnsightParts/ensightOutputFunctions.C b/applications/utilities/postProcessing/dataConversion/foamToEnsightParts/ensightOutputFunctions.C index 309b49a82d..2c92bee9cc 100644 --- a/applications/utilities/postProcessing/dataConversion/foamToEnsightParts/ensightOutputFunctions.C +++ b/applications/utilities/postProcessing/dataConversion/foamToEnsightParts/ensightOutputFunctions.C @@ -3,7 +3,7 @@ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | \\ / A nd | Copyright (C) 2011-2016 OpenFOAM Foundation - \\/ M anipulation | + \\/ M anipulation | Copyright (C) 2016 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -47,6 +47,8 @@ void Foam::ensightCaseEntry const label timeSet ) { + const ensight::VarName varName(fieldName); + caseFile.setf(ios_base::left); fileName dirName(dataMask); @@ -68,9 +70,9 @@ void Foam::ensightCaseEntry << ensightType.c_str() << " per measured node: " << ts << " " << setw(15) - << ("c" + Foam::name(cloudNo) + fieldName).c_str() + << ("c" + Foam::name(cloudNo) + varName).c_str() << " " - << (dirName/fieldName).c_str() + << (dirName/varName).c_str() << nl; } else @@ -78,9 +80,9 @@ void Foam::ensightCaseEntry caseFile << ensightType.c_str() << " per element: " - << setw(15) << fieldName + << setw(15) << varName << " " - << (dirName/fieldName).c_str() + << (dirName/varName).c_str() << nl; } } @@ -97,16 +99,16 @@ void Foam::ensightParticlePositions { Cloud parcels(mesh, cloudName, false); - fileName cloudDir = subDir/cloud::prefix/cloudName; - fileName postFileName = cloudDir/"positions"; + const fileName postFileName = + subDir/cloud::prefix/cloudName/"positions"; // the ITER/lagrangian subdirectory must exist - mkDir(dataDir/cloudDir); - ensightFile os(dataDir/postFileName, format); + mkDir(dataDir/postFileName.path()); + ensightFile os(dataDir, postFileName, format); // tag binary format (just like geometry files) os.writeBinaryHeader(); - os.write(postFileName); + os.write(postFileName); // description os.newline(); os.write("particle coordinates"); os.newline(); @@ -161,14 +163,19 @@ void Foam::ensightLagrangianField { Info<< " " << fieldObject.name() << flush; - fileName cloudDir = subDir/cloud::prefix/cloudName; - fileName postFileName = cloudDir/fieldObject.name(); + const fileName postFileName = + subDir/cloud::prefix/cloudName + /ensight::VarName(fieldObject.name()); - string title = - postFileName + " with " + pTraits::typeName + " values"; + // the ITER/lagrangian subdirectory was already created + // when writing positions - ensightFile os(dataDir/postFileName, format); - os.write(title); + ensightFile os(dataDir, postFileName, format); + os.write + ( + // description + string(postFileName + " with " + pTraits::typeName + " values") + ); os.newline(); IOField field(fieldObject); @@ -219,10 +226,10 @@ void Foam::ensightVolField { Info<< " " << fieldObject.name() << flush; - fileName postFileName = subDir/fieldObject.name(); + const fileName postFileName = subDir/ensight::VarName(fieldObject.name()); - ensightFile os(dataDir/postFileName, format); - os.write(postFileName); + ensightFile os(dataDir, postFileName, format); + os.write(postFileName); // description os.newline(); // ie, volField diff --git a/applications/utilities/postProcessing/dataConversion/foamToEnsightParts/foamToEnsightParts.C b/applications/utilities/postProcessing/dataConversion/foamToEnsightParts/foamToEnsightParts.C index 6aecc40fc2..0a3cdf6509 100644 --- a/applications/utilities/postProcessing/dataConversion/foamToEnsightParts/foamToEnsightParts.C +++ b/applications/utilities/postProcessing/dataConversion/foamToEnsightParts/foamToEnsightParts.C @@ -3,7 +3,7 @@ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | \\ / A nd | Copyright (C) 2011-2016 OpenFOAM Foundation - \\/ M anipulation | + \\/ M anipulation | Copyright (C) 2016 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -166,9 +166,9 @@ int main(int argc, char *argv[]) ensightDir = args.rootPath()/args.globalCaseName()/ensightDir; } - fileName dataDir = ensightDir/"data"; - fileName caseFileName = "Ensight.case"; - fileName dataMask = fileName("data")/ensightFile::mask(); + const fileName caseFileName = "Ensight.case"; + const fileName dataDir = ensightDir/"data"; + const fileName dataMask = dataDir.name()/ensightFile::mask(); // Ensight and Ensight/data directories must exist // do not remove old data - we might wish to convert new results @@ -178,6 +178,8 @@ int main(int argc, char *argv[]) Info<<"Warning: re-using existing directory" << nl << " " << ensightDir << endl; } + + // as per mkdir -p "Ensight/data" mkDir(ensightDir); mkDir(dataDir); @@ -216,6 +218,9 @@ int main(int argc, char *argv[]) // map times used Map timeIndices; + // TODO: Track the time indices used by the geometry + DynamicList