diff --git a/DG_Approximation.py b/DG_Approximation.py index 927b17f226f3b027af85cce04eb7c68d241a03a0..c83297ad0c6514b6079e900ca7e559c3d3e6d681 100644 --- a/DG_Approximation.py +++ b/DG_Approximation.py @@ -51,8 +51,8 @@ from Vectors_of_Polynomials import OrthonormalLegendre x = Symbol('x') xi = Symbol('z') -class DGScheme(object): +class DGScheme(object): def __init__(self, detector, **kwargs): # Unpack keyword arguments self.wave_speed = kwargs.pop('wave_speed', 1) @@ -128,9 +128,9 @@ class DGScheme(object): time_step = abs(self.cfl_number * self.cell_len / self.wave_speed) current_time = 0 - while (current_time < self.final_time): + while current_time < self.final_time: # Adjust for last cell - if (current_time+time_step > self.final_time): + if current_time+time_step > self.final_time: time_step = self.final_time-current_time self.cfl_number = self.wave_speed * time_step \ / self.num_grid_cells @@ -234,7 +234,7 @@ class DGScheme(object): new_row = [] for j in range(self.polynom_degree+1): new_entry = 0.0 - if (i == j): + if i == j: new_entry = 1.0 new_row.append(new_entry) mass_matrix.append(new_row) diff --git a/Initial_Condition.py b/Initial_Condition.py index c64a426aee963a37213fcf72b9fcaed8a3f14d21..3976d823b6dcb5d69b2adf44a688ccc6fe55bd84 100644 --- a/Initial_Condition.py +++ b/Initial_Condition.py @@ -6,21 +6,23 @@ TODO: Rename initial conditions -> Renamed returned names """ import numpy as np -class InitialCondition(object): +class InitialCondition(object): def __init__(self, left_bound, right_bound, config): self.left_bound = left_bound self.right_bound = right_bound self.interval_len = self.right_bound-self.left_bound + + self.function_name = 'None' pass def get_name(self): return self.function_name def calculate(self, x): - while (x < self.left_bound): + while x < self.left_bound: x = x + self.interval_len - while (x > self.right_bound): + while x > self.right_bound: x = x - self.interval_len return self._get_point(x) @@ -29,7 +31,6 @@ class InitialCondition(object): class InitialCondition1(InitialCondition): - def __init__(self, left_bound, right_bound, config): super().__init__(left_bound, right_bound, config) # Set name of function @@ -43,7 +44,6 @@ class InitialCondition1(InitialCondition): class InitialCondition2(InitialCondition): - def __init__(self, left_bound, right_bound, config): super().__init__(left_bound, right_bound, config) # Set name of function @@ -51,11 +51,11 @@ class InitialCondition2(InitialCondition): pass def _get_point(self, x): - if (x < -1): + if x < -1: x = x + 2 - if (x > 1): + if x > 1: x = x - 2 - if ((x >= -0.5) & (x <= 0.5)): + if (x >= -0.5) & (x <= 0.5): return 1 else: return 0 @@ -76,15 +76,15 @@ class InitialCondition3(InitialCondition): pass def _get_point(self, x): - if ((x >= -0.8) & (x <= -0.6)): + if (x >= -0.8) & (x <= -0.6): return 1/6 * (self._G(x, self.z-self.delta) + self._G(x, self.z+self.delta) + 4 * self._G(x, self.z)) - if ((x >= -0.4) & (x <= -0.2)): + if (x >= -0.4) & (x <= -0.2): return 1 - if ((x >= 0) & (x <= 0.2)): + if (x >= 0) & (x <= 0.2): return 1 - abs(10 * (x-0.1)) - if ((x >= 0.4) & (x <= 0.6)): + if (x >= 0.4) & (x <= 0.6): return 1/6 * (self._F(x, self.a-self.delta) + self._F(x, self.a+self.delta) + 4 * self._F(x, self.a)) diff --git a/Limiter.py b/Limiter.py index 8d653039886af4621795aa96a7da9305db94465f..44900e83e630fb955d042e87d0f5358cebafada0 100644 --- a/Limiter.py +++ b/Limiter.py @@ -5,30 +5,30 @@ """ import numpy as np -class Limiter(object): +class Limiter(object): def __init__(self, config): + self.function_name = 'None' pass def get_name(self): return self.function_name - def apply(): + def apply(self, projection, cell): pass -class NoLimiter(Limiter): +class NoLimiter(Limiter): def __init__(self, config): # Set name of function self.function_name = 'NoLimiter' pass def apply(self, projection, cell): - return np.transpose(self.projection)[cell] + return np.transpose(projection)[cell] class MinMod(Limiter): - def __init__(self, config): # Unpack necessary parameter self.erase_degree = config.pop('erase_degree', 0) @@ -36,7 +36,7 @@ class MinMod(Limiter): # Set name of function self.function_name = 'MinMod' + str(self.erase_degree) - # Pop unncessary parameter + # Pop unnecessary parameter config.pop('cell_len') self.projection = [] @@ -54,7 +54,7 @@ class MinMod(Limiter): adapted_projection = np.transpose(self.projection)[self.cell].copy() for i in range(len(adapted_projection)): - if (i > self.erase_degree): + if i > self.erase_degree: adapted_projection[i] = 0 return adapted_projection @@ -84,7 +84,6 @@ class MinMod(Limiter): class ModifiedMinMod(MinMod): - def __init__(self, config): # Unpack necessary configurations self.cell_len = config.pop('cell_len') @@ -99,7 +98,7 @@ class ModifiedMinMod(MinMod): pass def _determine_modification(self): - if (abs(self.cell_slope) <= self.mod_factor*self.cell_len**2): + if abs(self.cell_slope) <= self.mod_factor*self.cell_len**2: return True # former: u_tilde_entry return super()._determine_modification() diff --git a/Main.py b/Main.py index dac7651e3e481d998103fe45c255818f41fd9c0f..098ba84a6bb5f11e16be14a7a548f02dc2247570 100644 --- a/Main.py +++ b/Main.py @@ -2,8 +2,8 @@ """ @author: Laura C. Kühle -Code-Stil: E226, W503 -Docstring-Stil: D200, D400 +Code-Style: E226, W503 +Docstring-Style: D200, D400 """ from DG_Approximation import DGScheme diff --git a/Quadrature.py b/Quadrature.py index f682612f7a1ddb29f6ccce2bcee4b44dbdbc8a28..73424302882337388533ec25129b4e843c7c990e 100644 --- a/Quadrature.py +++ b/Quadrature.py @@ -5,9 +5,13 @@ """ import numpy.polynomial.legendre as leg -class Quadrature(object): +class Quadrature(object): def __init__(self, config): + self.function_name = 'None' + self.eval_points = None + self.weights = None + self.num_eval_points = None pass def get_name(self): @@ -24,7 +28,6 @@ class Quadrature(object): class Gauss(Quadrature): - def __init__(self, config): # Unpack necessary configurations self.num_eval_points = config.pop('num_eval_points', 6) diff --git a/Troubled_Cell_Detector.py b/Troubled_Cell_Detector.py index 725f3587dc30ce8ece034479ce32a01eab00d83e..cc7f7926720d58a8c5e64028b65e32447dd74b4c 100644 --- a/Troubled_Cell_Detector.py +++ b/Troubled_Cell_Detector.py @@ -8,8 +8,8 @@ TODO: Check cutoff_factor/whether _is_troubled_cell works correctly \ import numpy as np -class TroubledCellDetector(object): +class TroubledCellDetector(object): def __init__(self, config): pass @@ -21,7 +21,6 @@ class TroubledCellDetector(object): class NoDetection(TroubledCellDetector): - def __init__(self, config): # Set name of function self.function_name = 'NoDetection' @@ -32,7 +31,6 @@ class NoDetection(TroubledCellDetector): class Boxplot(TroubledCellDetector): - def __init__(self, config): # Set name of function self.function_name = 'Boxplot' @@ -65,17 +63,17 @@ class Boxplot(TroubledCellDetector): balance_factor = self.fold_len/4.0 % 1 # former: int(...) boundary_index = int(self.fold_len/4.0 - balance_factor) - first_quartil = (1-balance_factor)\ + first_quartile = (1-balance_factor)\ * sorted_fold[boundary_index-1][0]\ + balance_factor * sorted_fold[boundary_index][0] - third_quartil = (1-balance_factor)\ + third_quartile = (1-balance_factor)\ * sorted_fold[3*boundary_index-1][0]\ + balance_factor * sorted_fold[3*boundary_index][0] - lower_bound = first_quartil\ - - self.whisker_len * (third_quartil-first_quartil) - upper_bound = third_quartil\ - + self.whisker_len * (third_quartil-first_quartil) + lower_bound = first_quartile\ + - self.whisker_len * (third_quartile-first_quartile) + upper_bound = third_quartile\ + + self.whisker_len * (third_quartile-first_quartile) # Check for lower extreme outliers and add respective cells for cell in sorted_fold: @@ -95,7 +93,6 @@ class Boxplot(TroubledCellDetector): class Theoretical(TroubledCellDetector): - def __init__(self, config): # Set name of function self.function_name = 'Theoretical' diff --git a/Update_Scheme.py b/Update_Scheme.py index dfb518400b808255c0805bbfcab3a55e75ef7f04..097389ad3c1e30b3b7eda3785f8dbc1da6a44e01 100644 --- a/Update_Scheme.py +++ b/Update_Scheme.py @@ -10,13 +10,7 @@ phi = DG basis vector psi = wavelet vector TODO: Find better names for A, B, M1, and M2 -TODO: Contemplate how to handle closing brackets \ - (end of line?, indented in new line?) -> Done (in line seems to be fine) -TODO: Contemplate giving option to change number of iterations for history \ - -> Done (history_threshold) -TODO: Rename counter -> Done (iteration) -TODO: Investigate whether list comprehension improves runtime \ - (compared to for-loops) -> Done (not really but more compact) + """ import numpy as np from sympy import Symbol, integrate @@ -26,8 +20,8 @@ from Vectors_of_Polynomials import OrthonormalLegendre, AlpertsWavelet x = Symbol('x') xi = Symbol('z') -class UpdateScheme(object): +class UpdateScheme(object): def __init__(self, detector, limiter, init_cond, mesh, wave_speed, polynom_degree, num_grid_cells, final_time, history_threshold, left_bound, right_bound): @@ -78,6 +72,9 @@ class UpdateScheme(object): return self.current_projection, self.troubled_cells + def _apply_stability_method(self): + pass + def _reset(self): # Set fixed basis and wavelet vectors self.basis = OrthonormalLegendre( @@ -85,6 +82,7 @@ class UpdateScheme(object): self.wavelet = AlpertsWavelet(self.polynom_degree).get_vector(x) # Set additional necessary fixed instance variables + self.name = 'None' self.interval_len = self.right_bound-self.left_bound self.cell_len = self.interval_len / self.num_grid_cells self.M1 = self._set_multiwavelet_matrix(xi, -0.5*(xi-1), True) @@ -96,7 +94,7 @@ class UpdateScheme(object): new_row = [] for j in range(self.polynom_degree+1): new_entry = -1.0 - if ((j < i) & ((i+j) % 2 == 1)): + if (j < i) & ((i+j) % 2 == 1): new_entry = 1.0 new_row.append(new_entry*np.sqrt((i+0.5) * (j+0.5))) matrix.append(new_row) @@ -132,7 +130,7 @@ class UpdateScheme(object): entry = integrate(self.basis[i].subs(x, first_param) * self.wavelet[j].subs(x, second_param), (xi, -1, 1)) - if (is_M1): + if is_M1: entry = entry * (-1)**(j + self.polynom_degree + 1) row.append(np.float64(entry)) matrix.append(row) @@ -173,7 +171,6 @@ class UpdateScheme(object): class SSPRK3(UpdateScheme): - def __init__(self, detector, limiter, init_cond, mesh, wave_speed, polynom_degree, num_grid_cells, final_time, history_threshold, left_bound, right_bound): diff --git a/Vectors_of_Polynomials.py b/Vectors_of_Polynomials.py index 4301f4b4ac5be7952cdd4d37fae7d37dd30f1c40..43fa9608d13f6a1c4f7038f741bd426ba6202db6 100644 --- a/Vectors_of_Polynomials.py +++ b/Vectors_of_Polynomials.py @@ -5,18 +5,17 @@ """ import numpy as np -class Vector(object): +class Vector(object): def __init__(self, polynom_degree): self.polynom_degree = polynom_degree pass - def get_vector(self): + def get_vector(self, eval_points): pass class Legendre(Vector): - def get_vector(self, eval_point): vector = self._calculate_legendre_vector(eval_point) return vector @@ -24,10 +23,10 @@ class Legendre(Vector): def _calculate_legendre_vector(self, eval_point): vector = [] for degree in range(self.polynom_degree+1): - if (degree == 0): + if degree == 0: vector.append(1.0 + 0*eval_point) else: - if (degree == 1): + if degree == 1: vector.append(eval_point) else: poly = (2.0*degree - 1) / degree\ @@ -38,7 +37,6 @@ class Legendre(Vector): class OrthonormalLegendre(Legendre): - def get_vector(self, eval_point): leg_vector = self._calculate_legendre_vector(eval_point) vector = [leg_vector[degree] * np.sqrt(degree+0.5) @@ -51,19 +49,19 @@ class AlpertsWavelet(Vector): def get_vector(self, eval_point): degree = self.polynom_degree - if (degree == 0): + if degree == 0: return [np.sqrt(0.5) + eval_point*0] - if (degree == 1): + if degree == 1: return [np.sqrt(1.5) * (-1 + 2*eval_point), np.sqrt(0.5) * (-2 + 3*eval_point)] - if (degree == 2): + if degree == 2: return [1/3 * np.sqrt(0.5) * (1 - 24*eval_point + 30*(eval_point**2)), 1/2 * np.sqrt(1.5) * (3 - 16*eval_point + 15*(eval_point**2)), 1/3 * np.sqrt(2.5) * (4 - 15*eval_point + 12*(eval_point**2))] - if (degree == 3): + if degree == 3: return [np.sqrt(15/34) * (1 + 4*eval_point - 30*(eval_point**2) + 28*(eval_point**3)), np.sqrt(1/42) * (-4 + 105 * eval_point @@ -72,7 +70,7 @@ class AlpertsWavelet(Vector): - 105*(eval_point**2) + 64*(eval_point**3)), 1/2 * np.sqrt(5/34) * (-16 + 105*eval_point - 192*(eval_point**2) + 105*(eval_point**3))] - if (degree == 4): + if degree == 4: return [np.sqrt(1/186) * (1 + 30*eval_point + 210*(eval_point**2) - 840*(eval_point**3) + 630*(eval_point**4)), 0.5 * np.sqrt(1/38) * (-5 - 144*eval_point @@ -90,4 +88,3 @@ class AlpertsWavelet(Vector): raise ValueError('Invalid value: Albert\'s wavelet is only available \ up to degree 4 for this application') - return 0