add subdomain.write_pressure_to_interfaces(initial_values = True) to...

add subdomain.write_pressure_to_interfaces(initial_values = True) to LDDsimulation._init_initial_values and add initial_values flag to DomainPatch.write_pressure_to_interfaces to reflect this accordingly.

LDDsimulation/LDDsimulation.py

@@ -182,8 +182,8 @@ class LDDsimulation(object):
 
         ### actual iteration starts here
         # gobal stopping criterion for the iteration.
-        all_subdomains_done = False
-        while iteration < self._max_iter_num and not all_subdomains_done:
+        all_subdomains_are_done = False
+        while iteration < self._max_iter_num and not all_subdomains_are_done:
             # we need to loop over the subdomains and solve an L-scheme type step
             # on each subdomain. here it should be possible to parallelise in
             # the future.
@@ -213,12 +213,12 @@ class LDDsimulation(object):
                     # one subdomain is done. Check if all subdomains are done to
                     # stop the while loop if needed.
                     # For this, set
-                    all_subdomains_done = True
+                    all_subdomains_are_done = True
                     # then check if all domains are done. If not, reset
-                    # all_subdomains_done to False.
+                    # all_subdomains_are_done to False.
                     for subdomain, isdone in subdomain_calc_isdone.items():
                         if not isdone:
-                            all_subdomains_done = False
+                            all_subdomains_are_done = False
                     #TODO check if maybe I still need to do
                     # subdomain.iteration_number += 1
                     # or adapt the calc_gli_term method to reflect that one subdomain
@@ -421,6 +421,9 @@ class LDDsimulation(object):
 
     def _init_initial_values(self, interpolation_degree: int = 2):
         """ set initial values
+
+        This method populates the subdomain dictionaries as well as the interface
+        dictionaries with initial values needed for the LDDsolver iterations.
         """
         for num, subdomain in self.subdomain.items():
             mesh = subdomain.mesh
@@ -432,6 +435,7 @@ class LDDsimulation(object):
                 pw0 = df.Expression(p0['wetting'], domain = mesh, degree = interpolation_degree)
                 subdomain.pressure_prev_iter = {'wetting' : df.interpolate(pw0, V['wetting'])}
                 subdomain.pressure_prev_timestep = {'wetting' : df.interpolate(pw0, V['wetting'])}
+
             else:
                 pw0 = df.Expression(p0['wetting'], domain = mesh, degree = interpolation_degree)
                 pnw0 = df.Expression(p0['nonwetting'], domain = mesh, degree = interpolation_degree)
@@ -439,6 +443,8 @@ class LDDsimulation(object):
                                                 'nonwetting' : df.interpolate(pnw0, V['nonwetting'])}
                 subdomain.pressure_prev_timestep = {'wetting' : df.interpolate(pw0, V['wetting']),#
                                                     'nonwetting' : df.interpolate(pnw0, V['nonwetting'])}
+            # populate interface dictionaries with pressure values.
+            subdomain.write_pressure_to_interfaces(initial_values = True)
 
     def set_DirichletBC(self, interpolation_degree: int = 2, #
                         time: float = 0,#

LDDsimulation/domainPatch.py

@@ -179,14 +179,27 @@ class DomainPatch(df.SubDomain):
     # END constructor
 
     #### PUBLIC METHODS
-    def write_pressure_to_interfaces(self):
+    def write_pressure_to_interfaces(self, initial_values: bool = False):
         """ save the interface values of self.pressure to all neighbouring interfaces
 
         This method is used by the LDDsimulation.LDDsolver method.
+
+        # parameters
+        initial_values      bool        Flag to toggle wether to write self.pressure
+                                        or self.pressure_prev_iter to the interface.
+                                        This is needed for the first ever iteration
+                                        at time = 0.
         """
+        if initial_values:
+            # this assumes, that self.pressure_prev_iter has been set by
+            # LDDsimulation._init_initial_values().
+            from_func = self.pressure_prev_iter
+        else:
+            from_func = self.pressure
+
         if self.isRichards:
             for ind in self.has_interface:
-                interface[ind].read_pressure_dofs(from_function = self.pressure['wetting'], #
+                interface[ind].read_pressure_dofs(from_function = from_func['wetting'], #
                                         interface_dofs = self._dof_indices_of_interface[ind]['wetting'],#
                                         dof_to_vert_map = self.dof2vertex['wetting'],#
                                         local_to_parent_vertex_map = self.parent_mesh_index,#
@@ -195,7 +208,7 @@ class DomainPatch(df.SubDomain):
         else:
             for ind in self.has_interface:
                 for phase in ['wetting', 'nonwetting']:
-                    interface[ind].read_pressure_dofs(from_function = self.pressure[phase], #
+                    interface[ind].read_pressure_dofs(from_function = from_func[phase], #
                                             interface_dofs = self._dof_indices_of_interface[ind][phase],#
                                             dof_to_vert_map = self.dof2vertex[phase],#
                                             local_to_parent_vertex_map = self.parent_mesh_index,#