From 52ac10847f8501de2a52ff215a87ae474d92332c Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?K=C3=BChle=2C=20Laura=20Christine=20=28lakue103=29?= <laura.kuehle@uni-duesseldorf.de> Date: Fri, 10 Mar 2023 15:25:25 +0100 Subject: [PATCH] Added derivative basis to Basis class. --- Snakefile | 1 + scripts/tcd/Basis_Function.py | 94 +++++++++++++++++++++++++++++++++++ 2 files changed, 95 insertions(+) diff --git a/Snakefile b/Snakefile index d66ead6..b6035a9 100644 --- a/Snakefile +++ b/Snakefile @@ -29,6 +29,7 @@ TODO: ASk whether cfl number should be absolute value Urgent: TODO: Add Burgers class completely TODO: Add 'update_time_step()' to Burgers -> Done +TODO: Add derivative basis to Basis class -> Done TODO: Add 'solve_exactly()' to Burgers TODO: Add 'update_right_hand_side()' to Burgers TODO: Add initialization to Burgers diff --git a/scripts/tcd/Basis_Function.py b/scripts/tcd/Basis_Function.py index 4144fb9..69b1347 100644 --- a/scripts/tcd/Basis_Function.py +++ b/scripts/tcd/Basis_Function.py @@ -28,6 +28,8 @@ class Basis(ABC): Array of basis. wavelet : ndarray Array of wavelet. + derivative: ndarray + Array of derivative basis. inv_mass : ndarray Inverse mass matrix. basis_projection : Tuple[ndarray, ndarray] @@ -111,6 +113,29 @@ class Basis(ABC): """ pass + @property + @cache + def derivative(self) -> ndarray: + """Return derivative basis vector.""" + return self._build_derivative_vector(x) + + @abstractmethod + def _build_derivative_vector(self, eval_point: float) -> ndarray: + """Construct derivative basis vector. + + Parameters + ---------- + eval_point : float + Evaluation point. + + Returns + ------- + ndarray + Vector containing derivative basis evaluated at evaluation point. + + """ + pass + @property @cache def inverse_mass_matrix(self) -> ndarray: @@ -261,6 +286,55 @@ class Legendre(Basis): vector.append(poly) return np.array(vector) + def _build_derivative_vector(self, eval_point: float) -> ndarray: + """Construct derivative basis vector. + + Parameters + ---------- + eval_point : float + Evaluation point. + + Returns + ------- + ndarray + Vector containing derivative basis evaluated at evaluation point. + + """ + return self._calculate_derivative_legendre_vector(eval_point) + + def _calculate_derivative_legendre_vector(self, + eval_point: float) -> ndarray: + """Construct derivative Legendre vector. + + Parameters + ---------- + eval_point : float + Evaluation point. + + Returns + ------- + ndarray + Vector containing derivative Legendre polynomial evaluated at + evaluation point. + + """ + derivative_vector = [] + leg_vector = self._calculate_legendre_vector(eval_point) + for degree in range(self._polynomial_degree+1): + if degree == 0: + derivative_vector.append(0 * eval_point) + else: + if degree == 1: + derivative_vector.append(1.0+0 * eval_point) + else: + poly = ((2.0 * degree-1) / degree) * ( + leg_vector[degree-1]+( + eval_point * derivative_vector[-1]))-( + (degree-1) / degree) * \ + derivative_vector[-2] + derivative_vector.append(poly) + return np.array(derivative_vector) + class OrthonormalLegendre(Legendre): """Class for orthonormal Legendre basis.""" @@ -348,6 +422,26 @@ class OrthonormalLegendre(Legendre): raise ValueError('Invalid value: Alpert\'s wavelet is only available \ up to degree 4 for this application') + def _build_derivative_basis_vector(self, eval_point): + """Construct derivative basis vector. + + Parameters + ---------- + eval_point : float + Evaluation point. + + Returns + ------- + ndarray + Vector containing derivative Legendre polynomial evaluated at + evaluation point. + + """ + derivative_leg_vector = self._calculate_derivative_legendre_vector( + eval_point) + return np.array([derivative_leg_vector[degree] * np.sqrt(degree+0.5) + for degree in range(self._polynomial_degree+1)]) + def _build_inverse_mass_matrix(self) -> ndarray: """Construct inverse mass matrix. -- GitLab