From 9555fe9413b1a2010acce9bc37d16c4c99c2174e Mon Sep 17 00:00:00 2001 From: Antoine Marteau Date: Mon, 6 Jan 2025 13:00:10 +1100 Subject: [PATCH] fix P[Curl]Grad coverage and API for 1D + isHierarchical cov --- src/Polynomials/MonomialBases.jl | 5 +++++ src/Polynomials/NedelecPolyBases.jl | 4 ++-- test/PolynomialsTests/BernsteinBasesTests.jl | 6 ++++++ test/PolynomialsTests/LegendreBasesTests.jl | 2 ++ test/PolynomialsTests/MonomialBasesTests.jl | 2 ++ test/PolynomialsTests/PCurlGradBasesTests.jl | 20 ++++++++++++++++++-- test/PolynomialsTests/PGradBasesTests.jl | 11 +++++++++-- test/PolynomialsTests/QCurlGradBasesTests.jl | 1 - 8 files changed, 44 insertions(+), 7 deletions(-) diff --git a/src/Polynomials/MonomialBases.jl b/src/Polynomials/MonomialBases.jl index 32ce2fbf0..58a6282d7 100644 --- a/src/Polynomials/MonomialBases.jl +++ b/src/Polynomials/MonomialBases.jl @@ -26,6 +26,11 @@ MonomialBasis(args...) = UniformPolyBasis(Monomial, args...) function PGradBasis(::Type{Monomial},::Val{D},::Type{T},order::Int) where {D,T} NedelecPolyBasisOnSimplex{D}(Monomial,T,order) end +function PGradBasis(::Type{Monomial},::Val{1},::Type{T},order::Int) where T + @check T<:Real "T needs to be <:Real since represents the type of the components of the vector value" + V = VectorValue{1,T} + UniformPolyBasis(Monomial, Val(1), V, order+1) +end # 1D evaluation implementation diff --git a/src/Polynomials/NedelecPolyBases.jl b/src/Polynomials/NedelecPolyBases.jl index ebd65130b..807c2d4a8 100644 --- a/src/Polynomials/NedelecPolyBases.jl +++ b/src/Polynomials/NedelecPolyBases.jl @@ -2,6 +2,7 @@ struct NedelecPolyBasisOnSimplex{D,V,K,PT} <: PolynomialBasis{D,V,K,PT} order::Int function NedelecPolyBasisOnSimplex{D}(::Type{PT},::Type{T},order::Integer) where {D,PT,T} @check T<:Real "T needs to be <:Real since represents the type of the components of the vector value" + @notimplementedif !(D in (2,3)) K = Int(order)+1 V = VectorValue{D,T} new{D,V,K,PT}(Int(order)) @@ -242,7 +243,7 @@ Return a basis of ℕ𝔻ᴰₙ(△) = (ℙₙ)ᴰ ⊕ x × (ℙₙ \\ ℙₙ₋₁)ᴰ with n=`order`, the polynomial space for Nedelec elements on `D`-dimensional -simplices with scalar type `T`. +simplices with scalar type `T`. `D` must be 1, 2 or 3. The `order`=n argument has the following meaning: the curl of the functions in this basis is in ℙₙ. @@ -262,4 +263,3 @@ function PGradBasis(::Type{PT},::Val{D},::Type{T},order::Int) where {PT,D,T} @notimplemented "only implemented for monomials" end - diff --git a/test/PolynomialsTests/BernsteinBasesTests.jl b/test/PolynomialsTests/BernsteinBasesTests.jl index 33f007b4c..17ffc99ec 100644 --- a/test/PolynomialsTests/BernsteinBasesTests.jl +++ b/test/PolynomialsTests/BernsteinBasesTests.jl @@ -5,6 +5,8 @@ using Gridap.TensorValues using Gridap.Fields using Gridap.Polynomials +@test isHierarchical(Bernstein) == false + np = 3 x = [Point(0.),Point(1.),Point(.4)] xi = x[1] @@ -15,6 +17,7 @@ V = Float64 G = gradient_type(V,xi) H = gradient_type(G,xi) + # order 0 degenerated case order = 0 @@ -32,6 +35,7 @@ Hbx = repeat(permutedims(h),np) 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 @@ -52,6 +56,7 @@ Hbx = H[ 0. 0. 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 @@ -73,6 +78,7 @@ Hbx = H[ 2. -4. 2. 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 diff --git a/test/PolynomialsTests/LegendreBasesTests.jl b/test/PolynomialsTests/LegendreBasesTests.jl index 7072855af..404492b77 100644 --- a/test/PolynomialsTests/LegendreBasesTests.jl +++ b/test/PolynomialsTests/LegendreBasesTests.jl @@ -6,6 +6,8 @@ using Gridap.Fields using Gridap.Fields: Broadcasting using Gridap.Polynomials +@test isHierarchical(Legendre) == true + # Real-valued Q space with isotropic order x1 = Point(0.0) diff --git a/test/PolynomialsTests/MonomialBasesTests.jl b/test/PolynomialsTests/MonomialBasesTests.jl index 1114f5656..d3a734982 100644 --- a/test/PolynomialsTests/MonomialBasesTests.jl +++ b/test/PolynomialsTests/MonomialBasesTests.jl @@ -7,6 +7,8 @@ using Gridap.Polynomials using Gridap.Polynomials: _q_filter, _qs_filter, _p_filter, _ps_filter +@test isHierarchical(Monomial) == true + xi = Point(2,3) np = 5 x = fill(xi,np) diff --git a/test/PolynomialsTests/PCurlGradBasesTests.jl b/test/PolynomialsTests/PCurlGradBasesTests.jl index e446296d3..07d171796 100644 --- a/test/PolynomialsTests/PCurlGradBasesTests.jl +++ b/test/PolynomialsTests/PCurlGradBasesTests.jl @@ -10,12 +10,14 @@ xi = Point(4,2) np = 1 x = fill(xi,np) +PT = Monomial + order = 2 D = 2 T = Float64 V = VectorValue{D,T} G = gradient_type(V,xi) -b = PCurlGradBasis(Monomial, Val(D),T,order) +b = PCurlGradBasis(PT, Val(D),T,order) v = V[ (1.0, 0.0), (4.0, 0.0), (16.0, 0.0), (2.0, 0.0), (8.0, 0.0), (4.0, 0.0), # pterm ex @@ -57,9 +59,23 @@ D = 3 T = Float64 V = VectorValue{D,T} G = gradient_type(V,xi) -b = PCurlGradBasis(Monomial, Val(D),T,order) +b = PCurlGradBasis(PT, Val(D),T,order) @test length(b) == 15 @test get_order(b) == 2 + +# 1D + +order = 0 +D = 1 +T = Float64 +V = VectorValue{D,T} +b = PCurlGradBasis(PT,Val(D),T,order) + +@test b isa UniformPolyBasis{D,V,order+1,PT} + +@test_throws AssertionError PCurlGradBasis(PT,Val(D),V,order) + + end # module diff --git a/test/PolynomialsTests/PGradBasesTests.jl b/test/PolynomialsTests/PGradBasesTests.jl index d09cfb940..01138775d 100644 --- a/test/PolynomialsTests/PGradBasesTests.jl +++ b/test/PolynomialsTests/PGradBasesTests.jl @@ -11,6 +11,11 @@ np = 3 x = fill(xi,np) T = Float64 +D = 2 +for PT in (Legendre, Chebyshev, ModalC0, Bernstein) + @test_throws ErrorException PGradBasis(PT,Val(D),T,0) +end + PT = Monomial order = 0 @@ -46,20 +51,22 @@ bx = repeat(permutedims(v),np) test_field_array(b,x,bx,grad=∇bx) test_field_array(b,x[1],bx[1,:],grad=∇bx[1,:]) + # 1D order = 0 D = 1 T = Float64 V = VectorValue{D,T} -b = QGradBasis(PT,Val(D),T,order) +b = PGradBasis(PT,Val(D),T,order) @test b isa UniformPolyBasis{D,V,order+1,PT} -@test_throws AssertionError QGradBasis(PT,Val(D),V,order) +@test_throws AssertionError PGradBasis(PT,Val(D),V,order) # D > 3 not implemented @test_throws ErrorException NedelecPolyBasisOnSimplex{4}(PT, T, order) + end # module diff --git a/test/PolynomialsTests/QCurlGradBasesTests.jl b/test/PolynomialsTests/QCurlGradBasesTests.jl index 30384f53f..30a7d3207 100644 --- a/test/PolynomialsTests/QCurlGradBasesTests.jl +++ b/test/PolynomialsTests/QCurlGradBasesTests.jl @@ -92,7 +92,6 @@ order = 0 D = 1 T = Float64 V = VectorValue{D,T} -G = gradient_type(V,xi) b = QCurlGradBasis(PT,Val(D),T,order) @test b isa UniformPolyBasis{D,V,order+1,PT}