From 71dc3aca819d56fac1053c04170fb85c9075df2a Mon Sep 17 00:00:00 2001
From: David Seus <david.seus@ians.uni-stuttgart.de>
Date: Mon, 15 Apr 2019 16:56:24 +0200
Subject: [PATCH] fix _init_solution_files

---
 LDDsimulation/LDDsimulation.py | 54 +++++++++++++++++++---------------
 1 file changed, 30 insertions(+), 24 deletions(-)

diff --git a/LDDsimulation/LDDsimulation.py b/LDDsimulation/LDDsimulation.py
index 01ca7f6..fb412f3 100644
--- a/LDDsimulation/LDDsimulation.py
+++ b/LDDsimulation/LDDsimulation.py
@@ -125,9 +125,9 @@ class LDDsimulation(object):
         # These are interpolated and stored to each respective subdomain by
         # self._init_initial_values()
         self.initial_conditions = None
-        # dictionary of filenames for writing the solution to xdmf files for each
+        # dictionary of xdmf files for writing the solution for each
         # subdomain. This gets populated by self._init_solution_files()
-        self.solution_filename = None
+        self.solution_file = None
         ########## SOLVER DEFINITION AND PARAMETERS ########
         # print("\nLinear Algebra Backends:")
         # df.list_linear_algebra_backends()
@@ -214,6 +214,8 @@ class LDDsimulation(object):
         self._init_interfaces()
         print("initialise subdomains ...\n")
         self._init_subdomains()
+        print("creating xdmf files for each subdomain ...\n")
+        self._init_solution_files()
         # initial values get initialised here to be able to call the solver at different
         # timesteps without having to call self.run().
         self._init_initial_values()
@@ -246,8 +248,8 @@ class LDDsimulation(object):
                   "Initialising simulation now ... ")
             self.initialise()
 
+        # write initial values to file.
 
-        #
 
     def LDDsolver(self, time: float = None, #
                   debug: bool = False, #
@@ -486,35 +488,32 @@ class LDDsimulation(object):
         each subdomain.
         """
         # set up solution file names for each subdomain.
-        self.solution_filename = dict()
+        self.solution_file = dict()
         for subdom_ind, subdomain in self.subdomain.items():
             tau = subdomain.timestep_size
-            L = subodmain.L
+            L = subdomain.L
             h = subdomain.mesh.hmin()
             if subdomain.isRichards:
                 Lam = subdomain.lambda_param['wetting']
-                filename = self.output_dir+"solution_subdomain" + "{}_".format(subdom_ind)#
-                                           +"dt" + "{}".format(tau)#
-                                           # only show three digits of the mesh size.
-                                           +"_hmin" + "{number:.{digits}f}".format(number=h, digits=3)#
-                                           +"_lambda" + "{}".format(Lam)
-                                           + "_Lw" + "{}".format(L["wetting"])#
-                                           +".xdmf"
+                name1 ="solution_subdomain" + "{}_".format(subdom_ind)
+                # only show three digits of the mesh size.
+                name2 = "dt" + "{}".format(tau)+"_hmin" + "{number:.{digits}f}".format(number=h, digits=3)
+                name3 = "_lambda" + "{}".format(Lam) + "_Lw" + "{}".format(L["wetting"])
+                filename = self.output_dir + name1 + name2 + name3 +".xdmf"
             else:
                 Lam_w = subdomain.lambda_param['wetting']
                 Lam_nw = subdomain.lambda_param['nonwetting']
-                filename = self.output_dir+"solution_subdomain" + "{}_".format(subdom_ind)#
-                                           +"dt" + "{}".format(tau)#
-                                           # only show three digits of the mesh size.
-                                           +"_hmin" + "{number:.{digits}f}".format(number=h, digits=3)#
-                                           +"_lambda_w" + "{}".format(Lam_w)#
-                                           +"_lambda_nw" + "{}".format(Lam_nw)#
-                                           +"_Lw" + "{}".format(L["wetting"])#
-                                           +"_Lnw" + "{}".format(L["nonwetting"])
+                filename = self.output_dir+"solution_subdomain" + "{}_".format(subdom_ind)\
+                                           +"dt" + "{}".format(tau)\
+                                           +"_hmin" + "{number:.{digits}f}".format(number=h, digits=3)\
+                                           +"_lambda_w" + "{}".format(Lam_w)\
+                                           +"_lambda_nw" + "{}".format(Lam_nw)\
+                                           +"_Lw" + "{}".format(L["wetting"])\
+                                           +"_Lnw" + "{}".format(L["nonwetting"])\
                                            +".xdmf"
-
-            self.solution_filename.update( #
-                {ind: SolutionFile(self._mpi_communicator, filename)}
+            # create solution files and populate dictionary.
+            self.solution_file.update( #
+                {subdom_ind: SolutionFile(self._mpi_communicator, filename)}
                 )
 
 
@@ -771,6 +770,8 @@ class LDDsimulation(object):
 
         This method populates the subdomain dictionaries as well as the interface
         dictionaries with initial values needed for the LDDsolver iterations.
+        Additionally, the initial value gets written to the solution file of each
+        subdomain.
         """
         if interpolation_degree is None:
             interpolation_degree = self.interpolation_degree
@@ -782,7 +783,7 @@ class LDDsimulation(object):
             p0 = self.initial_conditions[num]
             for phase in subdomain.has_phases:
                 # note that the default interpolation degree is 2
-                pa0 = df.Expression(p0[phase], domain = mesh, degree = interpolation_degree)
+                pa0 = df.Expression(p0[phase], domain = mesh, degree = interpolation_degree, t=self.t)
                 pa0_interpolated = df.interpolate(pa0, V[phase])
                 # alocate memory
                 subdomain.pressure.update(
@@ -800,6 +801,11 @@ class LDDsimulation(object):
 
             # populate interface dictionaries with pressure values.
             subdomain.write_pressure_to_interfaces()
+            subdomain.write_solution_to_xdmf(#
+                file = self.solution_file[num], #
+                time = self.t,#
+                write_iter_for_fixed_time = False,
+                )
 
     def update_DirichletBC_dictionary(self, subdomain_index: int,#
                         interpolation_degree: int = None, #
-- 
GitLab