diff --git a/src/TypedPolynomials.jl b/src/TypedPolynomials.jl index 1b0f8b5..79eb06e 100644 --- a/src/TypedPolynomials.jl +++ b/src/TypedPolynomials.jl @@ -25,6 +25,7 @@ export @polyvar, variables, terms, differentiate, + antidifferentiate, subs include("types.jl") diff --git a/src/calculus.jl b/src/calculus.jl index 1f89b52..137a6b1 100644 --- a/src/calculus.jl +++ b/src/calculus.jl @@ -37,3 +37,38 @@ function _diff(m::Monomial{Vars}, end end, Val{N}())) end + +MP.antidifferentiate(v1::V, v2::V) where {V <: Variable} = v1*v2/2 +MP.antidifferentiate(v1::Variable, v2::Variable) = v1 * v2 + +function MP.antidifferentiate(m::Monomial{Vars}, v::V) where {Vars, V <: Variable} + if inmonomial(v, Vars...) + return _antidiff(m, v) + else + # Insert `v` in the monomial + return m * v + end +end + +_antidiff(m::Monomial, v::Variable) = _antidiff(m, exponents(m), v) + +function _antidiff(::Monomial{Vars}, + exponents::NTuple{N, Integer}, + v::Variable) where {Vars, N} + vi = find_variable_index(v, Vars) + new_m = Monomial{Vars,N}( + ntuple(i -> begin + if i == vi + (exponents[i] == 0) ? 1 : exponents[i] + 1 + else + exponents[i] + end + end, + Val{N}() + ) + ) + # Remark : as `exponents` are imposed to be `Integer` according to + # the method signature, we can use the Rational{Int} type here for the + # coefficient + return (1 // (exponents[vi] + 1)) * new_m +end diff --git a/test/calculus.jl b/test/calculus.jl index 79c8916..2d917c1 100644 --- a/test/calculus.jl +++ b/test/calculus.jl @@ -25,3 +25,31 @@ @test @inferred(differentiate(1, x)) == 0 end + +@testset "antiderivatives" begin + @polyvar x y z + + @test @inferred(antidifferentiate(x, x)) == 1//2*x^2 + @test @inferred(antidifferentiate(x, y)) == x*y + @test @inferred(antidifferentiate(y, x)) == x*y + @test @inferred(antidifferentiate(x^2, x)) == 1//3*x^3 + @test @inferred(antidifferentiate(x^2, y)) == x^2*y + @test @inferred(antidifferentiate(x^2 * y^3, y)) == x^2 * 1//4*y^4 + @test @inferred(antidifferentiate(x^2 * y^3, x)) == x^3 * 1//3*y^3 + @test @inferred(antidifferentiate(x^2 * y^3, z)) == x^2 * y^3 * z + @test @inferred(antidifferentiate(3x^2, x)) == x^3 + @test @inferred(antidifferentiate(3x^2 * y^0, y)) == 3x^2 * y + @test @inferred(antidifferentiate(3 * x^2 + 2 * x + 1, x)) == x^3 + x^2 + x + + @test @inferred(exponents(antidifferentiate(x^0,x))==(1,)) + @test @inferred(exponents(antidifferentiate(x^0*y*z^2,x))==(1,1,2)) + + m = x^2 + @test @wrappedallocs(antidifferentiate(m, x)) == 0 + @test @wrappedallocs(antidifferentiate(m, y)) == 0 + m = x^2 * y^3 + @test @wrappedallocs(antidifferentiate(m, x)) == 0 + @test @wrappedallocs(antidifferentiate(m, y)) == 0 + + @test @inferred(antidifferentiate(1, x)) == x +end