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

Moved plotting to new file 'Plotting'.

parent 3cb4a461
No related branches found
No related tags found
No related merge requests found
# -*- coding: utf-8 -*-
"""
@author: Laura C. Kühle
"""
import numpy as np
import matplotlib.pyplot as plt
import seaborn as sns
from sympy import Symbol
x = Symbol('x')
z = Symbol('z')
sns.set()
def plot_solution_and_approx(grid, exact, approx, color_exact, color_approx):
print(color_exact, color_approx)
plt.figure('exact_and_approx')
plt.plot(grid[0], exact[0], color_exact)
plt.plot(grid[0], approx[0], color_approx)
plt.xlabel('x')
plt.ylabel('u(x,t)')
plt.title('Solution and Approximation')
def plot_semilog_error(grid, pointwise_error):
plt.figure('semilog_error')
plt.semilogy(grid[0], pointwise_error[0])
plt.xlabel('x')
plt.ylabel('|u(x,t)-uh(x,t)|')
plt.title('Semilog Error plotted at Evaluation points')
def plot_error(grid, exact, approx):
plt.figure('error')
plt.plot(grid[0], exact[0]-approx[0])
plt.xlabel('X')
plt.ylabel('u(x,t)-uh(x,t)')
plt.title('Errors')
def plot_shock_tube(num_grid_cells, troubled_cell_history, time_history):
plt.figure('shock_tube')
for pos in range(len(time_history)):
current_cells = troubled_cell_history[pos]
for cell in current_cells:
plt.plot(cell, time_history[pos], 'k.')
plt.xlim((0, num_grid_cells // 2))
plt.xlabel('Cell')
plt.ylabel('Time')
plt.title('Shock Tubes')
def plot_details(fine_projection, fine_mesh, coarse_projection, basis, wavelet, multiwavelet_coeffs,
num_coarse_grid_cells, polynomial_degree):
averaged_projection = [[coarse_projection[degree][cell] * basis[degree].subs(x, value)
for cell in range(num_coarse_grid_cells)
for value in [-0.5, 0.5]]
for degree in range(polynomial_degree + 1)]
wavelet_projection = [[multiwavelet_coeffs[degree][cell] * wavelet[degree].subs(z, 0.5) * value
for cell in range(num_coarse_grid_cells)
for value in [(-1) ** (polynomial_degree + degree + 1), 1]]
for degree in range(polynomial_degree + 1)]
projected_coarse = np.sum(averaged_projection, axis=0)
projected_fine = np.sum([fine_projection[degree] * basis[degree].subs(x, 0)
for degree in range(polynomial_degree + 1)], axis=0)
projected_wavelet_coeffs = np.sum(wavelet_projection, axis=0)
plt.figure('coeff_details')
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 calculate_approximate_solution(projection, points, polynomial_degree, basis):
num_points = len(points)
basis_matrix = [[basis[degree].subs(x, points[point]) for point in range(num_points)]
for degree in range(polynomial_degree+1)]
approx = [[sum(projection[degree][cell] * basis_matrix[degree][point]
for degree in range(polynomial_degree+1))
for point in range(num_points)]
for cell in range(len(projection[0]))]
return np.reshape(np.array(approx), (1, len(approx) * num_points))
def calculate_exact_solution(mesh, cell_len, wave_speed, final_time, interval_len, quadrature, init_cond):
grid = []
exact = []
num_periods = np.floor(wave_speed * final_time / interval_len)
for cell in range(len(mesh)):
eval_points = mesh[cell]+cell_len / 2 * quadrature.get_eval_points()
eval_values = []
for point in range(len(eval_points)):
new_entry = init_cond.calculate(eval_points[point] - wave_speed * final_time + num_periods * interval_len)
eval_values.append(new_entry)
grid.append(eval_points)
exact.append(eval_values)
exact = np.reshape(np.array(exact), (1, len(exact) * len(exact[0])))
grid = np.reshape(np.array(grid), (1, len(grid) * len(grid[0])))
return grid, exact
......@@ -2,7 +2,7 @@
"""
@author: Laura C. Kühle, Soraya Terrab (sorayaterrab)
TODO: Move plotting to separate file (try to adjust for different equations)
TODO: Move plotting to separate file (try to adjust for different equations later) -> Done
TODO: Improve figure identifiers -> Done
"""
......@@ -13,6 +13,8 @@ import torch
from sympy import Symbol
import ANN_Model
from Plotting import plot_solution_and_approx, plot_semilog_error, plot_error, plot_shock_tube, plot_details, \
calculate_approximate_solution, calculate_exact_solution
x = Symbol('x')
z = Symbol('z')
......@@ -60,110 +62,40 @@ class TroubledCellDetector(object):
Here come some parameter.
"""
cell_averages = self._calculate_approximate_solution(projection, [0], 0)
left_reconstructions = self._calculate_approximate_solution(projection, [-1], self._polynomial_degree)
right_reconstructions = self._calculate_approximate_solution(projection, [1], self._polynomial_degree)
cell_averages = calculate_approximate_solution(projection, [0], 0, self._basis.get_basis_vector())
left_reconstructions = calculate_approximate_solution(projection, [-1], self._polynomial_degree,
self._basis.get_basis_vector())
right_reconstructions = calculate_approximate_solution(projection, [1], self._polynomial_degree,
self._basis.get_basis_vector())
middle_idx = stencil_length//2
return np.array(list(map(np.float64, zip(cell_averages[:, :middle_idx],
left_reconstructions[:, middle_idx], cell_averages[:, middle_idx],
right_reconstructions[:, middle_idx], cell_averages[:, middle_idx+1:]))))
def plot_results(self, projection, troubled_cell_history, time_history):
self._plot_shock_tube(troubled_cell_history, time_history)
plot_shock_tube(self._num_grid_cells, troubled_cell_history, time_history)
max_error = self._plot_mesh(projection)
print('p =', self._polynomial_degree)
print('N =', self._num_grid_cells)
print('maximum error =', max_error)
def _plot_shock_tube(self, troubled_cell_history, time_history):
plt.figure('shock_tube')
for pos in range(len(time_history)):
current_cells = troubled_cell_history[pos]
for cell in current_cells:
plt.plot(cell, time_history[pos], 'k.')
plt.xlim((0, self._num_grid_cells//2))
plt.xlabel('Cell')
plt.ylabel('Time')
plt.title('Shock Tubes')
def _plot_mesh(self, projection):
grid, exact = self._calculate_exact_solution(self._mesh[2:-2], self._cell_len)
approx = self._calculate_approximate_solution(projection[:, 1:-1], self._quadrature.get_eval_points(),
self._polynomial_degree)
grid, exact = calculate_exact_solution(self._mesh[2:-2], self._cell_len, self._wave_speed, self._final_time,
self._interval_len, self._quadrature, self._init_cond)
approx = calculate_approximate_solution(projection[:, 1:-1], self._quadrature.get_eval_points(),
self._polynomial_degree, self._basis.get_basis_vector())
pointwise_error = np.abs(exact-approx)
max_error = np.max(pointwise_error)
self._plot_solution_and_approx(grid, exact, approx, self._colors['exact'], self._colors['approx'])
plot_solution_and_approx(grid, exact, approx, self._colors['exact'], self._colors['approx'])
plt.legend(['Exact', 'Approx'])
self._plot_semilog_error(grid, pointwise_error)
self._plot_error(grid, exact, approx)
plot_semilog_error(grid, pointwise_error)
plot_error(grid, exact, approx)
return max_error
@staticmethod
def _plot_solution_and_approx(grid, exact, approx, color_exact, color_approx):
print(color_exact, color_approx)
plt.figure('exact_and_approx')
plt.plot(grid[0], exact[0], color_exact)
plt.plot(grid[0], approx[0], color_approx)
plt.xlabel('x')
plt.ylabel('u(x,t)')
plt.title('Solution and Approximation')
@staticmethod
def _plot_semilog_error(grid, pointwise_error):
plt.figure('semilog_error')
plt.semilogy(grid[0], pointwise_error[0])
plt.xlabel('x')
plt.ylabel('|u(x,t)-uh(x,t)|')
plt.title('Semilog Error plotted at Evaluation points')
@staticmethod
def _plot_error(grid, exact, approx):
plt.figure('error')
plt.plot(grid[0], exact[0]-approx[0])
plt.xlabel('X')
plt.ylabel('u(x,t)-uh(x,t)')
plt.title('Errors')
def _calculate_exact_solution(self, mesh, cell_len):
grid = []
exact = []
num_periods = np.floor(self._wave_speed * self._final_time / self._interval_len)
for cell in range(len(mesh)):
eval_points = mesh[cell] + cell_len/2 * self._quadrature.get_eval_points()
eval_values = []
for point in range(len(eval_points)):
new_entry = self._init_cond.calculate(eval_points[point] - self._wave_speed*self._final_time
+ num_periods*self._interval_len)
eval_values.append(new_entry)
grid.append(eval_points)
exact.append(eval_values)
exact = np.reshape(np.array(exact), (1, len(exact)*len(exact[0])))
grid = np.reshape(np.array(grid), (1, len(grid)*len(grid[0])))
return grid, exact
def _calculate_approximate_solution(self, projection, points, polynomial_degree):
num_points = len(points)
basis = self._basis.get_basis_vector()
basis_matrix = [[basis[degree].subs(x, points[point]) for point in range(num_points)]
for degree in range(polynomial_degree+1)]
approx = [[sum(projection[degree][cell] * basis_matrix[degree][point]
for degree in range(polynomial_degree+1))
for point in range(num_points)]
for cell in range(len(projection[0]))]
return np.reshape(np.array(approx), (1, len(approx) * num_points))
class NoDetection(TroubledCellDetector):
def get_cells(self, projection):
......@@ -237,40 +169,12 @@ class WaveletDetector(TroubledCellDetector):
return []
def plot_results(self, projection, troubled_cell_history, time_history):
self._plot_details(projection)
super().plot_results(projection, troubled_cell_history, time_history)
def _plot_details(self, projection):
fine_mesh = self._mesh[2:-2]
fine_projection = projection[:, 1:-1]
coarse_projection = self._calculate_coarse_projection(projection)
multiwavelet_coeffs = self._calculate_wavelet_coeffs(projection)
basis = self._basis.get_basis_vector()
wavelet = self._basis.get_wavelet_vector()
averaged_projection = [[coarse_projection[degree][cell] * basis[degree].subs(x, value)
for cell in range(self._num_coarse_grid_cells)
for value in [-0.5, 0.5]]
for degree in range(self._polynomial_degree + 1)]
wavelet_projection = [[multiwavelet_coeffs[degree][cell] * wavelet[degree].subs(z, 0.5) * value
for cell in range(self._num_coarse_grid_cells)
for value in [(-1) ** (self._polynomial_degree + degree + 1), 1]]
for degree in range(self._polynomial_degree + 1)]
projected_coarse = np.sum(averaged_projection, axis=0)
projected_fine = np.sum([fine_projection[degree] * basis[degree].subs(x, 0)
for degree in range(self._polynomial_degree + 1)], axis=0)
projected_wavelet_coeffs = np.sum(wavelet_projection, axis=0)
plt.figure('coeff_details')
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')
coarse_projection = self._calculate_coarse_projection(projection)
plot_details(projection[:, 1:-1], self._mesh[2:-2], coarse_projection, self._basis.get_basis_vector(),
self._basis.get_wavelet_vector(), multiwavelet_coeffs, self._num_coarse_grid_cells,
self._polynomial_degree)
super().plot_results(projection, troubled_cell_history, time_history)
def _calculate_coarse_projection(self, projection):
basis_projection_left, basis_projection_right = self._basis.get_basis_projections()
......@@ -289,18 +193,19 @@ class WaveletDetector(TroubledCellDetector):
return coarse_projection
def _plot_mesh(self, projection):
grid, exact = self._calculate_exact_solution(self._mesh[2:-2], self._cell_len)
approx = self._calculate_approximate_solution(projection[:, 1:-1], self._quadrature.get_eval_points(),
self._polynomial_degree)
grid, exact = calculate_exact_solution(self._mesh[2:-2], self._cell_len, self._wave_speed, self._final_time,
self._interval_len, self._quadrature, self._init_cond)
approx = calculate_approximate_solution(projection[:, 1:-1], self._quadrature.get_eval_points(),
self._polynomial_degree, self._basis.get_basis_vector())
pointwise_error = np.abs(exact-approx)
max_error = np.max(pointwise_error)
self._plot_coarse_mesh(projection)
self._plot_solution_and_approx(grid, exact, approx, self._colors['fine_exact'], self._colors['fine_approx'])
plot_solution_and_approx(grid, exact, approx, self._colors['fine_exact'], self._colors['fine_approx'])
plt.legend(['Exact (Coarse)', 'Approx (Coarse)', 'Exact (Fine)', 'Approx (Fine)'])
self._plot_semilog_error(grid, pointwise_error)
self._plot_error(grid, exact, approx)
plot_semilog_error(grid, pointwise_error)
plot_error(grid, exact, approx)
return max_error
......@@ -312,10 +217,11 @@ class WaveletDetector(TroubledCellDetector):
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._quadrature.get_eval_points(),
self._polynomial_degree)
self._plot_solution_and_approx(grid, exact, approx, self._colors['coarse_exact'], self._colors['coarse_approx'])
grid, exact = calculate_exact_solution(coarse_mesh[1:-1], coarse_cell_len, self._wave_speed, self._final_time,
self._interval_len, self._quadrature, self._init_cond)
approx = calculate_approximate_solution(coarse_projection, self._quadrature.get_eval_points(),
self._polynomial_degree, self._basis.get_basis_vector())
plot_solution_and_approx(grid, exact, approx, self._colors['coarse_exact'], self._colors['coarse_approx'])
class Boxplot(WaveletDetector):
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment