diff --git a/src/triangle.jl b/src/triangle.jl index e8f50ae..54c91e4 100644 --- a/src/triangle.jl +++ b/src/triangle.jl @@ -47,6 +47,13 @@ axes(P::TriangleWeight{T}) where T = (Inclusion(UnitTriangle{T}()),) Base.summary(io::IO, P::TriangleWeight) = print(io, "x^$(P.a)*y^$(P.b)*(1-x-y)^$(P.c) on the unit triangle") +function getindex(P::TriangleWeight{T}, xy::StaticVector{2}) where T + x,y = xy + convert(T, x^P.a * y^P.b * (1-x-y)^P.c) +end + +orthogonalityweight(P::JacobiTriangle) = TriangleWeight(P.a, P.b, P.c) + @simplify function *(Dx::PartialDerivative{1}, P::JacobiTriangle) a,b,c = P.a,P.b,P.c n = mortar(Fill.(oneto(∞),oneto(∞))) diff --git a/test/test_triangle.jl b/test/test_triangle.jl index 86aa461..0ae0646 100644 --- a/test/test_triangle.jl +++ b/test/test_triangle.jl @@ -333,4 +333,14 @@ import MultivariateOrthogonalPolynomials: tri_forwardrecurrence, grid, TriangleR C = P¹ \ (L * P) @test C[1:10,1:10] ≈ A[1:10,1:10] - B[1:10,1:10] end + + @testset "Weighted" begin + P¹ = JacobiTriangle(1,1,1) + W = Weighted(P¹) + @test W[SVector(0.1,0.2),1:5] ≈ 0.1 * 0.2 * (1 - 0.1 - 0.2) * P¹[SVector(0.1,0.2),1:5] + xy = axes(W,1) + ∂ˣ = PartialDerivative{1}(xy) + P = JacobiTriangle() + P \ W + end end