Skip to content
Snippets Groups Projects
Commit 4ca1263b authored by Laura Christine Kühle's avatar Laura Christine Kühle
Browse files

Split use of 'mesh' and 'grid' clearly.

parent ae80061f
No related branches found
No related tags found
No related merge requests found
......@@ -49,7 +49,7 @@ class TrainingDataGenerator:
self._quadrature_list = [Gauss({'num_nodes': pol_deg+1})
for pol_deg in range(5)]
self._mesh_list = [Mesh(left_bound=left_bound, right_bound=right_bound,
num_ghost_cells=0, num_grid_cells=2**exp)
num_ghost_cells=0, num_cells=2**exp)
for exp in range(3, 12)]
def build_training_data(self, init_cond_list, num_samples, balance=0.5,
......
......@@ -19,11 +19,12 @@ TODO: Discuss adding kwargs to attributes in documentation
TODO: Discuss descriptions (matrices, cfl number, right-hand side,
limiting slope, basis, wavelet, etc.)
TODO: Discuss referencing info on SSPRK3
TODO: Discuss name for quadrature mesh (now: grid)
Urgent:
TODO: Unify use of 'length' and 'len' in naming -> Done
TODO: Unify use of 'initial_condition' and 'init_cond' in naming -> Done
TODO: Unify use of 'mesh' and 'grid' in naming
TODO: Unify use of 'mesh' and 'grid' in naming -> Done
TODO: Check sign change in stiffness matrix to accommodate negative wave speed
TODO: Force input_size for each ANN model to be stencil length
TODO: Induce shift in IC class
......@@ -45,6 +46,7 @@ TODO: Add images to report
TODO: Add verbose output
Critical, but not urgent:
TODO: Restructure 'calculate_approximate_solution()'
TODO: Rework Theoretical TC for efficiency
TODO: Extract object initialization from DGScheme
TODO: Replace loops with list comprehension if feasible
......@@ -141,7 +143,7 @@ class DGScheme:
Polynomial degree. Default: 2.
cfl_number : float, optional
CFL number to ensure stability. Default: 0.2.
num_grid_cells : int, optional
num_mesh_cells : int, optional
Number of cells in the mesh. Usually exponential of 2. Default: 64.
final_time : float, optional
Final time for which approximation is calculated. Default: 1.
......@@ -191,7 +193,7 @@ class DGScheme:
self._basis = OrthonormalLegendre(kwargs.pop('polynomial_degree', 2))
# Initialize mesh with two ghost cells on each side
self._mesh = Mesh(num_grid_cells=kwargs.pop('num_grid_cells', 64),
self._mesh = Mesh(num_cells=kwargs.pop('num_mesh_cells', 64),
left_bound=kwargs.pop('left_bound', -1),
right_bound=kwargs.pop('right_bound', 1),
num_ghost_cells=2)
......@@ -231,7 +233,7 @@ class DGScheme:
config=self._detector_config, mesh=self._mesh, basis=self._basis)
self._update_scheme = getattr(Update_Scheme, self._update_scheme)(
polynomial_degree=self._basis.polynomial_degree,
num_grid_cells=self._mesh.num_grid_cells, detector=self._detector,
num_cells=self._mesh.num_cells, detector=self._detector,
limiter=self._limiter)
def approximate(self, data_file):
......@@ -326,7 +328,7 @@ def do_initial_projection(init_cond, mesh, basis, quadrature,
-------
ndarray
Matrix containing projection of size (N+2, p+1) with N being the
number of grid cells and p being the polynomial degree of the basis.
number of mesh cells and p being the polynomial degree of the basis.
"""
# Initialize matrix and set first entry to accommodate for ghost cell
......@@ -347,7 +349,7 @@ def do_initial_projection(init_cond, mesh, basis, quadrature,
output_matrix.append(basis.inverse_mass_matrix @ new_row)
# Set ghost cells to respective value
output_matrix[0] = output_matrix[mesh.num_grid_cells]
output_matrix[0] = output_matrix[mesh.num_cells]
output_matrix.append(output_matrix[1])
# print(np.array(output_matrix).shape)
......
......@@ -94,7 +94,7 @@ def plot_error(grid: ndarray, exact: ndarray, approx: ndarray) -> None:
plt.title('Errors')
def plot_shock_tube(num_grid_cells: int, troubled_cell_history: list,
def plot_shock_tube(num_cells: int, troubled_cell_history: list,
time_history: list) -> None:
"""Plots shock tube.
......@@ -103,7 +103,7 @@ def plot_shock_tube(num_grid_cells: int, troubled_cell_history: list,
Parameters
----------
num_grid_cells : int
num_cells : int
Number of cells in the mesh. Usually exponential of 2.
troubled_cell_history : list
List of detected troubled cells for each time step.
......@@ -116,7 +116,7 @@ def plot_shock_tube(num_grid_cells: int, troubled_cell_history: list,
current_cells = troubled_cell_history[pos]
for cell in current_cells:
plt.plot(cell, time_history[pos], 'k.')
plt.xlim((0, num_grid_cells))
plt.xlim((0, num_cells))
plt.xlabel('Cell')
plt.ylabel('Time')
plt.title('Shock Tubes')
......@@ -139,16 +139,16 @@ def plot_details(fine_projection: ndarray, fine_mesh: Mesh, basis: Basis,
Matrix of multiwavelet coefficients.
"""
num_coarse_grid_cells = len(coarse_projection[0])
num_coarse_mesh_cells = len(coarse_projection[0])
averaged_projection = [[coarse_projection[degree][cell]
* basis.basis[degree].subs(x, value)
for cell in range(num_coarse_grid_cells)
for cell in range(num_coarse_mesh_cells)
for value in [-0.5, 0.5]]
for degree in range(basis.polynomial_degree + 1)]
wavelet_projection = [[multiwavelet_coeffs[degree][cell]
* basis.wavelet[degree].subs(z, 0.5) * value
for cell in range(num_coarse_grid_cells)
for cell in range(num_coarse_mesh_cells)
for value in [(-1) ** (basis.polynomial_degree
+ degree + 1), 1]]
for degree in range(basis.polynomial_degree + 1)]
......@@ -375,8 +375,8 @@ def plot_results(projection: ndarray, troubled_cell_history: list,
Plots exact and approximate solution, errors, and troubled cells of a
projection given its evaluation history.
If coarse grid and projection are given, solutions are displayed for
both coarse and fine grid. Additionally, coefficient details are plotted.
If coarse mesh and projection are given, solutions are displayed for
both coarse and fine mesh. Additionally, coefficient details are plotted.
Parameters
----------
......@@ -401,7 +401,7 @@ def plot_results(projection: ndarray, troubled_cell_history: list,
colors: dict
Dictionary of colors used for plots.
coarse_projection: ndarray, optional
Matrix of projection on coarse grid for each polynomial degree.
Matrix of projection on coarse mesh for each polynomial degree.
Default: None.
multiwavelet_coeffs: ndarray, optional
Matrix of wavelet coefficients. Default: None.
......@@ -413,7 +413,7 @@ def plot_results(projection: ndarray, troubled_cell_history: list,
colors = _check_colors(colors)
# Plot troubled cells
plot_shock_tube(mesh.num_grid_cells, troubled_cell_history, time_history)
plot_shock_tube(mesh.num_cells, troubled_cell_history, time_history)
# Determine exact and approximate solution
grid, exact = calculate_exact_solution(
......@@ -422,9 +422,9 @@ def plot_results(projection: ndarray, troubled_cell_history: list,
projection[:, 1:-1], quadrature.nodes,
basis.polynomial_degree, basis.basis)
# Plot multiwavelet solution (fine and coarse grid)
# Plot multiwavelet solution (fine and coarse mesh)
if coarse_projection is not None:
coarse_mesh = Mesh(num_grid_cells=mesh.num_grid_cells//2,
coarse_mesh = Mesh(num_cells=mesh.num_cells//2,
num_ghost_cells=1, left_bound=mesh.bounds[0],
right_bound=mesh.bounds[1])
......@@ -462,7 +462,7 @@ def plot_results(projection: ndarray, troubled_cell_history: list,
plot_error(grid, exact, approx)
print('p =', basis.polynomial_degree)
print('N =', mesh.num_grid_cells)
print('N =', mesh.num_cells)
print('maximum error =', max_error)
......
......@@ -224,7 +224,7 @@ class WaveletDetector(TroubledCellDetector):
= self._basis.multiwavelet_projection
def _calculate_wavelet_coeffs(self, projection):
"""Calculates wavelet coefficients used for projection to coarser grid.
"""Calculates wavelet coefficients used for projection to coarser mesh.
Parameters
----------
......@@ -238,7 +238,7 @@ class WaveletDetector(TroubledCellDetector):
"""
output_matrix = []
for i in range(self._mesh.num_grid_cells):
for i in range(self._mesh.num_cells):
new_entry = 0.5*(
projection[:, i+1] @ self._wavelet_projection_left
+ projection[:, i+2] @ self._wavelet_projection_right)
......@@ -282,7 +282,7 @@ class WaveletDetector(TroubledCellDetector):
Returns
-------
ndarray
Matrix of projection on coarse grid for each polynomial degree.
Matrix of projection on coarse mesh for each polynomial degree.
"""
basis_projection_left, basis_projection_right\
......@@ -293,7 +293,7 @@ class WaveletDetector(TroubledCellDetector):
# Calculate projection on coarse mesh
output_matrix = []
for i in range(self._mesh.num_grid_cells//2):
for i in range(self._mesh.num_cells//2):
new_entry = 0.5 * (
projection[:, 2 * i] @ basis_projection_left
+ projection[:, 2 * i + 1] @ basis_projection_right)
......@@ -353,18 +353,18 @@ class Boxplot(WaveletDetector):
self._adjust_outer_fences = config.pop('adjust_outer_fences', True)
self._extreme_outlier_only = config.pop('extreme_outlier_only', True)
if self._mesh.num_grid_cells < self._fold_len:
self._fold_len = self._mesh.num_grid_cells
if self._mesh.num_cells < self._fold_len:
self._fold_len = self._mesh.num_cells
self._quantile_method = config.pop('quantile_method', 'weibull')
num_overlapping_cells = config.pop('num_overlapping_cells', 1)
num_folds = self._mesh.num_grid_cells//self._fold_len
num_folds = self._mesh.num_cells//self._fold_len
self._fold_indices = np.zeros([num_folds,
self._fold_len + 2 *
num_overlapping_cells]).astype(np.int32)
for fold in range(num_folds):
self._fold_indices[fold] = np.array(
[i % self._mesh.num_grid_cells for i in range(
[i % self._mesh.num_cells for i in range(
fold * self._fold_len - num_overlapping_cells,
(fold+1) * self._fold_len + num_overlapping_cells)])
......@@ -476,9 +476,9 @@ class Theoretical(WaveletDetector):
troubled_cells = []
max_avg = np.sqrt(0.5) \
* max(1, max(abs(projection[0][cell+1])
for cell in range(self._mesh.num_grid_cells)))
for cell in range(self._mesh.num_cells)))
for cell in range(self._mesh.num_grid_cells):
for cell in range(self._mesh.num_cells):
if self._is_troubled_cell(coeffs, cell, max_avg):
troubled_cells.append(cell)
......@@ -504,6 +504,6 @@ class Theoretical(WaveletDetector):
"""
max_value = abs(coeffs[cell])/max_avg
eps = self._cutoff_factor\
/ (self._mesh.cell_len*self._mesh.num_grid_cells)
/ (self._mesh.cell_len*self._mesh.num_cells)
return max_value > eps
......@@ -26,14 +26,14 @@ class UpdateScheme(ABC):
Performs time step.
"""
def __init__(self, polynomial_degree, num_grid_cells, detector, limiter):
def __init__(self, polynomial_degree, num_cells, detector, limiter):
"""Initializes UpdateScheme.
Parameters
----------
polynomial_degree : int
Polynomial degree.
num_grid_cells : int
num_cells : int
Number of cells in the mesh. Usually exponential of 2.
detector : TroubledCellDetector object
Troubled cell detector for evaluation.
......@@ -43,7 +43,7 @@ class UpdateScheme(ABC):
"""
# Unpack positional arguments
self._polynomial_degree = polynomial_degree
self._num_grid_cells = num_grid_cells
self._num_cells = num_cells
self._detector = detector
self._limiter = limiter
......@@ -169,8 +169,8 @@ class UpdateScheme(ABC):
degree.
"""
current_projection[:, 0] = current_projection[:, self._num_grid_cells]
current_projection[:, self._num_grid_cells+1] \
current_projection[:, 0] = current_projection[:, self._num_cells]
current_projection[:, self._num_cells+1] \
= current_projection[:, 1]
return current_projection
......@@ -313,14 +313,14 @@ class SSPRK3(UpdateScheme):
# Initialize vector and set first entry to accommodate for ghost cell
right_hand_side = [0]
for j in range(self._num_grid_cells):
for j in range(self._num_cells):
right_hand_side.append(2*(self._stiffness_matrix
@ current_projection[:, j+1]
+ self._boundary_matrix
@ current_projection[:, j]))
# Set ghost cells to respective value
right_hand_side[0] = right_hand_side[self._num_grid_cells]
right_hand_side[0] = right_hand_side[self._num_cells]
right_hand_side.append(right_hand_side[1])
return np.transpose(right_hand_side)
data_dir: 'Sep29'
data_dir: 'Oct03'
random_seed: 1234
# Parameter for Approximation with Troubled Cell Detection
......@@ -8,7 +8,7 @@ Approximation:
wave_speed: 1
polynomial_degree: 2
cfl_number: 0.2
num_grid_cells: 32 # 40 elements work well for Condition 3
num_mesh_cells: 32 # 40 elements work well for Condition 3
final_time: 1
left_bound: -1
right_bound: 1
......
......@@ -20,7 +20,7 @@ x = Symbol('x')
class Mesh:
"""Class for mesh/grid.
"""Class for mesh.
Each cell is characterized by its center.
......@@ -28,7 +28,7 @@ class Mesh:
----------
mode : str
Mode for mesh use. Either 'training' or 'evaluation'.
num_grid_cells : int
num_cells : int
Number of cells in the mesh (ghost cells notwithstanding). Usually
exponential of 2.
bounds : Tuple[float, float]
......@@ -42,14 +42,14 @@ class Mesh:
"""
def __init__(self, num_grid_cells: int, num_ghost_cells: int,
def __init__(self, num_cells: int, num_ghost_cells: int,
left_bound: float, right_bound: float,
training_data_mode: bool = False) -> None:
"""Initialize Mesh.
Parameters
----------
num_grid_cells : int
num_cells : int
Number of cells in the mesh (ghost cells notwithstanding). Has
to be an exponential of 2.
num_ghost_cells : int
......@@ -63,10 +63,10 @@ class Mesh:
generation. Default: False.
"""
self._num_grid_cells = num_grid_cells
self._num_cells = num_cells
self._mode = 'training' if training_data_mode else 'evaluation'
if not training_data_mode:
if not math.log(self._num_grid_cells, 2).is_integer():
if not math.log(self._num_cells, 2).is_integer():
raise ValueError('The number of cells in the mesh has to be '
'an exponential of 2')
self._num_ghost_cells = num_ghost_cells
......@@ -79,9 +79,9 @@ class Mesh:
return self._mode
@property
def num_grid_cells(self) -> int:
"""Return number of grid cells."""
return self._num_grid_cells
def num_cells(self) -> int:
"""Return number of mesh cells."""
return self._num_cells
@property
def bounds(self) -> Tuple[float, float]:
......@@ -96,7 +96,7 @@ class Mesh:
@property
def cell_len(self) -> float:
"""Return the length of a mesh cell."""
return self.interval_len/self.num_grid_cells
return self.interval_len/self.num_cells
@property
@cache
......@@ -114,7 +114,7 @@ class Mesh:
def create_data_dict(self) -> dict:
"""Return dictionary with data necessary to construct mesh."""
return {'num_grid_cells': self._num_grid_cells,
return {'num_cells': self._num_cells,
'left_bound': self._left_bound,
'right_bound': self._right_bound,
'num_ghost_cells': self._num_ghost_cells}
......@@ -134,17 +134,17 @@ class Mesh:
# Pick random point between left and right bound
point = np.random.uniform(self._left_bound, self._right_bound)
# Adjust grid spacing to be within interval if necessary
grid_spacing = self.cell_len
# Adjust mesh spacing to be within interval if necessary
mesh_spacing = self.cell_len
max_spacing = 2/stencil_len*min(point-self._left_bound,
self._right_bound-point)
while grid_spacing > max_spacing:
grid_spacing /= 2
while mesh_spacing > max_spacing:
mesh_spacing /= 2
# Return new mesh instance
return Mesh(left_bound=point - stencil_len/2 * grid_spacing,
right_bound=point + stencil_len/2 * grid_spacing,
num_grid_cells=stencil_len, num_ghost_cells=2,
return Mesh(left_bound=point - stencil_len/2 * mesh_spacing,
right_bound=point + stencil_len/2 * mesh_spacing,
num_cells=stencil_len, num_ghost_cells=2,
training_data_mode=True)
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment