From 26418ee9d3f7368d69ccb0b32ff7dfd2d47bc572 Mon Sep 17 00:00:00 2001 From: Will Bainbridge Date: Fri, 16 Feb 2018 11:01:48 +0000 Subject: [PATCH] BlendedInterfacialModel: Unified blending functions --- .../BlendedInterfacialModel.C | 482 ++++++------------ .../BlendedInterfacialModel.H | 42 +- 2 files changed, 180 insertions(+), 344 deletions(-) diff --git a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/BlendedInterfacialModel/BlendedInterfacialModel.C b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/BlendedInterfacialModel/BlendedInterfacialModel.C index 64345bf759..7ba04bfc56 100644 --- a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/BlendedInterfacialModel/BlendedInterfacialModel.C +++ b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/BlendedInterfacialModel/BlendedInterfacialModel.C @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | - \\ / A nd | Copyright (C) 2014-2016 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2014-2018 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -27,17 +27,42 @@ License #include "fixedValueFvsPatchFields.H" #include "surfaceInterpolate.H" +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + +// * * * * * * * * * * * * * Static Member Functions * * * * * * * * * * * * // + +template<> +inline tmp +blendedInterfacialModel::interpolate(tmp f) +{ + return f; +} + + +template<> +inline tmp +blendedInterfacialModel::interpolate(tmp f) +{ + return fvc::interpolate(f); +} + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace Foam + // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * // template -template +template void Foam::BlendedInterfacialModel::correctFixedFluxBCs ( - GeometricField& field + GeoField& field ) const { - typename GeometricField::Boundary& fieldBf = - field.boundaryFieldRef(); + typename GeoField::Boundary& fieldBf = field.boundaryFieldRef(); forAll(phase1_.phi()().boundaryField(), patchi) { @@ -55,6 +80,111 @@ void Foam::BlendedInterfacialModel::correctFixedFluxBCs } +template +template +< + class Type, + template class PatchField, + class GeoMesh, + class ... Args +> +Foam::tmp> +Foam::BlendedInterfacialModel::evaluate +( + tmp> + (ModelType::*method)(Args ...) const, + const word& name, + const dimensionSet& dims, + const bool subtract, + Args ... args +) const +{ + typedef GeometricField scalarGeoField; + typedef GeometricField typeGeoField; + + tmp f1, f2; + + if (model_.valid() || model1In2_.valid()) + { + f1 = + blendedInterfacialModel::interpolate + ( + blending_.f1(phase1_, phase2_) + ); + } + + if (model_.valid() || model2In1_.valid()) + { + f2 = + blendedInterfacialModel::interpolate + ( + blending_.f2(phase1_, phase2_) + ); + } + + tmp x + ( + new typeGeoField + ( + IOobject + ( + ModelType::typeName + ":" + name, + phase1_.mesh().time().timeName(), + phase1_.mesh(), + IOobject::NO_READ, + IOobject::NO_WRITE, + false + ), + phase1_.mesh(), + dimensioned("zero", dims, Zero) + ) + ); + + if (model_.valid()) + { + if (subtract) + { + FatalErrorInFunction + << "Cannot treat an interfacial model with no distinction " + << "between continuous and dispersed phases as signed" + << exit(FatalError); + } + + x.ref() += (model_().*method)(args ...)*(scalar(1) - f1() - f2()); + } + + if (model1In2_.valid()) + { + x.ref() += (model1In2_().*method)(args ...)*f1; + } + + if (model2In1_.valid()) + { + tmp dx = (model2In1_().*method)(args ...)*f2; + + if (subtract) + { + x.ref() -= dx; + } + else + { + x.ref() += dx; + } + } + + if + ( + correctFixedFluxBCs_ + && (model_.valid() || model1In2_.valid() || model2In1_.valid()) + ) + { + correctFixedFluxBCs(x.ref()); + } + + return x; +} + + // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // template @@ -169,59 +299,9 @@ template Foam::tmp Foam::BlendedInterfacialModel::K() const { - tmp f1, f2; + tmp (ModelType::*k)() const = &ModelType::K; - if (model_.valid() || model1In2_.valid()) - { - f1 = blending_.f1(phase1_, phase2_); - } - - if (model_.valid() || model2In1_.valid()) - { - f2 = blending_.f2(phase1_, phase2_); - } - - tmp x - ( - new volScalarField - ( - IOobject - ( - ModelType::typeName + ":K", - phase1_.mesh().time().timeName(), - phase1_.mesh(), - IOobject::NO_READ, - IOobject::NO_WRITE, - false - ), - phase1_.mesh(), - dimensionedScalar("zero", ModelType::dimK, 0) - ) - ); - - if (model_.valid()) - { - x.ref() += model_->K()*(scalar(1) - f1() - f2()); - } - if (model1In2_.valid()) - { - x.ref() += model1In2_->K()*f1; - } - if (model2In1_.valid()) - { - x.ref() += model2In1_->K()*f2; - } - - if - ( - correctFixedFluxBCs_ - && (model_.valid() || model1In2_.valid() || model2In1_.valid()) - ) - { - correctFixedFluxBCs(x.ref()); - } - - return x; + return evaluate(k, "K", ModelType::dimK, false); } @@ -229,59 +309,9 @@ template Foam::tmp Foam::BlendedInterfacialModel::K(const scalar residualAlpha) const { - tmp f1, f2; + tmp (ModelType::*k)(const scalar) const = &ModelType::K; - if (model_.valid() || model1In2_.valid()) - { - f1 = blending_.f1(phase1_, phase2_); - } - - if (model_.valid() || model2In1_.valid()) - { - f2 = blending_.f2(phase1_, phase2_); - } - - tmp x - ( - new volScalarField - ( - IOobject - ( - ModelType::typeName + ":K", - phase1_.mesh().time().timeName(), - phase1_.mesh(), - IOobject::NO_READ, - IOobject::NO_WRITE, - false - ), - phase1_.mesh(), - dimensionedScalar("zero", ModelType::dimK, 0) - ) - ); - - if (model_.valid()) - { - x.ref() += model_->K(residualAlpha)*(scalar(1) - f1() - f2()); - } - if (model1In2_.valid()) - { - x.ref() += model1In2_->K(residualAlpha)*f1; - } - if (model2In1_.valid()) - { - x.ref() += model2In1_->K(residualAlpha)*f2; - } - - if - ( - correctFixedFluxBCs_ - && (model_.valid() || model1In2_.valid() || model2In1_.valid()) - ) - { - correctFixedFluxBCs(x.ref()); - } - - return x; + return evaluate(k, "K", ModelType::dimK, false, residualAlpha); } @@ -289,67 +319,7 @@ template Foam::tmp Foam::BlendedInterfacialModel::Kf() const { - tmp f1, f2; - - if (model_.valid() || model1In2_.valid()) - { - f1 = fvc::interpolate - ( - blending_.f1(phase1_, phase2_) - ); - } - - if (model_.valid() || model2In1_.valid()) - { - f2 = fvc::interpolate - ( - blending_.f2(phase1_, phase2_) - ); - } - - tmp x - ( - new surfaceScalarField - ( - IOobject - ( - ModelType::typeName + ":Kf", - phase1_.mesh().time().timeName(), - phase1_.mesh(), - IOobject::NO_READ, - IOobject::NO_WRITE, - false - ), - phase1_.mesh(), - dimensionedScalar("zero", ModelType::dimK, 0) - ) - ); - - if (model_.valid()) - { - x.ref() += model_->Kf()*(scalar(1) - f1() - f2()); - } - - if (model1In2_.valid()) - { - x.ref() += model1In2_->Kf()*f1; - } - - if (model2In1_.valid()) - { - x.ref() += model2In1_->Kf()*f2; - } - - if - ( - correctFixedFluxBCs_ - && (model_.valid() || model1In2_.valid() || model2In1_.valid()) - ) - { - correctFixedFluxBCs(x.ref()); - } - - return x; + return evaluate(&ModelType::Kf, "Kf", ModelType::dimK, false); } @@ -358,61 +328,7 @@ template Foam::tmp> Foam::BlendedInterfacialModel::F() const { - tmp f1, f2; - - if (model_.valid() || model1In2_.valid()) - { - f1 = blending_.f1(phase1_, phase2_); - } - - if (model_.valid() || model2In1_.valid()) - { - f2 = blending_.f2(phase1_, phase2_); - } - - tmp> x - ( - new GeometricField - ( - IOobject - ( - ModelType::typeName + ":F", - phase1_.mesh().time().timeName(), - phase1_.mesh(), - IOobject::NO_READ, - IOobject::NO_WRITE, - false - ), - phase1_.mesh(), - dimensioned("zero", ModelType::dimF, Zero) - ) - ); - - if (model_.valid()) - { - x.ref() += model_->F()*(scalar(1) - f1() - f2()); - } - - if (model1In2_.valid()) - { - x.ref() += model1In2_->F()*f1; - } - - if (model2In1_.valid()) - { - x.ref() -= model2In1_->F()*f2; // note : subtraction - } - - if - ( - correctFixedFluxBCs_ - && (model_.valid() || model1In2_.valid() || model2In1_.valid()) - ) - { - correctFixedFluxBCs(x.ref()); - } - - return x; + return evaluate(&ModelType::F, "F", ModelType::dimF, true); } @@ -420,67 +336,7 @@ template Foam::tmp Foam::BlendedInterfacialModel::Ff() const { - tmp f1, f2; - - if (model_.valid() || model1In2_.valid()) - { - f1 = fvc::interpolate - ( - blending_.f1(phase1_, phase2_) - ); - } - - if (model_.valid() || model2In1_.valid()) - { - f2 = fvc::interpolate - ( - blending_.f2(phase1_, phase2_) - ); - } - - tmp x - ( - new surfaceScalarField - ( - IOobject - ( - ModelType::typeName + ":Ff", - phase1_.mesh().time().timeName(), - phase1_.mesh(), - IOobject::NO_READ, - IOobject::NO_WRITE, - false - ), - phase1_.mesh(), - dimensionedScalar("zero", ModelType::dimF*dimArea, 0) - ) - ); - - if (model_.valid()) - { - x.ref() += model_->Ff()*(scalar(1) - f1() - f2()); - } - - if (model1In2_.valid()) - { - x.ref() += model1In2_->Ff()*f1; - } - - if (model2In1_.valid()) - { - x.ref() -= model2In1_->Ff()*f2; // note : subtraction - } - - if - ( - correctFixedFluxBCs_ - && (model_.valid() || model1In2_.valid() || model2In1_.valid()) - ) - { - correctFixedFluxBCs(x.ref()); - } - - return x; + return evaluate(&ModelType::Ff, "Ff", ModelType::dimF*dimArea, true); } @@ -488,59 +344,7 @@ template Foam::tmp Foam::BlendedInterfacialModel::D() const { - tmp f1, f2; - - if (model_.valid() || model1In2_.valid()) - { - f1 = blending_.f1(phase1_, phase2_); - } - - if (model_.valid() || model2In1_.valid()) - { - f2 = blending_.f2(phase1_, phase2_); - } - - tmp x - ( - new volScalarField - ( - IOobject - ( - ModelType::typeName + ":D", - phase1_.mesh().time().timeName(), - phase1_.mesh(), - IOobject::NO_READ, - IOobject::NO_WRITE, - false - ), - phase1_.mesh(), - dimensionedScalar("zero", ModelType::dimD, 0) - ) - ); - - if (model_.valid()) - { - x.ref() += model_->D()*(scalar(1) - f1() - f2()); - } - if (model1In2_.valid()) - { - x.ref() += model1In2_->D()*f1; - } - if (model2In1_.valid()) - { - x.ref() += model2In1_->D()*f2; - } - - if - ( - correctFixedFluxBCs_ - && (model_.valid() || model1In2_.valid() || model2In1_.valid()) - ) - { - correctFixedFluxBCs(x.ref()); - } - - return x; + return evaluate(&ModelType::D, "D", ModelType::dimD, false); } diff --git a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/BlendedInterfacialModel/BlendedInterfacialModel.H b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/BlendedInterfacialModel/BlendedInterfacialModel.H index c51599d58c..4debb7db03 100644 --- a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/BlendedInterfacialModel/BlendedInterfacialModel.H +++ b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/BlendedInterfacialModel/BlendedInterfacialModel.H @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | - \\ / A nd | Copyright (C) 2014-2016 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2014-2018 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -37,7 +37,6 @@ SourceFiles #include "blendingMethod.H" #include "phasePair.H" #include "orderedPhasePair.H" - #include "geometricZeroField.H" // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // @@ -46,7 +45,22 @@ namespace Foam { /*---------------------------------------------------------------------------*\ - Class BlendedInterfacialModel Declaration + Class blendedInterfacialModel Declaration +\*---------------------------------------------------------------------------*/ + +class blendedInterfacialModel +{ + public: + + //- Convenience function to interpolate blending values. Needs to be + // specialised, so can't sit in the templated class. + template + static inline tmp interpolate(tmp f); +}; + + +/*---------------------------------------------------------------------------*\ + Class BlendedInterfacialModel Declaration \*---------------------------------------------------------------------------*/ template @@ -85,8 +99,26 @@ class BlendedInterfacialModel void operator=(const BlendedInterfacialModel&); //- Correct coeff/value on fixed flux boundary conditions - template - void correctFixedFluxBCs(GeometricField& field) const; + template + void correctFixedFluxBCs(GeoField& field) const; + + //- Return the blended coeff/value + template + < + class Type, + template class PatchField, + class GeoMesh, + class ... Args + > + tmp> evaluate + ( + tmp> + (ModelType::*method)(Args ...) const, + const word& name, + const dimensionSet& dims, + const bool subtract, + Args ... args + ) const; public: