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

Contained number of grid cells in mesh.

parent 61e02bff
Branches
No related tags found
No related merge requests found
Show whitespace changes
Inline Side-by-side
DG_Approximation.py
+ 16
− 23
View file @ df58c32b
......@@ -16,7 +16,7 @@ TODO: Put basis initialization for plots in function -> Done
TODO: Contain cell length in mesh -> Done
TODO: Contain bounds in mesh -> Done
TODO: Contain interval length in mesh -> Done
TODO: Contain number of grid cells in mesh
TODO: Contain number of grid cells in mesh -> Done
TODO: Create data dict for mesh separately
TODO: Create data dict for basis separately
TODO: Refactor eval_points in do_initial_projection()
......@@ -161,7 +161,6 @@ class DGScheme:
# Unpack keyword arguments
self._wave_speed = kwargs.pop('wave_speed', 1)
self._cfl_number = kwargs.pop('cfl_number', 0.2)
self._num_grid_cells = kwargs.pop('num_grid_cells', 64)
self._final_time = kwargs.pop('final_time', 1)
self._verbose = kwargs.pop('verbose', False)
self._history_threshold = kwargs.pop('history_threshold',
......@@ -175,16 +174,13 @@ class DGScheme:
self._quadrature = kwargs.pop('quadrature', 'Gauss')
self._quadrature_config = kwargs.pop('quadrature_config', {})
self._update_scheme = kwargs.pop('update_scheme', 'SSPRK3')
polynomial_degree = kwargs.pop('polynomial_degree', 2)
self._basis = OrthonormalLegendre(polynomial_degree)
self._basis = OrthonormalLegendre(kwargs.pop('polynomial_degree', 2))
# Initialize mesh with two ghost cells on each side
left_bound = kwargs.pop('left_bound', -1)
right_bound = kwargs.pop('right_bound', 1)
self._mesh = Mesh(num_grid_cells=self._num_grid_cells,
num_ghost_cells=2, left_bound=left_bound,
right_bound=right_bound)
self._mesh = Mesh(num_grid_cells=kwargs.pop('num_grid_cells', 64),
left_bound=kwargs.pop('left_bound', -1),
right_bound=kwargs.pop('right_bound', 1),
num_ghost_cells=2)
print(len(self._mesh.cells))
print(type(self._mesh.cells))
......@@ -219,12 +215,12 @@ class DGScheme:
config=self._quadrature_config)
self._detector = getattr(Troubled_Cell_Detector, self._detector)(
config=self._detector_config, mesh=self._mesh,
wave_speed=self._wave_speed, num_grid_cells=self._num_grid_cells,
final_time=self._final_time, basis=self._basis,
init_cond=self._init_cond, quadrature=self._quadrature)
wave_speed=self._wave_speed, final_time=self._final_time,
basis=self._basis, init_cond=self._init_cond,
quadrature=self._quadrature)
self._update_scheme = getattr(Update_Scheme, self._update_scheme)(
polynomial_degree=polynomial_degree,
num_grid_cells=self._num_grid_cells, detector=self._detector,
polynomial_degree=self._basis.polynomial_degree,
num_grid_cells=self._mesh.num_grid_cells, detector=self._detector,
limiter=self._limiter)
def approximate(self, data_file):
......@@ -244,8 +240,7 @@ class DGScheme:
"""
projection = do_initial_projection(
initial_condition=self._init_cond, mesh=self._mesh,
basis=self._basis, quadrature=self._quadrature,
num_grid_cells=self._num_grid_cells)
basis=self._basis, quadrature=self._quadrature)
time_step = abs(self._cfl_number * self._mesh.cell_len /
self._wave_speed)
......@@ -324,7 +319,7 @@ class DGScheme:
projection = do_initial_projection(
initial_condition=initial_condition, mesh=self._mesh,
basis=self._basis, quadrature=self._quadrature,
num_grid_cells=self._num_grid_cells, adjustment=adjustment)
adjustment=adjustment)
return self._basis.calculate_cell_average(
projection=projection[:, 1:-1], stencil_length=stencil_length,
......@@ -332,7 +327,7 @@ class DGScheme:
def do_initial_projection(initial_condition, mesh, basis, quadrature,
num_grid_cells, adjustment=0):
adjustment=0):
"""Calculates initial projection.
Calculates a projection at time step 0 and adds ghost cells on both
......@@ -348,8 +343,6 @@ def do_initial_projection(initial_condition, mesh, basis, quadrature,
Basis used for calculation.
quadrature: Quadrature object
Quadrature used for evaluation.
num_grid_cells : int
Number of cells in the mesh. Usually exponential of 2.
adjustment: float, optional
Extent of adjustment of each evaluation point in x-direction.
Default: 0.
......@@ -365,7 +358,7 @@ def do_initial_projection(initial_condition, mesh, basis, quadrature,
output_matrix = [0]
basis_vector = basis.basis
for cell in range(num_grid_cells):
for cell in range(mesh.num_grid_cells):
new_row = []
eval_point = mesh.bounds[0] + (cell+0.5)*mesh.cell_len
......@@ -382,7 +375,7 @@ def do_initial_projection(initial_condition, mesh, basis, quadrature,
output_matrix.append(basis.inverse_mass_matrix @ new_row)
# Set ghost cells to respective value
output_matrix[0] = output_matrix[num_grid_cells]
output_matrix[0] = output_matrix[mesh.num_grid_cells]
output_matrix.append(output_matrix[1])
# print(np.array(output_matrix).shape)
......
Plotting.py
+ 9
− 14
View file @ df58c32b
......@@ -125,8 +125,8 @@ def plot_shock_tube(num_grid_cells: int, troubled_cell_history: list,
def plot_details(fine_projection: ndarray, fine_mesh: Mesh, basis: Basis,
coarse_projection: ndarray, multiwavelet_coeffs: ndarray,
num_coarse_grid_cells: int) -> None:
coarse_projection: ndarray,
multiwavelet_coeffs: ndarray) -> None:
"""Plots details of projection to coarser mesh.
Parameters
......@@ -139,11 +139,9 @@ def plot_details(fine_projection: ndarray, fine_mesh: Mesh, basis: Basis,
Basis used for calculation.
multiwavelet_coeffs : ndarray
Matrix of multiwavelet coefficients.
num_coarse_grid_cells : int
Number of cells in the coarse mesh (half the cells of the fine mesh).
Usually exponential of 2.
"""
num_coarse_grid_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)
......@@ -369,8 +367,8 @@ def plot_approximation_results(data_file: str, directory: str, plot_name: str,
def plot_results(projection: ndarray, troubled_cell_history: list,
time_history: list, mesh: Mesh, num_grid_cells: int,
wave_speed: float, final_time: float, basis: Basis,
time_history: list, mesh: Mesh, wave_speed: float,
final_time: float, basis: Basis,
quadrature: Quadrature, init_cond: InitialCondition,
colors: dict = None, coarse_projection: ndarray = None,
multiwavelet_coeffs: ndarray = None) -> None:
......@@ -392,8 +390,6 @@ def plot_results(projection: ndarray, troubled_cell_history: list,
List of value of each time step.
mesh : Mesh
Mesh for calculation.
num_grid_cells : int
Number of cells in the mesh. Usually exponential of 2.
wave_speed : float
Speed of wave in rightward direction.
final_time : float
......@@ -419,7 +415,7 @@ def plot_results(projection: ndarray, troubled_cell_history: list,
colors = _check_colors(colors)
# Plot troubled cells
plot_shock_tube(num_grid_cells, troubled_cell_history, time_history)
plot_shock_tube(mesh.num_grid_cells, troubled_cell_history, time_history)
# Determine exact and approximate solution
grid, exact = calculate_exact_solution(
......@@ -430,7 +426,7 @@ def plot_results(projection: ndarray, troubled_cell_history: list,
# Plot multiwavelet solution (fine and coarse grid)
if coarse_projection is not None:
coarse_mesh = Mesh(num_grid_cells=num_grid_cells//2,
coarse_mesh = Mesh(num_grid_cells=mesh.num_grid_cells//2,
num_ghost_cells=1, left_bound=mesh.bounds[0],
right_bound=mesh.bounds[1])
......@@ -445,9 +441,8 @@ def plot_results(projection: ndarray, troubled_cell_history: list,
colors['coarse_approx'])
# Plot multiwavelet details
num_coarse_grid_cells = num_grid_cells//2
plot_details(projection[:, 1:-1], mesh, basis, coarse_projection,
multiwavelet_coeffs, num_coarse_grid_cells)
multiwavelet_coeffs)
plot_solution_and_approx(grid, exact, approx,
colors['fine_exact'],
......@@ -469,7 +464,7 @@ def plot_results(projection: ndarray, troubled_cell_history: list,
plot_error(grid, exact, approx)
print('p =', basis.polynomial_degree)
print('N =', num_grid_cells)
print('N =', mesh.num_grid_cells)
print('maximum error =', max_error)
......
Troubled_Cell_Detector.py
+ 5
− 10
View file @ df58c32b
......@@ -30,7 +30,7 @@ class TroubledCellDetector(ABC):
"""
def __init__(self, config, init_cond, quadrature, basis, mesh,
wave_speed=1, num_grid_cells=64, final_time=1):
wave_speed=1, final_time=1):
"""Initializes TroubledCellDetector.
Parameters
......@@ -47,15 +47,12 @@ class TroubledCellDetector(ABC):
Mesh for calculation.
wave_speed : float, optional
Speed of wave in rightward direction. Default: 1.
num_grid_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.
"""
self._mesh = mesh
self._wave_speed = wave_speed
self._num_grid_cells = num_grid_cells
self._final_time = final_time
self._basis = basis
self._init_cond = init_cond
......@@ -91,14 +88,12 @@ class TroubledCellDetector(ABC):
pass
def create_data_dict(self, projection):
left_bound, right_bound = self._mesh.bounds
return {'projection': projection, 'wave_speed': self._wave_speed,
'num_grid_cells': self._num_grid_cells,
'final_time': self._final_time,
'basis': {'polynomial_degree': self._basis.polynomial_degree},
'mesh': {'num_grid_cells': self._num_grid_cells,
'left_bound': left_bound,
'right_bound': right_bound,
'mesh': {'num_grid_cells': self._mesh.num_grid_cells,
'left_bound': self._mesh.bounds[0],
'right_bound': self._mesh.bounds[1],
'num_ghost_cells': 2}
}
......@@ -234,7 +229,7 @@ class WaveletDetector(TroubledCellDetector):
super()._reset(config)
# Set additional necessary parameter
self._num_coarse_grid_cells = self._num_grid_cells//2
self._num_coarse_grid_cells = self._mesh.num_grid_cells//2
self._wavelet_projection_left, self._wavelet_projection_right \
= self._basis.multiwavelet_projection
......
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment