From b2a50795cd43ba89c35f4a05ad7103af07110149 Mon Sep 17 00:00:00 2001
From: Niklas <niklas.vinnitchenko@tum.de>
Date: Tue, 11 Mar 2025 14:35:52 +0100
Subject: [PATCH 1/7] Initial comit

---
 fenicsxprecice/config.py                      |  53 ++--
 fenicsxprecice/fenicsxprecice.py              | 269 ++++++++++--------
 .../precice-config.xml                        |  38 +--
 .../solver-fenicsx/heat.py                    |  20 +-
 .../precice-adapter-config-D.json             |   5 +-
 .../precice-adapter-config-N.json             |   5 +-
 6 files changed, 218 insertions(+), 172 deletions(-)

diff --git a/fenicsxprecice/config.py b/fenicsxprecice/config.py
index ba7c938..b1fbe7c 100644
--- a/fenicsxprecice/config.py
+++ b/fenicsxprecice/config.py
@@ -18,9 +18,9 @@ def __init__(self, adapter_config_filename):
 
         self._config_file_name = None
         self._participant_name = None
-        self._coupling_mesh_name = None
-        self._read_data_name = None
-        self._write_data_name = None
+        self._coupling_mesh_names = None
+        self._read_data_names = {}
+        self._write_data_names = {}
 
         self.read_json(adapter_config_filename)
 
@@ -39,17 +39,18 @@ def read_json(self, adapter_config_filename):
         data = json.load(read_file)
         self._config_file_name = os.path.join(folder, data["config_file_name"])
         self._participant_name = data["participant_name"]
-        self._coupling_mesh_name = data["interface"]["coupling_mesh_name"]
+        self._coupling_mesh_names = list(data["interfaces"].keys())
 
-        try:
-            self._write_data_name = data["interface"]["write_data_name"]
-        except KeyError:
-            self._write_data_name = None  # not required for one-way coupling, if this participant reads data
+        for mesh_name in self._coupling_mesh_names:    
+            try:
+                self._write_data_names[mesh_name] = data["interfaces"][mesh_name]["write_data_name"]
+            except KeyError:
+                self._write_data_names[mesh_name] = None  # not required for one-way coupling, if this participant reads data
 
-        try:
-            self._read_data_name = data["interface"]["read_data_name"]
-        except KeyError:
-            self._read_data_name = None  # not required for one-way coupling, if this participant writes data
+            try:
+                self._read_data_names[mesh_name] = data["interfaces"][mesh_name]["read_data_name"]
+            except KeyError:
+                self._read_data_names[mesh_name] = None  # not required for one-way coupling, if this participant writes data
 
         read_file.close()
 
@@ -59,11 +60,27 @@ def get_config_file_name(self):
     def get_participant_name(self):
         return self._participant_name
 
-    def get_coupling_mesh_name(self):
-        return self._coupling_mesh_name
+    def get_coupling_mesh_names(self):
+        return self._coupling_mesh_names
 
-    def get_read_data_name(self):
-        return self._read_data_name
 
-    def get_write_data_name(self):
-        return self._write_data_name
+    def get_read_data_name(self, mesh_name):
+        """
+        Parameters
+        ----------
+        mesh_name : fenicsxprecice.Meshes member
+            Member of the enum fenicsxprecoce.Meshes
+    
+        """
+        return self._read_data_names[mesh_name.name]
+
+    def get_write_data_name(self, mesh_name):
+        """
+        Parameters
+        ----------
+        mesh_name : fenicsxprecice.Meshes member
+            Member of the enum fenicsxprecoce.Meshes
+    
+        """
+        return self._write_data_names[mesh_name.name]
+    
\ No newline at end of file
diff --git a/fenicsxprecice/fenicsxprecice.py b/fenicsxprecice/fenicsxprecice.py
index ef62bf9..0d98f38 100644
--- a/fenicsxprecice/fenicsxprecice.py
+++ b/fenicsxprecice/fenicsxprecice.py
@@ -11,6 +11,7 @@
 from .solverstate import SolverState
 from dolfinx import fem
 import copy
+from enum import Enum
 
 logger = logging.getLogger(__name__)
 logger.setLevel(level=logging.INFO)
@@ -48,6 +49,8 @@ def __init__(self, mpi_comm, adapter_config_filename='precice-adapter-config.jso
         """
 
         self._config = Config(adapter_config_filename)
+        
+        self.Meshes = Enum('Meshes', self._config.get_coupling_mesh_names())
 
         # Setup up MPI communicator
         self._comm = mpi_comm
@@ -60,17 +63,17 @@ def __init__(self, mpi_comm, adapter_config_filename='precice-adapter-config.jso
         )
 
         # FEniCSx related quantities
-        self._read_function_space = None  # initialized later
-        self._write_function_space = None  # initialized later
+        self._read_function_spaces = {}  # initialized later
+        self._write_function_spaces = {}  # initialized later
         self._dofmap = None  # initialized later using function space provided by user
 
         # coupling mesh related quantities
-        self._fenicsx_vertices = Vertices()
-        self._precice_vertex_ids = None  # initialized later
+        self._fenicsx_vertices = {}
+        self._precice_vertex_ids = {}  # initialized later
 
         # read data related quantities (read data is read from preCICE and applied in FEniCSx)
-        self._read_function_type = None  # stores whether read function is scalar or vector valued
-        self._write_function_type = None  # stores whether write function is scalar or vector valued
+        self._read_function_types = {}  # stores whether read function is scalar or vector valued
+        self._write_function_types = {}  # stores whether write function is scalar or vector valued
 
         # Interpolation strategy
         self._my_expression = SegregatedRBFInterpolationExpression
@@ -82,16 +85,21 @@ def __init__(self, mpi_comm, adapter_config_filename='precice-adapter-config.jso
         self._first_advance_done = False
 
         # Determine type of coupling in initialization
-        self._coupling_type = None
+        self._coupling_types = {}
 
         # Problem dimension in FEniCSx
         self._fenicsx_dims = None
         self._empty_rank = True
 
-    def create_coupling_expression(self):
+    def create_coupling_expression(self, mesh_name):
         """
         Creates a FEniCSx Expression in the form of an object of class GeneralInterpolationExpression or
         ExactInterpolationExpression. The adapter will hold this object till the coupling is on going.
+        
+        Parameters
+        ----------
+        mesh_name: member of fenicsxprecice.Meshes
+            Specifies for which field a coupling expression shall be created
 
         Returns
         -------
@@ -99,10 +107,10 @@ def create_coupling_expression(self):
             Reference to object of class GeneralInterpolationExpression or ExactInterpolationExpression.
         """
 
-        if not (self._read_function_type is FunctionType.SCALAR or self._read_function_type is FunctionType.VECTOR):
+        if not (self._read_function_types[mesh_name] is FunctionType.SCALAR or self._read_function_types[mesh_name] is FunctionType.VECTOR):
             raise Exception("No valid read_function is provided in initialization. Cannot create coupling expression")
 
-        coupling_expression = self._my_expression(self._read_function_space, self._read_function_type)
+        coupling_expression = self._my_expression(self._read_function_spaces[mesh_name], self._read_function_types[mesh_name])
 
         return coupling_expression
 
@@ -126,7 +134,7 @@ def update_coupling_expression(self, coupling_expression, data):
     def get_point_sources(self, data):
         raise Exception("PointSources are not implemented for the FEniCSx adapter.")
 
-    def read_data(self, dt):
+    def read_data(self, dt, mesh_name):
         """
         Read data from preCICE. Data is generated depending on the type of the read function (Scalar or Vector).
         For a scalar read function the data is a numpy array with shape (N) where N = number of coupling vertices
@@ -142,26 +150,26 @@ def read_data(self, dt):
             The coupling data. A dictionary containing nodal data with vertex coordinates as key and associated data as
             value.
         """
-        assert (self._coupling_type is CouplingMode.UNI_DIRECTIONAL_READ_COUPLING or
+        assert (self._coupling_types[mesh_name] is CouplingMode.UNI_DIRECTIONAL_READ_COUPLING or
                 CouplingMode.BI_DIRECTIONAL_COUPLING)
 
         read_data = None
 
         if not self._empty_rank:
             read_data = self._participant.read_data(
-                self._config.get_coupling_mesh_name(),
-                self._config.get_read_data_name(),
-                self._precice_vertex_ids,
+                mesh_name.name,
+                self._config.get_read_data_name(mesh_name),
+                self._precice_vertex_ids[mesh_name],
                 dt
             )
-            read_data = {tuple(key): value for key, value in zip(self._fenicsx_vertices.get_coordinates(), read_data)}
+            read_data = {tuple(key): value for key, value in zip(self._fenicsx_vertices[mesh_name].get_coordinates(), read_data)}
 
         else:
             pass
 
         return copy.deepcopy(read_data)
 
-    def write_data(self, write_function):
+    def write_data(self, write_function, mesh_name):
         """
         Writes data to preCICE. Depending on the dimensions of the simulation (2D-3D Coupling, 2D-2D coupling or
         Scalar/Vector write function) write_data is first converted into a format needed for preCICE.
@@ -172,40 +180,42 @@ def write_data(self, write_function):
             A FEniCSx function consisting of the data which this participant will write to preCICE in every time step.
         """
 
-        assert (self._coupling_type is CouplingMode.UNI_DIRECTIONAL_WRITE_COUPLING or
+        assert (self._coupling_types[mesh_name] is CouplingMode.UNI_DIRECTIONAL_WRITE_COUPLING or
                 CouplingMode.BI_DIRECTIONAL_COUPLING)
 
         w_func = write_function.copy()
 
         # Check that the function provided lives on the same function space provided during initialization
-        assert (self._write_function_type == determine_function_type(w_func))
-        assert (write_function.function_space == self._write_function_space)
+        assert (self._write_function_types[mesh_name] == determine_function_type(w_func))
+        assert (write_function.function_space == self._write_function_spaces[mesh_name])
 
         write_function_type = determine_function_type(write_function)
         assert (write_function_type in list(FunctionType))
-        write_data = convert_fenicsx_to_precice(write_function, self._fenicsx_vertices.get_coordinates())
+        write_data = convert_fenicsx_to_precice(write_function, self._fenicsx_vertices[mesh_name].get_coordinates())
         self._participant.write_data(
-            self._config.get_coupling_mesh_name(),
-            self._config.get_write_data_name(),
-            self._precice_vertex_ids,
+            mesh_name.name,
+            self._config.get_write_data_name(mesh_name),
+            self._precice_vertex_ids[mesh_name],
             write_data
         )
 
-    def initialize(self, coupling_subdomain, read_function_space=None, write_object=None):
+    def initialize(self, coupling_meshes):
         """
         Initializes the coupling and sets up the mesh where coupling happens in preCICE.
 
         Parameters
         ----------
-        coupling_subdomain : Object of class dolfinx.cpp.mesh.SubDomain
-            SubDomain of mesh which is the physical coupling boundary.
-        read_function_space : Object of class dolfinx.functions.functionspace.FunctionSpace
-            Function space on which the read function lives. If not provided then the adapter assumes that this
-            participant is a write-only participant.
-        write_object : Object of class dolfinx.functions.functionspace.FunctionSpace / dolfinx.functions.function.Function
-            Function space on which the write function lives or FEniCSx function related to the quantity to be written
-            by FEniCSx during each coupling iteration. If not provided then the adapter assumes that this participant is
-            a read-only participant.
+        coupling_meshes: A dictionary.
+            The key defines the name of the coupling mesh, the value is a list of the form [coupling_subdomain, read_function_space, write_object].\n
+            coupling_subdomain: Object of class dolfinx.cpp.mesh.SubDomain.
+                It is the subdomain of the mesh which is the physical coupling boundary.\n
+            read_function_space: Object of class dolfinx.functions.functionspace.FunctionSpace.
+                The function space on which the read function lives. If set to None then the adapter assumes that this
+                participant is a write-only participant.\n
+            write_object: Object of class dolfinx.functions.functionspace.FunctionSpace / dolfinx.functions.function.Function.
+                The Function space on which the write function lives or FEniCSx function related to the quantity to be written
+                by FEniCSx during each coupling iteration. If set to None then the adapter assumes that this participant is
+                a read-only participant.
 
         Returns
         -------
@@ -213,97 +223,106 @@ def initialize(self, coupling_subdomain, read_function_space=None, write_object=
             Recommended time step value from preCICE.
         """
 
-        write_function_space, write_function = None, None
-        if isinstance(write_object, fem.Function):  # precice.initialize_data() will be called using this Function
-            write_function_space = write_object.function_space
-            write_function = write_object
-        elif isinstance(write_object, fem.FunctionSpace):  # preCICE will use default zero values for initialization.
-            write_function_space = write_object
-            write_function = None
-        elif write_object is None:
-            pass
-        else:
-            raise Exception("Given write object is neither of type dolfinx.functions.function.Function or "
-                            "dolfinx.functions.functionspace.FunctionSpace")
-
-        if isinstance(read_function_space, fem.FunctionSpace):
-            pass
-        elif read_function_space is None:
-            pass
-        else:
-            raise Exception("Given read_function_space is not of type dolfinx.functions.functionspace.FunctionSpace")
-
-        if read_function_space is None and write_function_space:
-            self._coupling_type = CouplingMode.UNI_DIRECTIONAL_WRITE_COUPLING
-            assert (self._config.get_write_data_name())
-            print("Participant {} is write-only participant".format(self._config.get_participant_name()))
-            function_space = write_function_space
-        elif read_function_space and write_function_space is None:
-            self._coupling_type = CouplingMode.UNI_DIRECTIONAL_READ_COUPLING
-            assert (self._config.get_read_data_name())
-            print("Participant {} is read-only participant".format(self._config.get_participant_name()))
-            function_space = read_function_space
-        elif read_function_space and write_function_space:
-            self._coupling_type = CouplingMode.BI_DIRECTIONAL_COUPLING
-            assert (self._config.get_read_data_name() and self._config.get_write_data_name())
-            function_space = read_function_space
-        elif read_function_space is None and write_function_space is None:
-            raise Exception("Neither read_function_space nor write_function_space is provided. Please provide a "
-                            "write_object if this participant is used in one-way coupling and only writes data. "
-                            "Please provide a read_function_space if this participant is used in one-way coupling and "
-                            "only reads data. If two-way coupling is implemented then both read_function_space"
-                            " and write_object need to be provided.")
-        else:
-            raise Exception("Incorrect read and write function space combination provided. Please check input "
-                            "parameters in initialization")
-
-        if self._coupling_type is CouplingMode.UNI_DIRECTIONAL_READ_COUPLING or \
-                self._coupling_type is CouplingMode.BI_DIRECTIONAL_COUPLING:
-            self._read_function_type = determine_function_type(read_function_space)
-            self._read_function_space = read_function_space
-
-        if self._coupling_type is CouplingMode.UNI_DIRECTIONAL_WRITE_COUPLING or \
-                self._coupling_type is CouplingMode.BI_DIRECTIONAL_COUPLING:
-            # Ensure that function spaces of read and write functions are defined using the same mesh
-            self._write_function_type = determine_function_type(write_function_space)
-            self._write_function_space = write_function_space
-
-        # Set vertices on the coupling subdomain for this rank
-        self._fenicsx_dims = function_space.mesh.geometry.dim
-        # returns 3d coordinates (necessary later for writing the data!)
-        ids, coords = get_fenicsx_vertices(function_space, coupling_subdomain, self._fenicsx_dims)
-        # this isnt a problem in update_coupling_expression, because in this function
-        # , the two first dimensions are extracted. Exactly what we want!
-        self._fenicsx_vertices.set_ids(ids)
-        self._fenicsx_vertices.set_coordinates(coords)
-
-        # Set up mesh in preCICE
-        self._precice_vertex_ids = self._participant.set_mesh_vertices(
-            self._config.get_coupling_mesh_name(), self._fenicsx_vertices.get_coordinates()[
-                :, :2])  # give preCICE only 2D coordinates
-
-        if self._fenicsx_vertices.get_ids().size > 0:
-            self._empty_rank = False
-        else:
-            print("Rank {} has no part of coupling boundary.".format(self._comm.Get_rank()))
-
-        # Ensure that function spaces of read and write functions use the same mesh
-        if self._coupling_type is CouplingMode.BI_DIRECTIONAL_COUPLING:
-            assert (self._read_function_space.mesh is write_function_space.mesh
-                    ), "read_function_space and write_object need to be defined using the same mesh"
-
-        if self._fenicsx_dims != 2:
-            raise Exception("Currently the fenicsx-adapter only supports 2D cases")
-
-        if self._fenicsx_dims != self._participant.get_mesh_dimensions(self._config.get_coupling_mesh_name()):
-            raise Exception("Dimension of preCICE setup and FEniCSx do not match")
-
-        if self._participant.requires_initial_data():
-            if not write_function:
-                raise Exception(
-                    "preCICE requires you to write initial data. Please provide a write_function to initialize(...)")
-            self.write_data(write_function)
-
+        for c_mesh in coupling_meshes:
+            write_function_space, write_function = None, None
+            # [2]: write object!
+            write_object = coupling_meshes[c_mesh][2]
+            if isinstance(write_object, fem.Function):  # precice.initialize_data() will be called using this Function
+                write_function_space = write_object.function_space
+                write_function = write_object
+            elif isinstance(write_object, fem.FunctionSpace):  # preCICE will use default zero values for initialization.
+                write_function_space = write_object
+                write_function = None
+            elif write_object is None:
+                pass
+            else:
+                raise Exception(f"Given write object of {c_mesh} is neither of type dolfinx.functions.function.Function or "
+                                "dolfinx.functions.functionspace.FunctionSpace")
+
+            # [1]: read_function_space
+            read_function_space = coupling_meshes[c_mesh][1]
+            if isinstance(read_function_space, fem.FunctionSpace):
+                pass
+            elif read_function_space is None:
+                pass
+            else:
+                raise Exception(f"Given read_function_space of {c_mesh} is not of type dolfinx.functions.functionspace.FunctionSpace")
+
+            if read_function_space is None and write_function_space:
+                self._coupling_types[c_mesh] = CouplingMode.UNI_DIRECTIONAL_WRITE_COUPLING
+                assert (self._config.get_write_data_name(c_mesh))
+                print("Participant {} is write-only participant".format(self._config.get_participant_name()))
+                function_space = write_function_space
+            elif read_function_space and write_function_space is None:
+                self._coupling_types[c_mesh] = CouplingMode.UNI_DIRECTIONAL_READ_COUPLING
+                assert (self._config.get_read_data_name(c_mesh))
+                print("Participant {} is read-only participant".format(self._config.get_participant_name()))
+                function_space = read_function_space
+            elif read_function_space and write_function_space:
+                self._coupling_types[c_mesh] = CouplingMode.BI_DIRECTIONAL_COUPLING
+                assert (self._config.get_read_data_name(c_mesh) and self._config.get_write_data_name(c_mesh))
+                function_space = read_function_space
+            elif read_function_space is None and write_function_space is None:
+                raise Exception(f"Neither read_function_space nor write_function_space for {c_mesh} is provided. Please provide a "
+                                "write_object if this participant is used in one-way coupling and only writes data. "
+                                "Please provide a read_function_space if this participant is used in one-way coupling and "
+                                "only reads data. If two-way coupling is implemented then both read_function_space"
+                                " and write_object need to be provided.")
+            else:
+                raise Exception(f"Incorrect read and write function space combination provided for {c_mesh}. Please check input "
+                                "parameters in initialization")
+
+            coupling_type = self._coupling_types[c_mesh]
+            if coupling_type is CouplingMode.UNI_DIRECTIONAL_READ_COUPLING or \
+                    coupling_type is CouplingMode.BI_DIRECTIONAL_COUPLING:
+                self._read_function_types[c_mesh] = determine_function_type(read_function_space)
+                self._read_function_spaces[c_mesh] = read_function_space
+
+            if coupling_type is CouplingMode.UNI_DIRECTIONAL_WRITE_COUPLING or \
+                    coupling_type is CouplingMode.BI_DIRECTIONAL_COUPLING:
+                # Ensure that function spaces of read and write functions are defined using the same mesh
+                self._write_function_types[c_mesh] = determine_function_type(write_function_space)
+                self._write_function_spaces[c_mesh] = write_function_space
+
+            # Set vertices on the coupling subdomain for this rank (assumed to be equal across all fields that will be coupled)
+            self._fenicsx_dims = function_space.mesh.geometry.dim
+            # returns 3d coordinates (necessary later for writing the data!)
+            # coupling subdomain is at [0] !
+            ids, coords = get_fenicsx_vertices(function_space, coupling_meshes[c_mesh][0], self._fenicsx_dims)
+            # this isnt a problem in update_coupling_expression, because in this function
+            # , the two first dimensions are extracted. Exactly what we want!
+            self._fenicsx_vertices[c_mesh] = Vertices()
+            self._fenicsx_vertices[c_mesh].set_ids(ids)
+            self._fenicsx_vertices[c_mesh].set_coordinates(coords)
+
+            # Set up mesh in preCICE
+            self._precice_vertex_ids[c_mesh] = self._participant.set_mesh_vertices(
+                c_mesh.name, self._fenicsx_vertices[c_mesh].get_coordinates()[
+                    :, :2])  # give preCICE only 2D coordinates
+
+            if self._fenicsx_vertices[c_mesh].get_ids().size > 0:
+                self._empty_rank = False
+            else:
+                print("Rank {} has no part of coupling boundary.".format(self._comm.Get_rank()))
+
+            # Ensure that function spaces of read and write functions use the same mesh
+            if coupling_type is CouplingMode.BI_DIRECTIONAL_COUPLING:
+                assert (self._read_function_spaces[c_mesh].mesh is write_function_space.mesh
+                        ), "read_function_space and write_object need to be defined using the same mesh"
+
+            if self._fenicsx_dims != 2:
+                raise Exception("Currently the fenicsx-adapter only supports 2D cases")
+
+            if self._fenicsx_dims != self._participant.get_mesh_dimensions(c_mesh.name):
+                raise Exception("Dimension of preCICE setup and FEniCSx do not match")
+
+            if self._participant.requires_initial_data():
+                if not write_function:
+                    raise Exception(
+                        "preCICE requires you to write initial data. Please provide a write_function to initialize(...)")
+                self.write_data(write_function, c_mesh)
+
+        
         self._participant.initialize()
 
     def store_checkpoint(self, payload, t, n):
diff --git a/tutorials/partitioned-heat-conduction/precice-config.xml b/tutorials/partitioned-heat-conduction/precice-config.xml
index baab5d5..03d8261 100644
--- a/tutorials/partitioned-heat-conduction/precice-config.xml
+++ b/tutorials/partitioned-heat-conduction/precice-config.xml
@@ -10,37 +10,37 @@
   <data:scalar name="Temperature" waveform-degree="1" />
   <data:scalar name="Heat-Flux" waveform-degree="1" />
 
-  <mesh name="Dirichlet-Mesh" dimensions="2">
+  <mesh name="Dirichlet_Mesh" dimensions="2">
     <use-data name="Temperature" />
     <use-data name="Heat-Flux" />
   </mesh>
 
-  <mesh name="Neumann-Mesh" dimensions="2">
+  <mesh name="Neumann_Mesh" dimensions="2">
     <use-data name="Temperature" />
     <use-data name="Heat-Flux" />
   </mesh>
 
   <participant name="Dirichlet">
-    <provide-mesh name="Dirichlet-Mesh" />
-    <receive-mesh name="Neumann-Mesh" from="Neumann" />
-    <write-data name="Heat-Flux" mesh="Dirichlet-Mesh" />
-    <read-data name="Temperature" mesh="Dirichlet-Mesh" />
+    <provide-mesh name="Dirichlet_Mesh" />
+    <receive-mesh name="Neumann_Mesh" from="Neumann" />
+    <write-data name="Heat-Flux" mesh="Dirichlet_Mesh" />
+    <read-data name="Temperature" mesh="Dirichlet_Mesh" />
     <mapping:nearest-neighbor
       direction="read"
-      from="Neumann-Mesh"
-      to="Dirichlet-Mesh"
+      from="Neumann_Mesh"
+      to="Dirichlet_Mesh"
       constraint="consistent" />
   </participant>
 
   <participant name="Neumann">
-    <provide-mesh name="Neumann-Mesh" />
-    <receive-mesh name="Dirichlet-Mesh" from="Dirichlet" />
-    <write-data name="Temperature" mesh="Neumann-Mesh" />
-    <read-data name="Heat-Flux" mesh="Neumann-Mesh" />
+    <provide-mesh name="Neumann_Mesh" />
+    <receive-mesh name="Dirichlet_Mesh" from="Dirichlet" />
+    <write-data name="Temperature" mesh="Neumann_Mesh" />
+    <read-data name="Heat-Flux" mesh="Neumann_Mesh" />
     <mapping:nearest-neighbor
       direction="read"
-      from="Dirichlet-Mesh"
-      to="Neumann-Mesh"
+      from="Dirichlet_Mesh"
+      to="Neumann_Mesh"
       constraint="consistent" />
   </participant>
 
@@ -53,22 +53,22 @@
     <max-iterations value="100" />
     <exchange
       data="Heat-Flux"
-      mesh="Dirichlet-Mesh"
+      mesh="Dirichlet_Mesh"
       from="Dirichlet"
       to="Neumann"
       initialize="true"
       substeps="true" />
     <exchange
       data="Temperature"
-      mesh="Neumann-Mesh"
+      mesh="Neumann_Mesh"
       from="Neumann"
       to="Dirichlet"
       initialize="true"
       substeps="true" />
-    <relative-convergence-measure data="Heat-Flux" mesh="Dirichlet-Mesh" limit="1e-10" />
-    <relative-convergence-measure data="Temperature" mesh="Neumann-Mesh" limit="1e-10" />
+    <relative-convergence-measure data="Heat-Flux" mesh="Dirichlet_Mesh" limit="1e-10" />
+    <relative-convergence-measure data="Temperature" mesh="Neumann_Mesh" limit="1e-10" />
     <!-- <acceleration:IQN-ILS>
-      <data name="Temperature" mesh="Neumann-Mesh" />
+      <data name="Temperature" mesh="Neumann_Mesh" />
       <initial-relaxation value="0.1" />
       <max-used-iterations value="10" />
       <time-windows-reused value="5" />
diff --git a/tutorials/partitioned-heat-conduction/solver-fenicsx/heat.py b/tutorials/partitioned-heat-conduction/solver-fenicsx/heat.py
index 2ef955c..a7a08eb 100644
--- a/tutorials/partitioned-heat-conduction/solver-fenicsx/heat.py
+++ b/tutorials/partitioned-heat-conduction/solver-fenicsx/heat.py
@@ -128,9 +128,9 @@ def __call__(self, x):
     precice = Adapter(adapter_config_filename="precice-adapter-config-N.json", mpi_comm=MPI.COMM_WORLD)
 
 if problem is ProblemType.DIRICHLET:
-    precice.initialize(coupling_boundary, read_function_space=V, write_object=f_N)
+    precice.initialize({precice.Meshes.Dirichlet_Mesh: [coupling_boundary, V, f_N]})
 elif problem is ProblemType.NEUMANN:
-    precice.initialize(coupling_boundary, read_function_space=W, write_object=u_D)
+    precice.initialize({precice.Meshes.Neumann_Mesh: [coupling_boundary, W, u_D]})
 
 # get precice's dt
 precice_dt = precice.get_max_time_step_size()
@@ -144,7 +144,11 @@ def __call__(self, x):
 u, v = ufl.TrialFunction(V), ufl.TestFunction(V)
 F = u * v * ufl.dx + dt * ufl.dot(ufl.grad(u), ufl.grad(v)) * ufl.dx - (u_n + dt * f) * v * ufl.dx
 # create a coupling expression for the coupling_boundary and modify variational problem accordingly
-coupling_expression = precice.create_coupling_expression()
+if problem is ProblemType.DIRICHLET:
+    coupling_expression = precice.create_coupling_expression(precice.Meshes.Dirichlet_Mesh)
+else:
+    coupling_expression = precice.create_coupling_expression(precice.Meshes.Neumann_Mesh)
+
 if problem is ProblemType.DIRICHLET:
     # modify Dirichlet boundary condition on coupling interface
     bc_coup = fem.dirichletbc(coupling_expression, dofs_coupling)
@@ -195,7 +199,11 @@ def __call__(self, x):
     precice_dt = precice.get_max_time_step_size()
     dt = np.min([fenics_dt, precice_dt])
 
-    read_data = precice.read_data(dt)
+    if problem is ProblemType.DIRICHLET:
+        read_data = precice.read_data(dt, precice.Meshes.Dirichlet_Mesh)
+    else:
+        read_data = precice.read_data(dt, precice.Meshes.Neumann_Mesh)
+
     # Update the right hand side reusing the initial vector
     with b.localForm() as loc_b:
         loc_b.set(0)
@@ -218,10 +226,10 @@ def __call__(self, x):
         flux = determine_gradient(V_g, uh)
         flux_x = fem.Function(W)
         flux_x.interpolate(flux.sub(0))
-        precice.write_data(flux_x)
+        precice.write_data(flux_x, precice.Meshes.Dirichlet_Mesh)
     elif problem is ProblemType.NEUMANN:
         # Neumann problem reads flux and writes temperature on boundary to Dirichlet problem
-        precice.write_data(uh)
+        precice.write_data(uh, precice.Meshes.Neumann_Mesh)
 
     precice.advance(dt)
     precice_dt = precice.get_max_time_step_size()
diff --git a/tutorials/partitioned-heat-conduction/solver-fenicsx/precice-adapter-config-D.json b/tutorials/partitioned-heat-conduction/solver-fenicsx/precice-adapter-config-D.json
index 7989ba3..b78bec9 100644
--- a/tutorials/partitioned-heat-conduction/solver-fenicsx/precice-adapter-config-D.json
+++ b/tutorials/partitioned-heat-conduction/solver-fenicsx/precice-adapter-config-D.json
@@ -1,9 +1,10 @@
 {
   "participant_name": "Dirichlet",
   "config_file_name": "../precice-config.xml",
-  "interface": {
-      "coupling_mesh_name": "Dirichlet-Mesh",
+  "interfaces": {
+    "Dirichlet_Mesh": {
       "write_data_name": "Heat-Flux",
       "read_data_name": "Temperature"
+      }
     }
 }
diff --git a/tutorials/partitioned-heat-conduction/solver-fenicsx/precice-adapter-config-N.json b/tutorials/partitioned-heat-conduction/solver-fenicsx/precice-adapter-config-N.json
index 2afd52c..eec26f4 100644
--- a/tutorials/partitioned-heat-conduction/solver-fenicsx/precice-adapter-config-N.json
+++ b/tutorials/partitioned-heat-conduction/solver-fenicsx/precice-adapter-config-N.json
@@ -1,9 +1,10 @@
 {
   "participant_name": "Neumann",
   "config_file_name": "../precice-config.xml",
-  "interface": {
-      "coupling_mesh_name": "Neumann-Mesh",
+  "interfaces": {
+    "Neumann_Mesh": {
       "write_data_name": "Temperature",
       "read_data_name": "Heat-Flux"
+      }
     }
 }

From e219d9e902d6d4d5fac8aa53f92fe833cabac237 Mon Sep 17 00:00:00 2001
From: Niklas <niklas.vinnitchenko@tum.de>
Date: Tue, 11 Mar 2025 15:04:02 +0100
Subject: [PATCH 2/7] Preparation of proof of concept example

---
 tutorials/multiple-fields/left.py             |  47 +++++++
 .../precice-adapter-config-L.json             |  14 ++
 .../precice-adapter-config-R.json             |  14 ++
 tutorials/multiple-fields/precice-config.xml  | 127 ++++++++++++++++++
 tutorials/multiple-fields/right.py            |   0
 5 files changed, 202 insertions(+)
 create mode 100644 tutorials/multiple-fields/left.py
 create mode 100644 tutorials/multiple-fields/precice-adapter-config-L.json
 create mode 100644 tutorials/multiple-fields/precice-adapter-config-R.json
 create mode 100644 tutorials/multiple-fields/precice-config.xml
 create mode 100644 tutorials/multiple-fields/right.py

diff --git a/tutorials/multiple-fields/left.py b/tutorials/multiple-fields/left.py
new file mode 100644
index 0000000..1684ee1
--- /dev/null
+++ b/tutorials/multiple-fields/left.py
@@ -0,0 +1,47 @@
+from mpi4py import MPI
+from dolfinx import mesh
+from dolfinx.fem import functionspace
+from dolfinx import fem
+import numpy
+import ufl
+from dolfinx import default_scalar_type
+from dolfinx.fem.petsc import LinearProblem
+
+domain = mesh.create_unit_square(MPI.COMM_WORLD, 8, 8, mesh.CellType.quadrilateral)
+V = functionspace(domain, ("Lagrange", 1))
+uD = fem.Function(V)
+uD.interpolate(lambda x: 1 + x[0]**2 + 2 * x[1]**2)
+
+# Create facet to cell connectivity required to determine boundary facets
+tdim = domain.topology.dim
+fdim = tdim - 1
+domain.topology.create_connectivity(fdim, tdim)
+boundary_facets = mesh.exterior_facet_indices(domain.topology)
+
+boundary_dofs = fem.locate_dofs_topological(V, fdim, boundary_facets)
+bc = fem.dirichletbc(uD, boundary_dofs)
+
+u = ufl.TrialFunction(V)
+v = ufl.TestFunction(V)
+
+f = fem.Constant(domain, default_scalar_type(-6))
+
+a = ufl.dot(ufl.grad(u), ufl.grad(v)) * ufl.dx
+L = f * v * ufl.dx
+
+problem = LinearProblem(a, L, bcs=[bc], petsc_options={"ksp_type": "preonly", "pc_type": "lu"})
+uh = problem.solve()
+
+V2 = fem.functionspace(domain, ("Lagrange", 2))
+uex = fem.Function(V2)
+uex.interpolate(lambda x: 1 + x[0]**2 + 2 * x[1]**2)
+
+L2_error = fem.form(ufl.inner(uh - uex, uh - uex) * ufl.dx)
+error_local = fem.assemble_scalar(L2_error)
+error_L2 = numpy.sqrt(domain.comm.allreduce(error_local, op=MPI.SUM))
+
+error_max = numpy.max(numpy.abs(uD.x.array-uh.x.array))
+# Only print the error on one process
+if domain.comm.rank == 0:
+    print(f"Error_L2 : {error_L2:.2e}")
+    print(f"Error_max : {error_max:.2e}")
\ No newline at end of file
diff --git a/tutorials/multiple-fields/precice-adapter-config-L.json b/tutorials/multiple-fields/precice-adapter-config-L.json
new file mode 100644
index 0000000..cb7a629
--- /dev/null
+++ b/tutorials/multiple-fields/precice-adapter-config-L.json
@@ -0,0 +1,14 @@
+{
+  "participant_name": "Left",
+  "config_file_name": "precice-config.xml",
+  "interfaces": {
+    "Left1": {
+      "write_data_name": "LeftOut",
+      "read_data_name": "LeftIn"
+      },
+    "Left2": {
+      "write_data_name": "RightIn",
+      "read_data_name": "RightOut"
+    }
+    }
+}
diff --git a/tutorials/multiple-fields/precice-adapter-config-R.json b/tutorials/multiple-fields/precice-adapter-config-R.json
new file mode 100644
index 0000000..b1e1a35
--- /dev/null
+++ b/tutorials/multiple-fields/precice-adapter-config-R.json
@@ -0,0 +1,14 @@
+{
+  "participant_name": "Right",
+  "config_file_name": "precice-config.xml",
+  "interfaces": {
+    "Right1": {
+      "write_data_name": "LeftIn",
+      "read_data_name": "LeftOut"
+      },
+    "Right2": {
+      "write_data_name": "RightOut",
+      "read_data_name": "RightIn"
+    }
+    }
+}
diff --git a/tutorials/multiple-fields/precice-config.xml b/tutorials/multiple-fields/precice-config.xml
new file mode 100644
index 0000000..2955e79
--- /dev/null
+++ b/tutorials/multiple-fields/precice-config.xml
@@ -0,0 +1,127 @@
+<?xml version="1.0" encoding="UTF-8" ?>
+<precice-configuration>
+  <log>
+    <sink
+      filter="%Severity% > debug and %Rank% = 0"
+      format="---[precice] %ColorizedSeverity% %Message%"
+      enabled="true" />
+  </log>
+
+  <data:scalar name="LeftIn" waveform-degree="1" />
+  <data:scalar name="LeftOut" waveform-degree="1" />
+  <data:scalar name="RightIn" waveform-degree="1" />
+  <data:scalar name="RightOut" waveform-degree="1" />
+
+
+  <mesh name="Left1" dimensions="2">
+    <use-data name="LeftIn" />
+    <use-data name="LeftOut" />
+  </mesh>
+
+  <mesh name="Left2" dimensions="2">
+    <use-data name="RightIn" />
+    <use-data name="RightOut" />
+  </mesh>
+
+  <mesh name="Right1" dimensions="2">
+    <use-data name="LeftIn" />
+    <use-data name="LeftOut" />
+  </mesh>
+
+  <mesh name="Right2" dimensions="2">
+    <use-data name="RightIn" />
+    <use-data name="RightOut" />
+  </mesh>
+
+  <participant name="Left">
+    <provide-mesh name="Left1" />
+    <receive-mesh name="Right1" from="Right" />
+    <write-data name="LeftOut" mesh="Left1" />
+    <read-data name="LeftIn" mesh="Left1" />
+    <mapping:nearest-neighbor
+      direction="read"
+      from="Right1"
+      to="Left1"
+      constraint="consistent" />
+  </participant>
+
+    <participant name="Left">
+    <provide-mesh name="Left2" />
+    <receive-mesh name="Right2" from="Right" />
+    <write-data name="RightIn" mesh="Left2" />
+    <read-data name="RightOut" mesh="Left2" />
+    <mapping:nearest-neighbor
+      direction="read"
+      from="Right2"
+      to="Left2"
+      constraint="consistent" />
+  </participant>
+
+  <participant name="Right">
+    <provide-mesh name="Right1" />
+    <receive-mesh name="Left1" from="Left" />
+    <write-data name="LeftIn" mesh="Right1" />
+    <read-data name="LeftOut" mesh="Right1" />
+    <mapping:nearest-neighbor
+      direction="read"
+      from="Left1"
+      to="Right1"
+      constraint="consistent" />
+  </participant>
+
+    <participant name="Right">
+    <provide-mesh name="Right2" />
+    <receive-mesh name="Left2" from="Left" />
+    <write-data name="RightOut" mesh="Right2" />
+    <read-data name="RightIn" mesh="Right2" />
+    <mapping:nearest-neighbor
+      direction="read"
+      from="Left2"
+      to="Right2"
+      constraint="consistent" />
+  </participant>
+
+  <m2n:sockets acceptor="Left" connector="Right" exchange-directory=".." />
+
+  <coupling-scheme:serial-implicit>
+    <participants first="Left" second="Right" />
+    <max-time value="1.0" />
+    <time-window-size value="0.1" />
+    <max-iterations value="100" />
+    <exchange
+      data="LeftOut"
+      mesh="Left1"
+      from="Left"
+      to="Right"
+      initialize="true"
+      substeps="true" />
+    <exchange
+      data="LeftIn"
+      mesh="Right1"
+      from="Right"
+      to="Left"
+      initialize="true"
+      substeps="true" />
+    <exchange
+      data="RightOut"
+      mesh="Right2"
+      from="Right"
+      to="Left"
+      initialize="true"
+      substeps="true" />
+    <exchange
+      data="RightIn"
+      mesh="Left2"
+      from="Left"
+      to="Right"
+      initialize="true"
+      substeps="true" />
+    <relative-convergence-measure data="LeftOut" mesh="Left1" limit="1e-10" />
+    <relative-convergence-measure data="RightIn" mesh="Left2" limit="1e-10" />
+    <relative-convergence-measure data="LeftIn" mesh="Right1" limit="1e-10" />
+    <relative-convergence-measure data="RightOut" mesh="Right2" limit="1e-10" />
+    <acceleration:constant>
+      <relaxation value="0.5" />
+    </acceleration:constant>
+  </coupling-scheme:serial-implicit>
+</precice-configuration>
diff --git a/tutorials/multiple-fields/right.py b/tutorials/multiple-fields/right.py
new file mode 100644
index 0000000..e69de29

From 138cf3f6c361693903db8478000f883e7e980abf Mon Sep 17 00:00:00 2001
From: Niklas <niklas.vinnitchenko@tum.de>
Date: Tue, 11 Mar 2025 23:55:13 +0100
Subject: [PATCH 3/7] Completing trivial example with multiple data fields

---
 tutorials/multiple-fields/left.py            | 58 +++++++++++---------
 tutorials/multiple-fields/precice-config.xml | 15 ++---
 tutorials/multiple-fields/right.py           | 48 ++++++++++++++++
 3 files changed, 86 insertions(+), 35 deletions(-)

diff --git a/tutorials/multiple-fields/left.py b/tutorials/multiple-fields/left.py
index 1684ee1..88cd5a1 100644
--- a/tutorials/multiple-fields/left.py
+++ b/tutorials/multiple-fields/left.py
@@ -7,41 +7,47 @@
 from dolfinx import default_scalar_type
 from dolfinx.fem.petsc import LinearProblem
 
-domain = mesh.create_unit_square(MPI.COMM_WORLD, 8, 8, mesh.CellType.quadrilateral)
+import fenicsxprecice
+
+domain = mesh.create_rectangle(MPI.COMM_WORLD, [numpy.asarray((0,0)), numpy.asarray((1,1))], [10,10], mesh.CellType.quadrilateral)
+domain2 = mesh.create_rectangle(MPI.COMM_WORLD, [numpy.asarray((0,0)), numpy.asarray((1,1))], [11,11], mesh.CellType.quadrilateral)
 V = functionspace(domain, ("Lagrange", 1))
+V2 = functionspace(domain2, ("Lagrange", 1))
 uD = fem.Function(V)
 uD.interpolate(lambda x: 1 + x[0]**2 + 2 * x[1]**2)
+uD2 = fem.Function(V2)
+uD2.interpolate(lambda x: 1 + x[0]**2 + 3 * x[1]**2)
+
+def coupling_bc(x):
+    tol = 1E-14
+    return numpy.isclose(x[0], 1, tol)
 
-# Create facet to cell connectivity required to determine boundary facets
-tdim = domain.topology.dim
-fdim = tdim - 1
-domain.topology.create_connectivity(fdim, tdim)
-boundary_facets = mesh.exterior_facet_indices(domain.topology)
+precice = fenicsxprecice.Adapter(adapter_config_filename="precice-adapter-config-L.json", mpi_comm=MPI.COMM_SELF)
+precice.initialize({precice.Meshes.Left1:[coupling_bc, V, uD], 
+                    precice.Meshes.Left2:[coupling_bc, V2, uD2]})
 
-boundary_dofs = fem.locate_dofs_topological(V, fdim, boundary_facets)
-bc = fem.dirichletbc(uD, boundary_dofs)
+coupling_expression1 = precice.create_coupling_expression(precice.Meshes.Left1)
+coupling_expression2 = precice.create_coupling_expression(precice.Meshes.Left2)
 
-u = ufl.TrialFunction(V)
-v = ufl.TestFunction(V)
+while precice.is_coupling_ongoing():
 
-f = fem.Constant(domain, default_scalar_type(-6))
+    if precice.requires_writing_checkpoint():
+        precice.store_checkpoint(uD, 0, 0)
 
-a = ufl.dot(ufl.grad(u), ufl.grad(v)) * ufl.dx
-L = f * v * ufl.dx
+    read_data1 = precice.read_data(0, precice.Meshes.Left1)
+    read_data2 = precice.read_data(0, precice.Meshes.Left2)
 
-problem = LinearProblem(a, L, bcs=[bc], petsc_options={"ksp_type": "preonly", "pc_type": "lu"})
-uh = problem.solve()
+    precice.write_data(uD, precice.Meshes.Left1)
+    precice.write_data(uD2, precice.Meshes.Left2)
+    
+    precice.advance(1)
+    
+    if precice.requires_reading_checkpoint():
+        pass
+        
 
-V2 = fem.functionspace(domain, ("Lagrange", 2))
-uex = fem.Function(V2)
-uex.interpolate(lambda x: 1 + x[0]**2 + 2 * x[1]**2)
 
-L2_error = fem.form(ufl.inner(uh - uex, uh - uex) * ufl.dx)
-error_local = fem.assemble_scalar(L2_error)
-error_L2 = numpy.sqrt(domain.comm.allreduce(error_local, op=MPI.SUM))
+precice.finalize()
 
-error_max = numpy.max(numpy.abs(uD.x.array-uh.x.array))
-# Only print the error on one process
-if domain.comm.rank == 0:
-    print(f"Error_L2 : {error_L2:.2e}")
-    print(f"Error_max : {error_max:.2e}")
\ No newline at end of file
+print(read_data1)
+print(read_data2)
\ No newline at end of file
diff --git a/tutorials/multiple-fields/precice-config.xml b/tutorials/multiple-fields/precice-config.xml
index 2955e79..5da166a 100644
--- a/tutorials/multiple-fields/precice-config.xml
+++ b/tutorials/multiple-fields/precice-config.xml
@@ -36,19 +36,19 @@
   <participant name="Left">
     <provide-mesh name="Left1" />
     <receive-mesh name="Right1" from="Right" />
+    <provide-mesh name="Left2" />
+    <receive-mesh name="Right2" from="Right" />
+
     <write-data name="LeftOut" mesh="Left1" />
+    <write-data name="RightIn" mesh="Left2" />
+
     <read-data name="LeftIn" mesh="Left1" />
     <mapping:nearest-neighbor
       direction="read"
       from="Right1"
       to="Left1"
       constraint="consistent" />
-  </participant>
 
-    <participant name="Left">
-    <provide-mesh name="Left2" />
-    <receive-mesh name="Right2" from="Right" />
-    <write-data name="RightIn" mesh="Left2" />
     <read-data name="RightOut" mesh="Left2" />
     <mapping:nearest-neighbor
       direction="read"
@@ -67,9 +67,6 @@
       from="Left1"
       to="Right1"
       constraint="consistent" />
-  </participant>
-
-    <participant name="Right">
     <provide-mesh name="Right2" />
     <receive-mesh name="Left2" from="Left" />
     <write-data name="RightOut" mesh="Right2" />
@@ -86,7 +83,7 @@
   <coupling-scheme:serial-implicit>
     <participants first="Left" second="Right" />
     <max-time value="1.0" />
-    <time-window-size value="0.1" />
+    <time-window-size value="1.0" />
     <max-iterations value="100" />
     <exchange
       data="LeftOut"
diff --git a/tutorials/multiple-fields/right.py b/tutorials/multiple-fields/right.py
index e69de29..054e568 100644
--- a/tutorials/multiple-fields/right.py
+++ b/tutorials/multiple-fields/right.py
@@ -0,0 +1,48 @@
+from mpi4py import MPI
+from dolfinx import mesh
+from dolfinx.fem import functionspace
+from dolfinx import fem
+import numpy
+
+import fenicsxprecice
+
+domain = mesh.create_rectangle(MPI.COMM_WORLD, [numpy.asarray((1,0)), numpy.asarray((2,1))], [10,10], mesh.CellType.quadrilateral)
+domain2 = mesh.create_rectangle(MPI.COMM_WORLD, [numpy.asarray((1,0)), numpy.asarray((2,1))], [11,11], mesh.CellType.quadrilateral)
+V = functionspace(domain, ("Lagrange", 1))
+V2 = functionspace(domain2, ("Lagrange", 1))
+uD = fem.Function(V)
+uD.interpolate(lambda x: 1 + x[0]**2 + 2 * x[1]**2)
+uD2 = fem.Function(V2)
+uD2.interpolate(lambda x: 1 + x[0]**2 + 3 * x[1]**2)
+
+def coupling_bc(x):
+    tol = 1E-14
+    return numpy.isclose(x[0], 1, tol)
+
+precice = fenicsxprecice.Adapter(adapter_config_filename="precice-adapter-config-R.json", mpi_comm=MPI.COMM_WORLD)
+precice.initialize({precice.Meshes.Right1:[coupling_bc, V, uD], 
+                    precice.Meshes.Right2:[coupling_bc, V2, uD2]})
+
+coupling_expression1 = precice.create_coupling_expression(precice.Meshes.Right1)
+coupling_expression2 = precice.create_coupling_expression(precice.Meshes.Right2)
+
+while precice.is_coupling_ongoing():
+
+    if precice.requires_writing_checkpoint():
+        precice.store_checkpoint(uD, 0, 0)
+
+    read_data1 = precice.read_data(0, precice.Meshes.Right1)
+    read_data2 = precice.read_data(0, precice.Meshes.Right2)
+
+    precice.write_data(uD, precice.Meshes.Right1)
+    precice.write_data(uD2, precice.Meshes.Right2)
+
+    precice.advance(1)
+    
+    if precice.requires_reading_checkpoint():
+        pass
+
+precice.finalize()
+
+print(read_data1)
+print(read_data2)
\ No newline at end of file

From 4993a0aa6dc052263c30a0973dff07297acc17db Mon Sep 17 00:00:00 2001
From: Niklas <niklas.vinnitchenko@tum.de>
Date: Wed, 12 Mar 2025 12:28:51 +0100
Subject: [PATCH 4/7] Simplify example

---
 tutorials/multiple-fields/left.py  | 8 ++++----
 tutorials/multiple-fields/right.py | 8 ++++----
 2 files changed, 8 insertions(+), 8 deletions(-)

diff --git a/tutorials/multiple-fields/left.py b/tutorials/multiple-fields/left.py
index 88cd5a1..9ed774f 100644
--- a/tutorials/multiple-fields/left.py
+++ b/tutorials/multiple-fields/left.py
@@ -14,9 +14,9 @@
 V = functionspace(domain, ("Lagrange", 1))
 V2 = functionspace(domain2, ("Lagrange", 1))
 uD = fem.Function(V)
-uD.interpolate(lambda x: 1 + x[0]**2 + 2 * x[1]**2)
+uD.interpolate(lambda x: 1+x[0])
 uD2 = fem.Function(V2)
-uD2.interpolate(lambda x: 1 + x[0]**2 + 3 * x[1]**2)
+uD2.interpolate(lambda x: 2+x[0])
 
 def coupling_bc(x):
     tol = 1E-14
@@ -49,5 +49,5 @@ def coupling_bc(x):
 
 precice.finalize()
 
-print(read_data1)
-print(read_data2)
\ No newline at end of file
+print(read_data1) # expected: x[0]+3
+print(read_data2) # expected: x[0]+4
\ No newline at end of file
diff --git a/tutorials/multiple-fields/right.py b/tutorials/multiple-fields/right.py
index 054e568..b3183e4 100644
--- a/tutorials/multiple-fields/right.py
+++ b/tutorials/multiple-fields/right.py
@@ -11,9 +11,9 @@
 V = functionspace(domain, ("Lagrange", 1))
 V2 = functionspace(domain2, ("Lagrange", 1))
 uD = fem.Function(V)
-uD.interpolate(lambda x: 1 + x[0]**2 + 2 * x[1]**2)
+uD.interpolate(lambda x: 3+x[0])
 uD2 = fem.Function(V2)
-uD2.interpolate(lambda x: 1 + x[0]**2 + 3 * x[1]**2)
+uD2.interpolate(lambda x: 4+x[0])
 
 def coupling_bc(x):
     tol = 1E-14
@@ -44,5 +44,5 @@ def coupling_bc(x):
 
 precice.finalize()
 
-print(read_data1)
-print(read_data2)
\ No newline at end of file
+print(read_data1) # expected: x[0]+1
+print(read_data2) # expected: x[0]+2
\ No newline at end of file

From f0cd8fe22977ac0e976129c4463b112493bcacd2 Mon Sep 17 00:00:00 2001
From: Niklas <niklas.vinnitchenko@tum.de>
Date: Wed, 12 Mar 2025 12:56:29 +0100
Subject: [PATCH 5/7] Add changelog entry

---
 CHANGELOG.md | 1 +
 1 file changed, 1 insertion(+)

diff --git a/CHANGELOG.md b/CHANGELOG.md
index 4e4a2ba..24c787e 100644
--- a/CHANGELOG.md
+++ b/CHANGELOG.md
@@ -2,6 +2,7 @@
 
 ## current development
 
+* Support communication of multiple data fields. [#34](https://github.com/precice/fenicsx-adapter/pull/34)
 * Update to support dolfinx version 0.9.0 and preCICE v3. [#28](https://github.com/precice/fenicsx-adapter/pull/28)
 * Developed initial working version with example case `tutorials/partitioned-heat-conduction`. [#15](https://github.com/precice/fenicsx-adapter/pull/15)
 * Forked initial version of this adapter from [`precice/fenics-adapter@v1.2.0`](https://github.com/precice/fenics-adapter/releases/tag/v1.2.0).

From 46f1771e6082fe43d0f50e24b1df0209bfca56a7 Mon Sep 17 00:00:00 2001
From: Niklas <niklas.vinnitchenko@tum.de>
Date: Wed, 12 Mar 2025 13:02:32 +0100
Subject: [PATCH 6/7] Formatting

---
 fenicsxprecice/config.py           | 14 ++++++-------
 fenicsxprecice/fenicsxprecice.py   | 31 ++++++++++++++++++----------
 tutorials/multiple-fields/left.py  | 33 ++++++++++++++++++++----------
 tutorials/multiple-fields/right.py | 30 +++++++++++++++++++--------
 4 files changed, 70 insertions(+), 38 deletions(-)

diff --git a/fenicsxprecice/config.py b/fenicsxprecice/config.py
index b1fbe7c..8b6536b 100644
--- a/fenicsxprecice/config.py
+++ b/fenicsxprecice/config.py
@@ -41,16 +41,18 @@ def read_json(self, adapter_config_filename):
         self._participant_name = data["participant_name"]
         self._coupling_mesh_names = list(data["interfaces"].keys())
 
-        for mesh_name in self._coupling_mesh_names:    
+        for mesh_name in self._coupling_mesh_names:
             try:
                 self._write_data_names[mesh_name] = data["interfaces"][mesh_name]["write_data_name"]
             except KeyError:
-                self._write_data_names[mesh_name] = None  # not required for one-way coupling, if this participant reads data
+                # not required for one-way coupling, if this participant reads data
+                self._write_data_names[mesh_name] = None
 
             try:
                 self._read_data_names[mesh_name] = data["interfaces"][mesh_name]["read_data_name"]
             except KeyError:
-                self._read_data_names[mesh_name] = None  # not required for one-way coupling, if this participant writes data
+                # not required for one-way coupling, if this participant writes data
+                self._read_data_names[mesh_name] = None
 
         read_file.close()
 
@@ -63,14 +65,13 @@ def get_participant_name(self):
     def get_coupling_mesh_names(self):
         return self._coupling_mesh_names
 
-
     def get_read_data_name(self, mesh_name):
         """
         Parameters
         ----------
         mesh_name : fenicsxprecice.Meshes member
             Member of the enum fenicsxprecoce.Meshes
-    
+
         """
         return self._read_data_names[mesh_name.name]
 
@@ -80,7 +81,6 @@ def get_write_data_name(self, mesh_name):
         ----------
         mesh_name : fenicsxprecice.Meshes member
             Member of the enum fenicsxprecoce.Meshes
-    
+
         """
         return self._write_data_names[mesh_name.name]
-    
\ No newline at end of file
diff --git a/fenicsxprecice/fenicsxprecice.py b/fenicsxprecice/fenicsxprecice.py
index 0d98f38..14ac52c 100644
--- a/fenicsxprecice/fenicsxprecice.py
+++ b/fenicsxprecice/fenicsxprecice.py
@@ -49,7 +49,7 @@ def __init__(self, mpi_comm, adapter_config_filename='precice-adapter-config.jso
         """
 
         self._config = Config(adapter_config_filename)
-        
+
         self.Meshes = Enum('Meshes', self._config.get_coupling_mesh_names())
 
         # Setup up MPI communicator
@@ -95,7 +95,7 @@ def create_coupling_expression(self, mesh_name):
         """
         Creates a FEniCSx Expression in the form of an object of class GeneralInterpolationExpression or
         ExactInterpolationExpression. The adapter will hold this object till the coupling is on going.
-        
+
         Parameters
         ----------
         mesh_name: member of fenicsxprecice.Meshes
@@ -107,10 +107,13 @@ def create_coupling_expression(self, mesh_name):
             Reference to object of class GeneralInterpolationExpression or ExactInterpolationExpression.
         """
 
-        if not (self._read_function_types[mesh_name] is FunctionType.SCALAR or self._read_function_types[mesh_name] is FunctionType.VECTOR):
+        if not (self._read_function_types[mesh_name]
+                is FunctionType.SCALAR or self._read_function_types[mesh_name] is FunctionType.VECTOR):
             raise Exception("No valid read_function is provided in initialization. Cannot create coupling expression")
 
-        coupling_expression = self._my_expression(self._read_function_spaces[mesh_name], self._read_function_types[mesh_name])
+        coupling_expression = self._my_expression(
+            self._read_function_spaces[mesh_name],
+            self._read_function_types[mesh_name])
 
         return coupling_expression
 
@@ -162,7 +165,11 @@ def read_data(self, dt, mesh_name):
                 self._precice_vertex_ids[mesh_name],
                 dt
             )
-            read_data = {tuple(key): value for key, value in zip(self._fenicsx_vertices[mesh_name].get_coordinates(), read_data)}
+            read_data = {
+                tuple(key): value for key,
+                value in zip(
+                    self._fenicsx_vertices[mesh_name].get_coordinates(),
+                    read_data)}
 
         else:
             pass
@@ -230,7 +237,8 @@ def initialize(self, coupling_meshes):
             if isinstance(write_object, fem.Function):  # precice.initialize_data() will be called using this Function
                 write_function_space = write_object.function_space
                 write_function = write_object
-            elif isinstance(write_object, fem.FunctionSpace):  # preCICE will use default zero values for initialization.
+            # preCICE will use default zero values for initialization.
+            elif isinstance(write_object, fem.FunctionSpace):
                 write_function_space = write_object
                 write_function = None
             elif write_object is None:
@@ -246,7 +254,8 @@ def initialize(self, coupling_meshes):
             elif read_function_space is None:
                 pass
             else:
-                raise Exception(f"Given read_function_space of {c_mesh} is not of type dolfinx.functions.functionspace.FunctionSpace")
+                raise Exception(f"Given read_function_space of {
+                                c_mesh} is not of type dolfinx.functions.functionspace.FunctionSpace")
 
             if read_function_space is None and write_function_space:
                 self._coupling_types[c_mesh] = CouplingMode.UNI_DIRECTIONAL_WRITE_COUPLING
@@ -269,8 +278,8 @@ def initialize(self, coupling_meshes):
                                 "only reads data. If two-way coupling is implemented then both read_function_space"
                                 " and write_object need to be provided.")
             else:
-                raise Exception(f"Incorrect read and write function space combination provided for {c_mesh}. Please check input "
-                                "parameters in initialization")
+                raise Exception(f"Incorrect read and write function space combination provided for {
+                                c_mesh}. Please check input " "parameters in initialization")
 
             coupling_type = self._coupling_types[c_mesh]
             if coupling_type is CouplingMode.UNI_DIRECTIONAL_READ_COUPLING or \
@@ -284,7 +293,8 @@ def initialize(self, coupling_meshes):
                 self._write_function_types[c_mesh] = determine_function_type(write_function_space)
                 self._write_function_spaces[c_mesh] = write_function_space
 
-            # Set vertices on the coupling subdomain for this rank (assumed to be equal across all fields that will be coupled)
+            # Set vertices on the coupling subdomain for this rank (assumed to be
+            # equal across all fields that will be coupled)
             self._fenicsx_dims = function_space.mesh.geometry.dim
             # returns 3d coordinates (necessary later for writing the data!)
             # coupling subdomain is at [0] !
@@ -322,7 +332,6 @@ def initialize(self, coupling_meshes):
                         "preCICE requires you to write initial data. Please provide a write_function to initialize(...)")
                 self.write_data(write_function, c_mesh)
 
-        
         self._participant.initialize()
 
     def store_checkpoint(self, payload, t, n):
diff --git a/tutorials/multiple-fields/left.py b/tutorials/multiple-fields/left.py
index 9ed774f..ab62331 100644
--- a/tutorials/multiple-fields/left.py
+++ b/tutorials/multiple-fields/left.py
@@ -9,22 +9,34 @@
 
 import fenicsxprecice
 
-domain = mesh.create_rectangle(MPI.COMM_WORLD, [numpy.asarray((0,0)), numpy.asarray((1,1))], [10,10], mesh.CellType.quadrilateral)
-domain2 = mesh.create_rectangle(MPI.COMM_WORLD, [numpy.asarray((0,0)), numpy.asarray((1,1))], [11,11], mesh.CellType.quadrilateral)
+domain = mesh.create_rectangle(
+    MPI.COMM_WORLD, [
+        numpy.asarray(
+            (0, 0)), numpy.asarray(
+                (1, 1))], [
+                    10, 10], mesh.CellType.quadrilateral)
+domain2 = mesh.create_rectangle(
+    MPI.COMM_WORLD, [
+        numpy.asarray(
+            (0, 0)), numpy.asarray(
+                (1, 1))], [
+                    11, 11], mesh.CellType.quadrilateral)
 V = functionspace(domain, ("Lagrange", 1))
 V2 = functionspace(domain2, ("Lagrange", 1))
 uD = fem.Function(V)
-uD.interpolate(lambda x: 1+x[0])
+uD.interpolate(lambda x: 1 + x[0])
 uD2 = fem.Function(V2)
-uD2.interpolate(lambda x: 2+x[0])
+uD2.interpolate(lambda x: 2 + x[0])
+
 
 def coupling_bc(x):
     tol = 1E-14
     return numpy.isclose(x[0], 1, tol)
 
+
 precice = fenicsxprecice.Adapter(adapter_config_filename="precice-adapter-config-L.json", mpi_comm=MPI.COMM_SELF)
-precice.initialize({precice.Meshes.Left1:[coupling_bc, V, uD], 
-                    precice.Meshes.Left2:[coupling_bc, V2, uD2]})
+precice.initialize({precice.Meshes.Left1: [coupling_bc, V, uD],
+                    precice.Meshes.Left2: [coupling_bc, V2, uD2]})
 
 coupling_expression1 = precice.create_coupling_expression(precice.Meshes.Left1)
 coupling_expression2 = precice.create_coupling_expression(precice.Meshes.Left2)
@@ -39,15 +51,14 @@ def coupling_bc(x):
 
     precice.write_data(uD, precice.Meshes.Left1)
     precice.write_data(uD2, precice.Meshes.Left2)
-    
+
     precice.advance(1)
-    
+
     if precice.requires_reading_checkpoint():
         pass
-        
 
 
 precice.finalize()
 
-print(read_data1) # expected: x[0]+3
-print(read_data2) # expected: x[0]+4
\ No newline at end of file
+print(read_data1)  # expected: x[0]+3
+print(read_data2)  # expected: x[0]+4
diff --git a/tutorials/multiple-fields/right.py b/tutorials/multiple-fields/right.py
index b3183e4..8eee9bc 100644
--- a/tutorials/multiple-fields/right.py
+++ b/tutorials/multiple-fields/right.py
@@ -6,22 +6,34 @@
 
 import fenicsxprecice
 
-domain = mesh.create_rectangle(MPI.COMM_WORLD, [numpy.asarray((1,0)), numpy.asarray((2,1))], [10,10], mesh.CellType.quadrilateral)
-domain2 = mesh.create_rectangle(MPI.COMM_WORLD, [numpy.asarray((1,0)), numpy.asarray((2,1))], [11,11], mesh.CellType.quadrilateral)
+domain = mesh.create_rectangle(
+    MPI.COMM_WORLD, [
+        numpy.asarray(
+            (1, 0)), numpy.asarray(
+                (2, 1))], [
+                    10, 10], mesh.CellType.quadrilateral)
+domain2 = mesh.create_rectangle(
+    MPI.COMM_WORLD, [
+        numpy.asarray(
+            (1, 0)), numpy.asarray(
+                (2, 1))], [
+                    11, 11], mesh.CellType.quadrilateral)
 V = functionspace(domain, ("Lagrange", 1))
 V2 = functionspace(domain2, ("Lagrange", 1))
 uD = fem.Function(V)
-uD.interpolate(lambda x: 3+x[0])
+uD.interpolate(lambda x: 3 + x[0])
 uD2 = fem.Function(V2)
-uD2.interpolate(lambda x: 4+x[0])
+uD2.interpolate(lambda x: 4 + x[0])
+
 
 def coupling_bc(x):
     tol = 1E-14
     return numpy.isclose(x[0], 1, tol)
 
+
 precice = fenicsxprecice.Adapter(adapter_config_filename="precice-adapter-config-R.json", mpi_comm=MPI.COMM_WORLD)
-precice.initialize({precice.Meshes.Right1:[coupling_bc, V, uD], 
-                    precice.Meshes.Right2:[coupling_bc, V2, uD2]})
+precice.initialize({precice.Meshes.Right1: [coupling_bc, V, uD],
+                    precice.Meshes.Right2: [coupling_bc, V2, uD2]})
 
 coupling_expression1 = precice.create_coupling_expression(precice.Meshes.Right1)
 coupling_expression2 = precice.create_coupling_expression(precice.Meshes.Right2)
@@ -38,11 +50,11 @@ def coupling_bc(x):
     precice.write_data(uD2, precice.Meshes.Right2)
 
     precice.advance(1)
-    
+
     if precice.requires_reading_checkpoint():
         pass
 
 precice.finalize()
 
-print(read_data1) # expected: x[0]+1
-print(read_data2) # expected: x[0]+2
\ No newline at end of file
+print(read_data1)  # expected: x[0]+1
+print(read_data2)  # expected: x[0]+2

From 1faabaa2d4445a6540b575170b0870ca516fccbd Mon Sep 17 00:00:00 2001
From: Niklas <niklas.vinnitchenko@tum.de>
Date: Wed, 12 Mar 2025 13:09:44 +0100
Subject: [PATCH 7/7] Update documentation

---
 fenicsxprecice/fenicsxprecice.py | 8 ++++++++
 1 file changed, 8 insertions(+)

diff --git a/fenicsxprecice/fenicsxprecice.py b/fenicsxprecice/fenicsxprecice.py
index 14ac52c..88f112b 100644
--- a/fenicsxprecice/fenicsxprecice.py
+++ b/fenicsxprecice/fenicsxprecice.py
@@ -147,6 +147,12 @@ def read_data(self, dt, mesh_name):
         Note: For quasi 2D-3D coupled simulation (FEniCSx participant is 2D) the Z-component of the data and vertices
         is deleted.
 
+        Parameters
+        ----------
+        dt : offset within time window
+        mesh_name: member of fenicsxprecice.Meshes
+            Specifies for which mesh the data shall be read
+
         Returns
         -------
         data : dict_like
@@ -185,6 +191,8 @@ def write_data(self, write_function, mesh_name):
         ----------
         write_function : Object of class dolfinx.functions.function.Function
             A FEniCSx function consisting of the data which this participant will write to preCICE in every time step.
+        mesh_name : member of fenicsxprecice.Meshes
+            Specifies the mesh from which the data is written
         """
 
         assert (self._coupling_types[mesh_name] is CouplingMode.UNI_DIRECTIONAL_WRITE_COUPLING or