Adds functionality to allow FSI to work with solids4foam as the solid solver

.gitignore

@@ -37,6 +37,7 @@
 # OpenFOAM / WMake
 lnInclude/
 Make/linux64GccDPInt32Opt/
+Make/darwin64ClangDPInt32Opt/
 
 # Editors
 .cproject

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);

FSI/Displacement.C

@@ -7,8 +7,10 @@ preciceAdapter::FSI::Displacement::Displacement(
     const std::string namePointDisplacement,
     const std::string nameCellDisplacement)
 : pointDisplacement_(
-    const_cast<pointVectorField*>(
-        &mesh.lookupObject<pointVectorField>(namePointDisplacement))),
+    namePointDisplacement == "unused"
+        ? nullptr
+        : const_cast<pointVectorField*>(
+            &mesh.lookupObject<pointVectorField>(namePointDisplacement))),
   cellDisplacement_(
       const_cast<volVectorField*>(
           &mesh.lookupObject<volVectorField>(nameCellDisplacement))),
@@ -37,13 +39,54 @@ 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)
+    {
+        WarningIn(
+            "void preciceAdapter::FSI::Displacement::write"
+            "(double* buffer, bool meshConnectivity, const unsigned int dim))")
+            << "Please be aware of issues with using 'locationType faceNodes' "
+            << "in parallel. "
+            << "See https://github.com/precice/openfoam-adapter/issues/153."
+            << endl;
+
+        // 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];
+            }
+        }
+    }
 }
 
 

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);

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<word>("namePointDisplacement", "pointDisplacement");
@@ -65,6 +67,10 @@ bool preciceAdapter::FSI::FluidStructureInteraction::readConfig(const IOdictiona
     nameCellDisplacement_ = FSIdict.lookupOrDefault<word>("nameCellDisplacement", "cellDisplacement");
     DEBUG(adapterInfo("    cellDisplacement field name : " + nameCellDisplacement_));
 
+    // Read the name of the solidForce field (if different)
+    nameSolidForce_ = FSIdict.lookupOrDefault<word>("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."));
     }

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,

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<volVectorField>(nameSolidForce))
+    {
+        solidForce_ =
+            &const_cast<volVectorField&>(
+                mesh_.lookupObject<volVectorField>(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

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;

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<Foam::volSymmTensorField> 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<Foam::volSymmTensorField> preciceAdapter::FSI::ForceBase::devRhoReff()
     }
 }
 
-//lookup correct rho
+// lookup correct rho
 Foam::tmp<Foam::volScalarField> preciceAdapter::FSI::ForceBase::rho() const
 {
     // If volScalarField exists, read it from registry (for compressible cases)
@@ -88,10 +90,9 @@ Foam::tmp<Foam::volScalarField> preciceAdapter::FSI::ForceBase::rho() const
     }
 }
 
-//lookup correct mu
+// lookup correct mu
 Foam::tmp<Foam::volScalarField> 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);

changelog-entries/236.md

@@ -0,0 +1 @@
+- Added functionality to allow use of solids4foam with the OpenFOAM adapter