diff --git a/src/finiteVolume/cfdTools/general/MRF/MRFZone.C b/src/finiteVolume/cfdTools/general/MRF/MRFZone.C index e8ca39b089..ffa1d2bc6a 100644 --- a/src/finiteVolume/cfdTools/general/MRF/MRFZone.C +++ b/src/finiteVolume/cfdTools/general/MRF/MRFZone.C @@ -375,6 +375,45 @@ void Foam::MRFZone::relativeVelocity(volVectorField& U) const } +void Foam::MRFZone::absoluteVelocity(volVectorField& U) const +{ + const volVectorField& C = mesh_.C(); + + const vector& origin = origin_.value(); + const vector& Omega = Omega_.value(); + + const labelList& cells = mesh_.cellZones()[cellZoneID_]; + + forAll(cells, i) + { + label celli = cells[i]; + U[celli] += (Omega ^ (C[celli] - origin)); + } + + // Included patches + forAll(includedFaces_, patchi) + { + forAll(includedFaces_[patchi], i) + { + label patchFacei = includedFaces_[patchi][i]; + U.boundaryField()[patchi][patchFacei] = + (Omega ^ (C.boundaryField()[patchi][patchFacei] - origin)); + } + } + + // Excluded patches + forAll(excludedFaces_, patchi) + { + forAll(excludedFaces_[patchi], i) + { + label patchFacei = excludedFaces_[patchi][i]; + U.boundaryField()[patchi][patchFacei] += + (Omega ^ (C.boundaryField()[patchi][patchFacei] - origin)); + } + } +} + + void Foam::MRFZone::relativeFlux(surfaceScalarField& phi) const { relativeRhoFlux(geometricOneField(), phi); diff --git a/src/finiteVolume/cfdTools/general/MRF/MRFZone.H b/src/finiteVolume/cfdTools/general/MRF/MRFZone.H index d925f70392..c9b0429951 100644 --- a/src/finiteVolume/cfdTools/general/MRF/MRFZone.H +++ b/src/finiteVolume/cfdTools/general/MRF/MRFZone.H @@ -175,6 +175,9 @@ public: //- Make the given absolute velocity relative within the MRF region void relativeVelocity(volVectorField& U) const; + //- Make the given relative velocity absolute within the MRF region + void absoluteVelocity(volVectorField& U) const; + //- Make the given absolute flux relative within the MRF region void relativeFlux(surfaceScalarField& phi) const; diff --git a/src/finiteVolume/cfdTools/general/MRF/MRFZones.C b/src/finiteVolume/cfdTools/general/MRF/MRFZones.C index d95d20ab0e..5ceb48af51 100644 --- a/src/finiteVolume/cfdTools/general/MRF/MRFZones.C +++ b/src/finiteVolume/cfdTools/general/MRF/MRFZones.C @@ -86,6 +86,15 @@ void Foam::MRFZones::relativeVelocity(volVectorField& U) const } +void Foam::MRFZones::absoluteVelocity(volVectorField& U) const +{ + forAll(*this, i) + { + operator[](i).absoluteVelocity(U); + } +} + + void Foam::MRFZones::relativeFlux(surfaceScalarField& phi) const { forAll(*this, i) diff --git a/src/finiteVolume/cfdTools/general/MRF/MRFZones.H b/src/finiteVolume/cfdTools/general/MRF/MRFZones.H index 53a561a9e8..ae7f4c374f 100644 --- a/src/finiteVolume/cfdTools/general/MRF/MRFZones.H +++ b/src/finiteVolume/cfdTools/general/MRF/MRFZones.H @@ -81,6 +81,9 @@ public: //- Make the given absolute velocity relative within the MRF region void relativeVelocity(volVectorField& U) const; + //- Make the given relative velocity absolute within the MRF region + void absoluteVelocity(volVectorField& U) const; + //- Make the given absolute flux relative within the MRF region void relativeFlux(surfaceScalarField& phi) const;