diff --git a/LDDsimulation/LDDsimulation.py b/LDDsimulation/LDDsimulation.py index c83752a31382eb7a88d228a035a9a6f8b839d053..e15fb49679dd65380c9c14d9220d92c43c321f1f 100644 --- a/LDDsimulation/LDDsimulation.py +++ b/LDDsimulation/LDDsimulation.py @@ -938,12 +938,12 @@ class LDDsimulation(object): file.write(exact_capillary_pressure, timestep) saturation_w = df.project(S(exact_capillary_pressure), subdomain.function_space["pressure"]["wetting"]) # saturation_w.assign(Sat_w) - saturation_w.rename("Sw", "Sw") + saturation_w.rename("Sw_exact", "Sw_exact") file.write(saturation_w, timestep) # S_nw = 1-S(exact_capillary_pressure).vector().get_local() saturation_nw = df.project(1-S(exact_capillary_pressure), subdomain.function_space["pressure"]["wetting"]) # saturation_nw.assign(S_nw) - saturation_nw.rename("Snw", "Snw") + saturation_nw.rename("Snw_exact", "Snw_exact") file.write(saturation_nw, timestep) @@ -1459,63 +1459,63 @@ class LDDsimulation(object): ) - # def post_processing(self): - # """post processing of the simulation. - # calculate - # - pc_exact and pc_num - # - absolute differences to exact solution if present. - # - relative errors to exact solution if present - # """ - # for subdom_ind, subdomain in self.subdomain.items(): - # self._mesh = df.Mesh() - # solution_file = self.solution_file[subdom_ind] - # post_solution_file = self.postprocessed_solution_file[subdom_ind] - # pressure = dict() - # # the internal time series counter in df.XDMFFile starts at 0 and - # # refers to the first saved value. This corresponds to the initial - # # values. We then calculated self.number_of_timesteps many timesteps. - # t = self.starttime - # for internal_timestep in range(self.number_of_timesteps+1): - # for phase in subdomain.has_phases: - # pressure.update( - # {phase: df.Function(subdomain.function_space["pressure"][phase])} - # ) - # phase_pressure_string = "pressure_"+"{}".format(phase) - # solution_file.read_checkpoint( - # u=pressure[phase], - # name=phase_pressure_string, - # counter=internal_timestep - # ) - # pressure[phase].rename("pressure_"+"{}".format(phase), "pressure_"+"{}".format(phase)) - # post_solution_file.write(pressure[phase], t) - # t += self.timestep_size - # print(f"read pressure of {phase}phase:") - # print(pressure[phase].vector()[:]) - # solution_file.close() - # # # if we have an exact solution, write out the |u -uh|_L2 to see the - # # # absolute error. - # # if subdomain.pressure_exact is not None: - # # relative_L2_errornorm = dict() - # # for phase in subdomain.has_phases: - # # pa_exact = subdomain.pressure_exact[phase] - # # pa_exact.t = self.t - # # norm_pa_exact = df.norm(pa_exact, norm_type='L2', mesh=subdomain.mesh) - # # error_calculated = df.errornorm(pa_exact, subdomain.pressure[phase], 'L2', degree_rise=6) - # # if norm_pa_exact > self.tol: - # # relative_L2_errornorm.update( - # # {phase: error_calculated/norm_pa_exact} - # # ) - # # else: - # # relative_L2_errornorm.update( - # # {phase: error_calculated} - # # ) - # # errornorm_filename = self.output_dir+self.output_filename_parameter_part[subdom_ind]+\ - # # "_L2_errornorms_over_time" +".csv" - # # self.write_errornorms_to_csv( - # # filename = errornorm_filename, # - # # subdomain_index = subdom_ind, - # # errors = relative_L2_errornorm, - # # ) + def post_processing(self): + """post processing of the simulation. + calculate + - pc_exact and pc_num + - absolute differences to exact solution if present. + - relative errors to exact solution if present + """ + for subdom_ind, subdomain in self.subdomain.items(): + self._mesh = df.Mesh() + solution_file = self.solution_file[subdom_ind] + post_solution_file = self.postprocessed_solution_file[subdom_ind] + pressure = dict() + # the internal time series counter in df.XDMFFile starts at 0 and + # refers to the first saved value. This corresponds to the initial + # values. We then calculated self.number_of_timesteps many timesteps. + t = self.starttime + for internal_timestep in range(self.number_of_timesteps+1): + for phase in subdomain.has_phases: + pressure.update( + {phase: df.Function(subdomain.function_space["pressure"][phase])} + ) + phase_pressure_string = "pressure_"+"{}".format(phase) + solution_file.read_checkpoint( + u=pressure[phase], + name=phase_pressure_string, + counter=internal_timestep + ) + pressure[phase].rename("pressure_"+"{}".format(phase), "pressure_"+"{}".format(phase)) + post_solution_file.write(pressure[phase], t) + t += self.timestep_size + print(f"read pressure of {phase}phase:") + print(pressure[phase].vector()[:]) + solution_file.close() + # # if we have an exact solution, write out the |u -uh|_L2 to see the + # # absolute error. + # if subdomain.pressure_exact is not None: + # relative_L2_errornorm = dict() + # for phase in subdomain.has_phases: + # pa_exact = subdomain.pressure_exact[phase] + # pa_exact.t = self.t + # norm_pa_exact = df.norm(pa_exact, norm_type='L2', mesh=subdomain.mesh) + # error_calculated = df.errornorm(pa_exact, subdomain.pressure[phase], 'L2', degree_rise=6) + # if norm_pa_exact > self.tol: + # relative_L2_errornorm.update( + # {phase: error_calculated/norm_pa_exact} + # ) + # else: + # relative_L2_errornorm.update( + # {phase: error_calculated} + # ) + # errornorm_filename = self.output_dir+self.output_filename_parameter_part[subdom_ind]+\ + # "_L2_errornorms_over_time" +".csv" + # self.write_errornorms_to_csv( + # filename = errornorm_filename, # + # subdomain_index = subdom_ind, + # errors = relative_L2_errornorm, + # ) def _init_exact_solution_expression(self):