Skip to content

Commit

Permalink
Docstrings for some methods in basis Lobatto-Legendre (#1874)
Browse files Browse the repository at this point in the history
* Docstrings for some methods in basis LL

* double back slash
  • Loading branch information
DanielDoehring authored Mar 15, 2024
1 parent 38a9a52 commit 2dfde7f
Show file tree
Hide file tree
Showing 2 changed files with 72 additions and 7 deletions.
2 changes: 1 addition & 1 deletion docs/literate/src/files/scalar_linear_advection_1d.jl
Original file line number Diff line number Diff line change
Expand Up @@ -115,7 +115,7 @@ integral = sum(nodes.^3 .* weights)
# To approximate the solution, we need to get the polynomial coefficients $\{u_j^{Q_l}\}_{j=0}^N$
# for every element $Q_l$.

# After defining all nodes, we can implement the spatial coordinate $x$ and its initial value $u0$
# After defining all nodes, we can implement the spatial coordinate $x$ and its initial value $u0 = u(t_0)$
# for every node.
x = Matrix{Float64}(undef, length(nodes), n_elements)
for element in 1:n_elements
Expand Down
77 changes: 71 additions & 6 deletions src/solvers/dgsem/basis_lobatto_legendre.jl
Original file line number Diff line number Diff line change
Expand Up @@ -404,7 +404,8 @@ function calc_dsplit(nodes, weights)
return dsplit
end

# Calculate the polynomial derivative matrix D
# Calculate the polynomial derivative matrix D.
# This implements algorithm 37 "PolynomialDerivativeMatrix" from Kopriva's book.
function polynomial_derivative_matrix(nodes)
n_nodes = length(nodes)
d = zeros(n_nodes, n_nodes)
Expand All @@ -421,6 +422,7 @@ function polynomial_derivative_matrix(nodes)
end

# Calculate and interpolation matrix (Vandermonde matrix) between two given sets of nodes
# See algorithm 32 "PolynomialInterpolationMatrix" from Kopriva's book.
function polynomial_interpolation_matrix(nodes_in, nodes_out,
baryweights_in = barycentric_weights(nodes_in))
n_nodes_in = length(nodes_in)
Expand All @@ -433,6 +435,7 @@ function polynomial_interpolation_matrix(nodes_in, nodes_out,
return vandermonde
end

# This implements algorithm 32 "PolynomialInterpolationMatrix" from Kopriva's book.
function polynomial_interpolation_matrix!(vandermonde,
nodes_in, nodes_out,
baryweights_in)
Expand Down Expand Up @@ -463,7 +466,19 @@ function polynomial_interpolation_matrix!(vandermonde,
return vandermonde
end

# Calculate the barycentric weights for a given node distribution.
"""
barycentric_weights(nodes)
Calculate the barycentric weights for a given node distribution, i.e.,
```math
w_j = \\frac{1}{ \\prod_{k \\neq j} \\left( x_j - x_k \\right ) }
```
For details, see (especially Section 3)
- Jean-Paul Berrut and Lloyd N. Trefethen (2004).
Barycentric Lagrange Interpolation.
[DOI:10.1137/S0036144502417715](https://doi.org/10.1137/S0036144502417715)
"""
function barycentric_weights(nodes)
n_nodes = length(nodes)
weights = ones(n_nodes)
Expand Down Expand Up @@ -494,12 +509,31 @@ function calc_lhat(x, nodes, weights)
return lhat
end

# Calculate Lagrange polynomials for a given node distribution.
"""
lagrange_interpolating_polynomials(x, nodes, wbary)
Calculate Lagrange polynomials for a given node distribution with
associated barycentric weights `wbary` at a given point `x` on the
reference interval ``[-1, 1]``.
This returns all ``l_j(x)``, i.e., the Lagrange polynomials for each node ``x_j``.
Thus, to obtain the interpolating polynomial ``p(x)`` at ``x``, one has to
multiply the Lagrange polynomials with the nodal values ``u_j`` and sum them up:
``p(x) = \\sum_{j=1}^{n} u_j l_j(x)``.
For details, see e.g. Section 2 of
- Jean-Paul Berrut and Lloyd N. Trefethen (2004).
Barycentric Lagrange Interpolation.
[DOI:10.1137/S0036144502417715](https://doi.org/10.1137/S0036144502417715)
"""
function lagrange_interpolating_polynomials(x, nodes, wbary)
n_nodes = length(nodes)
polynomials = zeros(n_nodes)

for i in 1:n_nodes
# Avoid division by zero when `x` is close to node by using
# the Kronecker-delta property at nodes
# of the Lagrange interpolation polynomials.
if isapprox(x, nodes[i], rtol = eps(x))
polynomials[i] = 1
return polynomials
Expand All @@ -518,6 +552,17 @@ function lagrange_interpolating_polynomials(x, nodes, wbary)
return polynomials
end

"""
gauss_lobatto_nodes_weights(n_nodes::Integer)
Computes nodes ``x_j`` and weights ``w_j`` for the (Legendre-)Gauss-Lobatto quadrature.
This implements algorithm 25 "GaussLobattoNodesAndWeights" from the book
- David A. Kopriva, (2009).
Implementing spectral methods for partial differential equations:
Algorithms for scientists and engineers.
[DOI:10.1007/978-90-481-2261-5](https://doi.org/10.1007/978-90-481-2261-5)
"""
# From FLUXO (but really from blue book by Kopriva)
function gauss_lobatto_nodes_weights(n_nodes::Integer)
# From Kopriva's book
Expand Down Expand Up @@ -585,7 +630,7 @@ function gauss_lobatto_nodes_weights(n_nodes::Integer)
return nodes, weights
end

# From FLUXO (but really from blue book by Kopriva)
# From FLUXO (but really from blue book by Kopriva, algorithm 24)
function calc_q_and_l(N::Integer, x::Float64)
L_Nm2 = 1.0
L_Nm1 = x
Expand All @@ -609,7 +654,17 @@ function calc_q_and_l(N::Integer, x::Float64)
end
calc_q_and_l(N::Integer, x::Real) = calc_q_and_l(N, convert(Float64, x))

# From FLUXO (but really from blue book by Kopriva)
"""
gauss_nodes_weights(n_nodes::Integer)
Computes nodes ``x_j`` and weights ``w_j`` for the Gauss-Legendre quadrature.
This implements algorithm 23 "LegendreGaussNodesAndWeights" from the book
- David A. Kopriva, (2009).
Implementing spectral methods for partial differential equations:
Algorithms for scientists and engineers.
[DOI:10.1007/978-90-481-2261-5](https://doi.org/10.1007/978-90-481-2261-5)
"""
function gauss_nodes_weights(n_nodes::Integer)
# From Kopriva's book
n_iterations = 10
Expand Down Expand Up @@ -666,7 +721,17 @@ function gauss_nodes_weights(n_nodes::Integer)
end
end

# From FLUXO (but really from blue book by Kopriva)
"""
legendre_polynomial_and_derivative(N::Int, x::Real)
Computes the Legendre polynomial of degree `N` and its derivative at `x`.
This implements algorithm 22 "LegendrePolynomialAndDerivative" from the book
- David A. Kopriva, (2009).
Implementing spectral methods for partial differential equations:
Algorithms for scientists and engineers.
[DOI:10.1007/978-90-481-2261-5](https://doi.org/10.1007/978-90-481-2261-5)
"""
function legendre_polynomial_and_derivative(N::Int, x::Real)
if N == 0
poly = 1.0
Expand Down

0 comments on commit 2dfde7f

Please sign in to comment.