From f588a922345eb91519290217be6b16dc64b4ac88 Mon Sep 17 00:00:00 2001 From: Markus Muehlhaeusser Date: Fri, 17 Nov 2023 15:23:01 +0100 Subject: [PATCH] Add Alpha and AlphaGradient files --- FF/Alpha.C | 115 +++++++++++++++++++++++++++++++++++++++++++++ FF/Alpha.H | 43 +++++++++++++++++ FF/AlphaGradient.C | 73 ++++++++++++++++++++++++++++ FF/AlphaGradient.H | 42 +++++++++++++++++ 4 files changed, 273 insertions(+) create mode 100644 FF/Alpha.C create mode 100644 FF/Alpha.H create mode 100644 FF/AlphaGradient.C create mode 100644 FF/AlphaGradient.H diff --git a/FF/Alpha.C b/FF/Alpha.C new file mode 100644 index 00000000..070b2050 --- /dev/null +++ b/FF/Alpha.C @@ -0,0 +1,115 @@ +#include "Alpha.H" + +using namespace Foam; + +preciceAdapter::FF::Alpha::Alpha( + const Foam::fvMesh& mesh, + const std::string nameAlpha) +: Alpha_( + const_cast( + &mesh.lookupObject(nameAlpha))) +{ + dataType_ = scalar; +} + +std::size_t preciceAdapter::FF::Alpha::write(double* buffer, bool meshConnectivity, const unsigned int dim) +{ + int bufferIndex = 0; + + if (this->locationType_ == LocationType::volumeCenters) + { + if (cellSetNames_.empty()) + { + for (const auto& cell : Alpha_->internalField()) + { + buffer[bufferIndex++] = cell; + } + } + else + { + for (const auto& cellSetName : cellSetNames_) + { + cellSet overlapRegion(Alpha_->mesh(), cellSetName); + const labelList& cells = overlapRegion.toc(); + + for (const auto& currentCell : cells) + { + // Copy the alpha valus into the buffer + buffer[bufferIndex++] = Alpha_->internalField()[currentCell]; + } + } + } + } + + // For every boundary patch of the interface + for (uint j = 0; j < patchIDs_.size(); j++) + { + int patchID = patchIDs_.at(j); + + // For every cell of the patch + forAll(Alpha_->boundaryFieldRef()[patchID], i) + { + // Copy the Alpha into the buffer + buffer[bufferIndex++] = + Alpha_->boundaryFieldRef()[patchID][i]; + } + } + return bufferIndex; +} + +void preciceAdapter::FF::Alpha::read(double* buffer, const unsigned int dim) +{ + int bufferIndex = 0; + + if (this->locationType_ == LocationType::volumeCenters) + { + if (cellSetNames_.empty()) + { + for (auto& cell : Alpha_->ref()) + { + cell = buffer[bufferIndex++]; + } + } + else + { + for (const auto& cellSetName : cellSetNames_) + { + cellSet overlapRegion(Alpha_->mesh(), cellSetName); + const labelList& cells = overlapRegion.toc(); + + for (const auto& currentCell : cells) + { + // Copy the pressure into the buffer + Alpha_->ref()[currentCell] = buffer[bufferIndex++]; + } + } + } + } + + // For every boundary patch of the interface + for (uint j = 0; j < patchIDs_.size(); j++) + { + int patchID = patchIDs_.at(j); + // For every cell of the patch + forAll(Alpha_->boundaryFieldRef()[patchID], i) + { + Alpha_->boundaryFieldRef()[patchID][i] = buffer[bufferIndex++]; + } + } +} + +bool preciceAdapter::FF::Alpha::isLocationTypeSupported(const bool meshConnectivity) const +{ + if (meshConnectivity) + { + return (this->locationType_ == LocationType::faceCenters); + } + else + { + return (this->locationType_ == LocationType::faceCenters || this->locationType_ == LocationType::volumeCenters); + }} + +std::string preciceAdapter::FF::Alpha::getDataName() const +{ + return "Alpha"; +} diff --git a/FF/Alpha.H b/FF/Alpha.H new file mode 100644 index 00000000..7e5feb99 --- /dev/null +++ b/FF/Alpha.H @@ -0,0 +1,43 @@ +#ifndef FF_ALPHA_H +#define FF_ALPHA_H + +#include "CouplingDataUser.H" + +#include "fvCFD.H" +#include "cellSet.H" + +namespace preciceAdapter +{ +namespace FF +{ + +//- Class that writes and reads Alpha +class Alpha : public CouplingDataUser +{ + +private: + //- Alpha field + Foam::volScalarField* Alpha_; + +public: + //- Constructor + Alpha( + const Foam::fvMesh& mesh, + const std::string nameAlpha); + + //- Write the Alpha values into the buffer + std::size_t write(double* buffer, bool meshConnectivity, const unsigned int dim); + + //- Read the Alpha values from the buffer + void read(double* buffer, const unsigned int dim); + + bool isLocationTypeSupported(const bool meshConnectivity) const override; + + //- Get the name of the current data field + std::string getDataName() const override; +}; + +} +} + +#endif diff --git a/FF/AlphaGradient.C b/FF/AlphaGradient.C new file mode 100644 index 00000000..e80452be --- /dev/null +++ b/FF/AlphaGradient.C @@ -0,0 +1,73 @@ +#include "AlphaGradient.H" +#include "mixedFvPatchFields.H" + +using namespace Foam; + +preciceAdapter::FF::AlphaGradient::AlphaGradient( + const Foam::fvMesh& mesh, + const std::string nameAlpha) +: Alpha_( + const_cast( + &mesh.lookupObject(nameAlpha))) +{ + dataType_ = scalar; +} + +std::size_t preciceAdapter::FF::AlphaGradient::write(double* buffer, bool meshConnectivity, const unsigned int dim) +{ + int bufferIndex = 0; + + // For every boundary patch of the interface + for (uint j = 0; j < patchIDs_.size(); j++) + { + int patchID = patchIDs_.at(j); + + // Get the Alpha gradient boundary patch + const scalarField gradientPatch((Alpha_->boundaryFieldRef()[patchID]) + .snGrad()); + + // For every cell of the patch + forAll(gradientPatch, i) + { + // Copy the Alpha gradient into the buffer + buffer[bufferIndex++] = + -gradientPatch[i]; + } + } +} + +void preciceAdapter::FF::AlphaGradient::read(double* buffer, const unsigned int dim) +{ + int bufferIndex = 0; + + // For every boundary patch of the interface + for (uint j = 0; j < patchIDs_.size(); j++) + { + int patchID = patchIDs_.at(j); + + // Get the Alpha gradient boundary patch + scalarField& gradientPatch = + refCast( + Alpha_->boundaryFieldRef()[patchID]) + .gradient(); + + // For every cell of the patch + forAll(gradientPatch, i) + { + // Set the Alpha gradient as the buffer value + gradientPatch[i] = + buffer[bufferIndex++]; + } + } + return bufferIndex; +} + +bool preciceAdapter::FF::AlphaGradient::isLocationTypeSupported(const bool meshConnectivity) const +{ + return (this->locationType_ == LocationType::faceCenters); +} + +std::string preciceAdapter::FF::AlphaGradient::getDataName() const +{ + return "AlphaGradient"; +} diff --git a/FF/AlphaGradient.H b/FF/AlphaGradient.H new file mode 100644 index 00000000..b538f009 --- /dev/null +++ b/FF/AlphaGradient.H @@ -0,0 +1,42 @@ +#ifndef FF_ALPHA_GRADIENT_H +#define FF_ALPHA_GRADIENT_H + +#include "CouplingDataUser.H" + +#include "fvCFD.H" + +namespace preciceAdapter +{ +namespace FF +{ + +//- Class that writes and reads Alpha gradient +class AlphaGradient : public CouplingDataUser +{ + +private: + //- Alpha field + Foam::volScalarField* Alpha_; + +public: + //- Constructor + AlphaGradient( + const Foam::fvMesh& mesh, + const std::string nameAlpha); + + //- Write the Alpha gradient values into the buffer + std::size_t write(double* buffer, bool meshConnectivity, const unsigned int dim); + + //- Read the Alpha gradient values from the buffer + void read(double* buffer, const unsigned int dim); + + bool isLocationTypeSupported(const bool meshConnectivity) const override; + + //- Get the name of the current data field + std::string getDataName() const override; +}; + +} +} + +#endif