From ec5b93c902f5fa12b8b8a2319cd2528dd6b263c8 Mon Sep 17 00:00:00 2001 From: Andrew Heather <> Date: Wed, 24 Mar 2021 13:12:38 +0000 Subject: [PATCH] ENH: relative velocity usage for porosity in MRF regions - New optional keyword 'MRFRelativeVelocity' (default=true) usable by the porosity models. For example, in 'constant/porosityProperties': porosity1 { type DarcyForchheimer; ... MRFRelativeVelocity yes; ... } --- .../powerLawLopesdaCosta.C | 47 +++++----- .../powerLawLopesdaCosta.H | 25 ++++-- .../powerLawLopesdaCostaTemplates.C | 9 +- .../DarcyForchheimer/DarcyForchheimer.C | 70 +++++++-------- .../DarcyForchheimer/DarcyForchheimer.H | 25 ++++-- .../DarcyForchheimerTemplates.C | 11 +-- .../porosityModel/fixedCoeff/fixedCoeff.C | 49 +++++----- .../porosityModel/fixedCoeff/fixedCoeff.H | 25 ++++-- .../porosityModel/porosityModel.C | 89 ++++++++++++++++--- .../porosityModel/porosityModel.H | 24 +++-- .../porosityModel/porosityModelList.C | 25 +++--- .../general/porosityModel/powerLaw/powerLaw.C | 33 +++---- .../general/porosityModel/powerLaw/powerLaw.H | 23 +++-- .../powerLaw/powerLawTemplates.C | 10 +-- .../solidification/solidification.C | 39 ++++---- .../solidification/solidification.H | 33 ++++--- .../solidification/solidificationTemplates.C | 25 +++--- .../interRegionExplicitPorositySource.C | 4 +- 18 files changed, 336 insertions(+), 230 deletions(-) diff --git a/src/atmosphericModels/porosityModels/powerLawLopesdaCosta/powerLawLopesdaCosta.C b/src/atmosphericModels/porosityModels/powerLawLopesdaCosta/powerLawLopesdaCosta.C index c8e41752682..a16fe8e34a4 100644 --- a/src/atmosphericModels/porosityModels/powerLawLopesdaCosta/powerLawLopesdaCosta.C +++ b/src/atmosphericModels/porosityModels/powerLawLopesdaCosta/powerLawLopesdaCosta.C @@ -6,7 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2018 OpenFOAM Foundation - Copyright (C) 2020 OpenCFD Ltd. + Copyright (C) 2020-2021 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -314,12 +314,6 @@ Foam::porosityModels::powerLawLopesdaCosta::powerLawLopesdaCosta {} -// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * // - -Foam::porosityModels::powerLawLopesdaCosta::~powerLawLopesdaCosta() -{} - - // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // const Foam::scalarField& @@ -344,47 +338,47 @@ void Foam::porosityModels::powerLawLopesdaCosta::calcForce scalarField Udiag(U.size(), Zero); const scalarField& V = mesh_.V(); - apply(Udiag, V, rho, U); + apply(V, rho, U, Udiag); - force = Udiag*U; + force += Udiag*U; } void Foam::porosityModels::powerLawLopesdaCosta::correct ( - fvVectorMatrix& UEqn + const volVectorField& U, + scalarField& Udiag, + vectorField& Usource, + bool compressible ) const { - const vectorField& U = UEqn.psi(); const scalarField& V = mesh_.V(); - scalarField& Udiag = UEqn.diag(); - if (UEqn.dimensions() == dimForce) + if (compressible) { - const volScalarField& rho = - mesh_.lookupObject<volScalarField>(rhoName_); + const auto& rho = mesh_.lookupObject<volScalarField>(rhoName_); - apply(Udiag, V, rho, U); + apply(V, rho, U, Udiag); } else { - apply(Udiag, V, geometricOneField(), U); + apply(V, geometricOneField(), U, Udiag); } } void Foam::porosityModels::powerLawLopesdaCosta::correct ( - fvVectorMatrix& UEqn, + const volVectorField& U, const volScalarField& rho, - const volScalarField& mu + const volScalarField& mu, + scalarField& Udiag, + vectorField& Usource ) const { - const vectorField& U = UEqn.psi(); const scalarField& V = mesh_.V(); - scalarField& Udiag = UEqn.diag(); - apply(Udiag, V, rho, U); + apply(V, rho, U, Udiag); } @@ -394,18 +388,17 @@ void Foam::porosityModels::powerLawLopesdaCosta::correct volTensorField& AU ) const { - const vectorField& U = UEqn.psi(); + const volVectorField& U = UEqn.psi(); if (UEqn.dimensions() == dimForce) { - const volScalarField& rho = - mesh_.lookupObject<volScalarField>(rhoName_); + const auto& rho = mesh_.lookupObject<volScalarField>(rhoName_); - apply(AU, rho, U); + apply(rho, U, AU); } else { - apply(AU, geometricOneField(), U); + apply(geometricOneField(), U, AU); } } diff --git a/src/atmosphericModels/porosityModels/powerLawLopesdaCosta/powerLawLopesdaCosta.H b/src/atmosphericModels/porosityModels/powerLawLopesdaCosta/powerLawLopesdaCosta.H index 0fb1fd23f75..b52457bcc17 100644 --- a/src/atmosphericModels/porosityModels/powerLawLopesdaCosta/powerLawLopesdaCosta.H +++ b/src/atmosphericModels/porosityModels/powerLawLopesdaCosta/powerLawLopesdaCosta.H @@ -6,6 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2018 OpenFOAM Foundation + Copyright (C) 2021 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -133,19 +134,19 @@ class powerLawLopesdaCosta template<class RhoFieldType> void apply ( - scalarField& Udiag, const scalarField& V, const RhoFieldType& rho, - const vectorField& U + const vectorField& U, + scalarField& Udiag ) const; //- Apply resistance template<class RhoFieldType> void apply ( - tensorField& AU, const RhoFieldType& rho, - const vectorField& U + const vectorField& U, + tensorField& AU ) const; @@ -172,7 +173,7 @@ public: ); //- Destructor - virtual ~powerLawLopesdaCosta(); + virtual ~powerLawLopesdaCosta() = default; // Member Functions @@ -190,14 +191,22 @@ public: ) const; //- Add resistance - virtual void correct(fvVectorMatrix& UEqn) const; + virtual void correct + ( + const volVectorField& U, + scalarField& Udiag, + vectorField& Usource, + bool compressible + ) const; //- Add resistance virtual void correct ( - fvVectorMatrix& UEqn, + const volVectorField& U, const volScalarField& rho, - const volScalarField& mu + const volScalarField& mu, + scalarField& Udiag, + vectorField& Usource ) const; //- Add resistance diff --git a/src/atmosphericModels/porosityModels/powerLawLopesdaCosta/powerLawLopesdaCostaTemplates.C b/src/atmosphericModels/porosityModels/powerLawLopesdaCosta/powerLawLopesdaCostaTemplates.C index 7191843f0f7..637cfa0d4d6 100644 --- a/src/atmosphericModels/porosityModels/powerLawLopesdaCosta/powerLawLopesdaCostaTemplates.C +++ b/src/atmosphericModels/porosityModels/powerLawLopesdaCosta/powerLawLopesdaCostaTemplates.C @@ -6,6 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2018 OpenFOAM Foundation + Copyright (C) 2021 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -32,10 +33,10 @@ License template<class RhoFieldType> void Foam::porosityModels::powerLawLopesdaCosta::apply ( - scalarField& Udiag, const scalarField& V, const RhoFieldType& rho, - const vectorField& U + const vectorField& U, + scalarField& Udiag ) const { const scalar C1m1b2 = (C1_ - 1.0)/2.0; @@ -59,9 +60,9 @@ void Foam::porosityModels::powerLawLopesdaCosta::apply template<class RhoFieldType> void Foam::porosityModels::powerLawLopesdaCosta::apply ( - tensorField& AU, const RhoFieldType& rho, - const vectorField& U + const vectorField& U, + tensorField& AU ) const { const scalar C1m1b2 = (C1_ - 1.0)/2.0; diff --git a/src/finiteVolume/cfdTools/general/porosityModel/DarcyForchheimer/DarcyForchheimer.C b/src/finiteVolume/cfdTools/general/porosityModel/DarcyForchheimer/DarcyForchheimer.C index f2c9f856d5e..6f1fdf9f659 100644 --- a/src/finiteVolume/cfdTools/general/porosityModel/DarcyForchheimer/DarcyForchheimer.C +++ b/src/finiteVolume/cfdTools/general/porosityModel/DarcyForchheimer/DarcyForchheimer.C @@ -6,7 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2012-2016 OpenFOAM Foundation - Copyright (C) 2018-2020 OpenCFD Ltd. + Copyright (C) 2018-2021 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -175,79 +175,79 @@ void Foam::porosityModels::DarcyForchheimer::calcForce vectorField& force ) const { - scalarField Udiag(U.size(), Zero); - vectorField Usource(U.size(), Zero); const scalarField& V = mesh_.V(); - apply(Udiag, Usource, V, rho, mu, U); + scalarField Udiag(U.size(), Zero); + vectorField Usource(U.size(), Zero); + apply(V, rho, mu, U, Udiag, Usource); - force = Udiag*U - Usource; + force += Udiag*U - Usource; } void Foam::porosityModels::DarcyForchheimer::correct ( - fvVectorMatrix& UEqn + const volVectorField& U, + scalarField& Udiag, + vectorField& Usource, + bool compressible ) const { - const volVectorField& U = UEqn.psi(); const scalarField& V = mesh_.V(); - scalarField& Udiag = UEqn.diag(); - vectorField& Usource = UEqn.source(); - word rhoName(IOobject::groupName(rhoName_, U.group())); - word muName(IOobject::groupName(muName_, U.group())); - word nuName(IOobject::groupName(nuName_, U.group())); + const word rhoName(IOobject::groupName(rhoName_, U.group())); + const word muName(IOobject::groupName(muName_, U.group())); + const word nuName(IOobject::groupName(nuName_, U.group())); - if (UEqn.dimensions() == dimForce) + if (compressible) { const auto& rho = mesh_.lookupObject<volScalarField>(rhoName); + tmp<volScalarField> tmu; if (mesh_.foundObject<volScalarField>(muName)) { - const auto& mu = mesh_.lookupObject<volScalarField>(muName); - - apply(Udiag, Usource, V, rho, mu, U); + tmu = mesh_.lookupObject<volScalarField>(muName); } else { const auto& nu = mesh_.lookupObject<volScalarField>(nuName); - - apply(Udiag, Usource, V, rho, rho*nu, U); + tmu = rho*nu; } + + apply(V, rho, tmu(), U, Udiag, Usource); } else { + tmp<volScalarField> tnu; + if (mesh_.foundObject<volScalarField>(nuName)) { - const auto& nu = mesh_.lookupObject<volScalarField>(nuName); - - apply(Udiag, Usource, V, geometricOneField(), nu, U); + tnu = mesh_.lookupObject<volScalarField>(nuName); } else { const auto& rho = mesh_.lookupObject<volScalarField>(rhoName); const auto& mu = mesh_.lookupObject<volScalarField>(muName); - - apply(Udiag, Usource, V, geometricOneField(), mu/rho, U); + tnu = mu/rho; } + + apply(V, geometricOneField(), tnu(), U, Udiag, Usource); } } void Foam::porosityModels::DarcyForchheimer::correct ( - fvVectorMatrix& UEqn, + const volVectorField& U, const volScalarField& rho, - const volScalarField& mu + const volScalarField& mu, + scalarField& Udiag, + vectorField& Usource ) const { - const vectorField& U = UEqn.psi(); const scalarField& V = mesh_.V(); - scalarField& Udiag = UEqn.diag(); - vectorField& Usource = UEqn.source(); - apply(Udiag, Usource, V, rho, mu, U); + apply(V, rho, mu, U, Udiag, Usource); } @@ -259,16 +259,16 @@ void Foam::porosityModels::DarcyForchheimer::correct { const volVectorField& U = UEqn.psi(); - word rhoName(IOobject::groupName(rhoName_, U.group())); - word muName(IOobject::groupName(muName_, U.group())); - word nuName(IOobject::groupName(nuName_, U.group())); + const word rhoName(IOobject::groupName(rhoName_, U.group())); + const word muName(IOobject::groupName(muName_, U.group())); + const word nuName(IOobject::groupName(nuName_, U.group())); if (UEqn.dimensions() == dimForce) { const auto& rho = mesh_.lookupObject<volScalarField>(rhoName); const auto& mu = mesh_.lookupObject<volScalarField>(muName); - apply(AU, rho, mu, U); + apply(rho, mu, U, AU); } else { @@ -276,14 +276,14 @@ void Foam::porosityModels::DarcyForchheimer::correct { const auto& nu = mesh_.lookupObject<volScalarField>(nuName); - apply(AU, geometricOneField(), nu, U); + apply(geometricOneField(), nu, U, AU); } else { const auto& rho = mesh_.lookupObject<volScalarField>(rhoName); const auto& mu = mesh_.lookupObject<volScalarField>(muName); - apply(AU, geometricOneField(), mu/rho, U); + apply(geometricOneField(), mu/rho, U, AU); } } } diff --git a/src/finiteVolume/cfdTools/general/porosityModel/DarcyForchheimer/DarcyForchheimer.H b/src/finiteVolume/cfdTools/general/porosityModel/DarcyForchheimer/DarcyForchheimer.H index 98e57ba530f..90b7a74eeff 100644 --- a/src/finiteVolume/cfdTools/general/porosityModel/DarcyForchheimer/DarcyForchheimer.H +++ b/src/finiteVolume/cfdTools/general/porosityModel/DarcyForchheimer/DarcyForchheimer.H @@ -6,6 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2012-2017 OpenFOAM Foundation + Copyright (C) 2021 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -102,22 +103,22 @@ class DarcyForchheimer template<class RhoFieldType> void apply ( - scalarField& Udiag, - vectorField& Usource, const scalarField& V, const RhoFieldType& rho, const scalarField& mu, - const vectorField& U + const vectorField& U, + scalarField& Udiag, + vectorField& Usource ) const; //- Apply template<class RhoFieldType> void apply ( - tensorField& AU, const RhoFieldType& rho, const scalarField& mu, - const vectorField& U + const vectorField& U, + tensorField& AU ) const; //- No copy construct @@ -161,14 +162,22 @@ public: ) const; //- Add resistance - virtual void correct(fvVectorMatrix& UEqn) const; + virtual void correct + ( + const volVectorField& U, + scalarField& Udiag, + vectorField& Usource, + bool compressible + ) const; //- Add resistance virtual void correct ( - fvVectorMatrix& UEqn, + const volVectorField& U, const volScalarField& rho, - const volScalarField& mu + const volScalarField& mu, + scalarField& Udiag, + vectorField& Usource ) const; //- Add resistance diff --git a/src/finiteVolume/cfdTools/general/porosityModel/DarcyForchheimer/DarcyForchheimerTemplates.C b/src/finiteVolume/cfdTools/general/porosityModel/DarcyForchheimer/DarcyForchheimerTemplates.C index e3330d13d6c..982273a73d9 100644 --- a/src/finiteVolume/cfdTools/general/porosityModel/DarcyForchheimer/DarcyForchheimerTemplates.C +++ b/src/finiteVolume/cfdTools/general/porosityModel/DarcyForchheimer/DarcyForchheimerTemplates.C @@ -6,6 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2012-2016 OpenFOAM Foundation + Copyright (C) 2021 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -30,12 +31,12 @@ License template<class RhoFieldType> void Foam::porosityModels::DarcyForchheimer::apply ( - scalarField& Udiag, - vectorField& Usource, const scalarField& V, const RhoFieldType& rho, const scalarField& mu, - const vectorField& U + const vectorField& U, + scalarField& Udiag, + vectorField& Usource ) const { forAll(cellZoneIDs_, zoneI) @@ -64,10 +65,10 @@ void Foam::porosityModels::DarcyForchheimer::apply template<class RhoFieldType> void Foam::porosityModels::DarcyForchheimer::apply ( - tensorField& AU, const RhoFieldType& rho, const scalarField& mu, - const vectorField& U + const vectorField& U, + tensorField& AU ) const { forAll(cellZoneIDs_, zoneI) diff --git a/src/finiteVolume/cfdTools/general/porosityModel/fixedCoeff/fixedCoeff.C b/src/finiteVolume/cfdTools/general/porosityModel/fixedCoeff/fixedCoeff.C index 10da2270c2e..66bf2807196 100644 --- a/src/finiteVolume/cfdTools/general/porosityModel/fixedCoeff/fixedCoeff.C +++ b/src/finiteVolume/cfdTools/general/porosityModel/fixedCoeff/fixedCoeff.C @@ -6,7 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2012-2016 OpenFOAM Foundation - Copyright (C) 2018 OpenCFD Ltd. + Copyright (C) 2018-2021 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -47,11 +47,11 @@ namespace Foam void Foam::porosityModels::fixedCoeff::apply ( - scalarField& Udiag, - vectorField& Usource, const scalarField& V, const vectorField& U, - const scalar rho + const scalar rho, + scalarField& Udiag, + vectorField& Usource ) const { forAll(cellZoneIDs_, zoneI) @@ -77,9 +77,9 @@ void Foam::porosityModels::fixedCoeff::apply void Foam::porosityModels::fixedCoeff::apply ( - tensorField& AU, const vectorField& U, - const scalar rho + const scalar rho, + tensorField& AU ) const { forAll(cellZoneIDs_, zoneI) @@ -183,51 +183,46 @@ void Foam::porosityModels::fixedCoeff::calcForce const scalarField& V = mesh_.V(); const scalar rhoRef = coeffs_.get<scalar>("rhoRef"); - apply(Udiag, Usource, V, U, rhoRef); + apply(V, U, rhoRef, Udiag, Usource); - force = Udiag*U - Usource; + force += Udiag*U - Usource; } void Foam::porosityModels::fixedCoeff::correct ( - fvVectorMatrix& UEqn + const volVectorField& U, + scalarField& Udiag, + vectorField& Usource, + bool compressible ) const { - const vectorField& U = UEqn.psi(); const scalarField& V = mesh_.V(); - scalarField& Udiag = UEqn.diag(); - vectorField& Usource = UEqn.source(); scalar rho = 1.0; - if (UEqn.dimensions() == dimForce) + if (compressible) { coeffs_.readEntry("rhoRef", rho); } - apply(Udiag, Usource, V, U, rho); + apply(V, U, rho, Udiag, Usource); } void Foam::porosityModels::fixedCoeff::correct ( - fvVectorMatrix& UEqn, + const volVectorField& U, const volScalarField&, - const volScalarField& + const volScalarField&, + scalarField& Udiag, + vectorField& Usource ) const { - const vectorField& U = UEqn.psi(); const scalarField& V = mesh_.V(); - scalarField& Udiag = UEqn.diag(); - vectorField& Usource = UEqn.source(); - scalar rho = 1.0; - if (UEqn.dimensions() == dimForce) - { - coeffs_.readEntry("rhoRef", rho); - } + const scalar rho = coeffs_.get<scalar>("rhoRef"); - apply(Udiag, Usource, V, U, rho); + apply(V, U, rho, Udiag, Usource); } @@ -237,15 +232,13 @@ void Foam::porosityModels::fixedCoeff::correct volTensorField& AU ) const { - const vectorField& U = UEqn.psi(); - scalar rho = 1.0; if (UEqn.dimensions() == dimForce) { coeffs_.readEntry("rhoRef", rho); } - apply(AU, U, rho); + apply(UEqn.psi(), rho, AU); } diff --git a/src/finiteVolume/cfdTools/general/porosityModel/fixedCoeff/fixedCoeff.H b/src/finiteVolume/cfdTools/general/porosityModel/fixedCoeff/fixedCoeff.H index c198257d699..885c757183c 100644 --- a/src/finiteVolume/cfdTools/general/porosityModel/fixedCoeff/fixedCoeff.H +++ b/src/finiteVolume/cfdTools/general/porosityModel/fixedCoeff/fixedCoeff.H @@ -6,6 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2012-2017 OpenFOAM Foundation + Copyright (C) 2021 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -82,19 +83,19 @@ class fixedCoeff //- Apply void apply ( - scalarField& Udiag, - vectorField& Usource, const scalarField& V, const vectorField& U, - const scalar rho + const scalar rho, + scalarField& Udiag, + vectorField& Usource ) const; //- Apply void apply ( - tensorField& AU, const vectorField& U, - const scalar rho + const scalar rho, + tensorField& AU ) const; //- No copy construct @@ -138,14 +139,22 @@ public: ) const; //- Add resistance - virtual void correct(fvVectorMatrix& UEqn) const; + virtual void correct + ( + const volVectorField& U, + scalarField& Udiag, + vectorField& Usource, + bool compressible + ) const; //- Add resistance virtual void correct ( - fvVectorMatrix& UEqn, + const volVectorField& U, const volScalarField& rho, - const volScalarField& mu + const volScalarField& mu, + scalarField& Udiag, + vectorField& Usource ) const; //- Add resistance diff --git a/src/finiteVolume/cfdTools/general/porosityModel/porosityModel/porosityModel.C b/src/finiteVolume/cfdTools/general/porosityModel/porosityModel/porosityModel.C index e189d099e62..25583de4a44 100644 --- a/src/finiteVolume/cfdTools/general/porosityModel/porosityModel/porosityModel.C +++ b/src/finiteVolume/cfdTools/general/porosityModel/porosityModel/porosityModel.C @@ -28,6 +28,7 @@ License #include "porosityModel.H" #include "volFields.H" +#include "IOMRFZoneList.H" // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // @@ -97,7 +98,8 @@ Foam::porosityModel::porosityModel csysPtr_ ( coordinateSystem::New(mesh, coeffs_, coordinateSystem::typeName_()) - ) + ), + MRFRelativeVelocity_(true) { if (zoneName_.empty()) { @@ -146,12 +148,17 @@ Foam::porosityModel::porosityModel Info<< " local bounds: " << bb.span() << nl << endl; } + + if (dict.readIfPresent("MRFRelativeVelocity", MRFRelativeVelocity_)) + { + Info<< " MRFRelativeVelocity " << MRFRelativeVelocity_ << endl; + } } // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // -void Foam::porosityModel::transformModelData() +Foam::tmp<Foam::volVectorField> Foam::porosityModel::transformModelData() { if (!mesh_.upToDatePoints(*this)) { @@ -160,6 +167,32 @@ void Foam::porosityModel::transformModelData() // set model up-to-date wrt points mesh_.setUpToDatePoints(*this); } + + + auto* MRF = mesh_.findObject<IOMRFZoneList>("MRFProperties"); + if (MRF && MRFRelativeVelocity_) + { + auto tUref = + tmp<volVectorField>::New + ( + IOobject + ( + IOobject::scopedName(word("Uporous"), name_), + mesh_.time().timeName(), + mesh_, + IOobject::NO_READ, + IOobject::NO_WRITE + ), + mesh_, + dimensionedVector(dimVelocity, Zero) + ); + + MRF->makeRelative(tUref.ref()); + + return tUref; + } + + return nullptr; } @@ -170,13 +203,18 @@ Foam::tmp<Foam::vectorField> Foam::porosityModel::porosityModel::force const volScalarField& mu ) { - transformModelData(); + tmp<volVectorField> tUref = transformModelData(); - tmp<vectorField> tforce(new vectorField(U.size(), Zero)); + auto tforce = tmp<vectorField>::New(U.size(), Zero); if (!cellZoneIDs_.empty()) { this->calcForce(U, rho, mu, tforce.ref()); + + if (tUref.valid()) + { + this->calcForce(tUref(), rho, mu, tforce.ref()); + } } return tforce; @@ -190,16 +228,30 @@ void Foam::porosityModel::addResistance(fvVectorMatrix& UEqn) return; } - transformModelData(); - this->correct(UEqn); + bool compressible = UEqn.dimensions() == dimForce; + + tmp<volVectorField> tUref = transformModelData(); + + this->correct(UEqn.psi(), UEqn.diag(), UEqn.source(), compressible); + + if (tUref.valid()) + { + const auto& Uref = tUref(); + scalarField Udiag(Uref.size(), Zero); + vectorField Usource(Uref.size(), Zero); + + this->correct(Uref, Udiag, Usource, compressible); + + UEqn.source() -= Udiag*Uref - Usource; + } } void Foam::porosityModel::addResistance ( - fvVectorMatrix& UEqn, const volScalarField& rho, - const volScalarField& mu + const volScalarField& mu, + fvVectorMatrix& UEqn ) { if (cellZoneIDs_.empty()) @@ -207,8 +259,20 @@ void Foam::porosityModel::addResistance return; } - transformModelData(); - this->correct(UEqn, rho, mu); + tmp<volVectorField> tUref = transformModelData(); + + this->correct(UEqn.psi(), rho, mu, UEqn.diag(), UEqn.source()); + + if (tUref.valid()) + { + const auto& Uref = tUref(); + scalarField Udiag(Uref.size(), Zero); + vectorField Usource(Uref.size(), Zero); + + this->correct(Uref, rho, mu, Udiag, Usource); + + UEqn.source() -= Udiag*Uref - Usource; + } } @@ -224,7 +288,10 @@ void Foam::porosityModel::addResistance return; } - transformModelData(); + // Note: reference velocity contribution not handled! + // tmp<volVectorField> tUref = transformModelData(); + (void) transformModelData(); + this->correct(UEqn, AU); if (correctAUprocBC) diff --git a/src/finiteVolume/cfdTools/general/porosityModel/porosityModel/porosityModel.H b/src/finiteVolume/cfdTools/general/porosityModel/porosityModel/porosityModel.H index 229bbfc2c52..ce97f8b5b06 100644 --- a/src/finiteVolume/cfdTools/general/porosityModel/porosityModel/porosityModel.H +++ b/src/finiteVolume/cfdTools/general/porosityModel/porosityModel/porosityModel.H @@ -96,6 +96,9 @@ protected: //- Local coordinate system autoPtr<coordinateSystem> csysPtr_; + //- Flag to use relative velocity for MRF regions; default = yes + Switch MRFRelativeVelocity_; + // Protected Member Functions @@ -114,13 +117,21 @@ protected: vectorField& force ) const = 0; - virtual void correct(fvVectorMatrix& UEqn) const = 0; + virtual void correct + ( + const volVectorField& U, + scalarField& Udiag, + vectorField& Usource, + bool compressible + ) const = 0; virtual void correct ( - fvVectorMatrix& UEqn, + const volVectorField& U, const volScalarField& rho, - const volScalarField& mu + const volScalarField& mu, + scalarField& Udiag, + vectorField& Usource ) const = 0; virtual void correct @@ -232,7 +243,8 @@ public: const dictionary& dict() const; //- Transform the model data wrt mesh changes - virtual void transformModelData(); + // Returns the reference velocity (if applicable) + virtual tmp<volVectorField> transformModelData(); //- Return the force over the cell zone(s) virtual tmp<vectorField> force @@ -248,9 +260,9 @@ public: //- Add resistance virtual void addResistance ( - fvVectorMatrix& UEqn, const volScalarField& rho, - const volScalarField& mu + const volScalarField& mu, + fvVectorMatrix& UEqn ); //- Add resistance diff --git a/src/finiteVolume/cfdTools/general/porosityModel/porosityModel/porosityModelList.C b/src/finiteVolume/cfdTools/general/porosityModel/porosityModel/porosityModelList.C index 898dbae97e6..aac01ac808d 100644 --- a/src/finiteVolume/cfdTools/general/porosityModel/porosityModel/porosityModelList.C +++ b/src/finiteVolume/cfdTools/general/porosityModel/porosityModel/porosityModelList.C @@ -6,6 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2012-2014 OpenFOAM Foundation + Copyright (C) 2021 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -97,9 +98,8 @@ void Foam::porosityModelList::reset(const dictionary& dict) bool Foam::porosityModelList::read(const dictionary& dict) { bool allOk = true; - forAll(*this, i) + for (auto& pm : *this) { - porosityModel& pm = this->operator[](i); bool ok = pm.read(dict.subDict(pm.name())); allOk = (allOk && ok); } @@ -109,24 +109,21 @@ bool Foam::porosityModelList::read(const dictionary& dict) bool Foam::porosityModelList::writeData(Ostream& os) const { - forAll(*this, i) + for (const auto& pm : *this) { os << nl; - this->operator[](i).writeData(os); + pm.writeData(os); } return os.good(); } -void Foam::porosityModelList::addResistance -( - fvVectorMatrix& UEqn -) +void Foam::porosityModelList::addResistance(fvVectorMatrix& UEqn) { - forAll(*this, i) + for (auto& pm : *this) { - this->operator[](i).addResistance(UEqn); + pm.addResistance(UEqn); } } @@ -138,9 +135,9 @@ void Foam::porosityModelList::addResistance const volScalarField& mu ) { - forAll(*this, i) + for (auto& pm : *this) { - this->operator[](i).addResistance(UEqn, rho, mu); + pm.addResistance(rho, mu, UEqn); } } @@ -152,9 +149,9 @@ void Foam::porosityModelList::addResistance bool correctAUprocBC ) { - forAll(*this, i) + for (auto& pm : *this) { - this->operator[](i).addResistance(UEqn, AU, correctAUprocBC); + pm.addResistance(UEqn, AU, correctAUprocBC); } } diff --git a/src/finiteVolume/cfdTools/general/porosityModel/powerLaw/powerLaw.C b/src/finiteVolume/cfdTools/general/porosityModel/powerLaw/powerLaw.C index 417c292c649..5496ff54523 100644 --- a/src/finiteVolume/cfdTools/general/porosityModel/powerLaw/powerLaw.C +++ b/src/finiteVolume/cfdTools/general/porosityModel/powerLaw/powerLaw.C @@ -6,7 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2012-2017 OpenFOAM Foundation - Copyright (C) 2020 OpenCFD Ltd. + Copyright (C) 2020-2021 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -80,49 +80,50 @@ void Foam::porosityModels::powerLaw::calcForce scalarField Udiag(U.size(), Zero); const scalarField& V = mesh_.V(); - apply(Udiag, V, rho, U); + apply(V, rho, U, Udiag); - force = Udiag*U; + force += Udiag*U; } void Foam::porosityModels::powerLaw::correct ( - fvVectorMatrix& UEqn + const volVectorField& U, + scalarField& Udiag, + vectorField& Usource, + bool compressible ) const { - const volVectorField& U = UEqn.psi(); const scalarField& V = mesh_.V(); - scalarField& Udiag = UEqn.diag(); - if (UEqn.dimensions() == dimForce) + if (compressible) { const auto& rho = mesh_.lookupObject<volScalarField> ( IOobject::groupName(rhoName_, U.group()) ); - apply(Udiag, V, rho, U); + apply(V, rho, U, Udiag); } else { - apply(Udiag, V, geometricOneField(), U); + apply(V, geometricOneField(), U, Udiag); } } void Foam::porosityModels::powerLaw::correct ( - fvVectorMatrix& UEqn, + const volVectorField& U, const volScalarField& rho, - const volScalarField& mu + const volScalarField& mu, + scalarField& Udiag, + vectorField& Usource ) const { - const vectorField& U = UEqn.psi(); const scalarField& V = mesh_.V(); - scalarField& Udiag = UEqn.diag(); - apply(Udiag, V, rho, U); + apply(V, rho, U, Udiag); } @@ -141,11 +142,11 @@ void Foam::porosityModels::powerLaw::correct IOobject::groupName(rhoName_, U.group()) ); - apply(AU, rho, U); + apply(rho, U, AU); } else { - apply(AU, geometricOneField(), U); + apply(geometricOneField(), U, AU); } } diff --git a/src/finiteVolume/cfdTools/general/porosityModel/powerLaw/powerLaw.H b/src/finiteVolume/cfdTools/general/porosityModel/powerLaw/powerLaw.H index 0f3c483b744..097284fac2b 100644 --- a/src/finiteVolume/cfdTools/general/porosityModel/powerLaw/powerLaw.H +++ b/src/finiteVolume/cfdTools/general/porosityModel/powerLaw/powerLaw.H @@ -6,6 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2012-2017 OpenFOAM Foundation + Copyright (C) 2021 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -83,19 +84,19 @@ class powerLaw template<class RhoFieldType> void apply ( - scalarField& Udiag, const scalarField& V, const RhoFieldType& rho, - const vectorField& U + const vectorField& U, + scalarField& Udiag ) const; //- Apply resistance template<class RhoFieldType> void apply ( - tensorField& AU, const RhoFieldType& rho, - const vectorField& U + const vectorField& U, + tensorField& AU ) const; //- No copy construct @@ -139,14 +140,22 @@ public: ) const; //- Add resistance - virtual void correct(fvVectorMatrix& UEqn) const; + virtual void correct + ( + const volVectorField& U, + scalarField& Udiag, + vectorField& Usource, + bool compressible + ) const; //- Add resistance virtual void correct ( - fvVectorMatrix& UEqn, + const volVectorField& U, const volScalarField& rho, - const volScalarField& mu + const volScalarField& mu, + scalarField& Udiag, + vectorField& Usource ) const; //- Add resistance diff --git a/src/finiteVolume/cfdTools/general/porosityModel/powerLaw/powerLawTemplates.C b/src/finiteVolume/cfdTools/general/porosityModel/powerLaw/powerLawTemplates.C index 55701c70950..3706168198c 100644 --- a/src/finiteVolume/cfdTools/general/porosityModel/powerLaw/powerLawTemplates.C +++ b/src/finiteVolume/cfdTools/general/porosityModel/powerLaw/powerLawTemplates.C @@ -6,7 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2012-2016 OpenFOAM Foundation - Copyright (C) 2018 OpenCFD Ltd. + Copyright (C) 2018-2021 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -31,10 +31,10 @@ License template<class RhoFieldType> void Foam::porosityModels::powerLaw::apply ( - scalarField& Udiag, const scalarField& V, const RhoFieldType& rho, - const vectorField& U + const vectorField& U, + scalarField& Udiag ) const { const scalar C0 = C0_; @@ -56,9 +56,9 @@ void Foam::porosityModels::powerLaw::apply template<class RhoFieldType> void Foam::porosityModels::powerLaw::apply ( - tensorField& AU, const RhoFieldType& rho, - const vectorField& U + const vectorField& U, + tensorField& AU ) const { const scalar C0 = C0_; diff --git a/src/finiteVolume/cfdTools/general/porosityModel/solidification/solidification.C b/src/finiteVolume/cfdTools/general/porosityModel/solidification/solidification.C index 1281b2955fd..d1276d32d00 100644 --- a/src/finiteVolume/cfdTools/general/porosityModel/solidification/solidification.C +++ b/src/finiteVolume/cfdTools/general/porosityModel/solidification/solidification.C @@ -6,7 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2017 OpenFOAM Foundation - Copyright (C) 2020 OpenCFD Ltd. + Copyright (C) 2020-2021 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -62,12 +62,6 @@ Foam::porosityModels::solidification::solidification {} -// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * // - -Foam::porosityModels::solidification::~solidification() -{} - - // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // void Foam::porosityModels::solidification::calcTransformModelData() @@ -85,49 +79,50 @@ void Foam::porosityModels::solidification::calcForce scalarField Udiag(U.size(), Zero); const scalarField& V = mesh_.V(); - apply(Udiag, V, rho, U); + apply(V, rho, U, Udiag); - force = Udiag*U; + force += Udiag*U; } void Foam::porosityModels::solidification::correct ( - fvVectorMatrix& UEqn + const volVectorField& U, + scalarField& Udiag, + vectorField& Usource, + bool compressible ) const { - const volVectorField& U = UEqn.psi(); const scalarField& V = mesh_.V(); - scalarField& Udiag = UEqn.diag(); - if (UEqn.dimensions() == dimForce) + if (compressible) { const auto& rho = mesh_.lookupObject<volScalarField> ( IOobject::groupName(rhoName_, U.group()) ); - apply(Udiag, V, rho, U); + apply(V, rho, U, Udiag); } else { - apply(Udiag, V, geometricOneField(), U); + apply(V, geometricOneField(), U, Udiag); } } void Foam::porosityModels::solidification::correct ( - fvVectorMatrix& UEqn, + const volVectorField& U, const volScalarField& rho, - const volScalarField& mu + const volScalarField& mu, + scalarField& Udiag, + vectorField& Usource ) const { - const volVectorField& U = UEqn.psi(); const scalarField& V = mesh_.V(); - scalarField& Udiag = UEqn.diag(); - apply(Udiag, V, rho, U); + apply(V, rho, U, Udiag); } @@ -146,11 +141,11 @@ void Foam::porosityModels::solidification::correct IOobject::groupName(rhoName_, U.group()) ); - apply(AU, rho, U); + apply(rho, U, AU); } else { - apply(AU, geometricOneField(), U); + apply(geometricOneField(), U, AU); } } diff --git a/src/finiteVolume/cfdTools/general/porosityModel/solidification/solidification.H b/src/finiteVolume/cfdTools/general/porosityModel/solidification/solidification.H index 12a4d0de814..ec0eb904ebe 100644 --- a/src/finiteVolume/cfdTools/general/porosityModel/solidification/solidification.H +++ b/src/finiteVolume/cfdTools/general/porosityModel/solidification/solidification.H @@ -6,6 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2017 OpenFOAM Foundation + Copyright (C) 2021 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -120,40 +121,40 @@ class solidification template<class AlphaFieldType, class RhoFieldType> void apply ( - scalarField& Udiag, const scalarField& V, const AlphaFieldType& alpha, const RhoFieldType& rho, - const volVectorField& U + const volVectorField& U, + scalarField& Udiag ) const; //- Apply resistance template<class AlphaFieldType, class RhoFieldType> void apply ( - tensorField& AU, const AlphaFieldType& alpha, const RhoFieldType& rho, - const volVectorField& U + const volVectorField& U, + tensorField& AU ) const; //- Apply resistance template<class RhoFieldType> void apply ( - scalarField& Udiag, const scalarField& V, const RhoFieldType& rho, - const volVectorField& U + const volVectorField& U, + scalarField& Udiag ) const; //- Apply resistance template<class RhoFieldType> void apply ( - tensorField& AU, const RhoFieldType& rho, - const volVectorField& U + const volVectorField& U, + tensorField& AU ) const; //- No copy construct @@ -179,7 +180,7 @@ public: ); //- Destructor - virtual ~solidification(); + virtual ~solidification() = default; // Member Functions @@ -197,14 +198,22 @@ public: ) const; //- Add resistance - virtual void correct(fvVectorMatrix& UEqn) const; + virtual void correct + ( + const volVectorField& U, + scalarField& Udiag, + vectorField& Usource, + bool compressible + ) const; //- Add resistance virtual void correct ( - fvVectorMatrix& UEqn, + const volVectorField& U, const volScalarField& rho, - const volScalarField& mu + const volScalarField& mu, + scalarField& Udiag, + vectorField& Usource ) const; //- Add resistance diff --git a/src/finiteVolume/cfdTools/general/porosityModel/solidification/solidificationTemplates.C b/src/finiteVolume/cfdTools/general/porosityModel/solidification/solidificationTemplates.C index 08cf1eb3b59..ca0f43e595a 100644 --- a/src/finiteVolume/cfdTools/general/porosityModel/solidification/solidificationTemplates.C +++ b/src/finiteVolume/cfdTools/general/porosityModel/solidification/solidificationTemplates.C @@ -6,6 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2017 OpenFOAM Foundation + Copyright (C) 2021 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -33,11 +34,11 @@ License template<class AlphaFieldType, class RhoFieldType> void Foam::porosityModels::solidification::apply ( - scalarField& Udiag, const scalarField& V, const AlphaFieldType& alpha, const RhoFieldType& rho, - const volVectorField& U + const volVectorField& U, + scalarField& Udiag ) const { const auto& T = mesh_.lookupObject<volScalarField> @@ -61,10 +62,10 @@ void Foam::porosityModels::solidification::apply template<class AlphaFieldType, class RhoFieldType> void Foam::porosityModels::solidification::apply ( - tensorField& AU, const AlphaFieldType& alpha, const RhoFieldType& rho, - const volVectorField& U + const volVectorField& U, + tensorField& AU ) const { const auto& T = mesh_.lookupObject<volScalarField> @@ -88,15 +89,15 @@ void Foam::porosityModels::solidification::apply template<class RhoFieldType> void Foam::porosityModels::solidification::apply ( - scalarField& Udiag, const scalarField& V, const RhoFieldType& rho, - const volVectorField& U + const volVectorField& U, + scalarField& Udiag ) const { if (alphaName_ == "none") { - return apply(Udiag, V, geometricOneField(), rho, U); + return apply(V, geometricOneField(), rho, U, Udiag); } else { @@ -105,7 +106,7 @@ void Foam::porosityModels::solidification::apply IOobject::groupName(alphaName_, U.group()) ); - return apply(Udiag, V, alpha, rho, U); + return apply(V, alpha, rho, U, Udiag); } } @@ -113,14 +114,14 @@ void Foam::porosityModels::solidification::apply template<class RhoFieldType> void Foam::porosityModels::solidification::apply ( - tensorField& AU, const RhoFieldType& rho, - const volVectorField& U + const volVectorField& U, + tensorField& AU ) const { if (alphaName_ == "none") { - return apply(AU, geometricOneField(), rho, U); + return apply(geometricOneField(), rho, U, AU); } else { @@ -129,7 +130,7 @@ void Foam::porosityModels::solidification::apply IOobject::groupName(alphaName_, U.group()) ); - return apply(AU, alpha, rho, U); + return apply(alpha, rho, U, AU); } } diff --git a/src/fvOptions/sources/interRegion/interRegionExplicitPorositySource/interRegionExplicitPorositySource.C b/src/fvOptions/sources/interRegion/interRegionExplicitPorositySource/interRegionExplicitPorositySource.C index 4e00e9b2251..ab49de929bc 100644 --- a/src/fvOptions/sources/interRegion/interRegionExplicitPorositySource/interRegionExplicitPorositySource.C +++ b/src/fvOptions/sources/interRegion/interRegionExplicitPorositySource/interRegionExplicitPorositySource.C @@ -6,7 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2012-2016 OpenFOAM Foundation - Copyright (C) 2020 OpenCFD Ltd. + Copyright (C) 2020-2021 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -270,7 +270,7 @@ void Foam::fv::interRegionExplicitPorositySource::addSup muNbr.primitiveFieldRef() ); - porosityPtr_->addResistance(nbrEqn, rhoNbr, muNbr); + porosityPtr_->addResistance(rhoNbr, muNbr, nbrEqn); // Convert source from neighbour to local region fvMatrix<vector> porosityEqn(U, eqn.dimensions()); -- GitLab