From ea87ad0da288b24780915ee5f0eca240546ce336 Mon Sep 17 00:00:00 2001
From: David Seus <david.seus@ians.uni-stuttgart.de>
Date: Mon, 9 Sep 2019 11:43:21 +0200
Subject: [PATCH] fix meshstudy implementation
---
LDDsimulation/LDDsimulation.py | 58 +++++++++----------
LDDsimulation/domainPatch.py | 1 -
...TP-TP-2-patch-pure-dd-convergence-study.py | 58 +++++++++++++------
3 files changed, 66 insertions(+), 51 deletions(-)
diff --git a/LDDsimulation/LDDsimulation.py b/LDDsimulation/LDDsimulation.py
index d3c8c43..f9e7794 100644
--- a/LDDsimulation/LDDsimulation.py
+++ b/LDDsimulation/LDDsimulation.py
@@ -120,8 +120,6 @@ class LDDsimulation(object):
# self.calc_tol = self.tol
# list of timesteps that get analyed. Gets initiated by self._init_analyse_timesteps
self.analyse_timesteps = None
- # The number of subdomains are counted by self.init_meshes_and_markers()
- self._number_of_subdomains = 0
# variable to check if self.set_parameters() has been called
self._parameters_set = False
# flag holds true if self.initialise() has been executed.
@@ -257,7 +255,6 @@ class LDDsimulation(object):
self.saturation = saturation
# simulation start time and time variable
self.starttime = starttime
- self.t = starttime
# total number of timesteps.
self.number_of_timesteps = number_of_timesteps
self.number_of_timesteps_to_analyse = number_of_timesteps_to_analyse
@@ -287,6 +284,8 @@ class LDDsimulation(object):
"""
if not self._parameters_set:
raise RuntimeError('Parameters note set!\nYou forgott to run self.set_parameters(**kwds)')
+ # The number of subdomains are counted by self.init_meshes_and_markers()
+ self._number_of_subdomains = 0
self._add_parameters_to_output_dirname()
self._init_analyse_timesteps()
df.info(colored("initialise meshes and marker functions ...\n", "yellow"))
@@ -349,6 +348,7 @@ class LDDsimulation(object):
df.info(colored("Start time stepping","green"))
# write initial values to file.
self.timestep_num = int(1)
+ self.t = self.starttime
# note that at this point of the code self.t = self.starttime as set in
# the beginning.
@@ -554,7 +554,7 @@ class LDDsimulation(object):
# of the subdomain for communication.
subdomain.write_pressure_to_interfaces()
- if analyse_timestep:
+ if analyse_timestep and self.write_to_file['L_iterations_per_timestep']:
subdomain.write_solution_to_xdmf(#
file = solution_over_iteration_within_timestep[sd_index], #
time = time,#
@@ -1226,36 +1226,29 @@ class LDDsimulation(object):
subdomain.pressure_prev_iter[phase].assign(pa0_interpolated)
subdomain.pressure_prev_timestep[phase].assign(pa0_interpolated)
# write zeros to the abs_diff variables
+ if self.write2file['absolute_differences']:
+ absolute_difference = df.Function(subdomain.function_space["pressure"][phase])
+ # pressure_difference = pressure_exact[phase].vector()[:] - subdomain.pressure[phase].vector()[:]
+ # abs_diff_tmp = np.fabs(pressure_difference)
+ # absolute_difference.vector()[:] = abs_diff_tmp
+ if phase == "wetting":
+ data_set_name="abs_diff_w"
+ else:
+ data_set_name="abs_diff_nw"
+ self.write_function_to_xdmf(#
+ function=absolute_difference,
+ file=self.solution_file[num], #
+ time=self.starttime,#
+ data_set_name=data_set_name
+ )
- absolute_difference = df.Function(subdomain.function_space["pressure"][phase])
- # pressure_difference = pressure_exact[phase].vector()[:] - subdomain.pressure[phase].vector()[:]
- # abs_diff_tmp = np.fabs(pressure_difference)
- # absolute_difference.vector()[:] = abs_diff_tmp
- if phase == "wetting":
- data_set_name="abs_diff_w"
- else:
- data_set_name="abs_diff_nw"
- self.write_function_to_xdmf(#
- function=absolute_difference,
- file=self.solution_file[num], #
- time=self.starttime,#
- data_set_name=data_set_name
+ if self.write2file['solutions']:
+ subdomain.write_solution_to_xdmf(#
+ file = self.solution_file[num], #
+ time = self.starttime,#
+ write_iter_for_fixed_time = False,
)
- # this is also being done in the self.prepare_LDDsolver method!!!!
- # # populate interface dictionaries with pressure values. The calc_gli_term
- # # of each subdomain reads from the interface.pressure_values_prev
- # # and from interface.pressure_values dictionary so both need to be populated.
- # if self.interface_def_points is not None:
- # subdomain.write_pressure_to_interfaces()
- # subdomain.write_pressure_to_interfaces(previous_iter = True)
-
- subdomain.write_solution_to_xdmf(#
- file = self.solution_file[num], #
- time = self.starttime,#
- write_iter_for_fixed_time = False,
- )
-
def _init_DirichletBC_dictionary(self, interpolation_degree: int = None):
"""initialise dictionary that holds the dolfin.DirichletBC objects
@@ -1494,7 +1487,8 @@ class LDDsimulation(object):
def _init_simulation_output(self):
"""set up dictionary for simulation output"""
self.output = dict()
- for subdom_ind, subdomain in self.subdomain.items():
+ for subdom_ind in self.isRichards.keys():
+ subdomain = self.subdomain[subdom_ind]
self.output.update({subdom_ind: dict()})
self.output[subdom_ind].update({'mesh_size': subdomain.mesh_size})
self.output[subdom_ind].update({'errornorm': dict()})
diff --git a/LDDsimulation/domainPatch.py b/LDDsimulation/domainPatch.py
index b1c8b11..9b02d78 100644
--- a/LDDsimulation/domainPatch.py
+++ b/LDDsimulation/domainPatch.py
@@ -777,7 +777,6 @@ class DomainPatch(df.SubDomain):
This method sets the class variable
self.function_space["pressure"]
self.trialfunction
- self.pressure
self.testfunction
INPUT
diff --git a/TP-TP-2-patch-pure-dd-avoid-interface-at-origin/mesh_study_convergence/TP-TP-2-patch-pure-dd-convergence-study.py b/TP-TP-2-patch-pure-dd-avoid-interface-at-origin/mesh_study_convergence/TP-TP-2-patch-pure-dd-convergence-study.py
index 4e4aefc..2336cf7 100755
--- a/TP-TP-2-patch-pure-dd-avoid-interface-at-origin/mesh_study_convergence/TP-TP-2-patch-pure-dd-convergence-study.py
+++ b/TP-TP-2-patch-pure-dd-avoid-interface-at-origin/mesh_study_convergence/TP-TP-2-patch-pure-dd-convergence-study.py
@@ -24,16 +24,17 @@ solver_tol = 6E-7
max_iter_num = 1000
FEM_Lagrange_degree = 1
mesh_study = True
+resolutions = [5, 10, 15, 20, 25, 30, 35, 40, 45, 50]
############ GRID #######################
# mesh_resolution = 20
timestep_size = 0.0001
-number_of_timesteps = 100
+number_of_timesteps = 200
plot_timestep_every = 1
# decide how many timesteps you want analysed. Analysed means, that we write out
# subsequent errors of the L-iteration within the timestep.
number_of_timesteps_to_analyse = 0
-starttime = 0
+starttime = 0.0
Lw = 0.25 #/timestep_size
Lnw=Lw
@@ -47,6 +48,27 @@ analyse_condition = False
output_string = "./output/{}-{}_timesteps{}_P{}-solver_tol{}".format(datestr, use_case, number_of_timesteps, FEM_Lagrange_degree, solver_tol)
+# toggle what should be written to files
+if mesh_study:
+ write_to_file = {
+ 'meshes_and_markers': True,
+ 'L_iterations_per_timestep': False,
+ 'solutions': False,
+ 'absolute_differences': False,
+ 'condition_numbers': analyse_condition,
+ 'subsequent_errors': False
+ }
+else:
+ write_to_file = {
+ 'meshes_and_markers': True,
+ 'L_iterations_per_timestep': False,
+ 'solutions': True,
+ 'absolute_differences': True,
+ 'condition_numbers': analyse_condition,
+ 'subsequent_errors': True
+ }
+
+
##### Domain and Interface ####
# global simulation domain domain
sub_domain0_vertices = [df.Point(-1.0,-1.0), #
@@ -513,39 +535,39 @@ for subdomain in isRichards.keys():
#
# sa
+# toggle what should be written to files
if mesh_study:
write_to_file = {
'meshes_and_markers': True,
- 'L_iterations': True,
+ 'L_iterations_per_timestep': False,
'solutions': False,
'absolute_differences': False,
- 'condition_numbers': False,
+ 'condition_numbers': analyse_condition,
'subsequent_errors': False
}
else:
write_to_file = {
'meshes_and_markers': True,
- 'L_iterations': True,
+ 'L_iterations_per_timestep': False,
'solutions': True,
'absolute_differences': True,
- 'condition_numbers': True,
+ 'condition_numbers': analyse_condition,
'subsequent_errors': True
}
-# initialise LDD simulation class
-simulation = ldd.LDDsimulation(
- tol=1E-14,
- LDDsolver_tol=solver_tol,
- debug=debugflag,
- max_iter_num=max_iter_num,
- FEM_Lagrange_degree=FEM_Lagrange_degree,
- mesh_study=mesh_study
- )
-
-resolutions = [5, 10, 15, 20] #, 20, 25, 30, 35, 40, 45, 50]
-
for mesh_resolution in resolutions:
+ use_case = use_case + "-mesh-res_{}".format(mesh_resolution)
+ # initialise LDD simulation class
+ simulation = ldd.LDDsimulation(
+ tol=1E-14,
+ LDDsolver_tol=solver_tol,
+ debug=debugflag,
+ max_iter_num=max_iter_num,
+ FEM_Lagrange_degree=FEM_Lagrange_degree,
+ mesh_study=mesh_study
+ )
+
simulation.set_parameters(use_case=use_case,
output_dir=output_string,
subdomain_def_points=subdomain_def_points,
--
GitLab