diff --git a/channel-transport-reaction/chemical-fenics/chemical-reaction-advection-diffusion.py b/channel-transport-reaction/chemical-fenics/chemical-reaction-advection-diffusion.py
index c95795821..03d80968e 100644
--- a/channel-transport-reaction/chemical-fenics/chemical-reaction-advection-diffusion.py
+++ b/channel-transport-reaction/chemical-fenics/chemical-reaction-advection-diffusion.py
@@ -27,9 +27,8 @@ def inside(self, x, on_boundary):
 # Initialize preCICE
 precice = fenicsprecice.Adapter(
     adapter_config_filename="chemical-reaction-advection-diffusion.json")
-precice_dt = precice.initialize(
-    coupling_subdomain=CouplingDomain(),
-    read_function_space=W)
+precice.initialize(coupling_subdomain=CouplingDomain(), read_function_space=W)
+precice_dt = precice.get_max_time_step_size()
 
 flow_expr = precice.create_coupling_expression()
 
@@ -95,13 +94,14 @@ def inside(self, x, on_boundary):
   # No implicit coupling
 while precice.is_coupling_ongoing():
 
-    read_data = precice.read_data()
+    precice_dt = precice.get_max_time_step_size()
+    dt = np.min([default_dt, precice_dt])
+    read_data = precice.read_data(dt)
     precice.update_coupling_expression(flow_expr, read_data)
     flow_old.assign(flow)
     flow.interpolate(flow_expr)
     # If we add writing, do it here
 
-    dt = np.min([default_dt, precice_dt])
     k.assign(1. / dt)
 
     t += dt
@@ -128,6 +128,6 @@ def inside(self, x, on_boundary):
     vtkfileC << u_C, t
     vtkfileFlow << flow, t
 
-    precice_dt = precice.advance(dt)
+    precice.advance(dt)
 
 precice.finalize()
diff --git a/channel-transport-reaction/fluid-fenics/fluid.py b/channel-transport-reaction/fluid-fenics/fluid.py
index fafdfb677..436e4ecd9 100644
--- a/channel-transport-reaction/fluid-fenics/fluid.py
+++ b/channel-transport-reaction/fluid-fenics/fluid.py
@@ -118,10 +118,8 @@ def inside(self, x, on_boundary):
 
 
 precice = fenicsprecice.Adapter(adapter_config_filename="fluid-config.json")
-precice_dt = precice.initialize(
-    coupling_subdomain=CouplingDomain(),
-    write_object=u_)
-
+precice.initialize(coupling_subdomain=CouplingDomain(), write_object=u_)
+precice_dt = precice.get_max_time_step_size()
 
 dt = np.min([default_dt, precice_dt])
 k.assign(dt)
@@ -155,7 +153,8 @@ def inside(self, x, on_boundary):
 
     vtkfile << u_
 
-    precice_dt = precice.advance(dt)
+    precice.advance(dt)
+    precice_dt = precice.get_max_time_step_size()
     dt = np.min([default_dt, precice_dt])
     k.assign(dt)
 
diff --git a/elastic-tube-3d/solid-fenics/solid.py b/elastic-tube-3d/solid-fenics/solid.py
index bef8cfd4f..0278a2572 100644
--- a/elastic-tube-3d/solid-fenics/solid.py
+++ b/elastic-tube-3d/solid-fenics/solid.py
@@ -65,7 +65,8 @@ def neumann_boundary(x, on_boundary):
 precice = Adapter(adapter_config_filename="precice-adapter-config-fsi-s.json")
 
 # Initialize the coupling interface
-precice_dt = precice.initialize(coupling_boundary, read_function_space=V, write_object=V, fixed_boundary=fixed_boundary)
+precice.initialize(coupling_boundary, read_function_space=V, write_object=V, fixed_boundary=fixed_boundary)
+precice_dt = precice.get_max_time_step_size()
 
 fenics_dt = precice_dt  # if fenics_dt == precice_dt, no subcycling is applied
 dt = Constant(np.min([precice_dt, fenics_dt]))
@@ -168,15 +169,18 @@ def avg(x_old, x_new, alpha):
 
 u_n.rename("Displacement", "")
 u_np1.rename("Displacement", "")
-displacement_out << u_n
+displacement_out << (u_n, t)
 
 while precice.is_coupling_ongoing():
 
-    if precice.is_action_required(precice.action_write_iteration_checkpoint()):  # write checkpoint
+    if precice.requires_writing_checkpoint():  # write checkpoint
         precice.store_checkpoint(u_n, t, n)
 
+    precice_dt = precice.get_max_time_step_size()
+    dt = Constant(np.min([precice_dt, fenics_dt]))
+
     # read data from preCICE and get a new coupling expression
-    read_data = precice.read_data()
+    read_data = precice.read_data(dt)
 
     # Update the point sources on the coupling boundary with the new read data
     forces_x, forces_y, forces_z = precice.get_point_sources(read_data)
@@ -195,17 +199,15 @@ def avg(x_old, x_new, alpha):
     assert (b is not b_forces)
     solve(A, u_np1.vector(), b_forces)
 
-    dt = Constant(np.min([precice_dt, fenics_dt]))
-
     # Write relative displacements to preCICE
     u_delta.vector()[:] = u_np1.vector()[:] - u_n.vector()[:]
     precice.write_data(u_delta)
 
     # Call to advance coupling, also returns the optimum time step value
-    precice_dt = precice.advance(dt(0))
+    precice.advance(dt(0))
 
     # Either revert to old step if timestep has not converged or move to next timestep
-    if precice.is_action_required(precice.action_read_iteration_checkpoint()):  # roll back to checkpoint
+    if precice.requires_reading_checkpoint():  # roll back to checkpoint
         u_cp, t_cp, n_cp = precice.retrieve_checkpoint()
         u_n.assign(u_cp)
         t = t_cp
@@ -221,6 +223,6 @@ def avg(x_old, x_new, alpha):
             displacement_out << (u_n, t)
 
 # Plot tip displacement evolution
-displacement_out << u_n
+displacement_out << (u_n, t)
 
 precice.finalize()
diff --git a/flow-over-heated-plate/solid-fenics/solid.py b/flow-over-heated-plate/solid-fenics/solid.py
index 23ee3bbd4..f07cc29d5 100644
--- a/flow-over-heated-plate/solid-fenics/solid.py
+++ b/flow-over-heated-plate/solid-fenics/solid.py
@@ -115,13 +115,14 @@ def determine_heat_flux(V_g, u, k, flux):
 # Adapter definition and initialization
 precice = Adapter(adapter_config_filename="precice-adapter-config.json")
 
-precice_dt = precice.initialize(coupling_boundary, read_function_space=V, write_object=V_flux_y)
+precice.initialize(coupling_boundary, read_function_space=V, write_object=V_flux_y)
 
 # Create a FEniCS Expression to define and control the coupling boundary values
 coupling_expression = precice.create_coupling_expression()
 
 # Assigning appropriate dt
 dt = Constant(0)
+precice_dt = precice.get_max_time_step_size()
 dt.assign(np.min([fenics_dt, precice_dt]))
 
 # Define variational problem
@@ -149,7 +150,7 @@ def determine_heat_flux(V_g, u, k, flux):
 
 # Create output file
 file_out = File("output/%s.pvd" % precice.get_participant_name())
-file_out << u_n
+file_out << (u_n, t)
 
 print("output vtk for time = {}".format(float(t)))
 n = 0
@@ -159,16 +160,16 @@ def determine_heat_flux(V_g, u, k, flux):
 
 while precice.is_coupling_ongoing():
 
-    if precice.is_action_required(precice.action_write_iteration_checkpoint()):  # write checkpoint
+    if precice.requires_writing_checkpoint():  # write checkpoint
         precice.store_checkpoint(u_n, t, n)
 
-    read_data = precice.read_data()
+    precice_dt = precice.get_max_time_step_size()
+    dt.assign(np.min([fenics_dt, precice_dt]))
+    read_data = precice.read_data(dt)
 
     # Update the coupling expression with the new read data
     precice.update_coupling_expression(coupling_expression, read_data)
 
-    dt.assign(np.min([fenics_dt, precice_dt]))
-
     # Compute solution
     solve(a == L, u_np1, bcs)
 
@@ -177,9 +178,9 @@ def determine_heat_flux(V_g, u, k, flux):
     fluxes_y = fluxes.sub(1)  # only exchange y component of flux.
     precice.write_data(fluxes_y)
 
-    precice_dt = precice.advance(dt(0))
+    precice.advance(dt(0))
 
-    if precice.is_action_required(precice.action_read_iteration_checkpoint()):  # roll back to checkpoint
+    if precice.requires_reading_checkpoint():  # roll back to checkpoint
         u_cp, t_cp, n_cp = precice.retrieve_checkpoint()
         u_n.assign(u_cp)
         t = t_cp
@@ -193,7 +194,7 @@ def determine_heat_flux(V_g, u, k, flux):
         tol = 10e-5  # we need some tolerance, since otherwise output might be skipped.
         if abs((t + tol) % dt_out) < 2 * tol:  # output if t is a multiple of dt_out
             print("output vtk for time = {}".format(float(t)))
-            file_out << u_n
+            file_out << (u_n, t)
 
     # Update dirichlet BC
     u_D.t = t + float(dt)
diff --git a/partitioned-heat-conduction-complex/fenics/heat.py b/partitioned-heat-conduction-complex/fenics/heat.py
index 9ba1f519c..eb4b74e60 100644
--- a/partitioned-heat-conduction-complex/fenics/heat.py
+++ b/partitioned-heat-conduction-complex/fenics/heat.py
@@ -103,14 +103,15 @@ def determine_gradient(V_g, u, flux):
 # Initialize the adapter according to the specific participant
 if problem is ProblemType.DIRICHLET:
     precice = Adapter(adapter_config_filename="precice-adapter-config-D.json")
-    precice_dt = precice.initialize(coupling_boundary, read_function_space=V, write_object=f_N_function)
+    precice.initialize(coupling_boundary, read_function_space=V, write_object=f_N_function)
 elif problem is ProblemType.NEUMANN:
     precice = Adapter(adapter_config_filename="precice-adapter-config-N.json")
-    precice_dt = precice.initialize(coupling_boundary, read_function_space=V_g, write_object=u_D_function)
+    precice.initialize(coupling_boundary, read_function_space=V_g, write_object=u_D_function)
 
 boundary_marker = False
 
 dt = Constant(0)
+precice_dt = precice.get_max_time_step_size()
 dt.assign(np.min([fenics_dt, precice_dt]))
 
 # Define variational problem
@@ -173,7 +174,7 @@ def determine_gradient(V_g, u, flux):
 # output solution and reference solution at t=0, n=0
 n = 0
 print('output u^%d and u_ref^%d' % (n, n))
-temperature_out << u_n
+temperature_out << (u_n, t)
 ref_out << u_ref
 ranks << mesh_rank
 
@@ -191,10 +192,12 @@ def determine_gradient(V_g, u, flux):
 while precice.is_coupling_ongoing():
 
     # write checkpoint
-    if precice.is_action_required(precice.action_write_iteration_checkpoint()):
+    if precice.requires_writing_checkpoint():
         precice.store_checkpoint(u_n, t, n)
 
-    read_data = precice.read_data()
+    precice_dt = precice.get_max_time_step_size()
+    dt.assign(np.min([fenics_dt, precice_dt]))
+    read_data = precice.read_data(dt)
     if problem is ProblemType.DIRICHLET and (domain_part is DomainPart.CIRCULAR or domain_part is DomainPart.RECTANGLE):
         # We have to data for an arbitrary point that is not on the circle, to obtain exact solution.
         # See https://github.com/precice/fenics-adapter/issues/113 for details.
@@ -203,8 +206,6 @@ def determine_gradient(V_g, u, flux):
     # Update the coupling expression with the new read data
     precice.update_coupling_expression(coupling_expression, read_data)
 
-    dt.assign(np.min([fenics_dt, precice_dt]))
-
     # Compute solution u^n+1, use bcs u_D^n+1, u^n and coupling bcs
     solve(a == L, u_np1, bcs)
 
@@ -217,10 +218,10 @@ def determine_gradient(V_g, u, flux):
         # Neumann problem reads flux and writes temperature on boundary to Dirichlet problem
         precice.write_data(u_np1)
 
-    precice_dt = precice.advance(dt(0))
+    precice.advance(dt(0))
 
     # roll back to checkpoint
-    if precice.is_action_required(precice.action_read_iteration_checkpoint()):
+    if precice.requires_reading_checkpoint():
         u_cp, t_cp, n_cp = precice.retrieve_checkpoint()
         u_n.assign(u_cp)
         t = t_cp
@@ -237,7 +238,7 @@ def determine_gradient(V_g, u, flux):
         print('n = %d, t = %.2f: L2 error on domain = %.3g' % (n, t, error))
         # output solution and reference solution at t_n+1
         print('output u^%d and u_ref^%d' % (n, n))
-        temperature_out << u_n
+        temperature_out << (u_n, t)
         ref_out << u_ref
         error_out << error_pointwise
 
diff --git a/partitioned-heat-conduction-complex/fenics/precice-adapter-config-D.json b/partitioned-heat-conduction-complex/fenics/precice-adapter-config-D.json
index 9430d6e1e..2a31d960f 100644
--- a/partitioned-heat-conduction-complex/fenics/precice-adapter-config-D.json
+++ b/partitioned-heat-conduction-complex/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"
   }
 }
diff --git a/partitioned-heat-conduction-complex/fenics/precice-adapter-config-N.json b/partitioned-heat-conduction-complex/fenics/precice-adapter-config-N.json
index b21673b23..03d46a541 100644
--- a/partitioned-heat-conduction-complex/fenics/precice-adapter-config-N.json
+++ b/partitioned-heat-conduction-complex/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"
   }
 }
diff --git a/partitioned-heat-conduction-complex/precice-config.xml b/partitioned-heat-conduction-complex/precice-config.xml
index 3e38a7ee7..8fcca7d50 100644
--- a/partitioned-heat-conduction-complex/precice-config.xml
+++ b/partitioned-heat-conduction-complex/precice-config.xml
@@ -6,22 +6,22 @@
   </log>
 
   <data:scalar name="Temperature"/>
-  <data:vector name="Flux"/>
+  <data:vector name="Heat-Flux"/>
 
   <mesh name="Dirichlet-Mesh" dimensions="2">
     <use-data name="Temperature"/>
-    <use-data name="Flux"/>
+    <use-data name="Heat-Flux"/>
   </mesh>
 
   <mesh name="Neumann-Mesh" dimensions="2">
     <use-data name="Temperature"/>
-    <use-data name="Flux"/>
+    <use-data name="Heat-Flux"/>
   </mesh>
 
   <participant name="Dirichlet">
     <provide-mesh name="Dirichlet-Mesh" />
     <receive-mesh name="Neumann-Mesh" from="Neumann" />
-    <write-data name="Flux" mesh="Dirichlet-Mesh" />
+    <write-data name="Heat-Flux" mesh="Dirichlet-Mesh" />
     <read-data name="Temperature" mesh="Dirichlet-Mesh" />
     <mapping:nearest-projection direction="read" from="Neumann-Mesh" to="Dirichlet-Mesh" constraint="consistent" />
   </participant>
@@ -30,7 +30,7 @@
     <provide-mesh name="Neumann-Mesh" />
     <receive-mesh name="Dirichlet-Mesh" from="Dirichlet"/>
     <write-data name="Temperature" mesh="Neumann-Mesh"/>
-    <read-data name="Flux" mesh="Neumann-Mesh"/>
+    <read-data name="Heat-Flux" mesh="Neumann-Mesh"/>
     <mapping:nearest-projection direction="read" from="Dirichlet-Mesh" to="Neumann-Mesh" constraint="consistent" />
   </participant>
 
@@ -41,9 +41,9 @@
     <max-time value="1.0"/>
     <time-window-size value="0.1"/>
     <max-iterations value="100"/>
-    <exchange data="Flux" mesh="Dirichlet-Mesh" from="Dirichlet" to="Neumann" />
+    <exchange data="Heat-Flux" mesh="Dirichlet-Mesh" from="Dirichlet" to="Neumann" />
     <exchange data="Temperature" mesh="Neumann-Mesh" from="Neumann" to="Dirichlet" initialize="true"/>
-    <relative-convergence-measure data="Flux" mesh="Dirichlet-Mesh" limit="1e-5"/>
+    <relative-convergence-measure data="Heat-Flux" mesh="Dirichlet-Mesh" limit="1e-5"/>
     <relative-convergence-measure data="Temperature" mesh="Neumann-Mesh" limit="1e-5"/>
     <acceleration:IQN-ILS>
       <data name="Temperature" mesh="Neumann-Mesh"/>
diff --git a/partitioned-heat-conduction/fenics/heat.py b/partitioned-heat-conduction/fenics/heat.py
index e7b7f6529..3ca23990b 100644
--- a/partitioned-heat-conduction/fenics/heat.py
+++ b/partitioned-heat-conduction/fenics/heat.py
@@ -58,7 +58,7 @@ def determine_gradient(V_g, u, flux):
 command_group.add_argument("-d", "--dirichlet", help="create a dirichlet problem", dest="dirichlet",
                            action="store_true")
 command_group.add_argument("-n", "--neumann", help="create a neumann problem", dest="neumann", action="store_true")
-parser.add_argument("-e", "--error-tol", help="set error tolerance", type=float, default=10**-6,)
+parser.add_argument("-e", "--error-tol", help="set error tolerance", type=float, default=10**-8,)
 
 args = parser.parse_args()
 
@@ -101,10 +101,13 @@ def determine_gradient(V_g, u, flux):
 # Initialize the adapter according to the specific participant
 if problem is ProblemType.DIRICHLET:
     precice = Adapter(adapter_config_filename="precice-adapter-config-D.json")
-    precice_dt = precice.initialize(coupling_boundary, read_function_space=V, write_object=f_N_function)
+    precice.initialize(coupling_boundary, read_function_space=V, write_object=f_N_function)
+    precice_dt = precice.get_max_time_step_size()
 elif problem is ProblemType.NEUMANN:
     precice = Adapter(adapter_config_filename="precice-adapter-config-N.json")
-    precice_dt = precice.initialize(coupling_boundary, read_function_space=W, write_object=u_D_function)
+    precice.initialize(coupling_boundary, read_function_space=W, write_object=u_D_function)
+    precice_dt = precice.get_max_time_step_size()
+
 
 dt = Constant(0)
 dt.assign(np.min([fenics_dt, precice_dt]))
@@ -156,7 +159,7 @@ def determine_gradient(V_g, u, flux):
 # output solution and reference solution at t=0, n=0
 n = 0
 print('output u^%d and u_ref^%d' % (n, n))
-temperature_out << u_n
+temperature_out << (u_n, t)
 ref_out << u_ref
 ranks << mesh_rank
 
@@ -175,16 +178,15 @@ def determine_gradient(V_g, u, flux):
 while precice.is_coupling_ongoing():
 
     # write checkpoint
-    if precice.is_action_required(precice.action_write_iteration_checkpoint()):
+    if precice.requires_writing_checkpoint():
         precice.store_checkpoint(u_n, t, n)
 
-    read_data = precice.read_data()
+    dt.assign(np.min([fenics_dt, precice_dt]))
+    read_data = precice.read_data(dt)
 
     # Update the coupling expression with the new read data
     precice.update_coupling_expression(coupling_expression, read_data)
 
-    dt.assign(np.min([fenics_dt, precice_dt]))
-
     # Compute solution u^n+1, use bcs u_D^n+1, u^n and coupling bcs
     solve(a == L, u_np1, bcs)
 
@@ -198,10 +200,11 @@ def determine_gradient(V_g, u, flux):
         # Neumann problem reads flux and writes temperature on boundary to Dirichlet problem
         precice.write_data(u_np1)
 
-    precice_dt = precice.advance(dt(0))
+    precice.advance(dt)
+    precice_dt = precice.get_max_time_step_size()
 
     # roll back to checkpoint
-    if precice.is_action_required(precice.action_read_iteration_checkpoint()):
+    if precice.requires_reading_checkpoint():
         u_cp, t_cp, n_cp = precice.retrieve_checkpoint()
         u_n.assign(u_cp)
         t = t_cp
@@ -218,7 +221,7 @@ def determine_gradient(V_g, u, flux):
         print('n = %d, t = %.2f: L2 error on domain = %.3g' % (n, t, error))
         # output solution and reference solution at t_n+1
         print('output u^%d and u_ref^%d' % (n, n))
-        temperature_out << u_n
+        temperature_out << (u_n, t)
         ref_out << u_ref
         error_out << error_pointwise
 
diff --git a/partitioned-heat-conduction/fenics/problem_setup.py b/partitioned-heat-conduction/fenics/problem_setup.py
index b85f30358..ff2996fce 100644
--- a/partitioned-heat-conduction/fenics/problem_setup.py
+++ b/partitioned-heat-conduction/fenics/problem_setup.py
@@ -33,7 +33,7 @@ def inside(self, x, on_boundary):
 
 
 def get_geometry(domain_part):
-    nx = ny = 9
+    nx = ny = 11
 
     if domain_part is DomainPart.LEFT:
         p0 = Point(x_left, y_bottom)
diff --git a/partitioned-heat-conduction/fenics/run.sh b/partitioned-heat-conduction/fenics/run.sh
index e51fb750b..cf9534cf2 100755
--- a/partitioned-heat-conduction/fenics/run.sh
+++ b/partitioned-heat-conduction/fenics/run.sh
@@ -12,10 +12,10 @@ fi
 while getopts ":dn" opt; do
   case ${opt} in
   d)
-    python3 heat.py -d --error-tol 10e-3
+    python3 heat.py -d
     ;;
   n)
-    python3 heat.py -n --error-tol 10e-3
+    python3 heat.py -n
     ;;
   *)
     usage
diff --git a/perpendicular-flap/solid-fenics/solid.py b/perpendicular-flap/solid-fenics/solid.py
index 8c571c0b6..656901d14 100644
--- a/perpendicular-flap/solid-fenics/solid.py
+++ b/perpendicular-flap/solid-fenics/solid.py
@@ -69,8 +69,9 @@ def neumann_boundary(x, on_boundary):
 precice = Adapter(adapter_config_filename="precice-adapter-config-fsi-s.json")
 
 # Initialize the coupling interface
-precice_dt = precice.initialize(coupling_boundary, read_function_space=V, write_object=V, fixed_boundary=fixed_boundary)
+precice.initialize(coupling_boundary, read_function_space=V, write_object=V, fixed_boundary=fixed_boundary)
 
+precice_dt = precice.get_max_time_step_size()
 fenics_dt = precice_dt  # if fenics_dt == precice_dt, no subcycling is applied
 # fenics_dt = 0.02  # if fenics_dt < precice_dt, subcycling is applied
 dt = Constant(np.min([precice_dt, fenics_dt]))
@@ -174,15 +175,18 @@ def avg(x_old, x_new, alpha):
 
 u_n.rename("Displacement", "")
 u_np1.rename("Displacement", "")
-displacement_out << u_n
+displacement_out << (u_n, t)
 
 while precice.is_coupling_ongoing():
 
-    if precice.is_action_required(precice.action_write_iteration_checkpoint()):  # write checkpoint
+    if precice.requires_writing_checkpoint():  # write checkpoint
         precice.store_checkpoint(u_n, t, n)
 
+    precice_dt = precice.get_max_time_step_size()
+    dt = Constant(np.min([precice_dt, fenics_dt]))
+
     # read data from preCICE and get a new coupling expression
-    read_data = precice.read_data()
+    read_data = precice.read_data(dt)
 
     # Update the point sources on the coupling boundary with the new read data
     Forces_x, Forces_y = precice.get_point_sources(read_data)
@@ -199,16 +203,14 @@ def avg(x_old, x_new, alpha):
     assert (b is not b_forces)
     solve(A, u_np1.vector(), b_forces)
 
-    dt = Constant(np.min([precice_dt, fenics_dt]))
-
     # Write new displacements to preCICE
     precice.write_data(u_np1)
 
     # Call to advance coupling, also returns the optimum time step value
-    precice_dt = precice.advance(dt(0))
+    precice.advance(dt(0))
 
     # Either revert to old step if timestep has not converged or move to next timestep
-    if precice.is_action_required(precice.action_read_iteration_checkpoint()):  # roll back to checkpoint
+    if precice.requires_reading_checkpoint():  # roll back to checkpoint
         u_cp, t_cp, n_cp = precice.retrieve_checkpoint()
         u_n.assign(u_cp)
         t = t_cp
@@ -224,6 +226,6 @@ def avg(x_old, x_new, alpha):
             displacement_out << (u_n, t)
 
 # Plot tip displacement evolution
-displacement_out << u_n
+displacement_out << (u_n, t)
 
 precice.finalize()
diff --git a/volume-coupled-diffusion/fenics/volume-coupled-diffusion.py b/volume-coupled-diffusion/fenics/volume-coupled-diffusion.py
index eeffce490..d6694ede7 100644
--- a/volume-coupled-diffusion/fenics/volume-coupled-diffusion.py
+++ b/volume-coupled-diffusion/fenics/volume-coupled-diffusion.py
@@ -48,11 +48,11 @@ def inside(self, x, on_boundary):
 u_n = interpolate(u_ini, V)
 
 
-dt = precice.initialize(AllDomain(), read_function_space=V, write_object=u_n)
+precice.initialize(AllDomain(), read_function_space=V, write_object=u_n)
 volume_term = precice.create_coupling_expression()
 f = Function(V)
 
-
+dt = precice.get_max_time_step_size()
 dt_inv = Constant(1 / dt)
 
 diffusion_source = 1
@@ -87,16 +87,17 @@ def inside(self, x, on_boundary):
 # output solution and reference solution at t=0, n=0
 n = 0
 print('output u^%d and u_ref^%d' % (n, n))
-solution_out << u_n
+solution_out << (u_n, t)
 ranks << mesh_rank
 
 while precice.is_coupling_ongoing():
 
     # write checkpoint
-    if precice.is_action_required(precice.action_write_iteration_checkpoint()):
+    if precice.requires_writing_checkpoint():
         precice.store_checkpoint(u_n, t, n)
 
-    read_data = precice.read_data()
+    dt = precice.get_max_time_step_size()
+    read_data = precice.read_data(dt)
 
     # Update the coupling expression with the new read data
     precice.update_coupling_expression(volume_term, read_data)
@@ -111,10 +112,10 @@ def inside(self, x, on_boundary):
     # Write data to preCICE according to which problem is being solved
     precice.write_data(u_np1)
 
-    dt = precice.advance(dt)
+    precice.advance(dt)
 
     # roll back to checkpoint
-    if precice.is_action_required(precice.action_read_iteration_checkpoint()):
+    if precice.requires_reading_checkpoint():
         u_cp, t_cp, n_cp = precice.retrieve_checkpoint()
         u_n.assign(u_cp)
         t = t_cp
@@ -125,7 +126,7 @@ def inside(self, x, on_boundary):
         n += 1
 
     if precice.is_time_window_complete():
-        solution_out << u_n
+        solution_out << (u_n, t)
 
 # Hold plot
 precice.finalize()