From e5e1440020724508ebb3d1df6ef74aa01bb82646 Mon Sep 17 00:00:00 2001 From: Mark Olesen Date: Mon, 3 Apr 2023 14:23:14 +0200 Subject: [PATCH] ENH: align use of extrapolatedCalculated in faMatrix with fvMatrix STYLE: prefer GeometricField::New factory methods --- src/finiteArea/faMatrices/faMatrix/faMatrix.C | 61 +++------- .../faScalarMatrix/faScalarMatrix.C | 30 ++--- .../faScalarMatrix/faScalarMatrix.H | 5 +- .../fvMatrices/fvMatrix/fvMatrix.C | 113 +++++------------- .../fvScalarMatrix/fvScalarMatrix.C | 46 ++----- .../fvScalarMatrix/fvScalarMatrix.H | 4 +- 6 files changed, 75 insertions(+), 184 deletions(-) diff --git a/src/finiteArea/faMatrices/faMatrix/faMatrix.C b/src/finiteArea/faMatrices/faMatrix/faMatrix.C index b2d8169798..d4b8576638 100644 --- a/src/finiteArea/faMatrices/faMatrix/faMatrix.C +++ b/src/finiteArea/faMatrices/faMatrix/faMatrix.C @@ -29,6 +29,7 @@ License #include "areaFields.H" #include "edgeFields.H" #include "calculatedFaPatchFields.H" +#include "extrapolatedCalculatedFaPatchFields.H" #include "zeroGradientFaPatchFields.H" #include "IndirectList.H" #include "UniformList.H" @@ -610,7 +611,7 @@ void Foam::faMatrix::relax() template Foam::tmp Foam::faMatrix::D() const { - tmp tdiag(new scalarField(diag())); + auto tdiag = tmp::New(diag()); addCmptAvBoundaryDiag(tdiag.ref()); return tdiag; } @@ -619,20 +620,12 @@ Foam::tmp Foam::faMatrix::D() const template Foam::tmp Foam::faMatrix::A() const { - tmp tAphi + auto tAphi = areaScalarField::New ( - new areaScalarField - ( - IOobject - ( - "A("+psi_.name()+')', - psi_.instance(), - psi_.db() - ), - psi_.mesh(), - dimensions_/psi_.dimensions()/dimArea, - zeroGradientFaPatchScalarField::typeName - ) + "A(" + psi_.name() + ')', + psi_.mesh(), + dimensions_/psi_.dimensions()/dimArea, + faPatchFieldBase::extrapolatedCalculatedType() ); tAphi.ref().primitiveFieldRef() = D()/psi_.mesh().S(); @@ -646,22 +639,14 @@ template Foam::tmp> Foam::faMatrix::H() const { - tmp> tHphi + auto tHphi = GeometricField::New ( - new GeometricField - ( - IOobject - ( - "H("+psi_.name()+')', - psi_.instance(), - psi_.db() - ), - psi_.mesh(), - dimensions_/dimArea, - zeroGradientFaPatchScalarField::typeName - ) + "H(" + psi_.name() + ')', + psi_.mesh(), + dimensions_/dimArea, + faPatchFieldBase::extrapolatedCalculatedType() ); - GeometricField& Hphi = tHphi.ref(); + auto& Hphi = tHphi.ref(); // Loop over field components for (direction cmpt=0; cmpt::flux() const << abort(FatalError); } - // construct GeometricField - tmp> tfieldFlux + auto tfieldFlux = GeometricField::New ( - new GeometricField - ( - IOobject - ( - "flux("+psi_.name()+')', - psi_.instance(), - psi_.db() - ), - psi_.mesh(), - dimensions() - ) + "flux(" + psi_.name() + ')', + psi_.mesh(), + dimensions() ); - GeometricField& fieldFlux = - tfieldFlux.ref(); + auto& fieldFlux = tfieldFlux.ref(); for (direction cmpt=0; cmpt::nComponents; ++cmpt) { diff --git a/src/finiteArea/faMatrices/faScalarMatrix/faScalarMatrix.C b/src/finiteArea/faMatrices/faScalarMatrix/faScalarMatrix.C index 4fd605449d..d8a5ff4c7b 100644 --- a/src/finiteArea/faMatrices/faScalarMatrix/faScalarMatrix.C +++ b/src/finiteArea/faMatrices/faScalarMatrix/faScalarMatrix.C @@ -6,7 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2016-2017 Wikki Ltd - Copyright (C) 2019-2021 OpenCFD Ltd. + Copyright (C) 2019-2023 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -24,13 +24,11 @@ License You should have received a copy of the GNU General Public License along with OpenFOAM. If not, see . -Description - Finite-Area scalar matrix member functions and operators - \*---------------------------------------------------------------------------*/ #include "faScalarMatrix.H" -#include "zeroGradientFaPatchFields.H" +#include "extrapolatedCalculatedFaPatchFields.H" +#include "profiling.H" #include "PrecisionAdaptor.H" // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // @@ -136,24 +134,14 @@ Foam::tmp Foam::faMatrix::residual() const template<> Foam::tmp Foam::faMatrix::H() const { - tmp tHphi + auto tHphi = areaScalarField::New ( - new areaScalarField - ( - IOobject - ( - "H("+psi_.name()+')', - psi_.instance(), - psi_.db(), - IOobject::NO_READ, - IOobject::NO_WRITE - ), - psi_.mesh(), - dimensions_/dimArea, - zeroGradientFaPatchScalarField::typeName - ) + "H(" + psi_.name() + ')', + psi_.mesh(), + dimensions_/dimArea, + faPatchFieldBase::extrapolatedCalculatedType() ); - areaScalarField& Hphi = tHphi.ref(); + auto& Hphi = tHphi.ref(); Hphi.primitiveFieldRef() = (lduMatrix::H(psi_.primitiveField()) + source_); addBoundarySource(Hphi.primitiveFieldRef()); diff --git a/src/finiteArea/faMatrices/faScalarMatrix/faScalarMatrix.H b/src/finiteArea/faMatrices/faScalarMatrix/faScalarMatrix.H index e8652d8475..272cbd1b25 100644 --- a/src/finiteArea/faMatrices/faScalarMatrix/faScalarMatrix.H +++ b/src/finiteArea/faMatrices/faScalarMatrix/faScalarMatrix.H @@ -38,10 +38,11 @@ Author \*---------------------------------------------------------------------------*/ -#ifndef faScalarMatrix_H -#define faScalarMatrix_H +#ifndef Foam_faScalarMatrix_H +#define Foam_faScalarMatrix_H #include "faMatrix.H" +#include "faMatricesFwd.H" // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // diff --git a/src/finiteVolume/fvMatrices/fvMatrix/fvMatrix.C b/src/finiteVolume/fvMatrices/fvMatrix/fvMatrix.C index 38dd6520a8..14bacc7f0e 100644 --- a/src/finiteVolume/fvMatrices/fvMatrix/fvMatrix.C +++ b/src/finiteVolume/fvMatrices/fvMatrix/fvMatrix.C @@ -6,7 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2011-2017 OpenFOAM Foundation - Copyright (C) 2016-2022 OpenCFD Ltd. + Copyright (C) 2016-2023 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -1285,7 +1285,7 @@ void Foam::fvMatrix::boundaryManipulate template Foam::tmp Foam::fvMatrix::D() const { - tmp tdiag(new scalarField(diag())); + auto tdiag = tmp::New(diag()); addCmptAvBoundaryDiag(tdiag.ref()); return tdiag; } @@ -1318,22 +1318,12 @@ Foam::tmp> Foam::fvMatrix::DD() const template Foam::tmp Foam::fvMatrix::A() const { - tmp tAphi + auto tAphi = volScalarField::New ( - new volScalarField - ( - IOobject - ( - "A("+psi_.name()+')', - psi_.instance(), - psi_.mesh(), - IOobject::NO_READ, - IOobject::NO_WRITE - ), - psi_.mesh(), - dimensions_/psi_.dimensions()/dimVol, - extrapolatedCalculatedFvPatchScalarField::typeName - ) + "A(" + psi_.name() + ')', + psi_.mesh(), + dimensions_/psi_.dimensions()/dimVol, + fvPatchFieldBase::extrapolatedCalculatedType() ); tAphi.ref().primitiveFieldRef() = D()/psi_.mesh().V(); @@ -1347,24 +1337,14 @@ template Foam::tmp> Foam::fvMatrix::H() const { - tmp> tHphi + auto tHphi = GeometricField::New ( - new GeometricField - ( - IOobject - ( - "H("+psi_.name()+')', - psi_.instance(), - psi_.mesh(), - IOobject::NO_READ, - IOobject::NO_WRITE - ), - psi_.mesh(), - dimensions_/dimVol, - extrapolatedCalculatedFvPatchScalarField::typeName - ) + "H(" + psi_.name() + ')', + psi_.mesh(), + dimensions_/dimVol, + fvPatchFieldBase::extrapolatedCalculatedType() ); - GeometricField& Hphi = tHphi.ref(); + auto& Hphi = tHphi.ref(); // Loop over field components for (direction cmpt=0; cmpt::H() const template Foam::tmp Foam::fvMatrix::H1() const { - tmp tH1 + auto tH1 = volScalarField::New ( - new volScalarField - ( - IOobject - ( - "H(1)", - psi_.instance(), - psi_.mesh(), - IOobject::NO_READ, - IOobject::NO_WRITE - ), - psi_.mesh(), - dimensions_/(dimVol*psi_.dimensions()), - extrapolatedCalculatedFvPatchScalarField::typeName - ) + "H(1)", + psi_.mesh(), + dimensions_/(dimVol*psi_.dimensions()), + fvPatchFieldBase::extrapolatedCalculatedType() ); - volScalarField& H1_ = tH1.ref(); + auto& H1_ = tH1.ref(); H1_.primitiveFieldRef() = lduMatrix::H1(); @@ -1452,7 +1422,6 @@ Foam::tmp Foam::fvMatrix::H1() const } - template Foam::tmp> Foam::fvMatrix:: @@ -1475,25 +1444,13 @@ flux() const << abort(FatalError); } - // construct GeometricField - tmp> tfieldFlux + auto tfieldFlux = GeometricField::New ( - new GeometricField - ( - IOobject - ( - "flux("+psi_.name()+')', - psi_.instance(), - psi_.mesh(), - IOobject::NO_READ, - IOobject::NO_WRITE - ), - psi_.mesh(), - dimensions() - ) + "flux(" + psi_.name() + ')', + psi_.mesh(), + dimensions() ); - GeometricField& fieldFlux = - tfieldFlux.ref(); + auto& fieldFlux = tfieldFlux.ref(); fieldFlux.setOriented(); @@ -2898,24 +2855,14 @@ Foam::operator& const DimensionedField& psi ) { - tmp> tMphi + auto tMphi = GeometricField::New ( - new GeometricField - ( - IOobject - ( - "M&" + psi.name(), - psi.instance(), - psi.mesh(), - IOobject::NO_READ, - IOobject::NO_WRITE - ), - psi.mesh(), - M.dimensions()/dimVol, - extrapolatedCalculatedFvPatchScalarField::typeName - ) + "M&" + psi.name(), + psi.mesh(), + M.dimensions()/dimVol, + fvPatchFieldBase::extrapolatedCalculatedType() ); - GeometricField& Mphi = tMphi.ref(); + auto& Mphi = tMphi.ref(); // Loop over field components if (M.hasDiag()) diff --git a/src/finiteVolume/fvMatrices/fvScalarMatrix/fvScalarMatrix.C b/src/finiteVolume/fvMatrices/fvScalarMatrix/fvScalarMatrix.C index 8d8fb32705..591450e206 100644 --- a/src/finiteVolume/fvMatrices/fvScalarMatrix/fvScalarMatrix.C +++ b/src/finiteVolume/fvMatrices/fvScalarMatrix/fvScalarMatrix.C @@ -6,7 +6,7 @@ \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2011-2017 OpenFOAM Foundation - Copyright (C) 2016-2021 OpenCFD Ltd. + Copyright (C) 2016-2023 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. @@ -356,24 +356,14 @@ Foam::tmp Foam::fvMatrix::residual() const template<> Foam::tmp Foam::fvMatrix::H() const { - tmp tHphi + auto tHphi = volScalarField::New ( - new volScalarField - ( - IOobject - ( - "H("+psi_.name()+')', - psi_.instance(), - psi_.mesh(), - IOobject::NO_READ, - IOobject::NO_WRITE - ), - psi_.mesh(), - dimensions_/dimVol, - extrapolatedCalculatedFvPatchScalarField::typeName - ) + "H(" + psi_.name() + ')', + psi_.mesh(), + dimensions_/dimVol, + fvPatchFieldBase::extrapolatedCalculatedType() ); - volScalarField& Hphi = tHphi.ref(); + auto& Hphi = tHphi.ref(); Hphi.primitiveFieldRef() = (lduMatrix::H(psi_.primitiveField()) + source_); addBoundarySource(Hphi.primitiveFieldRef()); @@ -388,24 +378,14 @@ Foam::tmp Foam::fvMatrix::H() const template<> Foam::tmp Foam::fvMatrix::H1() const { - tmp tH1 + auto tH1 = volScalarField::New ( - new volScalarField - ( - IOobject - ( - "H(1)", - psi_.instance(), - psi_.mesh(), - IOobject::NO_READ, - IOobject::NO_WRITE - ), - psi_.mesh(), - dimensions_/(dimVol*psi_.dimensions()), - extrapolatedCalculatedFvPatchScalarField::typeName - ) + "H(1)", + psi_.mesh(), + dimensions_/(dimVol*psi_.dimensions()), + fvPatchFieldBase::extrapolatedCalculatedType() ); - volScalarField& H1_ = tH1.ref(); + auto& H1_ = tH1.ref(); H1_.primitiveFieldRef() = lduMatrix::H1(); //addBoundarySource(Hphi.primitiveField()); diff --git a/src/finiteVolume/fvMatrices/fvScalarMatrix/fvScalarMatrix.H b/src/finiteVolume/fvMatrices/fvScalarMatrix/fvScalarMatrix.H index 82545cec64..075dce5df4 100644 --- a/src/finiteVolume/fvMatrices/fvScalarMatrix/fvScalarMatrix.H +++ b/src/finiteVolume/fvMatrices/fvScalarMatrix/fvScalarMatrix.H @@ -34,8 +34,8 @@ SourceFiles \*---------------------------------------------------------------------------*/ -#ifndef fvScalarMatrix_H -#define fvScalarMatrix_H +#ifndef Foam_fvScalarMatrix_H +#define Foam_fvScalarMatrix_H #include "fvMatrix.H" #include "fvMatricesFwd.H"