Add FMU version to oscillator example

oscillator/README.md

@@ -1,7 +1,7 @@
 ---
 title: Oscillator
 permalink: tutorials-oscillator.html
-keywords: Python, ODE
+keywords: Python, ODE, FMI
 summary: We solve an oscillator with two masses in a partitioned fashion. Each mass is solved by an independent ODE.
 ---
 
@@ -19,15 +19,14 @@ Note that this case applies a Schwarz-type coupling method and not (like most ot
 
 ## Available solvers
 
-This tutorial is only available in Python. You need to have preCICE and the Python bindings installed on your system.
+There are two different implementations:
 
 - *Python*: An example solver using the preCICE [Python bindings](https://www.precice.org/installation-bindings-python.html). This solver also depends on the Python libraries `numpy`, which you can get from your system package manager or with `pip3 install --user <package>`.
+- *FMI*: An example solver using the [preCICE-FMI Runner](https://github.com/precice/fmi-runner). The Runner executes the FMU model `Oscillator.fmu` for computation. The compiled FMU model for Linux is part of this repository. For other systems, please recompile the model from the provided [C-files](https://github.com/precice/tutorials/tree/master/oscillator/fmi/fmu). If you want to change the model parameters or the initial conditions of the simulation, have a look inside the setting files for [MassLeft](https://github.com/precice/tutorials/tree/master/oscillator/fmi/MassLeft) and [MassRight](https://github.com/precice/tutorials/tree/master/oscillator/fmi/MassRight). For more information, please refer to [2]
 
 ## Running the Simulation
 
-### Python
-
-Open two separate terminals and start both participants by calling:
+Open two separate terminals and start both participants. For example, you can run a simulation where the left participant is computed in Python and the right participant is computed with FMI with these commands:
 
 ```bash
 cd python
@@ -37,19 +36,21 @@ cd python
 and
 
 ```bash
-cd python
+cd fmi
 ./run.sh -r
 ```
 
+Of course, you can also use the same solver for both sides.
+
 ## Post-processing
 
-Each simulation run creates two files containing position and velocity of the two masses over time. These files are called `trajectory-Mass-Left.csv` and `trajectory-Mass-Right.csv`. You can use the script `plot-trajectory.py` for post-processing. Type `python3 plot-trajectory --help` to see available options. You can, for example plot the trajectory by running
+Each simulation run creates two files containing position and velocity of the two masses over time. These files are called `trajectory-Mass-Left.csv` and `trajectory-Mass-Right.csv`. You can use the script `plot-trajectory.py` for post-processing. Type `python3 plot-trajectory --help` to see available options. You can, for example, plot the trajectory of the Python solver by running
 
 ```bash
 python3 plot-trajectory.py python/output/trajectory-Mass-Left.csv TRAJECTORY
 ```
 
-This allows you to study the effect of different time stepping schemes on energy conservation. Newmark beta conserves energy:
+The solvers allow you to study the effect of different time stepping schemes on energy conservation. Newmark beta conserves energy:
 
 ![Trajectory for Newmark beta scheme](images/tutorials-oscillator-trajectory-newmark-beta.png)
 
@@ -62,3 +63,5 @@ For details, refer to [1].
 ## References
 
 [1] V. Schüller, B. Rodenberg, B. Uekermann and H. Bungartz, A Simple Test Case for Convergence Order in Time and Energy Conservation of Black-Box Coupling Schemes, in: WCCM-APCOM2022. [URL](https://www.scipedia.com/public/Rodenberg_2022a)
+
+[2] L. Willeke, [A preCICE-FMI Runner to couple controller models to PDEs](http://dx.doi.org/10.18419/opus-13130), Master Thesis, University of Stuttgart, 2023

oscillator/fmi/MassLeft/fmi-settings.json

@@ -0,0 +1,20 @@
+{
+    "model_params": {
+        "mass.m": 1,
+        "spring_fixed.c": 39.478,
+        "spring_middle.c": 157.914
+    },
+    "initial_conditions": {
+        "mass.u": 1,
+        "mass.v": 0, 
+        "mass.a": -197.392
+    },
+    "simulation_params": {
+        "fmu_file_name": "./Oscillator.fmu",
+        "output_file_name": "./output/trajectory-Mass-Left.csv",
+        "output": ["mass.u", "mass.v"],
+        "fmu_read_data_names": ["force_in"],
+        "fmu_write_data_names": ["force_out"],
+        "fmu_instance_name": "instance_left"
+    }
+}

oscillator/fmi/MassLeft/precice-settings.json

@@ -0,0 +1,9 @@
+{
+    "coupling_params": {
+        "participant_name": "Mass-Left",
+        "config_file_name": "../precice-config.xml",
+        "mesh_name": "Mass-Left-Mesh",
+        "write_data": {"name": "Force-Left", "type": "scalar"},
+        "read_data": {"name": "Force-Right", "type": "scalar"}
+    }
+}

oscillator/fmi/MassRight/fmi-settings.json

@@ -0,0 +1,20 @@
+{
+    "model_params": {
+        "mass.m": 1,
+        "spring_fixed.c": 39.478,
+        "spring_middle.c": 157.914
+    },
+    "initial_conditions": {
+        "mass.u": 0,
+        "mass.v": 0, 
+        "mass.a": 157.914
+    },
+    "simulation_params": {
+        "fmu_file_name": "./Oscillator.fmu",
+        "output_file_name": "./output/trajectory-Mass-Right.csv",
+        "output": ["mass.u", "mass.v"],
+        "fmu_read_data_names": ["force_in"],
+        "fmu_write_data_names": ["force_out"],
+        "fmu_instance_name": "instance_right"
+    }
+}

oscillator/fmi/MassRight/precice-settings.json

@@ -0,0 +1,9 @@
+{
+    "coupling_params": {
+        "participant_name": "Mass-Right",
+        "config_file_name": "../precice-config.xml",
+        "mesh_name": "Mass-Right-Mesh",
+        "write_data": {"name": "Force-Right", "type": "scalar"},
+        "read_data": {"name": "Force-Left", "type": "scalar"}
+    }
+}

oscillator/fmi/calculate-error.py

@@ -0,0 +1,91 @@
+import numpy as np
+import pandas as pd
+import json
+import os
+import argparse
+from enum import Enum
+import sys
+
+class Participant(Enum):
+	MASS_LEFT = "Mass-Left"
+	MASS_RIGHT = "Mass-Right"
+
+parser = argparse.ArgumentParser()
+parser.add_argument("fmi_setting_file_left", help="Path to the fmi setting file for MassLeft.", type=str)
+parser.add_argument("precice_setting_file_left", help="Path to the precice setting file for MassLeft.", type=str)
+parser.add_argument("fmi_setting_file_right", help="Path to the fmi setting file for MassRight.", type=str)
+parser.add_argument("precice_setting_file_right", help="Path to the precice setting file for MassRight.", type=str)
+parser.add_argument("participant_name", help="Participant for which the error should be calculated", type=str, 
+    choices=[p.value for p in Participant])
+args = parser.parse_args()
+
+# Get input files
+fmi_file_left       = args.fmi_setting_file_left
+precice_file_left   = args.precice_setting_file_left
+fmi_file_right      = args.fmi_setting_file_right
+precice_file_right  = args.precice_setting_file_right
+participant_name    = args.participant_name
+
+# Read json files
+folder = os.path.dirname(os.path.join(os.getcwd(), os.path.dirname(sys.argv[0]), fmi_file_left))
+path = os.path.join(folder, os.path.basename(fmi_file_left))
+read_file = open(path, "r")
+fmi_data_left = json.load(read_file)
+
+folder = os.path.dirname(os.path.join(os.getcwd(), os.path.dirname(sys.argv[0]), precice_file_left))
+path = os.path.join(folder, os.path.basename(precice_file_left))
+read_file = open(path, "r")
+precice_data_left = json.load(read_file)
+
+folder = os.path.dirname(os.path.join(os.getcwd(), os.path.dirname(sys.argv[0]), fmi_file_right))
+path = os.path.join(folder, os.path.basename(fmi_file_right))
+read_file = open(path, "r")
+fmi_data_right = json.load(read_file)
+
+folder = os.path.dirname(os.path.join(os.getcwd(), os.path.dirname(sys.argv[0]), precice_file_right))
+path = os.path.join(folder, os.path.basename(precice_file_right))
+read_file = open(path, "r")
+precice_data_right = json.load(read_file)
+
+
+# Define variables
+k_1         = fmi_data_left["model_params"]["spring_fixed.c"]
+k_2         = fmi_data_right["model_params"]["spring_fixed.c"]
+u0_1        = fmi_data_left["initial_conditions"]["mass.u"]
+u0_2        = fmi_data_right["initial_conditions"]["mass.u"]
+k_12_left   = fmi_data_left["model_params"]["spring_middle.c"]
+k_12_right  = fmi_data_right["model_params"]["spring_middle.c"]
+
+if k_12_left == k_12_right:
+	k_12 = k_12_left
+else:
+	raise Exception("k_12 has to be equal in both participants. Please adjust input values.")
+
+
+# Define analytical solution and read computed results
+K = np.array([[k_1 + k_12, -k_12], [-k_12, k_2 + k_12]])
+eigenvalues, eigenvectors = np.linalg.eig(K)
+omega = np.sqrt(eigenvalues)
+A, B = eigenvectors
+c = np.linalg.solve(eigenvectors, [u0_1, u0_2])
+
+if participant_name == Participant.MASS_LEFT.value:
+    filename = fmi_data_left["simulation_params"]["output_file_name"]
+    df = pd.read_csv(filename, delimiter=',')	
+
+    def u_analytical(t): return c[0] * A[0] * np.cos(omega[0] * t) + c[1] * A[1] * np.cos(omega[1] * t)
+
+elif participant_name == Participant.MASS_RIGHT.value:
+    filename = fmi_data_right["simulation_params"]["output_file_name"]
+    df = pd.read_csv(filename, delimiter=',')
+
+    def u_analytical(t): return c[0] * B[0] * np.cos(omega[0] * t) + c[1] * B[1] * np.cos(omega[1] * t)
+
+times       = df.iloc[:,0]
+positions   = df.iloc[:,1]
+
+
+# Calculate error
+error = np.max(abs(u_analytical(np.array(times)) - np.array(positions)))
+print("Error w.r.t analytical solution:")
+print(f"{error}")

oscillator/fmi/clean.sh

@@ -0,0 +1,8 @@
+#!/bin/sh
+set -e -u
+
+. ../../tools/cleaning-tools.sh
+
+rm -rfv ./output/
+
+clean_precice_logs .