diff --git a/src/bernstein.jl b/src/bernstein.jl index 6483c0bdd2c2e4b48aa640fb4b2c39623cb0885a..969676db9a14c10070d1e46c4634ada39de0bc57 100644 --- a/src/bernstein.jl +++ b/src/bernstein.jl @@ -5,32 +5,6 @@ using dg1d using LinearAlgebra -struct BernsteinElement - Npts::Int64 - N::Int64 - x::Vector{Float64} - M::Matrix{Float64} - invM::Matrix{Float64} - S::Matrix{Float64} - MDM::Matrix{Float64} # = (M^-1 D^T M) = M^-1 S^T - MB_lhs::Vector{Float64} # = M^(-1) Bn(-1) - MB_rhs::Vector{Float64} # = M^(-1) Bn(-1) -end - - -function BernsteinElement(Npts) - N = Npts-1 - x = collect(range(-1.0,1.0,length=Npts)) - M = mass_matrix(Npts) - invM = inv(M) - S = stiffness_matrix(Npts) - MDM = invM * transpose(S) - MB_lhs = invM[:,1] # because Bn(-1) = δ(n,1) - MB_rhs = invM[:,end] # because Bn(1) = δ(n,Npts) - return BernsteinElement(Npts, N, x, M, invM, S, MDM, MB_lhs, MB_rhs) -end - - # Normalized so that bernstein_polynomial(-1,0,N) == bernstein_polynomial(1,N,N) == 1 function bernstein_polynomial(x, n::Integer, N::Integer) return binomial(N, n) / 2^N * (1 + x)^n * (1 - x)^(N-n)