diff --git a/.gitignore b/.gitignore index d273fd25..b7f69be9 100644 --- a/.gitignore +++ b/.gitignore @@ -37,6 +37,7 @@ # OpenFOAM / WMake lnInclude/ Make/linux64GccDPInt32Opt/ +Make/darwin64ClangDPInt32Opt/ # Editors .cproject diff --git a/Adapter.C b/Adapter.C index e4da2dca..faaab7b2 100644 --- a/Adapter.C +++ b/Adapter.C @@ -757,6 +757,12 @@ void preciceAdapter::Adapter::storeMeshPoints() void preciceAdapter::Adapter::reloadMeshPoints() { + if (!mesh_.moving()) + { + DEBUG(adapterInfo("Mesh points not moved as the mesh is not moving")); + return; + } + // In Foam::polyMesh::movePoints. // TODO: The function movePoints overwrites the pointer to the old mesh. // Therefore, if you revert the mesh, the oldpointer will be set to the points, which are the new values. @@ -1301,8 +1307,8 @@ void preciceAdapter::Adapter::readMeshCheckpoint() { DEBUG(adapterInfo("Reading a mesh checkpoint...")); - //TODO only the meshPhi field is here, which is a surfaceScalarField. The other fields can be removed. - // Reload all the fields of type mesh surfaceScalarField + // TODO only the meshPhi field is here, which is a surfaceScalarField. The other fields can be removed. + // Reload all the fields of type mesh surfaceScalarField for (uint i = 0; i < meshSurfaceScalarFields_.size(); i++) { // Load the volume field @@ -1529,8 +1535,8 @@ void preciceAdapter::Adapter::teardown() } meshVolVectorFieldCopies_.clear(); - //TODO for the internal volume - // volScalarInternal + // TODO for the internal volume + // volScalarInternal for (uint i = 0; i < volScalarInternalFieldCopies_.size(); i++) { delete volScalarInternalFieldCopies_.at(i); diff --git a/FSI/Displacement.C b/FSI/Displacement.C index 90a6e2de..398a65ed 100644 --- a/FSI/Displacement.C +++ b/FSI/Displacement.C @@ -7,8 +7,10 @@ preciceAdapter::FSI::Displacement::Displacement( const std::string namePointDisplacement, const std::string nameCellDisplacement) : pointDisplacement_( - const_cast( - &mesh.lookupObject(namePointDisplacement))), + namePointDisplacement == "unused" + ? nullptr + : const_cast( + &mesh.lookupObject(namePointDisplacement))), cellDisplacement_( const_cast( &mesh.lookupObject(nameCellDisplacement))), @@ -37,13 +39,52 @@ void preciceAdapter::FSI::Displacement::initialize() void preciceAdapter::FSI::Displacement::write(double* buffer, bool meshConnectivity, const unsigned int dim) { /* TODO: Implement - * We need two nested for-loops for each patch, - * the outer for the locations and the inner for the dimensions. - * See the preCICE writeBlockVectorData() implementation. - */ - FatalErrorInFunction - << "Writing displacements is not supported." - << exit(FatalError); + * We need two nested for-loops for each patch, + * the outer for the locations and the inner for the dimensions. + * See the preCICE writeBlockVectorData() implementation. + */ + + // Copy the displacement field from OpenFOAM to the buffer + + if (this->locationType_ == LocationType::faceCenters) + { + // For every boundary patch of the interface + for (const label patchID : patchIDs_) + { + // Write the displacement to the preCICE buffer + // For every cell of the patch + forAll(cellDisplacement_->boundaryField()[patchID], i) + { + for (unsigned int d = 0; d < dim; ++d) + buffer[i * dim + d] = + cellDisplacement_->boundaryField()[patchID][i][d]; + } + } + } + else if (this->locationType_ == LocationType::faceNodes) + { + DEBUG(adapterInfo( + "Please be aware of issues with using 'locationType faceNodes' " + "in parallel. \n" + "See https://github.com/precice/openfoam-adapter/issues/153.", + "warning")); + + // For every boundary patch of the interface + for (const label patchID : patchIDs_) + { + // Write the displacement to the preCICE buffer + // For every cell of the patch + forAll(pointDisplacement_->boundaryField()[patchID], i) + { + const labelList& meshPoints = + mesh_.boundaryMesh()[patchID].meshPoints(); + + for (unsigned int d = 0; d < dim; ++d) + buffer[i * dim + d] = + pointDisplacement_->internalField()[meshPoints[i]][d]; + } + } + } } diff --git a/FSI/DisplacementDelta.C b/FSI/DisplacementDelta.C index ca4e4c2a..53b56240 100644 --- a/FSI/DisplacementDelta.C +++ b/FSI/DisplacementDelta.C @@ -37,10 +37,10 @@ void preciceAdapter::FSI::DisplacementDelta::initialize() void preciceAdapter::FSI::DisplacementDelta::write(double* buffer, bool meshConnectivity, const unsigned int dim) { /* TODO: Implement - * We need two nested for-loops for each patch, - * the outer for the locations and the inner for the dimensions. - * See the preCICE writeBlockVectorData() implementation. - */ + * We need two nested for-loops for each patch, + * the outer for the locations and the inner for the dimensions. + * See the preCICE writeBlockVectorData() implementation. + */ FatalErrorInFunction << "Writing displacementDeltas is not supported." << exit(FatalError); diff --git a/FSI/FSI.C b/FSI/FSI.C index a10e0326..5050bb31 100644 --- a/FSI/FSI.C +++ b/FSI/FSI.C @@ -27,7 +27,9 @@ bool preciceAdapter::FSI::FluidStructureInteraction::configure(const IOdictionar // addWriters() and addReaders(). // Check the solver type and determine it if needed if ( - solverType_.compare("compressible") == 0 || solverType_.compare("incompressible") == 0) + solverType_.compare("compressible") == 0 + || solverType_.compare("incompressible") == 0 + || solverType_.compare("solid") == 0) { DEBUG(adapterInfo("Known solver type: " + solverType_)); } @@ -38,7 +40,7 @@ bool preciceAdapter::FSI::FluidStructureInteraction::configure(const IOdictionar } else { - DEBUG(adapterInfo("Determining the solver type for the FSI module... (override by setting solverType to one of {compressible, incompressible})")); + DEBUG(adapterInfo("Determining the solver type for the FSI module... (override by setting solverType to one of {compressible, incompressible, solid})")); solverType_ = determineSolverType(); } @@ -54,8 +56,8 @@ bool preciceAdapter::FSI::FluidStructureInteraction::readConfig(const IOdictiona DEBUG(adapterInfo(" user-defined solver type : " + solverType_)); /* TODO: Read the names of any needed fields and parameters. - * Include the force here? - */ + * Include the force here? + */ // Read the name of the pointDisplacement field (if different) namePointDisplacement_ = FSIdict.lookupOrDefault("namePointDisplacement", "pointDisplacement"); @@ -65,6 +67,10 @@ bool preciceAdapter::FSI::FluidStructureInteraction::readConfig(const IOdictiona nameCellDisplacement_ = FSIdict.lookupOrDefault("nameCellDisplacement", "cellDisplacement"); DEBUG(adapterInfo(" cellDisplacement field name : " + nameCellDisplacement_)); + // Read the name of the solidForce field (if different) + nameSolidForce_ = FSIdict.lookupOrDefault("nameSolidForce", "solidForce"); + DEBUG(adapterInfo(" solidForce field name : " + nameSolidForce_)); + return true; } @@ -117,7 +123,7 @@ bool preciceAdapter::FSI::FluidStructureInteraction::addWriters(std::string data { interface->addCouplingDataWriter( dataName, - new Force(mesh_, solverType_) /* TODO: Add any other arguments here */ + new Force(mesh_, solverType_, nameSolidForce_) /* TODO: Add any other arguments here */ ); DEBUG(adapterInfo("Added writer: Force.")); } @@ -165,7 +171,7 @@ bool preciceAdapter::FSI::FluidStructureInteraction::addReaders(std::string data { interface->addCouplingDataReader( dataName, - new Force(mesh_, solverType_) /* TODO: Add any other arguments here */ + new Force(mesh_, solverType_, nameSolidForce_) /* TODO: Add any other arguments here */ ); DEBUG(adapterInfo("Added reader: Force.")); } diff --git a/FSI/FSI.H b/FSI/FSI.H index 3e93d77e..008da0a6 100644 --- a/FSI/FSI.H +++ b/FSI/FSI.H @@ -34,6 +34,9 @@ protected: //- Name of the pointDisplacement field std::string nameCellDisplacement_ = "cellDisplacement"; + //- Name of the force field used by the solid + std::string nameSolidForce_ = "solidForce"; + /* TODO: Declare here any parameters that should be read from / the configuration file. See CHT/CHT.H for reference. / We want to support in-house solvers with different field names, diff --git a/FSI/Force.C b/FSI/Force.C index 46c298fb..d2867f30 100644 --- a/FSI/Force.C +++ b/FSI/Force.C @@ -4,7 +4,8 @@ using namespace Foam; preciceAdapter::FSI::Force::Force( const Foam::fvMesh& mesh, - const std::string solverType) + const std::string solverType, + const std::string nameSolidForce) : ForceBase(mesh, solverType) { Force_ = new volVectorField( @@ -19,6 +20,17 @@ preciceAdapter::FSI::Force::Force( "fdim", dimensionSet(1, 1, -2, 0, 0, 0, 0), Foam::vector::zero)); + + if (mesh_.foundObject(nameSolidForce)) + { + solidForce_ = + &const_cast( + mesh_.lookupObject(nameSolidForce)); + } + else + { + solidForce_ = nullptr; + } } void preciceAdapter::FSI::Force::write(double* buffer, bool meshConnectivity, const unsigned int dim) @@ -28,7 +40,37 @@ void preciceAdapter::FSI::Force::write(double* buffer, bool meshConnectivity, co void preciceAdapter::FSI::Force::read(double* buffer, const unsigned int dim) { - this->readFromBuffer(buffer); + // Copy the force field from the buffer to OpenFOAM + + // Here we assume that a force volVectorField exists, which is used by + // the OpenFOAM solver + + // Set boundary forces + for (unsigned int j = 0; j < patchIDs_.size(); j++) + { + // Get the ID of the current patch + const unsigned int patchID = patchIDs_.at(j); + + if (this->locationType_ == LocationType::faceCenters) + { + // Make a force field + vectorField& force = solidForce_->boundaryFieldRef()[patchID]; + + // Copy the forces from the buffer into the force field + forAll(force, i) + { + for (unsigned int d = 0; d < dim; ++d) + force[i][d] = buffer[i * dim + d]; + } + } + else if (this->locationType_ == LocationType::faceNodes) + { + // Here we could easily interpolate the face values to point values + // and assign them to some field, but I guess there is no need + // unless it will be used + notImplemented("Read forces not implemented for faceNodes!"); + } + } } bool preciceAdapter::FSI::Force::isLocationTypeSupported(const bool meshConnectivity) const diff --git a/FSI/Force.H b/FSI/Force.H index 062e42c0..851863f3 100644 --- a/FSI/Force.H +++ b/FSI/Force.H @@ -15,11 +15,15 @@ private: //- Force field Foam::volVectorField* Force_; + //- Solid force field + Foam::volVectorField* solidForce_; + public: //- Constructor Force( const Foam::fvMesh& mesh, - const std::string solverType); + const std::string solverType, + const std::string nameSolidForce); //- Write the forces values into the buffer void write(double* buffer, bool meshConnectivity, const unsigned int dim) override; diff --git a/FSI/ForceBase.C b/FSI/ForceBase.C index 37e1c9e8..b9aaafba 100644 --- a/FSI/ForceBase.C +++ b/FSI/ForceBase.C @@ -9,22 +9,24 @@ preciceAdapter::FSI::ForceBase::ForceBase( : mesh_(mesh), solverType_(solverType) { - //What about type "basic"? - if (solverType_.compare("incompressible") != 0 && solverType_.compare("compressible") != 0) + // What about type "basic"? + if (solverType_.compare("incompressible") != 0 + && solverType_.compare("compressible") != 0 + && solverType_.compare("solid") != 0) { FatalErrorInFunction << "Force based calculations only support " - << "compressible or incompressible solver types." + << "compressible, incompressible, or solid solver types." << exit(FatalError); } dataType_ = vector; } -//Calculate viscous force +// Calculate viscous force Foam::tmp preciceAdapter::FSI::ForceBase::devRhoReff() const { - //For turbulent flows + // For turbulent flows typedef compressible::turbulenceModel cmpTurbModel; typedef incompressible::turbulenceModel icoTurbModel; @@ -52,7 +54,7 @@ Foam::tmp preciceAdapter::FSI::ForceBase::devRhoReff() } } -//lookup correct rho +// lookup correct rho Foam::tmp preciceAdapter::FSI::ForceBase::rho() const { // If volScalarField exists, read it from registry (for compressible cases) @@ -88,10 +90,9 @@ Foam::tmp preciceAdapter::FSI::ForceBase::rho() const } } -//lookup correct mu +// lookup correct mu Foam::tmp preciceAdapter::FSI::ForceBase::mu() const { - if (solverType_.compare("incompressible") == 0) { typedef immiscibleIncompressibleTwoPhaseMixture iitpMixture; @@ -191,10 +192,10 @@ void preciceAdapter::FSI::ForceBase::writeToBuffer(double* buffer, void preciceAdapter::FSI::ForceBase::readFromBuffer(double* buffer) const { /* TODO: Implement - * We need two nested for-loops for each patch, - * the outer for the locations and the inner for the dimensions. - * See the preCICE readBlockVectorData() implementation. - */ + * We need two nested for-loops for each patch, + * the outer for the locations and the inner for the dimensions. + * See the preCICE readBlockVectorData() implementation. + */ FatalErrorInFunction << "Reading forces is not supported." << exit(FatalError); diff --git a/changelog-entries/236.md b/changelog-entries/236.md new file mode 100644 index 00000000..e80f533f --- /dev/null +++ b/changelog-entries/236.md @@ -0,0 +1 @@ +- Added functionality to allow use of solids4foam with the OpenFOAM adapter diff --git a/docs/README.md b/docs/README.md index ab1265f1..a1a437c5 100644 --- a/docs/README.md +++ b/docs/README.md @@ -18,9 +18,9 @@ This adapter can read/write the following fields: - Heat flux (read + write) - Sink temperature (read + write) - Heat transfer coefficient (read + write) -- Force (write) +- Force (read + write) - Stress (write) -- Displacement (read) +- Displacement (read + write) - Displacement delta (read) - Pressure (read + write) - Pressure gradient (read + write) diff --git a/docs/config.md b/docs/config.md index 6371bac7..b40a8ab1 100644 --- a/docs/config.md +++ b/docs/config.md @@ -120,7 +120,7 @@ writeData ); ``` -For fluid-structure interaction, `writeData` can be `Force` or `Stress`, where `Stress` is essentially a force vector scaled by the cell face in spatial coordinates (with any postfix), thus, a conservative quantity as well.`readData` can be `Displacement` and `DisplacementDelta` (with any postfix). `DisplacementDelta` refers to the last coupling time step, which needs to considered in the case of subcycling. +For fluid-structure interaction, `writeData` can be `Force` or `Stress`, where `Stress` is essentially a force vector scaled by the cell face in spatial coordinates (with any postfix), thus, a conservative quantity as well. `readData` can be `Displacement` and `DisplacementDelta` (with any postfix). `DisplacementDelta` refers to the last coupling time step, which needs to considered in the case of subcycling. Structure solvers (such as solids4Foam) can also write `Displacement` and read `Force`. {% warning %} You will run into problems when you use `Displacement(Delta)` as write data set and execute RBF mappings in parallel. This would affect users who use OpenFOAM and the adapter as the Solid participant in order to compute solid mechanics with OpenFOAM (currently not officially supported at all). Have a look [at this issue on GitHub](https://github.com/precice/openfoam-adapter/issues/153) for details. @@ -312,7 +312,7 @@ and depends on the density (`rho [ 1 -3 0 0 0 0 0 ]`) and heat capacity (`Cp #### Fluid-structure interaction -The adapter's FSI functionality supports both compressible and incompressible solvers. +The adapter's FSI functionality supports both compressible and incompressible solvers, as well as solid (e.g., solids4Foam) solvers. For incompressible solvers, it tries to read uniform values for the density and kinematic viscosity (if it is not already available) from the `FSI` subdictionary of `preciceDict`: diff --git a/docs/openfoam-support.md b/docs/openfoam-support.md index 29beb897..0eba1fb8 100644 --- a/docs/openfoam-support.md +++ b/docs/openfoam-support.md @@ -40,7 +40,7 @@ Known not supported versions: OpenFOAM v1606+ or older, OpenFOAM 3 or older, foa ## Supported OpenFOAM solvers -We support mainstream OpenFOAM solvers such as pimpleFoam (for FSI) or buoyantPimpleFoam, buoyantSimpleFoam, laplacianFoam (for CHT). Our community has tried the adapter with multiple different solvers that support function objects. +We support mainstream OpenFOAM solvers such as pimpleFoam and solids4Foam for FSI, buoyantPimpleFoam, buoyantSimpleFoam, and laplacianFoam for CHT, or pimpleFoam and sonicLiquidFoam for FF. Our community has, additionally, tried the adapter with multiple different solvers that support function objects. ## Notes on OpenFOAM features