diff --git a/Two-phase-Richards/multi-patch/layered_soil_with_inner_patch/TP-R-layered_soil_with_inner_patch-realistic-different-permeabilties.py b/Two-phase-Richards/multi-patch/layered_soil_with_inner_patch/TP-R-layered_soil_with_inner_patch-realistic-different-permeabilties.py
index 395fd982591fe37d39670961fd08379635572c17..60bb3414de7063adb6cb99f5f222ff88f8f54500 100755
--- a/Two-phase-Richards/multi-patch/layered_soil_with_inner_patch/TP-R-layered_soil_with_inner_patch-realistic-different-permeabilties.py
+++ b/Two-phase-Richards/multi-patch/layered_soil_with_inner_patch/TP-R-layered_soil_with_inner_patch-realistic-different-permeabilties.py
@@ -43,10 +43,10 @@ mesh_study = False
 resolutions = {
                 # 1: 5e-5,  # h=2
                 # 2: 5e-5,  # h=1.1180
-                # 4: 5e-5,  # h=0.5590
+                # 4: 5e-5,  # h=0.51590
                 # 8: 5e-5,  # h=0.2814
                 # 16: 3e-5, # h=0.1412
-                32: 5e-6,
+                32: 2e-6,
                 # 64: 3e-6,
                 # 128: 2e-6
                 }
@@ -55,55 +55,55 @@ resolutions = {
 # The list is looped over and a simulation is run with t_0 as initial time
 #  for each element t_0 in starttimes.
 starttimes = {0: 0.0}
-timestep_size = 0.0025
-number_of_timesteps = 400
+# starttimes = {0: 0.0, 1:0.3, 2:0.6, 3:0.9}
+timestep_size = 0.001
+number_of_timesteps = 1500
 
 # LDD scheme parameters  ######################################################
+Lw1 = 0.01  # /timestep_size
+Lnw1 = Lw1
 
-Lw1 = 0.25  # /timestep_size
-Lnw1 = 0.025
+Lw2 = 0.01 # /timestep_size
+Lnw2 = Lw2
 
-Lw2 = 0.25  # /timestep_size
-Lnw2 = 0.025
+Lw3 = 0.002 # /timestep_size
+Lnw3 = 0.002
 
-Lw3 = 0.25  # /timestep_size
-Lnw3 = 0.025
+Lw4 = 0.002 # /timestep_size
+Lnw4 = 0.002
 
-Lw4 = 0.25  # /timestep_size
-Lnw4 = 0.025
+Lw5 = 0.002 # /timestep_size
+Lnw5 = 0.002
 
-Lw5 = 0.25  # /timestep_size
-Lnw5 = 0.025
+Lw6 = 0.002  # /timestep_size
+Lnw6 = 0.002
 
-Lw6 = 0.25  # /timestep_size
-Lnw6 = 0.025
+lambda12_w = 1
+lambda12_nw = 0.4
 
-lambda12_w = 40
-lambda12_nw = 40
+lambda23_w = 1
+lambda23_nw = 0.4
 
-lambda23_w = 40
-lambda23_nw = 40
+lambda24_w = 1
+lambda24_nw= 0.4
 
-lambda24_w = 40
-lambda24_nw= 40
+lambda25_w= 1
+lambda25_nw= 0.4
 
-lambda25_w= 40
-lambda25_nw= 40
+lambda34_w = 1
+lambda34_nw = 0.4
 
-lambda34_w = 40
-lambda34_nw = 40
+lambda36_w = 1
+lambda36_nw = 0.4
 
-lambda36_w = 40
-lambda36_nw = 40
+lambda45_w = 1
+lambda45_nw = 0.4
 
-lambda45_w = 40
-lambda45_nw = 40
+lambda46_w = 1
+lambda46_nw = 0.4
 
-lambda46_w = 40
-lambda46_nw = 40
-
-lambda56_w = 40
-lambda56_nw = 40
+lambda56_w = 1
+lambda56_nw = 0.4
 
 include_gravity = False
 debugflag = False
@@ -113,7 +113,7 @@ analyse_condition = False
 # when number_of_timesteps is high, it might take a long time to write all
 # timesteps to disk. Therefore, you can choose to only write data of every
 # plot_timestep_every timestep to disk.
-plot_timestep_every = 1
+plot_timestep_every = 3
 # Decide how many timesteps you want analysed. Analysed means, that
 # subsequent errors of the L-iteration within the timestep are written out.
 number_of_timesteps_to_analyse = 5
@@ -222,11 +222,11 @@ gravity_acceleration = 9.81
 # https://www.geotechdata.info/parameter/soil-porosity.html
 # Dict of the form: { subdom_num : porosity }
 porosity = {
-    1: 0.37,  #0.2,  # Clayey gravels, clayey sandy gravels
-    2: 0.37,  #0.22, # Silty gravels, silty sandy gravels
-    3: 0.2,  #0.37, # Clayey sands
-    4: 0.01,  #0.2 # Silty or sandy clay
-    5: 0.025,  #
+    1: 0.57,  #0.2,  # Clayey gravels, clayey sandy gravels
+    2: 0.57,  #0.22, # Silty gravels, silty sandy gravels
+    3: 0.005,  #0.57, # Clayey sands
+    4: 0.005,  #0.2 # Silty or sandy clay
+    5: 0.005,  #
     6: 0.005,  #
 }
 
@@ -291,10 +291,10 @@ lambda_param = {
 intrinsic_permeability = {
     1: 0.1,  # sand
     2: 0.1,  # sand, there is a range
-    3: 0.05,  #10e-2,  # clay has a range
-    4: 0.025,  #10e-3
-    5: 0.05,  #10e-2,  # clay has a range
-    6: 0.025,  #10e-3
+    3: 0.001,  #10e-2,  # clay has a range
+    4: 0.001,  #10e-3
+    5: 0.001,  #10e-2,  # clay has a range
+    6: 0.001,  #10e-3
 }
 
 # relative permeabilties
@@ -453,10 +453,12 @@ if __name__ == '__main__':
         for mesh_resolution, solver_tol in resolutions.items():
             simulation_parameter.update({"solver_tol": solver_tol})
             hlp.info(simulation_parameter["use_case"])
+            processQueue = mp.Queue()
             LDDsim = mp.Process(
                         target=hlp.run_simulation,
                         args=(
                             simulation_parameter,
+                            processQueue,
                             starttime,
                             mesh_resolution
                             )
@@ -468,3 +470,10 @@ if __name__ == '__main__':
             #     starttime=starttime,
             #     parameter=simulation_parameter
             #     )
+
+        # LDDsim.join()
+        if mesh_study:
+            simulation_output_dir = processQueue.get()
+            hlp.merge_spacetime_errornorms(isRichards=isRichards,
+                                           resolutions=resolutions,
+                                           output_dir=simulation_output_dir)