Skip to content

Commit

Permalink
Chebyshev test coverage and method ambiguities
Browse files Browse the repository at this point in the history
  • Loading branch information
Antoinemarteau committed Jan 6, 2025
1 parent 6b07a81 commit 439ea2e
Show file tree
Hide file tree
Showing 3 changed files with 135 additions and 12 deletions.
26 changes: 14 additions & 12 deletions src/Polynomials/ChebyshevBases.jl
Original file line number Diff line number Diff line change
Expand Up @@ -24,34 +24,37 @@ High level constructors of [`ChebyshevBasis`](@ref).
"""
ChebyshevBasis(args...; kind=:T) = UniformPolyBasis(Chebyshev{kind}, args...)

UniformPolyBasis{D}(::Type{Chebyshev{:U}}, args...) where D = @notimplemented "1D evaluation for second kind needed here"
function UniformPolyBasis(
::Type{Chebyshev{:U}}, ::Val{D}, ::Type{V}, ::Int64) where {D, V}

@notimplemented "1D evaluation for second kind need to be implemented here"
end


# 1D evaluation implementation

function _evaluate_1d!(
::Type{<:Chebyshev},::Val{0},c::AbstractMatrix{T},x,d) where T<:Number
::Type{Chebyshev{kind}},::Val{0},c::AbstractMatrix{T},x,d) where {kind,T<:Number}

@inbounds c[d,1] = one(T)
end

function _evaluate_1d!(
::Type{Chebyshev{:T}},::Val{K},c::AbstractMatrix{T},x,d) where {K,T<:Number}
::Type{Chebyshev{kind}},::Val{K},c::AbstractMatrix{T},x,d) where {kind,K,T<:Number}

n = K + 1 # n > 1
ξ = (2*x[d] - 1) # ξ ∈ [-1,1]
ξ2 = 2*ξ

@inbounds c[d,1] = one(T)
@inbounds c[d,2] = ξ
@inbounds c[d,2] = (kind == :T) ? ξ : ξ2
for i in 3:n
@inbounds c[d,i] = c[d,i-1]*ξ2 - c[d,i-2]
end
end


function _gradient_1d!(
::Type{<:Chebyshev},::Val{0},g::AbstractMatrix{T},x,d) where T<:Number
::Type{Chebyshev{:T}},::Val{0},g::AbstractMatrix{T},x,d) where T<:Number

@inbounds g[d,1] = zero(T)
end
Expand All @@ -76,13 +79,12 @@ function _gradient_1d!(
end


function _hessian_1d!(
::Type{Chebyshev{:T}},::Val{K},h::AbstractMatrix{T},x,d) where {K,T<:Number}
function _gradient_1d!(
::Type{Chebyshev{:U}},::Val{0},g::AbstractMatrix{T},x,d) where T<:Number

@notimplemented
@inbounds g[d,1] = zero(T)
end


_evaluate_1d!(::Type{Chebyshev{:U}},::Val{K},c::AbstractMatrix{T},x,d) where {K,T<:Number} = @notimplemented
_gradient_1d!(::Type{Chebyshev{:U}},::Val{K},g::AbstractMatrix{T},x,d) where {K,T<:Number} = @notimplemented
_gradient_1d!(::Type{Chebyshev{:U}},::Val{K},h::AbstractMatrix{T},x,d) where {K,T<:Number} = @notimplemented
_hessian_1d!( ::Type{Chebyshev{:U}},::Val{K},h::AbstractMatrix{T},x,d) where {K,T<:Number} = @notimplemented

119 changes: 119 additions & 0 deletions test/PolynomialsTests/ChebyshevBasesTests.jl
Original file line number Diff line number Diff line change
@@ -0,0 +1,119 @@
module ChebyshevBasisTests

using Test
using Gridap.TensorValues
using Gridap.Fields
using Gridap.Polynomials
using Gridap.Polynomials: binoms
using ForwardDiff

@test isHierarchical(Chebyshev) == true

np = 3
x = [Point(2.),Point(3.),Point(4.)]
x1 = x[1]


# Only test 1D evaluations as tensor product structure is tested in monomial tests

V = Float64
G = gradient_type(V,x1)
H = gradient_type(G,x1)

chebyshev_T(N) = t -> begin
ξ = 2t-1
sq = sqrt*ξ-1)
.5*( (ξ - sq)^N ++ sq)^N )
end
chebyshev_U(N) = t -> begin
ξ = 2t-1
sq = sqrt*ξ-1)
.5*( (ξ + sq)^(N+1) -- sq)^(N+1) )/sq
end
_∇(b) = t -> ForwardDiff.derivative(b, t)
_H(b) = t -> ForwardDiff.derivative(y -> ForwardDiff.derivative(b, y), t)

######################################
# First Kind Chebyshev ( Hessian TODO)
######################################

# order 0 degenerated case
order = 0
b = ChebyshevBasis(Val(1),V,order)
@test get_order(b) == 0
@test get_orders(b) == (0,)

bx = [ chebyshev_T(n)( xi[1]) for xi in x, n in 0:order]
∇bx = [ G(_∇(chebyshev_T(n))(xi[1])) for xi in x, n in 0:order]
Hbx = [ H(_H(chebyshev_T(n))(xi[1])) for xi in x, n in 0:order]

test_field_array(b,x,bx,, grad=∇bx)#,gradgrad=Hbx)
test_field_array(b,x[1],bx[1,:],,grad=∇bx[1,:])#,gradgrad=Hbx[1,:])


# Order 1
order = 1
b = ChebyshevBasis(Val(1),V,order)

bx = [ chebyshev_T(n)( xi[1]) for xi in x, n in 0:order]
∇bx = [ G(_∇(chebyshev_T(n))(xi[1])) for xi in x, n in 0:order]
Hbx = [ H(_H(chebyshev_T(n))(xi[1])) for xi in x, n in 0:order]

test_field_array(b,x,bx,,grad=∇bx)#,gradgrad=Hbx)
test_field_array(b,x[1],bx[1,:],,grad=∇bx[1,:])#,gradgrad=Hbx[1,:])


# Order 2
order = 2
b = ChebyshevBasis(Val(1),V,order)


bx = [ chebyshev_T(n)( xi[1]) for xi in x, n in 0:order]
∇bx = [ G(_∇(chebyshev_T(n))(xi[1])) for xi in x, n in 0:order]
Hbx = [ H(_H(chebyshev_T(n))(xi[1])) for xi in x, n in 0:order]

test_field_array(b,x,bx,,grad=∇bx)#,gradgrad=Hbx)
test_field_array(b,x[1],bx[1,:],,grad=∇bx[1,:])#,gradgrad=Hbx[1,:])


# Order 3
order = 3
b = ChebyshevBasis(Val(1),V,order)

bx = [ chebyshev_T(n)( xi[1]) for xi in x, n in 0:order]
∇bx = [ G(_∇(chebyshev_T(n))(xi[1])) for xi in x, n in 0:order]
Hbx = [ H(_H(chebyshev_T(n))(xi[1])) for xi in x, n in 0:order]

test_field_array(b,x,bx,, grad=∇bx)#,gradgrad=Hbx)
test_field_array(b,x[1],bx[1,:],,grad=∇bx[1,:])#,gradgrad=Hbx[1,:])


# Order 4
order = 4
b = ChebyshevBasis(Val(1),V,order)

bx = [ chebyshev_T(n)( xi[1]) for xi in x, n in 0:order]
∇bx = [ G(_∇(chebyshev_T(n))(xi[1])) for xi in x, n in 0:order]
Hbx = [ H(_H(chebyshev_T(n))(xi[1])) for xi in x, n in 0:order]

test_field_array(b,x,bx,,grad=∇bx)#,gradgrad=Hbx)
test_field_array(b,x[1],bx[1,:],,grad=∇bx[1,:])#,gradgrad=Hbx[1,:])


############################
# Second Kind Chebyshev TODO
############################

@test_throws ErrorException ChebyshevBasis(Val(1),Float64,0;kind=:U)

order = 1
d = 1
b = ChebyshevBasis(Val(1),V,order)
Hb = FieldGradientArray{2}(b)
r, c, g, h = return_cache(Hb,x)

@test_throws ErrorException Polynomials._gradient_1d!(Chebyshev{:U},Val(order),g,x1,d)
@test_throws ErrorException Polynomials._hessian_1d!( Chebyshev{:U},Val(order),h,x1,d)


end # module
2 changes: 2 additions & 0 deletions test/PolynomialsTests/runtests.jl
Original file line number Diff line number Diff line change
Expand Up @@ -18,6 +18,8 @@ using Test

@testset "LegendreBases" begin include("LegendreBasesTests.jl") end

@testset "ChebyshevBases" begin include("ChebyshevBasesTests.jl") end

@testset "BernsteinBases" begin include("BernsteinBasesTests.jl") end

#@testset "ChangeBasis" begin include("ChangeBasisTests.jl") end
Expand Down

0 comments on commit 439ea2e

Please sign in to comment.