Add partitioned-heat OpenFOAM participant with solver

changelog-entries/223.md

@@ -0,0 +1 @@
+- Added OpenFOAM solver as well as the participants (Dirichlet and Neumann) to our partiitoned-heat tutorial [#223](https://github.com/precice/tutorials/pull/223)

partitioned-heat-conduction/README.md

@@ -15,50 +15,52 @@ We solve a partitioned heat equation. For information on the non-partitioned cas
 
 Case setup from [3]. `D` denotes the Dirichlet participant and `N` denotes the Neumann participant.
 
-The heat equation is solved on a rectangular domain `Omega = [0,2] x [0,1]` with given Dirichlet boundary conditions. We split the domain at `x_c = 1` using a straight vertical line, the coupling interface. The left part of the domain will be referred to as the Dirichlet partition and the right part as the Neumann partition. To couple the two participants we use Dirichlet-Neumann coupling. Here, the Dirichlet participant receives Dirichlet boundary conditions (`Temperature`) at the coupling interface and solves the heat equation using these boundary conditions on the left part of the domain. Then the Dirichlet participant computes the resulting heat flux (`Flux`) from the solution and sends it to the Neumann participant. The Neumann participant uses the flux as a Neumann boundary condition to solve the heat equation on the right part of the domain. We then extract the temperature from the solution and send it back to the Dirichlet participant. This establishes the coupling between the two participants.
+The heat equation is solved on a rectangular domain `Omega = [0,2] x [0,1]` with given Dirichlet boundary conditions. We split the domain at `x_c = 1` using a straight vertical line, the coupling interface. The left part of the domain will be referred to as the Dirichlet partition and the right part as the Neumann partition. To couple the two participants we use Dirichlet-Neumann coupling. Here, the Dirichlet participant receives Dirichlet boundary conditions (`Temperature`) at the coupling interface and solves the heat equation using these boundary conditions on the left part of the domain. Then the Dirichlet participant computes the resulting heat flux (`Heat-Flux`) from the solution and sends it to the Neumann participant. The Neumann participant uses the flux as a Neumann boundary condition to solve the heat equation on the right part of the domain. We then extract the temperature from the solution and send it back to the Dirichlet participant. This establishes the coupling between the two participants.
 
 This simple case allows us to compare the solution for the partitioned case to a known analytical solution (method of manufactures solutions, see [1, p.37ff]). For more usage examples and details, please refer to [3, sect. 4.1].
 
 ## Available solvers and dependencies
 
-You can either couple a solver with itself or different solvers with each other. In any case you will need to have preCICE and the python bindings installed on your system.
+You can either couple a solver with itself or different solvers with each other. In any case you will need to have preCICE
 
 * FEniCS. Install [FEniCS](https://fenicsproject.org/download/) and the [FEniCS-adapter](https://github.com/precice/fenics-adapter). The code is largely based on this [fenics-tutorial](https://github.com/hplgit/fenics-tutorial/blob/master/pub/python/vol1/ft03_heat.py) from [1].
 
 * Nutils. Install [Nutils](http://www.nutils.org/en/latest/).
 
-## Running the simulation
+* OpenFOAM. Install OpenFOAM and the [OpenFOAM adapter](https://www.precice.org/adapter-openfoam-overview.html). This tutorial uses a custom solver, which you can find in `tutorials/partitioned-heat-conduction/openfoam-solver` and build using `cd tutorials/partitioned-heat-conduction/openfoam-solver && wmake`. Have a look at the section below (Notes on the OpenFOAM case) for further information.
 
-This tutorial is for FEniCS and Nutils. You can find the corresponding `run.sh` script in the folders `fenics` and `nutils`.
+### Notes on the OpenFOAM case
 
-For choosing whether you want to run the Dirichlet-kind and a Neumann-kind participant, please provide the following commandline input:
+Running this tutorial with OpenFOAM is a bit of a challenge and requires some special considerations:
 
-* `-d` flag will enforce Dirichlet boundary conditions on the coupling interface.
-* `-n` flag will enforce Neumann boundary conditions on the coupling interface.
+* First of all, OpenFOAM does not provide a Laplace solver with a non-zero right-hand side. Therefore, we provide a modified Laplace solver together with the tutorial, which needs to be compiled before running the tutorial. The solver can be compiled by executing `wmake` in the solver directory `./openfoam-solver/`. The generated executable will be stored in the `FOAM_USER_APPBIN` by default. Afterwards, the custom solver `heatTransfer` can be started from the respective OpenFOAM case directory, as usual.
 
-For running the case, open two terminals run:
+* The second challenge is given by the time- and space-dependent Dirichlet boundary conditions required for domain boundaries not belonging to the interface. For this purpose, a valid installation of `groovyBC` (part of `swak4Foam`) is required.
 
-```bash
-cd fenics
-./run.sh -d
-```
+* The third challenge is given by the space-dependent initial conditions. We use `funkySetFields` (installed with OpenFOAM) to evaluate the initial condition. You can directly execute the `./run.sh` script, which calls the `setInitialField.sh` in order to evaluate the required initial condition and store it in the `0` directory. Note that `run.sh`  deletes the `0` time directory and copies it again from `0.orig`. If you start modifying the initial or boundary conditions, make sure you modify the files located in the `0.orig` directory in combination with the default `run.sh` scripts.
 
-and
+## Running the simulation
 
-```bash
-cd fenics
-./run.sh -n
-```
+You can find the corresponding `run.sh` script in each participant solver.
+
+In case of `fenics` and `nutils` the Dirichlet-kind and a Neumann-kind participant are currently merged into a single participant directory. Therefore, please provide the following command line input argument:
 
-If you want to use Nutils for one or both sides of the setup, just `cd nutils`. The FEniCS case also supports parallel runs. Here, you cannot use the `run.sh` script, but must simply execute
+* `-d` flag will enforce Dirichlet boundary conditions on the coupling interface.
+* `-n` flag will enforce Neumann boundary conditions on the coupling interface.
+
+For running the case, a Dirichlet and a Neumann participant need to be executed, e.g., `./run.sh -d` and `./run.sh. -n`
+
+The FEniCS case also supports parallel runs. Here, you cannot use the `run.sh` script, but must simply execute
 
 ```bash
 mpirun -n <N_PROC> python3 heat.py -d
 ```
 
+OpenFOAM supports parallel runs as usual. However, you need to execute the command manually by running: `mpirun -np <N_PROC> ./heatTransfer`.
+
 ### Note on the combination of Nutils & FEniCS
 
-You can mix the Nutils and FEniCS solver, if you like. Note that the error for a pure FEniCS simulation is lower than for a mixed one. We did not yet study the origin of this error, but assume that this is due to the fact that Nutils uses Gauss points as coupling mesh and therefore entails extrapolation in the data mapping at the top and bottom corners.
+You can mix the Nutils and FEniCS solver, if you like. Note that the error for a pure FEniCS simulation is lower than for a mixed one, because the FEniCS participants use the same coupling mesh, i.e., the mapping error becomes significantly smaller.
 
 ## Visualization
 

partitioned-heat-conduction/fenics/heat.py

@@ -170,7 +170,7 @@ def determine_gradient(V_g, u, flux):
 
 if problem is ProblemType.DIRICHLET:
     flux = Function(V_g)
-    flux.rename("Flux", "")
+    flux.rename("Heat-Flux", "")
 
 while precice.is_coupling_ongoing():
 

partitioned-heat-conduction/fenics/precice-adapter-config-D.json

@@ -3,7 +3,7 @@
   "config_file_name": "../precice-config.xml",
   "interface": {
       "coupling_mesh_name": "Dirichlet-Mesh",
-      "write_data_name": "Flux",
+      "write_data_name": "Heat-Flux",
       "read_data_name": "Temperature"
     }
 }

partitioned-heat-conduction/fenics/precice-adapter-config-N.json

@@ -4,6 +4,6 @@
   "interface": {
       "coupling_mesh_name": "Neumann-Mesh",
       "write_data_name": "Temperature",
-      "read_data_name": "Flux"
+      "read_data_name": "Heat-Flux"
     }
 }

partitioned-heat-conduction/nutils/heat.py

@@ -68,8 +68,8 @@ def main(side='Dirichlet'):
     vertex_ids = interface.set_mesh_vertices(mesh_id, vertices)
 
     # coupling data
-    write_data = "Temperature" if side == "Neumann" else "Flux"
-    read_data = "Flux" if side == "Neumann" else "Temperature"
+    write_data = "Temperature" if side == "Neumann" else "Heat-Flux"
+    read_data = "Heat-Flux" if side == "Neumann" else "Temperature"
     write_data_id = interface.get_data_id(write_data, mesh_id)
     read_data_id = interface.get_data_id(read_data, mesh_id)
 

partitioned-heat-conduction/openfoam-dirichlet/0.orig/T

@@ -0,0 +1,39 @@
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       volScalarField;
+    location    "0";
+    object      T;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+dimensions      [0 0 0 1 0 0 0];
+
+
+internalField   uniform 0;
+
+boundaryField
+{
+    interface
+    {
+        type            fixedValue;
+        value           uniform 2;
+    }
+
+    DirichletBoundary
+    {
+        type            groovyBC;
+        variables       "val=1+pow(pos().x,2)+(3*pow(pos().y,2))+1.3*time();";
+        valueExpression "val";
+        value           uniform 0;
+        evaluateDuringConstruction 1;
+    }
+
+    defaultFaces
+    {
+        type            empty;
+    }
+}
+
+// ************************************************************************* //

partitioned-heat-conduction/openfoam-dirichlet/clean.sh

@@ -0,0 +1,6 @@
+#!/bin/sh
+set -e -u
+
+. ../../tools/cleaning-tools.sh
+
+clean_openfoam .

partitioned-heat-conduction/openfoam-dirichlet/constant/transportProperties

@@ -0,0 +1,10 @@
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       dictionary;
+    location    "constant";
+    object      transportProperties;
+}
+
+DT               DT  [ 0  2 -1  0 0 0 0 ] 1;

partitioned-heat-conduction/openfoam-dirichlet/run.sh

@@ -0,0 +1,9 @@
+#!/bin/sh
+set -e -u
+
+blockMesh
+touch openfoam-dirichlet.foam
+./setInitialField.sh
+
+../../tools/run-openfoam.sh "$@"
+. ../../tools/openfoam-remove-empty-dirs.sh && openfoam_remove_empty_dirs

partitioned-heat-conduction/openfoam-dirichlet/setInitialField.sh

@@ -0,0 +1,8 @@
+#!/bin/sh
+set -e -u
+
+# Remove the old directory and copy the uninitialized field
+rm -rf ./0
+cp -r ./0.orig 0
+# Initialize the new field
+funkySetFields -keepPatches -field T -expression '1+pow(pos().x,2)+(3*pow(pos().y,2))+1.3*time()' -time '0'

partitioned-heat-conduction/openfoam-dirichlet/system/blockMeshDict

@@ -0,0 +1,58 @@
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       dictionary;
+    object      blockMeshDict;
+}
+
+vertices
+(
+
+    (0 0 0)
+    (1 0 0)
+    (1 1 0)
+    (0 1 0)
+
+    (0 0 .1)
+    (1 0 .1)
+    (1 1 .1)
+    (0 1 .1)
+);
+
+blocks
+(
+    hex (0 1 2 3 4 5 6 7) (100 100 1) simpleGrading (1 1 1)
+);
+
+edges
+(
+);
+
+boundary
+(
+
+    interface
+    {
+        type patch;
+        faces
+        (
+            (1 2 6 5)
+        );
+    }
+
+    DirichletBoundary
+    {
+        type patch;
+        faces
+        (
+            (4 7 3 0)
+            (7 6 2 3)
+            (4 0 1 5)
+        );
+    }
+);
+
+mergePatchPairs
+(
+);

partitioned-heat-conduction/openfoam-dirichlet/system/controlDict

@@ -0,0 +1,52 @@
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       dictionary;
+    location    "system";
+    object      controlDict;
+}
+
+// Make sure you build the solver located
+// in ./openfoam-solver before running
+// the case
+application     heatTransfer;
+
+libs ( "libgroovyBC.so" ) ;
+
+startFrom       startTime;
+
+startTime       0;
+
+stopAt          endTime;
+
+endTime         1;
+
+deltaT          0.1;
+
+writeControl    runTime;
+
+writeInterval   0.1;
+
+purgeWrite      0;
+
+writeFormat     ascii;
+
+writePrecision  6;
+
+writeCompression off;
+
+timeFormat      general;
+
+timePrecision   6;
+
+runTimeModifiable false;
+
+functions
+{
+    preCICE_Adapter
+    {
+        type preciceAdapterFunctionObject;
+        libs ("libpreciceAdapterFunctionObject.so");
+    }
+}

partitioned-heat-conduction/openfoam-dirichlet/system/decomposeParDict

@@ -0,0 +1,16 @@
+FoamFile {
+    version 2.0;
+    class dictionary;
+    object decomposeParDict;
+    format ascii;
+}
+
+numberOfSubdomains 2;
+
+method          simple;
+
+simpleCoeffs
+{
+    n               (2 1 1);
+    delta           0.001;
+}

partitioned-heat-conduction/openfoam-dirichlet/system/fvSchemes

@@ -0,0 +1,40 @@
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       dictionary;
+    location    "system";
+    object      fvSchemes;
+}
+
+ddtSchemes
+{
+    default         Euler;
+}
+
+gradSchemes
+{
+    default         Gauss linear;
+    grad(T)         Gauss linear;
+}
+
+divSchemes
+{
+    default         none;
+}
+
+laplacianSchemes
+{
+    default         none;
+    laplacian(DT,T) Gauss linear corrected;
+}
+
+interpolationSchemes
+{
+    default         linear;
+}
+
+snGradSchemes
+{
+    default         corrected;
+}

partitioned-heat-conduction/openfoam-dirichlet/system/fvSolution

@@ -0,0 +1,24 @@
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       dictionary;
+    location    "system";
+    object      fvSolution;
+}
+
+solvers
+{
+    T
+    {
+        solver          PCG;
+        preconditioner  DIC;
+        tolerance       1e-06;
+        relTol          0;
+    }
+}
+
+SIMPLE
+{
+    nNonOrthogonalCorrectors 2;
+}