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

Reworked details plot. Rewrote reshaping in plot methods.

parent c6eb35c7
No related branches found
No related tags found
No related merge requests found
Show whitespace changes
Inline Side-by-side
Showing
with 90 additions and 32 deletions
DG_Approximation.py
+ 90
− 32
View file @ 2fe6898c
......@@ -4,8 +4,8 @@
Plotter:
TODO: Contemplate using Seaborn instead of matplotlib
TODO: Check whether details work/why not -> Rework Figure 4
TODO: Rewrite reshaping in plot methods
TODO: Check whether details work/why not -> Rework Figure 4 -> Done
TODO: Rewrite reshaping in plot methods -> Done
TODO: Double-check everything!
TODO: Replace loops with list comprehension if feasible
......@@ -24,10 +24,15 @@ TODO: Write documentation for all methods
TODO: Add a verbose option
TODO: Check whether consistency is given/possible for each class instance
TODO: Remove unnecessary code in plot methods
TODO: Make sure time is changed to current_time everywhere -> Done
TODO: Make sure all instance variables are actually necessary
TODO: Make sure instance variables are only set in __init__()
TODO: Label plot for shock tubes -> Done
TODO: Extract calculation of coarse projection (_calculate_coarse_projection()) -> Done
TODO: Contemplate moving plots to pertaining files
TODO: Discuss details plot!
TODO: Fix typo in Vectors_of_Polynomials
"""
import numpy as np
......@@ -43,6 +48,7 @@ import Limiter
import Quadrature
import Update_Scheme
from Vectors_of_Polynomials import OrthonormalLegendre
from Vectors_of_Polynomials import AlpertsWavelet
x = Symbol('x')
xi = Symbol('z')
......@@ -130,6 +136,8 @@ class DGScheme(object):
self._plot_shock_tube()
self._plot_coarse_mesh(projection, 'k-', 'y')
self._plot_details(projection)
approx = self._plot_fine_mesh(projection, 'k-.', 'b')
# What is that??? What is it for?
......@@ -298,14 +306,13 @@ class DGScheme(object):
return grid, exact
def _plot_coarse_mesh(self, projection, color_exact, color_approx):
def _calculate_coarse_projection(self, projection):
basis_matrix_left = self._set_basis_matrix(xi, 0.5*(xi-1))
basis_matrix_right = self._set_basis_matrix(xi, 0.5*(xi+1))
coarse_cell_len = 2*self.cell_len
coarse_mesh = np.arange(self.left_bound - (0.5*coarse_cell_len), self.right_bound + (1.5*coarse_cell_len),
coarse_cell_len)
# Calculate projection on coarse mesh
################################################################################################################
# Remove later, but keep now for better time comparison
tic = timeit.default_timer()
transposed_vector = np.transpose(projection[:, 1:-1])
output_matrix = []
......@@ -318,29 +325,36 @@ class DGScheme(object):
toc = timeit.default_timer()
print('Trans:', toc-tic)
# !!!
# print()
tic = timeit.default_timer()
################################################################################################################
transposed_vector = projection[:, 1:-1]
output_matrix = []
for i in range(int(len(transposed_vector[0])/2)):
# print(transposed_vector[:, 2*i].shape)
new_entry = 0.5*(transposed_vector[:, 2*i] @ basis_matrix_left
+ transposed_vector[:, 2*i+1] @ basis_matrix_right)
output_matrix.append(new_entry)
# print(new_entry)
coarse_projection1 = np.transpose(np.array(output_matrix))
# Remove later, but keep now for better time comparison
toc = timeit.default_timer()
print('Vecto:', toc-tic)
# print(coarse_projection == coarse_projection1)
return coarse_projection
def _plot_coarse_mesh(self, projection, color_exact, color_approx):
coarse_cell_len = 2*self.cell_len
coarse_mesh = np.arange(self.left_bound - (0.5*coarse_cell_len), self.right_bound + (1.5*coarse_cell_len),
coarse_cell_len)
coarse_projection = self._calculate_coarse_projection(projection)
# Plot exact and approximate solutions for coarse mesh
grid, exact = self._calculate_exact_solution(coarse_mesh[1:-1], coarse_cell_len)
approx = self._calculate_approximate_solution(coarse_projection)
self._plot_solution_and_approx(grid, exact, approx, color_exact, color_approx)
# plt.legend(['Exact-Coarse', 'Approx-Coarse'])
def _plot_fine_mesh(self, projection, color_exact, color_approx):
grid, exact = self._calculate_exact_solution(self.mesh[2:-2], self.cell_len)
......@@ -383,24 +397,70 @@ class DGScheme(object):
plt.title('Errors')
# =============================================================================
# def _plot_details(self, fine_projection, coarse, projection):
# projected_coarse = np.sum(self.coarse_mesh, axis=0)
def _plot_details(self, projection):
fine_mesh = self.mesh[2:-2]
fine_projection = projection[:, 1:-1]
coarse_projection = self._calculate_coarse_projection(projection)
averaged_projection = [[coarse_projection[degree][cell]*self.basis[degree].subs(x, value)
for cell in range(len(coarse_projection[0]))
for value in [-0.5, 0.5]]
for degree in range(self.polynom_degree+1)]
multiwavelet_coeffs = self.detector._calculate_wavelet_coeffs(projection)
wavelet = AlpertsWavelet(self.polynom_degree).get_vector(x)
wavelet_projection = [[multiwavelet_coeffs[degree][cell]*wavelet[degree].subs(x, 0.5) * value
for cell in range(len(coarse_projection[0]))
for value in [(-1)**(self.polynom_degree+degree+1), 1]]
for degree in range(self.polynom_degree+1)]
avgMatrix = []
testMatrix = []
for degree in range(self.polynom_degree + 1):
leftMesh = coarse_projection[degree] * self.basis[degree].subs(x, -1 / 2)
rightMesh = coarse_projection[degree] * self.basis[degree].subs(x, 1 / 2)
leftTest = multiwavelet_coeffs[degree] * wavelet[degree].subs(x, 1 / 2) * (-1) ** (self.polynom_degree + 1 + degree)
rightTest = multiwavelet_coeffs[degree] * wavelet[degree].subs(x, 1 / 2)
newRowMesh = []
newRowTest = []
for i in range(len(coarse_projection[0])):
newRowMesh.append(leftMesh[i])
newRowMesh.append(rightMesh[i])
newRowTest.append(leftTest[i])
newRowTest.append(rightTest[i])
avgMatrix.append(newRowMesh)
testMatrix.append(newRowTest)
print(avgMatrix == averaged_projection)
print(testMatrix == wavelet_projection)
# uAvg = np.sum(avgMatrix, axis=0)
# uH = np.sum([newMatrix[degree] * phiVector[degree].subs(x, 0)
# for degree in range(p + 1)], axis=0)
# testD = np.sum(testMatrix, axis=0)
projected_coarse = np.sum(averaged_projection, axis=0)
# print(type(projected_coarse))
# projected_fine = np.sum([self.fine_mesh[degree]
# * self.basis[degree].subs(x, xi)
# for degree in range(self.polynom_degree+1)],
# axis=0)
# print(projected_coarse)
projected_fine = np.sum([fine_projection[degree] * self.basis[degree].subs(x, 0)
for degree in range(self.polynom_degree+1)], axis=0)
# print(type(projected_fine))
# print(projected_fine)
# print('Proj: ', projected_coarse.shape, projected_fine.shape)
# projected_wavelet_coeffs = np.sum(self.multiwavelet_coeffs, axis=0)
#
# plt.figure(4)
# plt.plot(self.fine_mesh, projected_fine-projected_coarse, 'm-.')
# plt.plot(self.fine_mesh, projected_wavelet_coeffs, 'y')
# plt.legend(['Fine-Coarse', 'Wavelet Coeff'])
# plt.xlabel('X')
# plt.ylabel('Detail Coefficients')
# plt.title('Wavelet Coefficients')
projected_wavelet_coeffs = np.sum(wavelet_projection, axis=0)
# print(type(projected_wavelet_coeffs))
# print(projected_wavelet_coeffs)
# print('Wave: ', projected_wavelet_coeffs.shape)
# print(fine_mesh.shape)
plt.figure(4)
plt.plot(fine_mesh, projected_fine-projected_coarse, 'm-.')
plt.plot(fine_mesh, projected_wavelet_coeffs, 'y')
plt.legend(['Fine-Coarse', 'Wavelet Coeff'])
plt.xlabel('X')
plt.ylabel('Detail Coefficients')
plt.title('Wavelet Coefficients')
# =============================================================================
def _plot_shock_tube(self):
......@@ -412,10 +472,8 @@ class DGScheme(object):
current_cells = troubled_cell_history[pos]
for cell in current_cells:
plt.plot(cell, time_history[pos], 'k.')
# =============================================================================
# print(pos, time_history[pos])
# print(current_cells)
# =============================================================================
plt.xlabel('Cell')
plt.ylabel('Time')
plt.title('Shock Tubes')
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment