Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

[WIP] Add symbolic integral #240

Merged
merged 22 commits into from
Jun 15, 2021
Merged
Show file tree
Hide file tree
Changes from 10 commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Jump to
Jump to file
Failed to load files.
Diff view
Diff view
1 change: 0 additions & 1 deletion Project.toml
Original file line number Diff line number Diff line change
Expand Up @@ -41,7 +41,6 @@ Reexport = "0.2, 1"
ReferenceTests = "0.9"
Requires = "1.1"
RuntimeGeneratedFunctions = "0.4.3, 0.5"
SafeTestsets = "0"
SciMLBase = "1.8"
Setfield = "0.7"
SpecialFunctions = "0.7, 0.8, 0.9, 0.10, 1.0"
Expand Down
6 changes: 6 additions & 0 deletions src/Symbolics.jl
Original file line number Diff line number Diff line change
Expand Up @@ -6,6 +6,9 @@ using LinearAlgebra

using Reexport

using DomainSets

import DomainSets: Domain
@reexport using SymbolicUtils

import SymbolicUtils: Term, Add, Mul, Pow, Sym, symtype,
Expand Down Expand Up @@ -52,6 +55,9 @@ export Differential, expand_derivatives

include("diff.jl")

export Integral
include("integral.jl")

include("linear_algebra.jl")

import Libdl
Expand Down
21 changes: 21 additions & 0 deletions src/diff.jl
Original file line number Diff line number Diff line change
Expand Up @@ -113,6 +113,25 @@ function expand_derivatives(O::Symbolic, simplify=false; occurances=nothing)
return expand_derivatives(operation(arg)(inner), simplify)
end
end
elseif isa(operation(arg) , Integral)
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Whitespace.

if isa(operation(arg).domain , AbstractInterval)
domain = operation(arg).domain
a , b = DomainSets.endpoints(domain)
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Whitespace.

c = 0
inner_function = expand_derivatives_(arguments(arg)[1])
if isa(value(a), Term)
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Why isa Term? Shouldn't it be istree? Check SymbolicUtils' docs.

t1 = replaceSym(operation(arg).x ,value(a) , inner_function )
t2 = D(a)
c -= t1*t2
elseif isa(value(b) , Term)
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I think you mean if instead of elseif? So if both a and b are terms, it would work out...

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I resolved this

t1 = replaceSym(operation(arg).x ,value(b) , inner_function )
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

formatting needs to be cleaned up, we also never use camelCase function names. Good to stick to the package's style when you're contributing to it. Where is the expand_derivatives_ function...?

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

you could use substitute here.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Yes, I infact removed replaceSym

t2 = D(b)
c += t1*t2
end
inner = expand_derivatives_(D(arguments(arg)[1]))
c += operation(arg)(inner)
return c
end
end

l = length(arguments(arg))
Expand Down Expand Up @@ -151,6 +170,8 @@ function expand_derivatives(O::Symbolic, simplify=false; occurances=nothing)
x = +((!_iszero(c) ? vcat(c, exprs) : exprs)...)
return simplify ? SymbolicUtils.simplify(x) : x
end
elseif isa(operation(O) , Integral )
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Could you be more careful about whitespace?

return operation(O)(expand_derivatives_(arguments(O)[1]))
elseif !hasderiv(O)
return O
else
Expand Down
34 changes: 34 additions & 0 deletions src/integral.jl
Original file line number Diff line number Diff line change
@@ -0,0 +1,34 @@
# basic integral struct with upper bound and lower bound.
struct Integral{X, T <: Domain} <: Function
x
domain::T
Integral(x,domain) = new{typeof(x),typeof(domain)}(Symbolics.value(x), domain)
end

(I::Integral)(x) = Term{SymbolicUtils.symtype(x)}(I, [x])
(I::Integral)(x::Num) = Num(I(Symbolics.value(x)))
SymbolicUtils.promote_symtype(::Integral, x) = x

function Base.show(io::IO, I::Integral)
print(io, "Integral(", I.x, ", ", I.domain, ")")
end

Base.:(==)(I1::Integral, I2::Integral) = (isequal(I1.x, I2.x) && isequal(I1.domain, I2.domain))
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

When comparing numbers, we should use == instead of isequal. Consider -0.0 vs 0.0.

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Ah, I guess domains could be symbolic.

(D::Differential)(I::Integral{X,T}) where{X,T} = I∘D
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Indentation.


function replaceSym(a::Sym, b, O)
if isa(O , Sym)
if isequal(O , a)
return b
else
return O
end
else
args_replace = Vector{Symbolic{Real}}()
args_ = arguments(O)
for i in 1:length(args_)
push!(args_replace, replaceSym(a , b , args_[i]))
end
return operation(O)(args_replace...)
end
end
7 changes: 7 additions & 0 deletions test/integral.jl
Original file line number Diff line number Diff line change
@@ -0,0 +1,7 @@
using Symbolics, Test
using DomainSets
@variables x y

I1 = Integral( x , ClosedInterval(1, 5))
I2 = Integral( x , ClosedInterval(1, 5))
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Also here.

@test I1 == I2