From 1d9f39761ea4f833a70c3bebc9f2fa2ac0dc0d74 Mon Sep 17 00:00:00 2001 From: Will Bainbridge Date: Wed, 16 Sep 2020 16:15:35 +0100 Subject: [PATCH] interpolation: Added field evaluation The field evaluations have been optimised using an additional fieldInterpolation base class so that a virtual call happens only once per field. This is the same pattern as that used to optimise Function1. --- .../interpolation/interpolation.C | 123 ++++++++++++++++-- .../interpolation/interpolation.H | 81 +++++++++--- .../interpolation/interpolationNew.C | 66 ---------- .../interpolationCell/interpolationCell.C | 4 +- .../interpolationCell/interpolationCell.H | 2 +- .../interpolationCellPatchConstrained.C | 4 +- .../interpolationCellPatchConstrained.H | 2 +- .../interpolationCellPoint.C | 6 +- .../interpolationCellPoint.H | 4 +- .../interpolationCellPointI.H | 6 +- .../findCellPointFaceTet.H | 10 +- .../findCellPointFaceTriangle.H | 8 +- .../interpolationCellPointFace.C | 30 ++--- .../interpolationCellPointFace.H | 5 +- .../interpolationPointMVC.C | 4 +- .../interpolationPointMVC.H | 2 +- .../interpolationPointMVCI.H | 4 +- 17 files changed, 225 insertions(+), 136 deletions(-) delete mode 100644 src/finiteVolume/interpolation/interpolation/interpolation/interpolationNew.C diff --git a/src/finiteVolume/interpolation/interpolation/interpolation/interpolation.C b/src/finiteVolume/interpolation/interpolation/interpolation/interpolation.C index e6ee5791b3..d9ac7742b0 100644 --- a/src/finiteVolume/interpolation/interpolation/interpolation/interpolation.C +++ b/src/finiteVolume/interpolation/interpolation/interpolation/interpolation.C @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | Website: https://openfoam.org - \\ / A nd | Copyright (C) 2011-2018 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2011-2020 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -37,16 +37,123 @@ Foam::interpolation::interpolation ) : psi_(psi), - pMesh_(psi.mesh()), - pMeshPoints_(pMesh_.points()), - pMeshFaces_(pMesh_.faces()), - pMeshFaceCentres_(pMesh_.faceCentres()), - pMeshFaceAreas_(pMesh_.faceAreas()) + mesh_(psi.mesh()) {} -// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // +// * * * * * * * * * * * * * * * * Selectors * * * * * * * * * * * * * * * * // + +template +Foam::autoPtr> Foam::interpolation::New +( + const word& interpolationType, + const GeometricField& psi +) +{ + typename dictionaryConstructorTable::iterator cstrIter = + dictionaryConstructorTablePtr_->find(interpolationType); + + if (cstrIter == dictionaryConstructorTablePtr_->end()) + { + FatalErrorInFunction + << "Unknown interpolation type " << interpolationType + << " for field " << psi.name() << nl << nl + << "Valid interpolation types : " << endl + << dictionaryConstructorTablePtr_->sortedToc() + << exit(FatalError); + } + + return autoPtr>(cstrIter()(psi)); +} + + +template +Foam::autoPtr> Foam::interpolation::New +( + const dictionary& interpolationSchemes, + const GeometricField& psi +) +{ + return New(word(interpolationSchemes.lookup(psi.name())), psi); +} + + +// * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * // + +template +Type Foam::interpolation::interpolate +( + const barycentric& coordinates, + const tetIndices& tetIs, + const label facei +) const +{ + return + interpolate + ( + tetIs.tet(mesh_).barycentricToPoint(coordinates), + tetIs.cell(), + facei + ); +} + + +template +Foam::tmp> +Foam::fieldInterpolation::interpolate +( + const vectorField& position, + const labelField& celli, + const labelField& facei +) const +{ + tmp> tField(new Field(position.size())); + + Field& field = tField.ref(); + + forAll(field, i) + { + field[i] = + static_cast(*this).interpolate + ( + position[i], + celli[i], + isNull(facei) ? -1 : facei[i] + ); + } + + return tField; +} + + +template +Foam::tmp> +Foam::fieldInterpolation::interpolate +( + const Field& coordinates, + const labelField& celli, + const labelField& tetFacei, + const labelField& tetPti, + const labelField& facei +) const +{ + tmp> tField(new Field(coordinates.size())); + + Field& field = tField.ref(); + + forAll(field, i) + { + field[i] = + static_cast(*this).interpolate + ( + coordinates[i], + tetIndices(celli[i], tetFacei[i], tetPti[i]), + isNull(facei) ? -1 : facei[i] + ); + } + + return tField; +} -#include "interpolationNew.C" // ************************************************************************* // diff --git a/src/finiteVolume/interpolation/interpolation/interpolation/interpolation.H b/src/finiteVolume/interpolation/interpolation/interpolation/interpolation.H index 2cc4f2371d..56557a4f47 100644 --- a/src/finiteVolume/interpolation/interpolation/interpolation/interpolation.H +++ b/src/finiteVolume/interpolation/interpolation/interpolation/interpolation.H @@ -54,18 +54,15 @@ class polyMesh; template class interpolation { - protected: // Protected data + //- The vol field to interpolate const GeometricField& psi_; - const polyMesh& pMesh_; - const vectorField& pMeshPoints_; - const faceList& pMeshFaces_; - const vectorField& pMeshFaceCentres_; - const vectorField& pMeshFaceAreas_; + //- Reference to the mesh + const polyMesh& mesh_; public: @@ -135,25 +132,73 @@ public: const label facei = -1 ) const = 0; + //- As above, but for a field + virtual tmp> interpolate + ( + const vectorField& position, + const labelField& celli, + const labelField& facei = NullObjectRef() + ) const = 0; + //- Interpolate field to the given coordinates in the tetrahedron - // defined by the given indices. Calls interpolate function - // above here execpt where overridden by derived + // defined by the given indices. Calls interpolate function + // above here except where overridden by derived // interpolation types. virtual Type interpolate ( const barycentric& coordinates, const tetIndices& tetIs, const label facei = -1 - ) const - { - return - interpolate - ( - tetIs.tet(pMesh_).barycentricToPoint(coordinates), - tetIs.cell(), - facei - ); - } + ) const; + + //- As above, but for a field + virtual tmp> interpolate + ( + const Field& coordinates, + const labelField& celli, + const labelField& tetFacei, + const labelField& tetPti, + const labelField& facei = NullObjectRef() + ) const = 0; +}; + + +/*---------------------------------------------------------------------------*\ + Class fieldInterpolation Declaration +\*---------------------------------------------------------------------------*/ + +template +class fieldInterpolation +: + public interpolation +{ +public: + + // Constructors + + //- Inherit base class constructors + using interpolation::interpolation; + + + // Member Functions + + //- Interpolate field to the given points in the given cells + virtual tmp> interpolate + ( + const vectorField& position, + const labelField& celli, + const labelField& facei = NullObjectRef() + ) const; + + //- Interpolate field to the given coordinates in the given tetrahedra + virtual tmp> interpolate + ( + const Field& coordinates, + const labelField& celli, + const labelField& tetFacei, + const labelField& tetPti, + const labelField& facei = NullObjectRef() + ) const; }; diff --git a/src/finiteVolume/interpolation/interpolation/interpolation/interpolationNew.C b/src/finiteVolume/interpolation/interpolation/interpolation/interpolationNew.C deleted file mode 100644 index e59527eb35..0000000000 --- a/src/finiteVolume/interpolation/interpolation/interpolation/interpolationNew.C +++ /dev/null @@ -1,66 +0,0 @@ -/*---------------------------------------------------------------------------*\ - ========= | - \\ / F ield | OpenFOAM: The Open Source CFD Toolbox - \\ / O peration | Website: https://openfoam.org - \\ / A nd | Copyright (C) 2011-2018 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 . - -\*---------------------------------------------------------------------------*/ - -#include "interpolation.H" -#include "volFields.H" - -// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // - -template -Foam::autoPtr> Foam::interpolation::New -( - const word& interpolationType, - const GeometricField& psi -) -{ - typename dictionaryConstructorTable::iterator cstrIter = - dictionaryConstructorTablePtr_->find(interpolationType); - - if (cstrIter == dictionaryConstructorTablePtr_->end()) - { - FatalErrorInFunction - << "Unknown interpolation type " << interpolationType - << " for field " << psi.name() << nl << nl - << "Valid interpolation types : " << endl - << dictionaryConstructorTablePtr_->sortedToc() - << exit(FatalError); - } - - return autoPtr>(cstrIter()(psi)); -} - - -template -Foam::autoPtr> Foam::interpolation::New -( - const dictionary& interpolationSchemes, - const GeometricField& psi -) -{ - return New(word(interpolationSchemes.lookup(psi.name())), psi); -} - - -// ************************************************************************* // diff --git a/src/finiteVolume/interpolation/interpolation/interpolationCell/interpolationCell.C b/src/finiteVolume/interpolation/interpolation/interpolationCell/interpolationCell.C index ba4a59f407..4704f6f642 100644 --- a/src/finiteVolume/interpolation/interpolation/interpolationCell/interpolationCell.C +++ b/src/finiteVolume/interpolation/interpolation/interpolationCell/interpolationCell.C @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | Website: https://openfoam.org - \\ / A nd | Copyright (C) 2011-2018 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2011-2020 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -34,7 +34,7 @@ Foam::interpolationCell::interpolationCell const GeometricField& psi ) : - interpolation(psi) + fieldInterpolation>(psi) {} diff --git a/src/finiteVolume/interpolation/interpolation/interpolationCell/interpolationCell.H b/src/finiteVolume/interpolation/interpolation/interpolationCell/interpolationCell.H index 75aba212b8..3c751721d9 100644 --- a/src/finiteVolume/interpolation/interpolation/interpolationCell/interpolationCell.H +++ b/src/finiteVolume/interpolation/interpolation/interpolationCell/interpolationCell.H @@ -48,7 +48,7 @@ class fvMesh; template class interpolationCell : - public interpolation + public fieldInterpolation> { public: diff --git a/src/finiteVolume/interpolation/interpolation/interpolationCellPatchConstrained/interpolationCellPatchConstrained.C b/src/finiteVolume/interpolation/interpolation/interpolationCellPatchConstrained/interpolationCellPatchConstrained.C index 1c7612f798..9bdb31a513 100644 --- a/src/finiteVolume/interpolation/interpolation/interpolationCellPatchConstrained/interpolationCellPatchConstrained.C +++ b/src/finiteVolume/interpolation/interpolation/interpolationCellPatchConstrained/interpolationCellPatchConstrained.C @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | Website: https://openfoam.org - \\ / A nd | Copyright (C) 2011-2018 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2011-2020 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -34,7 +34,7 @@ Foam::interpolationCellPatchConstrained::interpolationCellPatchConstrained const GeometricField& psi ) : - interpolation(psi) + fieldInterpolation>(psi) {} diff --git a/src/finiteVolume/interpolation/interpolation/interpolationCellPatchConstrained/interpolationCellPatchConstrained.H b/src/finiteVolume/interpolation/interpolation/interpolationCellPatchConstrained/interpolationCellPatchConstrained.H index 61116c074a..301e94d160 100644 --- a/src/finiteVolume/interpolation/interpolation/interpolationCellPatchConstrained/interpolationCellPatchConstrained.H +++ b/src/finiteVolume/interpolation/interpolation/interpolationCellPatchConstrained/interpolationCellPatchConstrained.H @@ -50,7 +50,7 @@ class fvMesh; template class interpolationCellPatchConstrained : - public interpolation + public fieldInterpolation> { public: diff --git a/src/finiteVolume/interpolation/interpolation/interpolationCellPoint/interpolationCellPoint.C b/src/finiteVolume/interpolation/interpolation/interpolationCellPoint/interpolationCellPoint.C index 233400987b..1b4d15e791 100644 --- a/src/finiteVolume/interpolation/interpolation/interpolationCellPoint/interpolationCellPoint.C +++ b/src/finiteVolume/interpolation/interpolation/interpolationCellPoint/interpolationCellPoint.C @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | Website: https://openfoam.org - \\ / A nd | Copyright (C) 2011-2018 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2011-2020 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -34,7 +34,7 @@ Foam::interpolationCellPoint::interpolationCellPoint const GeometricField& psi ) : - interpolation(psi), + fieldInterpolation>(psi), psip_ ( volPointInterpolation::New(psi.mesh()).interpolate @@ -57,7 +57,7 @@ Foam::interpolationCellPoint::interpolationCellPoint tmp> psip ) : - interpolation(psi), + fieldInterpolation>(psi), psip_(psip) { // Uses cellPointWeight to do interpolation which needs tet decomposition diff --git a/src/finiteVolume/interpolation/interpolation/interpolationCellPoint/interpolationCellPoint.H b/src/finiteVolume/interpolation/interpolation/interpolationCellPoint/interpolationCellPoint.H index 1937f799e0..dad6e858b1 100644 --- a/src/finiteVolume/interpolation/interpolation/interpolationCellPoint/interpolationCellPoint.H +++ b/src/finiteVolume/interpolation/interpolation/interpolationCellPoint/interpolationCellPoint.H @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | Website: https://openfoam.org - \\ / A nd | Copyright (C) 2011-2018 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2011-2020 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -48,7 +48,7 @@ namespace Foam template class interpolationCellPoint : - public interpolation + public fieldInterpolation> { protected: diff --git a/src/finiteVolume/interpolation/interpolation/interpolationCellPoint/interpolationCellPointI.H b/src/finiteVolume/interpolation/interpolation/interpolationCellPoint/interpolationCellPointI.H index d16e191909..b23235e6d6 100644 --- a/src/finiteVolume/interpolation/interpolation/interpolationCellPoint/interpolationCellPointI.H +++ b/src/finiteVolume/interpolation/interpolation/interpolationCellPoint/interpolationCellPointI.H @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | Website: https://openfoam.org - \\ / A nd | Copyright (C) 2011-2018 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2011-2020 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -51,7 +51,7 @@ inline Type Foam::interpolationCellPoint::interpolate const label facei ) const { - return interpolate(cellPointWeight(this->pMesh_, position, celli, facei)); + return interpolate(cellPointWeight(this->mesh_, position, celli, facei)); } @@ -79,7 +79,7 @@ inline Type Foam::interpolationCellPoint::interpolate } } - const triFace triIs = tetIs.faceTriIs(this->pMesh_); + const triFace triIs = tetIs.faceTriIs(this->mesh_); return this->psi_[tetIs.cell()]*coordinates[0] diff --git a/src/finiteVolume/interpolation/interpolation/interpolationCellPointFace/findCellPointFaceTet.H b/src/finiteVolume/interpolation/interpolation/interpolationCellPointFace/findCellPointFaceTet.H index 722f837b44..cc9abe2ba8 100644 --- a/src/finiteVolume/interpolation/interpolation/interpolationCellPointFace/findCellPointFaceTet.H +++ b/src/finiteVolume/interpolation/interpolation/interpolationCellPointFace/findCellPointFaceTet.H @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | Website: https://openfoam.org - \\ / A nd | Copyright (C) 2011-2018 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2011-2020 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -53,8 +53,8 @@ bool interpolationCellPointFace::findTet { bool foundTet = false; - const labelList& thisFacePoints = this->pMeshFaces_[nFace]; - tetPoints[2] = this->pMeshFaceCentres_[nFace]; + const labelList& thisFacePoints = this->mesh_.faces()[nFace]; + tetPoints[2] = this->mesh_.faceCentres()[nFace]; label pointi = 0; @@ -65,8 +65,8 @@ bool interpolationCellPointFace::findTet tetPointLabels[0] = thisFacePoints[pointi]; tetPointLabels[1] = thisFacePoints[nextPointLabel]; - tetPoints[0] = this->pMeshPoints_[tetPointLabels[0]]; - tetPoints[1] = this->pMeshPoints_[tetPointLabels[1]]; + tetPoints[0] = this->mesh_.points()[tetPointLabels[0]]; + tetPoints[1] = this->mesh_.points()[tetPointLabels[1]]; bool inside = true; scalar dist = 0.0; diff --git a/src/finiteVolume/interpolation/interpolation/interpolationCellPointFace/findCellPointFaceTriangle.H b/src/finiteVolume/interpolation/interpolation/interpolationCellPointFace/findCellPointFaceTriangle.H index 10dffdcaf7..efcda512f6 100644 --- a/src/finiteVolume/interpolation/interpolation/interpolationCellPointFace/findCellPointFaceTriangle.H +++ b/src/finiteVolume/interpolation/interpolation/interpolationCellPointFace/findCellPointFaceTriangle.H @@ -21,8 +21,8 @@ bool interpolationCellPointFace::findTriangle { bool foundTriangle = false; vector tetPoints[3]; - const labelList& facePoints = this->pMeshFaces_[nFace]; - tetPoints[2] = this->pMeshFaceCentres_[nFace]; + const labelList& facePoints = this->mesh_.faces()[nFace]; + tetPoints[2] = this->mesh_.faceCentres()[nFace]; label pointi = 0; @@ -35,8 +35,8 @@ bool interpolationCellPointFace::findTriangle tetPointLabels[0] = facePoints[pointi]; tetPointLabels[1] = facePoints[nextPointLabel]; - tetPoints[0] = this->pMeshPoints_[tetPointLabels[0]]; - tetPoints[1] = this->pMeshPoints_[tetPointLabels[1]]; + tetPoints[0] = this->mesh_.points()[tetPointLabels[0]]; + tetPoints[1] = this->mesh_.points()[tetPointLabels[1]]; vector fc = (tetPoints[0] + tetPoints[1] + tetPoints[2])/3.0; diff --git a/src/finiteVolume/interpolation/interpolation/interpolationCellPointFace/interpolationCellPointFace.C b/src/finiteVolume/interpolation/interpolation/interpolationCellPointFace/interpolationCellPointFace.C index a1059a39c4..690e27bd65 100644 --- a/src/finiteVolume/interpolation/interpolation/interpolationCellPointFace/interpolationCellPointFace.C +++ b/src/finiteVolume/interpolation/interpolation/interpolationCellPointFace/interpolationCellPointFace.C @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | Website: https://openfoam.org - \\ / A nd | Copyright (C) 2011-2018 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2011-2020 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -39,7 +39,7 @@ Foam::interpolationCellPointFace::interpolationCellPointFace const GeometricField& psi ) : - interpolation(psi), + fieldInterpolation>(psi), psip_ ( volPointInterpolation::New(psi.mesh()).interpolate @@ -73,8 +73,8 @@ Type Foam::interpolationCellPointFace::interpolate // only use face information when the position is on a face if (facei < 0) { - const vector& cellCentre = this->pMesh_.cellCentres()[celli]; - const labelList& cellFaces = this->pMesh_.cells()[celli]; + const vector& cellCentre = this->mesh_.cellCentres()[celli]; + const labelList& cellFaces = this->mesh_.cells()[celli]; vector projection = position - cellCentre; tetPoints[3] = cellCentre; @@ -92,10 +92,10 @@ Type Foam::interpolationCellPointFace::interpolate { label nFace = cellFaces[facei]; - vector normal = this->pMeshFaceAreas_[nFace]; + vector normal = this->mesh_.faceAreas()[nFace]; normal /= mag(normal); - const vector& faceCentreTmp = this->pMeshFaceCentres_[nFace]; + const vector& faceCentreTmp = this->mesh_.faceCentres()[nFace]; scalar multiplierNumerator = (faceCentreTmp - cellCentre) & normal; scalar multiplierDenominator = projection & normal; @@ -167,7 +167,7 @@ Type Foam::interpolationCellPointFace::interpolate while (facei < cellFaces.size() && !foundTet) { label nFace = cellFaces[facei]; - if (nFace < this->pMeshFaceAreas_.size()) + if (nFace < this->mesh_.faceAreas().size()) { foundTet = findTet ( @@ -221,7 +221,7 @@ Type Foam::interpolationCellPointFace::interpolate else { label patchi = - this->pMesh_.boundaryMesh().whichPatch(closestFace); + this->mesh_.boundaryMesh().whichPatch(closestFace); // If the boundary patch is not empty use the face value // else use the cell value @@ -229,7 +229,7 @@ Type Foam::interpolationCellPointFace::interpolate { ts[2] = this->psi_.boundaryField()[patchi] [ - this->pMesh_.boundaryMesh()[patchi].whichFace + this->mesh_.boundaryMesh()[patchi].whichFace ( closestFace ) @@ -265,7 +265,7 @@ Type Foam::interpolationCellPointFace::interpolate else { label patchi = - this->pMesh_.boundaryMesh().whichPatch(closestFace); + this->mesh_.boundaryMesh().whichPatch(closestFace); // If the boundary patch is not empty use the face value // else use the cell value @@ -273,7 +273,7 @@ Type Foam::interpolationCellPointFace::interpolate { t = this->psi_.boundaryField()[patchi] [ - this->pMesh_.boundaryMesh()[patchi].whichFace + this->mesh_.boundaryMesh()[patchi].whichFace ( closestFace ) @@ -312,14 +312,14 @@ Type Foam::interpolationCellPointFace::interpolate } else { - label patchi = this->pMesh_.boundaryMesh().whichPatch(facei); + label patchi = this->mesh_.boundaryMesh().whichPatch(facei); // If the boundary patch is not empty use the face value // else use the cell value if (this->psi_.boundaryField()[patchi].size()) { t += phi[2]*this->psi_.boundaryField()[patchi] - [this->pMesh_.boundaryMesh()[patchi].whichFace(facei)]; + [this->mesh_.boundaryMesh()[patchi].whichFace(facei)]; } else { @@ -336,14 +336,14 @@ Type Foam::interpolationCellPointFace::interpolate } else { - label patchi = this->pMesh_.boundaryMesh().whichPatch(facei); + label patchi = this->mesh_.boundaryMesh().whichPatch(facei); // If the boundary patch is not empty use the face value // else use the cell value if (this->psi_.boundaryField()[patchi].size()) { t = this->psi_.boundaryField()[patchi] - [this->pMesh_.boundaryMesh()[patchi].whichFace(facei)]; + [this->mesh_.boundaryMesh()[patchi].whichFace(facei)]; } else { diff --git a/src/finiteVolume/interpolation/interpolation/interpolationCellPointFace/interpolationCellPointFace.H b/src/finiteVolume/interpolation/interpolation/interpolationCellPointFace/interpolationCellPointFace.H index 8f6ca40ed6..8d5893149a 100644 --- a/src/finiteVolume/interpolation/interpolation/interpolationCellPointFace/interpolationCellPointFace.H +++ b/src/finiteVolume/interpolation/interpolation/interpolationCellPointFace/interpolationCellPointFace.H @@ -46,7 +46,7 @@ namespace Foam template class interpolationCellPointFace : - public interpolation + public fieldInterpolation> { // Private Data @@ -95,6 +95,9 @@ public: // Member Functions + //- Inherit interpolate from interpolation + using interpolation::interpolate; + //- Interpolate field to the given point in the given cell Type interpolate ( diff --git a/src/finiteVolume/interpolation/interpolation/interpolationPointMVC/interpolationPointMVC.C b/src/finiteVolume/interpolation/interpolation/interpolationPointMVC/interpolationPointMVC.C index cf87b23973..b926726f48 100644 --- a/src/finiteVolume/interpolation/interpolation/interpolationPointMVC/interpolationPointMVC.C +++ b/src/finiteVolume/interpolation/interpolation/interpolationPointMVC/interpolationPointMVC.C @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | Website: https://openfoam.org - \\ / A nd | Copyright (C) 2011-2018 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2011-2020 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -34,7 +34,7 @@ Foam::interpolationPointMVC::interpolationPointMVC const GeometricField& psi ) : - interpolation(psi), + fieldInterpolation>(psi), psip_ ( volPointInterpolation::New(psi.mesh()).interpolate diff --git a/src/finiteVolume/interpolation/interpolation/interpolationPointMVC/interpolationPointMVC.H b/src/finiteVolume/interpolation/interpolation/interpolationPointMVC/interpolationPointMVC.H index f655ab41bf..c777011ada 100644 --- a/src/finiteVolume/interpolation/interpolation/interpolationPointMVC/interpolationPointMVC.H +++ b/src/finiteVolume/interpolation/interpolation/interpolationPointMVC/interpolationPointMVC.H @@ -48,7 +48,7 @@ namespace Foam template class interpolationPointMVC : - public interpolation + public fieldInterpolation> { protected: diff --git a/src/finiteVolume/interpolation/interpolation/interpolationPointMVC/interpolationPointMVCI.H b/src/finiteVolume/interpolation/interpolation/interpolationPointMVC/interpolationPointMVCI.H index 872106a707..bf1fe8980a 100644 --- a/src/finiteVolume/interpolation/interpolation/interpolationPointMVC/interpolationPointMVCI.H +++ b/src/finiteVolume/interpolation/interpolation/interpolationPointMVC/interpolationPointMVCI.H @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | Website: https://openfoam.org - \\ / A nd | Copyright (C) 2011-2018 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2011-2020 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -45,7 +45,7 @@ inline Type Foam::interpolationPointMVC::interpolate { return interpolate ( - pointMVCWeight(this->pMesh_, position, celli, facei) + pointMVCWeight(this->mesh_, position, celli, facei) ); }