set up a gravity exapmle

Usecases/Two-phase-Richards/multi-patch/layered_soil/TP-R-layered_soil-g-but-same-perm-coarse-dt-longterm.py

@@ -49,13 +49,13 @@ datestr = date.strftime("%Y-%m-%d")
 
 
 # Name of the usecase that will be printed during simulation.
-use_case = "TP-R-layered-soil-realistic-g-same-intrinsic-perm"
+use_case = "TP-R-layered-soil-realistic-g-same-intrinsic-vanG-Mualem"
 # The name of this very file. Needed for creating log output.
 thisfile = "TP-R-layered_soil-g-but-same-perm-coarse-dt-longterm.py"
 
 # GENERAL SOLVER CONFIG  ######################################################
 # maximal iteration per timestep
-max_iter_num = 1000
+max_iter_num = 2000
 FEM_Lagrange_degree = 1
 
 # GRID AND MESH STUDY SPECIFICATIONS  #########################################
@@ -66,7 +66,7 @@ resolutions = {
                 # 4: 1e-6,
                 # 8: 1e-5,
                 # 16: 5e-6,
-                32: 7e-6,
+                32: 5e-5,
                 # 64: 2e-6,
                 # 128: 1e-6,
                 # 256: 1e-6,
@@ -75,27 +75,34 @@ resolutions = {
 # starttimes gives a list of starttimes to run the simulation from.
 # 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.01
-number_of_timesteps = 400
+# starttimes = {0: 0.0}
+starttimes = {
+    0: 0.0,
+    250: 0.25,
+    500: 0.5,
+    750: 0.75,
+    1000: 1.0
+    }
+timestep_size = 0.001
+number_of_timesteps = 300
 
 
 # LDD scheme parameters  ######################################################
-Lw1 = 0.05  # /timestep_size
+Lw1 = 0.025  # /timestep_size
 Lnw1 = Lw1
-Lw2 = 0.05  # /timestep_size
+Lw2 = 0.25  # /timestep_size
 Lnw2 = Lw2
-Lw3 = 0.05  # /timestep_size
-Lnw3 = Lw3
-Lw4 = 0.05  # /timestep_size
-Lnw4 = Lw4
+Lw3 = 0.25  # /timestep_size
+Lnw3 = 0.025
+Lw4 = 0.075  # /timestep_size
+Lnw4 = 0.025
 
-lambda12_w = 4
-lambda12_nw = 4
+lambda12_w = 1
+lambda12_nw = 1
 lambda23_w = 4
-lambda23_nw = 4
-lambda34_w = 4
-lambda34_nw = 4
+lambda23_nw = 1
+lambda34_w = 2
+lambda34_nw = 2
 
 include_gravity = True
 debugflag = False
@@ -108,7 +115,7 @@ analyse_condition = False
 plot_timestep_every = 1
 # 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
+number_of_timesteps_to_analyse = 7
 
 # fine grained control over data to be written to disk in the mesh study case
 # as well as for a regular simuation for a fixed grid.
@@ -267,15 +274,17 @@ intrinsic_permeability = {
 #     4: {"wetting": {"monomial": {"power": 3}},
 #         "nonwetting": {"monomial": {"power": 3}} },
 # }
+
+
 rel_perm_definition = {
     1: {"wetting": {"vanGenuchtenMualem": {"n": 3}},
         "nonwetting": {"vanGenuchtenMualem": {"n": 3}} },
     2: {"wetting": {"vanGenuchtenMualem": {"n": 3}},
         "nonwetting": {"vanGenuchtenMualem": {"n": 3}} },
-    3: {"wetting": {"vanGenuchtenMualem": {"n": 6}},
-        "nonwetting": {"vanGenuchtenMualem": {"n": 6}} },
-    4: {"wetting": {"vanGenuchtenMualem": {"n": 6}},
-        "nonwetting": {"vanGenuchtenMualem": {"n": 6}} },
+    3: {"wetting": {"vanGenuchtenMualem": {"n": 5}},
+        "nonwetting": {"vanGenuchtenMualem": {"n": 5}} },
+    4: {"wetting": {"vanGenuchtenMualem": {"n": 5}},
+        "nonwetting": {"vanGenuchtenMualem": {"n": 5}} },
 }
 
 rel_perm_dict = fts.generate_relative_permeability_dicts(rel_perm_definition)
@@ -286,8 +295,8 @@ ka_prime = rel_perm_dict["ka_prime"]
 Spc_on_subdomains = {
     1: {"vanGenuchten": {"n": 3, "alpha": 0.001}},
     2: {"vanGenuchten": {"n": 3, "alpha": 0.001}},
-    3: {"vanGenuchten": {"n": 6, "alpha": 0.001}},
-    4: {"vanGenuchten": {"n": 6, "alpha": 0.001}},
+    3: {"vanGenuchten": {"n": 5, "alpha": 0.001}},
+    4: {"vanGenuchten": {"n": 5, "alpha": 0.001}},
 }
 Spc = fts.generate_Spc_dicts(Spc_on_subdomains)
 S_pc_sym = Spc["symbolic"]
@@ -300,11 +309,18 @@ sat_pressure_relationship = Spc["dolfin"]
 x, y = sym.symbols('x[0], x[1]')  # needed by UFL
 t = sym.symbols('t', positive=True)
 
+# p_e_sym_2patch = {
+#     1: {'wetting': -6 - (1+t*t)*(1 + x*x + y*y),
+#         'nonwetting': 0.0*t},
+#     2: {'wetting': -6.0 - (1.0 + t*t)*(1.0 + x*x),
+#         'nonwetting': (-1-t*(1.1+y + x**2))*y**2},
+# }
+
 p_e_sym_2patch = {
-    1: {'wetting': -6 - (1+t*t)*(1 + x*x + y*y),
+    1: {'wetting': -7.0 - (1.0 + t*t)*(1.0 + x*y + y*y),
         'nonwetting': 0.0*t},
-    2: {'wetting': -6.0 - (1.0 + t*t)*(1.0 + x*x),
-        'nonwetting': (-1-t*(1.1+y + x**2))*y**2},
+    2: {'wetting': -7.0 - (1.0 + t*t)*(1.0 + x*y),
+        'nonwetting': (-1-t*(1.1 + y + x**2))*y**2},
 }
 
 p_e_sym = {
@@ -435,7 +451,7 @@ if __name__ == '__main__':
             #     parameter=simulation_parameter
             #     )
 
-        LDDsim.join()
+        # LDDsim.join()
         if mesh_study:
             simulation_output_dir = processQueue.get()
             hlp.merge_spacetime_errornorms(isRichards=isRichards,

Usecases/Two-phase-Richards/multi-patch/layered_soil/TP-R-layered_soil-g-but-same-perm.py

@@ -49,7 +49,7 @@ datestr = date.strftime("%Y-%m-%d")
 
 
 # Name of the usecase that will be printed during simulation.
-use_case = "TP-R-layered-soil-realistic-same-perm-vanG-Mualem-analyise-timesteps"
+use_case = "TP-R-layered-soil-realistic-same-perm-vanG-Mualem-g-test"
 # The name of this very file. Needed for creating log output.
 thisfile = "TP-R-layered_soil-g-but-same-perm.py"
 
@@ -65,8 +65,8 @@ resolutions = {
                 # 2: 1e-6,
                 # 4: 1e-6,
                 # 8: 1e-5,
-                # 16: 5e-6,
-                32: 4e-6,
+                16: 5e-6,
+                # 32: 4e-6,
                 # 64: 2e-6,
                 # 128: 1e-6,
                 # 256: 1e-6,
@@ -80,34 +80,32 @@ resolutions = {
 # timestep_num: initial_time
 # starttimes = {0: 0.0}
 starttimes = {
-    190: 0.19,
-    440: 0.44,
-    690: 0.69,
-    940: 0.94,
-    1140: 1.14,
-    # 1440: 1.440
+    0: 0.0,
+    1: 0.3,
+    2: 0.6,
+    3: 0.9,
     }
-timestep_size = 0.001
-number_of_timesteps = 60
+timestep_size = 0.0001
+number_of_timesteps = 1
 
 # LDD scheme parameters  ######################################################
 Lw1 = 0.025  # /timestep_size
 Lnw1 = Lw1
 Lw2 = 0.005  # /timestep_size
 Lnw2 = Lw2
-Lw3 = 0.005  # /timestep_size
+Lw3 = 0.025  # /timestep_size
 Lnw3 = 0.025
-Lw4 = 0.005  # /timestep_size
+Lw4 = 0.025  # /timestep_size
 Lnw4 = 0.025
 
-lambda12_w = 0.5
-lambda12_nw = 4
-lambda23_w = 0.5
-lambda23_nw = 4
-lambda34_w = 0.5
-lambda34_nw = 4
+lambda12_w = 1
+lambda12_nw = 1
+lambda23_w = 1
+lambda23_nw = 1
+lambda34_w = 1
+lambda34_nw = 1
 
-include_gravity = False
+include_gravity = True
 debugflag = False
 analyse_condition = False
 
@@ -126,14 +124,14 @@ if mesh_study:
     write_to_file = {
         # output the relative errornorm (integration in space) w.r.t. an exact
         # solution for each timestep into a csv file.
-        'space_errornorms': False,
+        'space_errornorms': True,
         # save the mesh and marker functions to disk
-        'meshes_and_markers': False,
+        'meshes_and_markers': True,
         # save xdmf/h5 data for each LDD iteration for timesteps determined by
         # number_of_timesteps_to_analyse. I/O intensive!
-        'L_iterations_per_timestep': False,
+        'L_iterations_per_timestep': True,
         # save solution to xdmf/h5.
-        'solutions': False,
+        'solutions': True,
         # save absolute differences w.r.t an exact solution to xdmf/h5 file
         # to monitor where on the domains errors happen
         'absolute_differences': True,
@@ -299,11 +297,18 @@ sat_pressure_relationship = Spc["dolfin"]
 x, y = sym.symbols('x[0], x[1]')  # needed by UFL
 t = sym.symbols('t', positive=True)
 
+# p_e_sym_2patch = {
+#     1: {'wetting': -6 - (1+t*t)*(1 + x*x + y*y),
+#         'nonwetting': 0.0*t},  # -1-t*(1.1 + y + x**2)**2},
+#     2: {'wetting': -6.0 - (1.0 + t*t)*(1.0 + x*x),
+#         'nonwetting': (-1-t*(1.1+y + x**2))*y**2},
+#         # 'nonwetting': (-1-t*(1.1 + x**2)**2 - sym.sqrt(5+t**2))*y**2},
+# }
 p_e_sym_2patch = {
-    1: {'wetting': -6 - (1+t*t)*(1 + x*x + y*y),
+    1: {'wetting': (-7.0 - (1.0 + t*t)*(1.0 + x*x + y*y)),
         'nonwetting': 0.0*t},  # -1-t*(1.1 + y + x**2)**2},
-    2: {'wetting': -6.0 - (1.0 + t*t)*(1.0 + x*x),
-        'nonwetting': (-1-t*(1.1+y + x**2))*y**2},
+    2: {'wetting': -7.0 - (1.0 + t*t)*(1.0 + x*x),
+        'nonwetting': (-2-t*(1.1+y + x**2))*y**2},
         # 'nonwetting': (-1-t*(1.1 + x**2)**2 - sym.sqrt(5+t**2))*y**2},
 }