Skip to content

Implementing domains as an interface #141

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

Merged
merged 13 commits into from
Oct 9, 2023
8 changes: 4 additions & 4 deletions Project.toml
Original file line number Diff line number Diff line change
@@ -1,6 +1,6 @@
name = "DomainSets"
uuid = "5b8099bc-c8ec-5219-889f-1d9e522a28bf"
version = "0.6.7"
version = "0.7.0"

[deps]
CompositeTypes = "b152e2b5-7a66-4b01-a709-34e65c35f657"
Expand All @@ -11,15 +11,15 @@ StaticArrays = "90137ffa-7385-5640-81b9-e52037218182"
Statistics = "10745b16-79ce-11e8-11f9-7d13ad32a3b2"

[compat]
julia = "1.6"
CompositeTypes = "0.1.2"
IntervalSets = "0.7.4"
StaticArrays = "0.12.2, 1"
StableRNGs = "1"
julia = "1.6"
StaticArrays = "0.12.2, 1"

[extras]
Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40"
StableRNGs = "860ef19b-820b-49d6-a774-d7a799459cd3"
Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40"

[targets]
test = ["Test", "StableRNGs"]
Empty file added docs/src/interface.md
Empty file.
2 changes: 1 addition & 1 deletion src/DomainSets.jl
Original file line number Diff line number Diff line change
Expand Up @@ -43,7 +43,6 @@ export ..

import CompositeTypes: component, components


################################
## Exhaustive list of exports
################################
Expand Down Expand Up @@ -192,6 +191,7 @@ include("maps/basic.jl")
include("maps/affine.jl")
include("maps/arithmetics.jl")

include("generic/core.jl")
include("generic/domain.jl")
include("generic/geometry.jl")
include("generic/canonical.jl")
Expand Down
4 changes: 2 additions & 2 deletions src/domains/ball.jl
Original file line number Diff line number Diff line change
Expand Up @@ -56,7 +56,7 @@ isequaldomain(d1::Ball, d2::Ball) = isclosedset(d1)==isclosedset(d2) &&
radius(d1)==radius(d2) && center(d1)==center(d2)
hash(d::Ball, h::UInt) = hashrec("Ball", isclosedset(d), radius(d), center(d), h)

function issubset(d1::Ball, d2::Ball)
function issubset_domain(d1::Ball, d2::Ball)
if dimension(d1) == dimension(d2)
if center(d1) == center(d2)
if radius(d1) < radius(d2)
Expand Down Expand Up @@ -119,7 +119,7 @@ approx_indomain(x, d::ClosedUnitBall, tolerance) = norm(x) <= 1+tolerance
isequaldomain(d1::UnitBall, d2::UnitBall) = isclosedset(d1)==isclosedset(d2) &&
dimension(d1) == dimension(d2)

issubset(d1::UnitBall, d2::UnitBall) =
issubset_domain(d1::UnitBall, d2::UnitBall) =
dimension(d1) == dimension(d2) && (isclosedset(d2) || isopenset(d1))

convert(::Type{SublevelSet}, d::UnitBall{T,C}) where {T,C} =
Expand Down
34 changes: 18 additions & 16 deletions src/domains/boundingbox.jl
Original file line number Diff line number Diff line change
Expand Up @@ -3,52 +3,54 @@ boundingbox(d::Vector{T}) where {T <: Number} = minimum(d)..maximum(d)
boundingbox(d::Set{T}) where {T<:Number} = minimum(d)..maximum(d)

"Return the bounding box of the union of two or more bounding boxes."
unionbox(d::Domain) = d
unionbox(d1::Domain, d2::Domain) = unionbox(promote_domains(d1, d2)...)
unionbox(d1::Domain, d2::Domain, domains...) =
unionbox(d) = d
unionbox(d1, d2) = unionbox1(promote_domains(d1, d2)...)
unionbox(d1, d2, domains...) =
unionbox(unionbox(d1,d2), domains...)

unionbox(d1::Domain{T}, d2::Domain{T}) where {T} = unionbox1(d1, d2)
unionbox1(d1, d2) = unionbox2(d1, d2)
unionbox2(d1, d2) = fullspace(d1)
unionbox1(d1::EmptySpace, d2) = d2
unionbox1(d1::FullSpace, d2) = d1
unionbox2(d1, d2::EmptySpace) = d1
unionbox2(d1, d2::FullSpace) = d2

unionbox(d1::D, d2::D) where {D<:FixedInterval} = d1
unionbox1(d1::D, d2::D) where {D<:FixedInterval} = d1
Copy link
Member

Choose a reason for hiding this comment

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

What's unionbox1 and below intersectbox1?

Copy link
Member Author

@daanhb daanhb Sep 5, 2023

Choose a reason for hiding this comment

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

unionbox(d1,d2) is a function of two arguments, which sometimes simplifies for certain d1 (independently of d2 or for some possibilities of d2). Similarly, it may simplify for certain d2. The call chain goes:
unionbox(d1,d2) -> unionbox1(d1,d2) -> unionbox2(d1,d2) -> some_default_algorithm

The idea is that unionbox1 can dispatch safely (without ambiguity) solely on d1, and unionbox2 can dispatch on d2. I've used this pattern in many places and it works rather well. If a combination of d1 and d2 is specific enough, then one can dispatch on it using just the unionbox function.


function unionbox(d1::AbstractInterval{T}, d2::AbstractInterval{T}) where {T}
function unionbox1(d1::AbstractInterval, d2::AbstractInterval)
a, b = endpoints(d1)
c, d = endpoints(d2)
A = min(a, c)
B = max(b, d)
isinf(A) && isinf(B) ? fullspace(T) : A..B
isinf(A) && isinf(B) ? fullspace(d1) : A..B
end

unionbox(d1::HyperRectangle{T}, d2::HyperRectangle{T}) where {T} =
function unionbox(d1::HyperRectangle, d2::HyperRectangle)
@assert dimension(d1) == dimension(d2)
T = promote_type(eltype(d1),eltype(d2))
Rectangle{T}(map(unionbox, components(d1), components(d2)))
end

"Return the bounding box of the intersection of two or more bounding boxes."
intersectbox(d::Domain) = d
intersectbox(d1::Domain, d2::Domain) = intersectbox(promote_domains(d1, d2)...)
intersectbox(d1::Domain, d2::Domain, domains...) =
intersectbox(d) = d
intersectbox(d1, d2) = intersectbox1(promote_domains(d1, d2)...)
intersectbox(d1, d2, domains...) =
intersectbox(intersectbox(d1,d2), domains...)

intersectbox(d1::Domain{T}, d2::Domain{T}) where {T} = intersectbox1(d1, d2)
intersectbox1(d1, d2) = intersectbox2(d1, d2)
intersectbox2(d1, d2) = fullspace(d1)
intersectbox1(d1::EmptySpace, d2) = d1
intersectbox1(d1::FullSpace, d2) = d2
intersectbox2(d1, d2::EmptySpace) = d2
intersectbox2(d1, d2::FullSpace) = d1

intersectbox(d1::D, d2::D) where {D<:FixedInterval} = d1
intersectbox1(d1::D, d2::D) where {D<:FixedInterval} = d1

intersectbox(d1::AbstractInterval{T}, d2::AbstractInterval{T}) where {T} =
intersectbox1(d1::AbstractInterval, d2::AbstractInterval) =
intersectdomain(d1, d2)

function intersectbox(d1::HyperRectangle{T}, d2::HyperRectangle{T}) where {T}
function intersectbox1(d1::HyperRectangle, d2::HyperRectangle)
@assert dimension(d1) == dimension(d2)
T = eltype(d1)
d = Rectangle{T}(map(intersectbox, components(d1), components(d2)))
isempty(d) ? emptyspace(d) : d
end
Expand Down
4 changes: 2 additions & 2 deletions src/domains/cube.jl
Original file line number Diff line number Diff line change
Expand Up @@ -351,8 +351,8 @@ ProductDomain(domains::NTuple{N,D}) where {N,D <: FixedInterval} =
ProductDomain(domains::SVector{N,<:FixedInterval}) where {N} =
FixedIntervalProduct(domains)
ProductDomain{T}(domains::D...) where {N,S,T<:SVector{N,S},D <: FixedInterval} =
FixedIntervalProduct(convert.(Ref(Domain{S}), domains))
FixedIntervalProduct(convert_eltype.(Ref(S), domains))
ProductDomain{T}(domains::NTuple{N,D}) where {N,S,T<:SVector{N,S},D <: FixedInterval} =
FixedIntervalProduct(convert.(Ref(Domain{S}), domains))
FixedIntervalProduct(convert_eltype.(Ref(S), domains))
ProductDomain{T}(domains::SVector{N,<:FixedInterval}) where {N,T<:SVector{N}} =
FixedIntervalProduct(domains)
25 changes: 13 additions & 12 deletions src/domains/indicator.jl
Original file line number Diff line number Diff line change
Expand Up @@ -13,7 +13,7 @@ implementing `in` directly.
abstract type AbstractIndicatorFunction{T} <: Domain{T} end

"The indicator function of a domain is the function `f(x) = x ∈ D`."
indicatorfunction(d::Domain) = x -> x ∈ d
indicatorfunction(d) = x -> x ∈ d

indomain(x, d::AbstractIndicatorFunction) = _indomain(x, d, indicatorfunction(d))
_indomain(x, d::AbstractIndicatorFunction, f) = f(x)
Expand All @@ -35,8 +35,8 @@ indicatorfunction(d::IndicatorFunction) = d.f
similardomain(d::IndicatorFunction, ::Type{T}) where {T} = IndicatorFunction{T}(d.f)

convert(::Type{IndicatorFunction}, d::AbstractIndicatorFunction) = d
convert(::Type{IndicatorFunction}, d::Domain{T}) where {T} =
IndicatorFunction{T}(indicatorfunction(d))
convert(::Type{IndicatorFunction}, d) =
IndicatorFunction{deltype(d)}(indicatorfunction(checkdomain(d)))

isequaldomain(d1::IndicatorFunction, d2::IndicatorFunction) =
indicatorfunction(d1)==indicatorfunction(d2)
Expand All @@ -45,18 +45,19 @@ intersectdomain1(d1::IndicatorFunction, d2) = BoundedIndicatorFunction(d1.f, d2)
intersectdomain2(d1, d2::IndicatorFunction) = BoundedIndicatorFunction(d2.f, d1)

"An indicator function with a known bounding domain."
struct BoundedIndicatorFunction{F,D,T} <: AbstractIndicatorFunction{T}
struct BoundedIndicatorFunction{T,F,D} <: AbstractIndicatorFunction{T}
f :: F
domain :: D
end

BoundedIndicatorFunction(f, domain::Domain{T}) where {T} =
BoundedIndicatorFunction{typeof(f),typeof(domain),T}(f, domain)

function BoundedIndicatorFunction(f, domain)
T = eltype(domain)
BoundedIndicatorFunction{typeof(f),typeof(domain),T}(f, domain)
end
BoundedIndicatorFunction(f, domain) =
BoundedIndicatorFunction{deltype(domain)}(f, domain)
BoundedIndicatorFunction{T}(f, domain::Domain{T}) where {T} =
BoundedIndicatorFunction{T,typeof(f),typeof(domain)}(f, domain)
BoundedIndicatorFunction{T}(f, domain) where {T} =
_BoundedIndicatorFunction(T, f, checkdomain(domain))
_BoundedIndicatorFunction(::Type{T}, f, domain) where {T} =
BoundedIndicatorFunction{T,typeof(f),typeof(domain)}(f, domain)

indicatorfunction(d::BoundedIndicatorFunction) = d.f

Expand All @@ -70,7 +71,7 @@ hash(d::BoundedIndicatorFunction, h::UInt) =
hashrec(indicatorfunction(d), boundingdomain(d), h)

similardomain(d::BoundedIndicatorFunction, ::Type{T}) where {T} =
BoundedIndicatorFunction(d.f, convert(Domain{T}, d.domain))
BoundedIndicatorFunction(d.f, convert_eltype(T, d.domain))

boundingbox(d::BoundedIndicatorFunction) = boundingbox(boundingdomain(d))

Expand Down
24 changes: 14 additions & 10 deletions src/domains/interval.jl
Original file line number Diff line number Diff line change
Expand Up @@ -2,7 +2,7 @@
iscompact(d::TypedEndpointsInterval{:closed,:closed}) = true
iscompact(d::TypedEndpointsInterval) = false

isinterval(d::Domain) = false
isinterval(d) = false
isinterval(d::AbstractInterval) = true

hash(d::AbstractInterval, h::UInt) =
Expand Down Expand Up @@ -54,9 +54,9 @@ end
boundary(d::AbstractInterval) = Point(leftendpoint(d)) ∪ Point(rightendpoint(d))
corners(d::AbstractInterval) = [leftendpoint(d), rightendpoint(d)]

normal(d::AbstractInterval, x) = (abs(minimum(d)-x) < abs(maximum(d)-x)) ? -one(eltype(d)) : one(eltype(d))
normal(d::AbstractInterval, x) = (abs(minimum(d)-x) < abs(maximum(d)-x)) ? -one(deltype(d)) : one(deltype(d))

distance_to(d::AbstractInterval, x) = x ∈ d ? zero(eltype(d)) : min(abs(x-supremum(d)), abs(x-infimum(d)))
distance_to(d::AbstractInterval, x) = x ∈ d ? zero(deltype(d)) : min(abs(x-supremum(d)), abs(x-infimum(d)))

boundingbox(d::AbstractInterval) = d

Expand All @@ -79,6 +79,10 @@ and `ChebyshevInterval`.
abstract type FixedInterval{L,R,T} <: TypedEndpointsInterval{L,R,T} end
const ClosedFixedInterval{T} = FixedInterval{:closed,:closed,T}

convert(::Type{AbstractInterval{T}}, d::FixedInterval{L,R,T}) where {L,R,T} = d
convert(::Type{AbstractInterval{T}}, d::FixedInterval{L,R,S}) where {L,R,S,T} =
similardomain(d, T)

closure(d::AbstractInterval) = ClosedInterval(endpoints(d)...)
closure(d::ClosedInterval) = d

Expand Down Expand Up @@ -244,8 +248,8 @@ function similar_interval(d::HalfLine{T,C}, a::S, b::S) where {T,S,C}
HalfLine{promote_type(float(T),S),C}()
end

point_in_domain(d::NonnegativeRealLine) = zero(eltype(d))
point_in_domain(d::PositiveRealLine) = one(eltype(d))
point_in_domain(d::NonnegativeRealLine) = zero(deltype(d))
Copy link
Member

Choose a reason for hiding this comment

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

Probably this name is clearer, but this is called a "choice function" in the Axiom of Choice so we could just call this choice(d)

Copy link
Member Author

Choose a reason for hiding this comment

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

Hmm, it would be good if choice(s) was defined for a Set too, for example, but it isn't. This function really is just for testing purposes, but it could have any name.

Copy link
Member

Choose a reason for hiding this comment

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

Just add choice(d::AbstractSet) = first(d)?

Copy link
Member

Choose a reason for hiding this comment

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

iterate(d)[1] is kind of the same...

Copy link
Member Author

Choose a reason for hiding this comment

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

the point_in_domain function roughly speaking picks some point in the interior of a domain, if possible. It may be good to think about how to "sample" points from a domain more generally (like rand), but that is a separate issue.

Copy link
Member Author

Choose a reason for hiding this comment

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

Okay, I'm renaming to choice and adding a deprecation for point_in_domain

point_in_domain(d::PositiveRealLine) = one(deltype(d))


"The negative halfline `(-∞,0]` or `(-∞,0)`, right-closed or right-open."
Expand Down Expand Up @@ -276,8 +280,8 @@ function similar_interval(d::NegativeHalfLine{T,C}, a::S, b::S) where {S,T,C}
NegativeHalfLine{promote_type(S,float(T)),C}()
end

point_in_domain(d::NegativeRealLine) = -one(eltype(d))
point_in_domain(d::NonpositiveRealLine) = zero(eltype(d))
point_in_domain(d::NegativeRealLine) = -one(deltype(d))
point_in_domain(d::NonpositiveRealLine) = zero(deltype(d))


"The real line `(-∞,∞)`."
Expand Down Expand Up @@ -318,11 +322,11 @@ intersect(d1::FixedInterval, d2::FixedInterval) = intersectdomain(d1, d2)

# Promotion to joint type T
uniondomain(d1::TypedEndpointsInterval, d2::TypedEndpointsInterval) =
uniondomain(promote(d1,d2)...)
uniondomain(promote_domains(d1,d2)...)
intersectdomain(d1::TypedEndpointsInterval, d2::TypedEndpointsInterval) =
intersectdomain(promote(d1,d2)...)
intersectdomain(promote_domains(d1,d2)...)
setdiffdomain(d1::AbstractInterval, d2::AbstractInterval) =
setdiffdomain(promote(d1,d2)...)
setdiffdomain(promote_domains(d1,d2)...)



Expand Down
11 changes: 7 additions & 4 deletions src/domains/trivial.jl
Original file line number Diff line number Diff line change
Expand Up @@ -13,7 +13,7 @@ EmptySpace(::Type{T}) where {T} = EmptySpace{T}()
similardomain(::EmptySpace, ::Type{T}) where {T} = EmptySpace{T}()

"Return the empty space with the same element type as the given domain."
emptyspace(d) = emptyspace(eltype(d))
emptyspace(d) = emptyspace(deltype(d))
emptyspace(::Type{T}) where {T} = EmptySpace{T}()

indomain(x::T, d::EmptySpace{T}) where {T} = false
Expand Down Expand Up @@ -64,14 +64,17 @@ struct FullSpace{T} <: Domain{T} end
const AnyFullSpace = FullSpace{Any}

FullSpace() = FullSpace{Float64}()
FullSpace(d) = FullSpace{eltype(d)}()
FullSpace(d) = FullSpace{deltype(d)}()

"Return the full space with the same element type as the given domain."
fullspace(d) = fullspace(eltype(d))
fullspace(d) = fullspace(deltype(d))
fullspace(::Type{T}) where {T} = FullSpace{T}()

isfullspace(d::FullSpace) = true
isfullspace(d::Domain) = false
isfullspace(d) = _isfullspace(d, DomainStyle(d))
_isfullspace(d, ::IsDomain) = false
_isfullspace(d, ::NotDomain) = error("isfullspace invoked on non-domain type")

similardomain(::FullSpace, ::Type{T}) where {T} = FullSpace{T}()

Expand Down Expand Up @@ -125,7 +128,7 @@ type `T`.
struct TypeDomain{T} <: Domain{T} end

"Return the domain for the element type of the given domain."
typedomain(d) = typedomain(eltype(d))
typedomain(d) = typedomain(deltype(d))
typedomain(::Type{T}) where {T} = TypeDomain{T}()

iscompatiblepair(x::T, d::TypeDomain{T}) where {T} = true
Expand Down
Loading