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

Moved inverse mass matrix to 'Basis_Function'.

parent f05ef836
No related branches found
No related tags found
No related merge requests found
...@@ -22,6 +22,8 @@ class Vector: ...@@ -22,6 +22,8 @@ class Vector:
Array of basis. Array of basis.
wavelet : ndarray wavelet : ndarray
Array of wavelet. Array of wavelet.
inv_mass: ndarray
Inverse mass matrix.
Methods Methods
------- -------
...@@ -29,6 +31,8 @@ class Vector: ...@@ -29,6 +31,8 @@ class Vector:
Returns basis vector. Returns basis vector.
get_wavelet_vector get_wavelet_vector
Returns wavelet vector. Returns wavelet vector.
get_inverse_mass_matrix()
Returns inverse mass matrix.
get_basis_projections() get_basis_projections()
Returns basis projections. Returns basis projections.
get_wavelet_projections() get_wavelet_projections()
...@@ -47,6 +51,7 @@ class Vector: ...@@ -47,6 +51,7 @@ class Vector:
self._polynomial_degree = polynomial_degree self._polynomial_degree = polynomial_degree
self._basis = self._build_basis_vector(x) self._basis = self._build_basis_vector(x)
self._wavelet = self._build_wavelet_vector(z) self._wavelet = self._build_wavelet_vector(z)
self._inv_mass = self._build_inverse_mass_matrix()
def get_basis_vector(self): def get_basis_vector(self):
"""Returns basis vector.""" """Returns basis vector."""
...@@ -88,6 +93,21 @@ class Vector: ...@@ -88,6 +93,21 @@ class Vector:
""" """
return [] return []
def get_inverse_mass_matrix(self):
"""Returns inverse mass matrix."""
return self._inv_mass
def _build_inverse_mass_matrix(self):
"""Constructs inverse mass matrix.
Returns
-------
ndarray
Inverse mass matrix.
"""
pass
def get_basis_projections(self): def get_basis_projections(self):
"""Returns basis projection.""" """Returns basis projection."""
pass pass
...@@ -238,6 +258,18 @@ class OrthonormalLegendre(Legendre): ...@@ -238,6 +258,18 @@ class OrthonormalLegendre(Legendre):
raise ValueError('Invalid value: Alpert\'s wavelet is only available \ raise ValueError('Invalid value: Alpert\'s wavelet is only available \
up to degree 4 for this application') up to degree 4 for this application')
def _build_inverse_mass_matrix(self):
mass_matrix = []
for i in range(self._polynomial_degree+1):
new_row = []
for j in range(self._polynomial_degree+1):
new_entry = 0.0
if i == j:
new_entry = 1.0
new_row.append(new_entry)
mass_matrix.append(new_row)
return np.array(mass_matrix)
def get_basis_projections(self): def get_basis_projections(self):
"""Returns basis projection. """Returns basis projection.
......
...@@ -4,7 +4,7 @@ ...@@ -4,7 +4,7 @@
Urgent: Urgent:
TODO: Extract do_initial_projection() from DGScheme -> Done TODO: Extract do_initial_projection() from DGScheme -> Done
TODO: Move inverse mass matrix to basis TODO: Move inverse mass matrix to basis -> Done
TODO: Extract calculate_cell_average() from TCD TODO: Extract calculate_cell_average() from TCD
TODO: Extract calculate_[...]_solution() from Plotting TODO: Extract calculate_[...]_solution() from Plotting
TODO: Extract plotting from TCD completely TODO: Extract plotting from TCD completely
...@@ -364,18 +364,6 @@ def do_initial_projection(initial_condition, basis, quadrature, ...@@ -364,18 +364,6 @@ def do_initial_projection(initial_condition, basis, quadrature,
number of grid cells and p being the polynomial degree. number of grid cells and p being the polynomial degree.
""" """
# Set inverse mass matrix
mass_matrix = []
for i in range(polynomial_degree+1):
new_row = []
for j in range(polynomial_degree+1):
new_entry = 0.0
if i == j:
new_entry = 1.0
new_row.append(new_entry)
mass_matrix.append(new_row)
inv_mass = np.array(mass_matrix)
# Initialize matrix and set first entry to accommodate for ghost cell # Initialize matrix and set first entry to accommodate for ghost cell
output_matrix = [0] output_matrix = [0]
basis_vector = basis.get_basis_vector() basis_vector = basis.get_basis_vector()
...@@ -386,19 +374,16 @@ def do_initial_projection(initial_condition, basis, quadrature, ...@@ -386,19 +374,16 @@ def do_initial_projection(initial_condition, basis, quadrature,
eval_point = left_bound + (cell+0.5)*cell_len eval_point = left_bound + (cell+0.5)*cell_len
for degree in range(polynomial_degree + 1): for degree in range(polynomial_degree + 1):
new_entry = sum( new_row.append(np.float64(sum(initial_condition.calculate(
initial_condition.calculate(
eval_point + cell_len/2 eval_point + cell_len/2
* quadrature.get_eval_points()[point] * quadrature.get_eval_points()[point] - adjustment)
- adjustment)
* basis_vector[degree].subs( * basis_vector[degree].subs(
x, quadrature.get_eval_points()[point]) x, quadrature.get_eval_points()[point])
* quadrature.get_weights()[point] * quadrature.get_weights()[point]
for point in range(quadrature.get_num_points())) for point in range(quadrature.get_num_points()))))
new_row.append(np.float64(new_entry))
new_row = np.array(new_row) new_row = np.array(new_row)
output_matrix.append(inv_mass @ new_row) output_matrix.append(basis.get_inverse_mass_matrix() @ new_row)
# Set ghost cells to respective value # Set ghost cells to respective value
output_matrix[0] = output_matrix[num_grid_cells] output_matrix[0] = output_matrix[num_grid_cells]
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment