diff --git a/FunctionMaps/LICENSE b/FunctionMaps/LICENSE deleted file mode 100644 index f52eae8..0000000 --- a/FunctionMaps/LICENSE +++ /dev/null @@ -1,21 +0,0 @@ -MIT License - -Copyright (c) 2024 Daan Huybrechs and contributors - -Permission is hereby granted, free of charge, to any person obtaining a copy -of this software and associated documentation files (the "Software"), to deal -in the Software without restriction, including without limitation the rights -to use, copy, modify, merge, publish, distribute, sublicense, and/or sell -copies of the Software, and to permit persons to whom the Software is -furnished to do so, subject to the following conditions: - -The above copyright notice and this permission notice shall be included in all -copies or substantial portions of the Software. - -THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR -IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, -FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE -AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER -LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, -OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE -SOFTWARE. diff --git a/FunctionMaps/Project.toml b/FunctionMaps/Project.toml deleted file mode 100644 index f4bfe38..0000000 --- a/FunctionMaps/Project.toml +++ /dev/null @@ -1,24 +0,0 @@ -name = "FunctionMaps" -uuid = "a75834df-64dc-405c-a530-e22a155ff74e" -authors = ["Daan Huybrechs and contributors"] -version = "0.1.0" - -[deps] -CompositeTypes = "b152e2b5-7a66-4b01-a709-34e65c35f657" -LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e" -StaticArrays = "90137ffa-7385-5640-81b9-e52037218182" - -[compat] -julia = "1.6" -Aqua = "0.8" -CompositeTypes = "0.1.3" -LinearAlgebra = "1" -StaticArrays = "1" -Test = "1" - -[extras] -Aqua = "4c88cf16-eb10-579e-8560-4a9242c79595" -Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40" - -[targets] -test = ["Test", "Aqua"] diff --git a/FunctionMaps/README.md b/FunctionMaps/README.md deleted file mode 100644 index e347aa0..0000000 --- a/FunctionMaps/README.md +++ /dev/null @@ -1,9 +0,0 @@ -# FunctionMaps - -[![Stable](https://img.shields.io/badge/docs-stable-blue.svg)](https://daanhb.github.io/FunctionMaps.jl/stable/) -[![Dev](https://img.shields.io/badge/docs-dev-blue.svg)](https://daanhb.github.io/FunctionMaps.jl/dev/) -[![Build Status](https://github.com/daanhb/FunctionMaps.jl/actions/workflows/CI.yml/badge.svg?branch=main)](https://github.com/daanhb/FunctionMaps.jl/actions/workflows/CI.yml?query=branch%3Amain) -[![Coverage](https://codecov.io/gh/daanhb/FunctionMaps.jl/branch/main/graph/badge.svg)](https://codecov.io/gh/daanhb/FunctionMaps.jl) - -FunctionMaps.jl is a package designed to represent scalar- and vector-valued functions. It makes it easy to make new maps from existing ones, via composition -or by taking products. diff --git a/FunctionMaps/src/FunctionMaps.jl b/FunctionMaps/src/FunctionMaps.jl deleted file mode 100644 index 9105721..0000000 --- a/FunctionMaps/src/FunctionMaps.jl +++ /dev/null @@ -1,62 +0,0 @@ -module FunctionMaps - -using CompositeTypes, CompositeTypes.Display, CompositeTypes.Indexing -using LinearAlgebra -using StaticArrays - -import Base: - convert, show, - ==, hash, # for the equality of maps - ∘, # used to compose maps - isreal # to check whether maps are real - -import CompositeTypes: component, components - -# Exhaustive list of exports: - -# from util/common.jl -export prectype, numtype - -# from generic/map.jl -export Map, MapRef, - applymap, isequalmap, - domaintype, codomaintype, - inverse, leftinverse, rightinverse, - mapsize, jacobian, jacdet, diffvolume, - isreal -# from generic/canonical.jl -export canonicalmap -# from generic/composite.jl -export composedmap -# from generic/product.jl -export productmap - -# from types/basic.jl -export IdentityMap, - StaticIdentityMap, VectorIdentityMap, - ZeroMap, UnityMap, ConstantMap, - isconstantmap, mapconstant -# from types/affine.jl -export AffineMap, Translation, LinearMap, - affinematrix, affinevector, - islinearmap, isaffinemap - -include("util/common.jl") - -include("generic/map.jl") -include("generic/interface.jl") -include("generic/canonical.jl") -include("generic/lazy.jl") -include("generic/inverse.jl") -include("generic/jacobian.jl") -include("generic/composite.jl") -include("generic/product.jl") -include("generic/isomorphism.jl") -include("types/basic.jl") -include("types/affine.jl") -include("types/coordinates.jl") -include("types/arithmetics.jl") - -include("deprecated.jl") - -end diff --git a/FunctionMaps/src/deprecated.jl b/FunctionMaps/src/deprecated.jl deleted file mode 100644 index f575f76..0000000 --- a/FunctionMaps/src/deprecated.jl +++ /dev/null @@ -1,16 +0,0 @@ - -# some renames, to be deprecated in a later breaking release -islinear(m::Map) = islinearmap(m) -isaffine(m::Map) = isaffinemap(m) -isconstant(m::Map) = isconstantmap(m) -constant(m::Map) = mapconstant(m) -isidentity(m::Map) = isidentitymap(m) - -isreal(::Type{T}) where T = isrealtype(T) - -matrix(m::Map) = affinematrix(m) -vector(m::Map) = affinevector(m) - -@deprecate convert_domaintype(map::Map, ::Type{T}) where {T} convert_domaintype(T, map) -@deprecate convert_numtype(map::Map{T}, ::Type{U}) where {T,U} convert_numtype(U, map) -@deprecate convert_prectype(map::Map{T}, ::Type{U}) where {T,U} convert_prectype(U, map) diff --git a/FunctionMaps/src/generic/canonical.jl b/FunctionMaps/src/generic/canonical.jl deleted file mode 100644 index b951982..0000000 --- a/FunctionMaps/src/generic/canonical.jl +++ /dev/null @@ -1,94 +0,0 @@ -""" - canonicalmap([ctype::CanonicalType, ]map) - -Return an associated canonical map, if any, of the given map. - -Optionally, a canonical type argument may specify an alternative canonical map. -Canonical maps help with converting between equal maps of different types. -""" -canonicalmap(m) = m - -"Does the map have a canonical map?" -hascanonicalmap(m) = !(m === canonicalmap(m)) - - -"Supertype of different kinds of canonical objects." -abstract type CanonicalType end - -canonicalmap(ctype::CanonicalType, m) = m -hascanonicalmap(ctype::CanonicalType, m) = !(m === canonicalmap(ctype, m)) - - -""" - Equal <: CanonicalType - -A canonical object that is equal but simpler. -""" -struct Equal <: CanonicalType end - -canonicalmap(::Equal, m) = m -equalmap(m) = canonicalmap(Equal(), m) -hasequalmap(m) = hascanonicalmap(Equal(), m) - -""" - simplify(m) - -Simplify the given map to an equal map. -""" -simplify(m) = equalmap(m) -""" - simplifies(m) - -Does the map simplify? -""" -simplifies(m) = hasequalmap(m) - -"Convert the given map to a map defined in FunctionMaps.jl." -tofunctionmap(m::Map) = m -tofunctionmap(m::MapRef) = tofunctionmap(functionmap(m)) - - -""" - Equivalent <: CanonicalType - -A canonical object that is equivalent but may have different type. -""" -struct Equivalent <: CanonicalType end - -canonicalmap(::Equivalent, m) = canonicalmap(Equal(), m) - -canonicalmap(::Equivalent, m::Map{<:StaticVector{1,T}}) where {T} = convert(Map{T}, m) -canonicalmap(::Equivalent, m::Map{NTuple{N,T}}) where {N,T} = convert(Map{SVector{N,T}}, m) - -equivalentmap(m) = canonicalmap(Equivalent(), m) -hasequivalentmap(m) = hascanonicalmap(Equivalent(), m) - -""" - CanonicalExtensionType <: CanonicalType - -Canonical types used to translate between packages. -""" -abstract type CanonicalExtensionType <: CanonicalType -end - -"Return the extension type associated with the given object." -canonicalextensiontype(m) = canonicalextensiontype(typeof(m)) - - -==(m1::Map, m2::Map) = isequalmap(m1, m2) # Method from Base - -""" - isequalmap(map1, map2) - -Are the two given maps equal? -""" -isequalmap(m1, m2) = isequalmap1(m1, m2) -isequalmap1(m1, m2) = simplifies(m1) ? isequalmap(simplify(m1), m2) : isequalmap2(m1, m2) -isequalmap2(m1, m2) = simplifies(m2) ? isequalmap(m1, simplify(m2)) : default_isequalmap(m1, m2) -default_isequalmap(m1, m2) = m1===m2 - -# Associated with == we have to define the hashes of maps -Base.hash(m::Map, h::UInt) = map_hash(m, h) - -map_hash(m) = map_hash(m, zero(UInt)) -map_hash(m, h::UInt) = invoke(hash, Tuple{Any,UInt}, m, h) diff --git a/FunctionMaps/src/generic/composite.jl b/FunctionMaps/src/generic/composite.jl deleted file mode 100644 index d6923ac..0000000 --- a/FunctionMaps/src/generic/composite.jl +++ /dev/null @@ -1,228 +0,0 @@ - -"The composition of several maps." -struct ComposedMap{T,MAPS} <: CompositeLazyMap{T} - maps :: MAPS -end - -ComposedMap(maps...) = ComposedMap{domaintype(maps[1])}(maps...) -ComposedMap{T}(maps...) where {T} = ComposedMap{T,typeof(maps)}(maps) - -# TODO: make proper conversion -similarmap(m::ComposedMap, ::Type{T}) where {T} = ComposedMap{T}(m.maps...) - -codomaintype(m::ComposedMap) = codomaintype(m.maps[end]) - -# Maps are applied in the order that they appear in m.maps -applymap(m::ComposedMap, x) = applymap_rec(x, m.maps...) -applymap_rec(x) = x -applymap_rec(x, map1, maps...) = applymap_rec(applymap(map1, x), maps...) - -# The size of a composite map depends on the first and the last map to be applied -# We check whether they are scalar_to_vector, vector_to_vector, etcetera -mapsize(m::ComposedMap) = _composed_mapsize(m, m.maps[end], m.maps[1], mapsize(m.maps[end]), mapsize(m.maps[1])) -_composed_mapsize(m, m_end, m1, S_end::Tuple{Int,Int}, S1::Tuple{Int,Int}) = (S_end[1],S1[2]) -_composed_mapsize(m, m_end, m1, S_end::Tuple{Int,Int}, S1::Tuple{Int}) = - is_vector_to_scalar(m_end) ? () : (S_end[1],) -_composed_mapsize(m, m_end, m1, S_end::Tuple{Int,Int}, S1::Tuple{}) = - is_vector_to_scalar(m_end) ? () : (S_end[1],) -_composed_mapsize(m, m_end, m1, S_end::Tuple{Int}, S1::Tuple{Int,Int}) = (S_end[1],S1[2]) -_composed_mapsize(m, m_end, m1, S_end::Tuple{Int}, S1::Tuple{Int}) = (S_end[1],) -_composed_mapsize(m, m_end, m1, S_end::Tuple{Int}, S1::Tuple{}) = (S_end[1],) -_composed_mapsize(m, m_end, m1, S_end::Tuple{}, S1::Tuple{Int,Int}) = (1,S1[2]) -_composed_mapsize(m, m_end, m1, S_end::Tuple{}, S1::Tuple{Int}) = () -_composed_mapsize(m, m_end, m1, S_end::Tuple{}, S1::Tuple{}) = () - -function jacobian(m::ComposedMap, x) - f, fd = backpropagate(x, reverse(components(m))...) - fd -end -backpropagate(x, m1) = (applymap(m1, x), jacobian(m1, x)) -function backpropagate(x, m2, ms...) - f, fd = backpropagate(x, ms...) - applymap(m2, f), jacobian(m2, f) * fd -end - -for op in (:inverse, :leftinverse, :rightinverse) - @eval $op(cmap::ComposedMap) = ComposedMap(reverse(map($op, components(cmap)))...) -end - -inverse(m::ComposedMap, x) = inverse_rec(x, reverse(components(m))...) -inverse_rec(x) = x -inverse_rec(x, map1, maps...) = inverse_rec(inverse(map1, x), maps...) - -leftinverse(m::ComposedMap, x) = leftinverse_rec(x, reverse(components(m))...) -leftinverse_rec(x) = x -leftinverse_rec(x, map1, maps...) = leftinverse_rec(leftinverse(map1, x), maps...) - -rightinverse(m::ComposedMap, x) = rightinverse_rec(x, reverse(components(m))...) -rightinverse_rec(x) = x -rightinverse_rec(x, map1, maps...) = rightinverse_rec(rightinverse(map1, x), maps...) - - -promote_composing_maps(m1, m2) = - _promote_composing_maps(codomaintype(m1), domaintype(m2), m1, m2) -function _promote_composing_maps(::Type{S}, ::Type{T}, m1, m2) where {S,T} - U = promote_type(S, T) - convert_codomaintype(U, m1), convert_domaintype(U, m2) -end - -composedmap() = () -composedmap(m) = m -composedmap(m1, m2) = composedmap1(promote_composing_maps(m1, m2)...) -composedmap1(m1, m2) = composedmap2(m1, m2) -composedmap2(m1, m2) = default_composedmap(m1, m2) -default_composedmap(m1, m2) = ComposedMap(m1, m2) - -composedmap(m1, m2, maps...) = composedmap(composedmap(m1, m2), maps...) - -composedmap(m1::ComposedMap, m2::ComposedMap) = - ComposedMap(components(m1)..., components(m2)...) -function composedmap1(m1::ComposedMap, m2) - T = promote_type(codomaintype(m1), domaintype(m2)) - ComposedMap(components(m1)..., convert_domaintype(T, m2)) -end -function composedmap2(m1, m2::ComposedMap) - T = promote_type(codomaintype(m1), domaintype(m2)) - ComposedMap(convert_codomaintype(T, m1), components(m2)...) -end - -# Arguments to ∘ should be reversed before passing on to mapcompose -Base.:∘(map1::AbstractMap, map2::AbstractMap) = composedmap(map2, map1) - - -isequalmap(m1::ComposedMap, m2::ComposedMap) = - ncomponents(m1) == ncomponents(m2) && all(map(isequalmap, components(m1), components(m2))) -map_hash(m::ComposedMap, h::UInt) = hashrec("ComposedMap", collect(components(m)), h) - -Display.combinationsymbol(m::ComposedMap) = Display.Symbol('∘') -Display.displaystencil(m::ComposedMap) = - composite_displaystencil(m; reversecomponents=true) -map_object_parentheses(m::ComposedMap) = true -map_stencil_parentheses(m::ComposedMap) = true -show(io::IO, mime::MIME"text/plain", m::ComposedMap) = composite_show(io, mime, m) -show(io::IO, m::ComposedMap) = composite_show_compact(io, m) - -## Lazy multiplication - -"The lazy multiplication of one or more maps." -struct MulMap{T,MAPS} <: CompositeLazyMap{T} - maps :: MAPS -end - -MulMap(map1::Map{T}, maps::Map{T}...) where {T} = MulMap{T}(map1, maps...) -MulMap{T}(maps::Map{T}...) where {T} = MulMap{T,typeof(maps)}(maps) -MulMap{T}(maps...) where {T} = _mulmap(T, convert.(Map{T}, maps)...) -_mulmap(::Type{T}, maps...) where {T} = MulMap{T,typeof(maps)}(maps) - -similarmap(m::MulMap, ::Type{T}) where {T} = MulMap{T}(m.maps...) - -applymap(m::MulMap, x) = reduce(*, applymap.(components(m), Ref(x))) - -multiply_map() = () -multiply_map(m) = m -multiply_map(m1, m2) = multiply_map1(m1, m2) -multiply_map1(m1, m2) = multiply_map2(m1, m2) -multiply_map2(m1, m2) = default_multiply_map(m1, m2) -default_multiply_map(m1, m2) = MulMap(m1, m2) - -multiply_map(m1, m2, maps...) = multiply_map(multiply_map(m1, m2), maps...) - -multiply_map(m1::MulMap, m2::MulMap) = - MulMap(components(m1)..., components(m2)...) -multiply_map1(m1::MulMap, m2) = MulMap(components(m1)..., m2) -multiply_map2(m1, m2::MulMap) = MulMap(m1, components(m2)...) - -function Display.displaystencil(m::MulMap) - A = Any[] - list = components(m) - push!(A, "x -> ") - push!(A, Display.SymbolObject(list[1])) - push!(A, "(x)") - for i in 2:length(list) - push!(A, " * ") - push!(A, Display.SymbolObject(list[i])) - push!(A, "(x)") - end - A -end -show(io::IO, mime::MIME"text/plain", m::MulMap) = composite_show(io, mime, m) - - -## Lazy sum - -"The lazy sum of one or more maps." -struct SumMap{T,MAPS} <: CompositeLazyMap{T} - maps :: MAPS -end - -SumMap(map1::Map{T}, maps::Map{T}...) where {T} = SumMap{T}(map1, maps...) -SumMap{T}(maps::Map{T}...) where {T} = SumMap{T,typeof(maps)}(maps) -SumMap{T}(maps...) where {T} = _summap(T, convert.(Map{T}, maps)...) -_summap(::Type{T}, maps...) where {T} = SumMap{T,typeof(maps)}(maps) - -similarmap(m::SumMap, ::Type{T}) where {T} = SumMap{T}(m.maps...) - -applymap(m::SumMap, x) = reduce(+, applymap.(components(m), Ref(x))) - -sum_map() = () -sum_map(m) = m -sum_map(m1, m2) = sum_map1(m1, m2) -sum_map1(m1, m2) = sum_map2(m1, m2) -sum_map2(m1, m2) = default_sum_map(m1, m2) -default_sum_map(m1, m2) = SumMap(m1, m2) - -sum_map(m1, m2, maps...) = sum_map(sum_map(m1, m2), maps...) - -sum_map(m1::SumMap, m2::SumMap) = - SumMap(components(m1)..., components(m2)...) -sum_map1(m1::SumMap, m2) = SumMap(components(m1)..., m2) -sum_map2(m1, m2::SumMap) = SumMap(m1, components(m2)...) - -function Display.displaystencil(m::SumMap) - A = Any[] - list = components(m) - push!(A, "x -> ") - push!(A, Display.SymbolObject(list[1])) - push!(A, "(x)") - for i in 2:length(list) - push!(A, " + ") - push!(A, Display.SymbolObject(list[i])) - push!(A, "(x)") - end - A -end -show(io::IO, mime::MIME"text/plain", m::SumMap) = composite_show(io, mime, m) - - -# Define the jacobian of a composite map -jacobian(m::ComposedMap) = composite_jacobian(reverse(components(m))...) -composite_jacobian(map1) = jacobian(map1) -composite_jacobian(map1, map2) = multiply_map(composedmap(map2, jacobian(map1)), jacobian(map2)) -function composite_jacobian(map1, map2, maps...) - rest = ComposedMap(reverse(maps)..., map2) - f1 = composedmap(rest, jacobian(map1)) - f2 = composite_jacobian(map2, maps...) - multiply_map(f1, f2) -end - -jacobian(m::MulMap) = mul_jacobian(components(m)...) -mul_jacobian() = () -mul_jacobian(map1) = jacobian(map1) -mul_jacobian(map1, map2) = sum_map(multiply_map(jacobian(map1), map2), multiply_map(map1, jacobian(map2))) -function mul_jacobian(map1, map2, maps...) - rest = multiply_map(map2, maps...) - mul_jacobian(map1, rest) -end -function jacobian(m::MulMap, x) - z = map(t -> applymap(t,x), components(m)) - zd = map(t -> jacobian(t, x), components(m)) - sum(prod(z[1:i-1]) * zd[i] * prod(z[i+1:end]) for i in 1:ncomponents(m)) -end - - -jacobian(m::SumMap) = sum_jacobian(components(m)...) -sum_jacobian() = () -sum_jacobian(map1) = jacobian(map1) -sum_jacobian(maps...) = sum_map(map(jacobian, maps)...) - -jacobian(m::SumMap, x) = sum(jacobian(mc, x) for mc in components(m)) diff --git a/FunctionMaps/src/generic/interface.jl b/FunctionMaps/src/generic/interface.jl deleted file mode 100644 index 24ac3cc..0000000 --- a/FunctionMaps/src/generic/interface.jl +++ /dev/null @@ -1,73 +0,0 @@ - -""" - MapStyle(m) - -Trait to indicate whether or not `m` implements the map interface. -""" -abstract type MapStyle end - -""" - IsMap() - -indicates an object implements the map interface. -""" -struct IsMap <: MapStyle end - -""" - NotMap() - -indicates an object does not implement the Map interface. -""" -struct NotMap <: MapStyle end - - -MapStyle(m) = MapStyle(typeof(m)) -# the default is no map -MapStyle(::Type) = NotMap() -# subtypes of Map are maps -MapStyle(::Type{<:Map}) = IsMap() - -""" - functionmap(m) - -Return a map associated with the object `m`. -""" -functionmap(m::Map) = m - -""" - MapRef(m) - -A reference to a map. - -In a function call, `MapRef(x)` can be used to indicate that `x` should be -treated as a map, e.g., `foo(x, MapRef(m))`. -""" -struct MapRef{M} - map :: M -end - -functionmap(m::MapRef) = m.map -domaintype(m::MapRef) = domaintype(functionmap(m)) - - -""" -`AnyMap` is the union of `Map` and `MapRef`. - -In both cases `map(m::AnyMap)` returns the map itself. -""" -const AnyMap = Union{Map,MapRef} - -""" - checkmap(m) - -Checks that `m` is a map or refers to a map and if so returns that map, -throws an error otherwise. -""" -checkmap(m::Map) = m -# we trust the explicit intention of a user providing a map reference -checkmap(m::MapRef) = functionmap(m) -# for other objects we check MapStyle -checkmap(m) = _checkmap(m, MapStyle(m)) -_checkmap(m, ::IsMap) = m -_checkmap(m, ::NotMap) = - throw(ArgumentError("Map does not implement map interface as indicated by MapStyle.")) diff --git a/FunctionMaps/src/generic/inverse.jl b/FunctionMaps/src/generic/inverse.jl deleted file mode 100644 index 3af081e..0000000 --- a/FunctionMaps/src/generic/inverse.jl +++ /dev/null @@ -1,56 +0,0 @@ - - -"A lazy inverse stores a map `m` and returns `inverse(m, x)`." -struct LazyInverse{T,M} <: SimpleLazyMap{T} - map :: M -end - -LazyInverse(m) = LazyInverse{codomaintype(m)}(m) -LazyInverse{T}(m) where {T} = LazyInverse{T,typeof(m)}(m) - -applymap(m::LazyInverse, x) = inverse(supermap(m), x) - -Display.displaystencil(m::LazyInverse) = ["LazyInverse(", supermap(m), ")"] -show(io::IO, mime::MIME"text/plain", m::LazyInverse) = composite_show(io, mime, m) - - -""" - inverse(m[, x]) - -Return the inverse of `m`. The two-argument function evaluates the inverse -at the point `x`. -""" -inverse(m) = LazyInverse(m) -inverse(m::LazyInverse) = supermap(m) -# Concrete maps should implement inverse(m, x) - -Base.:\(m::AbstractMap, x) = inverse(m, x) - -implements_inverse(m) = !(inverse(m) isa LazyInverse) - -""" - leftinverse(m[, x]) - -Return a left inverse of the given map. This left inverse `mli` is not unique, -but in any case it is such that `(mli ∘ m) * x = x` for each `x` in the domain -of `m`. - -The two-argument function applies the left inverse to the point `x`. -""" -leftinverse(m) = inverse(m) -leftinverse(m, x) = inverse(m, x) - -""" - rightinverse(m[, x]) - -Return a right inverse of the given map. This right inverse `mri` is not unique, -but in any case it is such that `(m ∘ mri) * y = y` for each `y` in the range -of `m`. - -The two-argument function applies the right inverse to the point `x`. -""" -rightinverse(m) = inverse(m) -rightinverse(m, x) = inverse(m, x) - -isequalmap(m1::LazyInverse, m2::LazyInverse) = isequalmap(supermap(m1), supermap(m2)) -map_hash(m::LazyInverse, h::UInt) = hashrec("LazyInverse", supermap(m), h) diff --git a/FunctionMaps/src/generic/isomorphism.jl b/FunctionMaps/src/generic/isomorphism.jl deleted file mode 100644 index 5155f03..0000000 --- a/FunctionMaps/src/generic/isomorphism.jl +++ /dev/null @@ -1,124 +0,0 @@ - -""" - Isomorphism{T,U} <: TypedMap{T,U} - -An isomorphism is a bijection between types that preserves norms. -""" -abstract type Isomorphism{T,U} <: TypedMap{T,U} end - -show(io::IO, m::Isomorphism{T,U}) where {T,U} = print(io, "x : $(T) -> x : $(U)") -map_object_parentheses(m::Isomorphism) = true - -struct VectorToNumber{T} <: Isomorphism{SVector{1,T},T} -end -struct NumberToVector{T} <: Isomorphism{T,SVector{1,T}} -end - -""" - VectorToNumber() - VectorToNumber{T}() - -Map a length 1 vector `x` to the number `x[1]`. - -See also: [`NumberToVector`](@ref). -""" -VectorToNumber() = VectorToNumber{Float64}() - -""" - NumberToVector() - NumberToVector{T}() - -Map a number `x` to the length 1 vector `[x]`. - -See also: [`VectorToNumber`](@ref). -""" -NumberToVector() = NumberToVector{Float64}() - -mapsize(::VectorToNumber) = (1,1) -mapsize(::NumberToVector) = (1,) - -inverse(::VectorToNumber{T}) where {T} = NumberToVector{T}() -inverse(::NumberToVector{T}) where {T} = VectorToNumber{T}() -inverse(m::VectorToNumber, x) = inverse(m)(x) -inverse(m::NumberToVector, x) = inverse(m)(x) - -applymap(::VectorToNumber, x) = x[1] -applymap(::NumberToVector, x) = SVector(x) - -jacobian(::VectorToNumber{T}, x) where {T} = transpose(SVector(one(T))) -jacobian(::VectorToNumber{T}) where {T} = ConstantMap{SVector{1,T}}(transpose(SVector(one(T)))) - -jacobian(::NumberToVector{T}, x) where {T} = SVector(one(T)) -jacobian(::NumberToVector{T}) where {T} = ConstantMap{T}(SVector(one(T))) - - -struct VectorToComplex{T} <: Isomorphism{SVector{2,T},Complex{T}} -end -struct ComplexToVector{T} <: Isomorphism{Complex{T},SVector{2,T}} -end - -""" - VectorToComplex() - VectorToComplex{T}() - -Map a length 2 vector ``[a;b]`` to the complex number ``a+bi``. - -See also: [`ComplexToVector`](@ref). -""" -VectorToComplex() = VectorToComplex{Float64}() - -""" - ComplexToVector() - ComplexToVector{T}() - -Map a complex number ``a+bi`` to the length 2 vector ``[a; b]``. - -See also: [`VectorToComplex`](@ref). -""" -ComplexToVector() = ComplexToVector{Float64}() - -mapsize(::VectorToComplex) = (1,2) -mapsize(::ComplexToVector) = (2,) - -applymap(::VectorToComplex, x) = x[1] + im*x[2] -applymap(::ComplexToVector, x) = SVector(real(x), imag(x)) - -inverse(::VectorToComplex{T}) where {T} = ComplexToVector{T}() -inverse(::ComplexToVector{T}) where {T} = VectorToComplex{T}() -inverse(m::VectorToComplex, x) = inverse(m)(x) -inverse(m::ComplexToVector, x) = inverse(m)(x) - -jacobian(::VectorToComplex{T}, x) where {T} = transpose(SVector(one(T),one(T)*im)) -jacobian(::VectorToComplex{T}) where {T} = ConstantMap{SVector{2,T}}(transpose(SVector(one(T),one(T)*im))) - - -"Map a static vector to a tuple." -struct VectorToTuple{N,T} <: Isomorphism{SVector{N,T},NTuple{N,T}} -end -"Map a tuple to a static vector." -struct TupleToVector{N,T} <: Isomorphism{NTuple{N,T},SVector{N,T}} -end - -inverse(::VectorToTuple{N,T}) where {N,T} = TupleToVector{N,T}() -inverse(::TupleToVector{N,T}) where {N,T} = VectorToTuple{N,T}() -inverse(m::VectorToTuple, x) = inverse(m)(x) -inverse(m::TupleToVector, x) = inverse(m)(x) - -applymap(::VectorToTuple, x) = tuple(x...) -applymap(::TupleToVector, x) = SVector(x) - - -"Map a nested vector or tuple to a flat vector." -struct NestedToFlat{N,T,U,DIM} <: Isomorphism{U,SVector{N,T}} -end -"Map a flattened vector to a nested one." -struct FlatToNested{N,T,U,DIM} <: Isomorphism{SVector{N,T},U} -end - -inverse(::NestedToFlat{N,T,U,DIM}) where {N,T,U,DIM} = FlatToNested{N,T,U,DIM}() -inverse(::FlatToNested{N,T,U,DIM}) where {N,T,U,DIM} = NestedToFlat{N,T,U,DIM}() -inverse(m::NestedToFlat, x) = inverse(m)(x) -inverse(m::FlatToNested, x) = inverse(m)(x) - -applymap(::NestedToFlat{N,T,U,DIM}, x) where {N,T,U,DIM} = convert_tocartesian(x, Val{DIM}()) -applymap(::FlatToNested{N,T,U,DIM}, x) where {N,T,U,DIM} = convert_fromcartesian(x, Val{DIM}()) diff --git a/FunctionMaps/src/generic/jacobian.jl b/FunctionMaps/src/generic/jacobian.jl deleted file mode 100644 index a6b30dd..0000000 --- a/FunctionMaps/src/generic/jacobian.jl +++ /dev/null @@ -1,179 +0,0 @@ - -"A lazy Jacobian `J` stores a map `m` and returns `J(x) = jacobian(m, x)`." -struct LazyJacobian{T,M} <: SimpleLazyMap{T} - map :: M -end - -LazyJacobian(m) = LazyJacobian{domaintype(m),typeof(m)}(m) -applymap(m::LazyJacobian, x) = jacobian(supermap(m), x) - -mapsize(m::LazyJacobian) = mapsize(supermap(m)) - -Display.displaystencil(m::LazyJacobian) = ["LazyJacobian(", supermap(m), ")"] -show(io::IO, mime::MIME"text/plain", m::LazyJacobian) = composite_show(io, mime, m) - - -""" - jacobian(m[, x]) - -Return the jacobian map. The two-argument version evaluates the jacobian -at a point `x`. -""" -jacobian(m) = LazyJacobian(m) -jacobian!(y, m, x) = y .= jacobian(m, x) - - -"A `DeterminantMap` returns the determinant of the result of a given map." -struct DeterminantMap{T,M} <: SimpleLazyMap{T} - map :: M -end - -DeterminantMap(m) = DeterminantMap{domaintype(m)}(m) -DeterminantMap{T}(m) where {T} = DeterminantMap{T,typeof(m)}(m) -DeterminantMap{T}(m::Map{T}) where {T} = DeterminantMap{T,typeof(m)}(m) -DeterminantMap{T}(m::Map{S}) where {S,T} = DeterminantMap{T}(convert(Map{T}, m)) - -determinantmap(m) = DeterminantMap(m) - -applymap(m::DeterminantMap, x) = det(supermap(m)(x)) - - -""" - jacdet(m[, x]) - -Return the determinant of the jacobian as a map. The two-argument version -evaluates the jacobian determinant at a point `x`. -""" -jacdet(m) = determinantmap(jacobian(m)) -jacdet(m, x) = det(jacobian(m, x)) - - -"An `AbsMap` returns the absolute value of the result of a given map." -struct AbsMap{T,M} <: SimpleLazyMap{T} - map :: M -end - -AbsMap(m) = AbsMap{domaintype(m)}(m) -AbsMap{T}(m) where {T} = AbsMap{T,typeof(m)}(m) -AbsMap{T}(m::Map{T}) where {T} = AbsMap{T,typeof(m)}(m) -AbsMap{T}(m::Map{S}) where {S,T} = AbsMap{T}(convert(Map{T}, m)) - -absmap(m) = AbsMap(m) - -applymap(m::AbsMap, x) = abs(supermap(m)(x)) - - -"A lazy volume element evaluates to `diffvolume(m, x)` on the fly." -struct LazyDiffVolume{T,M} <: SimpleLazyMap{T} - map :: M -end - -LazyDiffVolume(m) = LazyDiffVolume{domaintype(m)}(m) -LazyDiffVolume{T}(m) where {T} = LazyDiffVolume{T,typeof(m)}(m) -LazyDiffVolume{T}(m::Map{T}) where {T} = LazyDiffVolume{T,typeof(m)}(m) -LazyDiffVolume{T}(m::Map{S}) where {S,T} = LazyDiffVolume{T}(convert(Map{T}, m)) - -applymap(m::LazyDiffVolume, x) = diffvolume(supermap(m), x) - -Display.displaystencil(m::LazyDiffVolume) = ["LazyDiffVolume(", supermap(m), ")"] -show(io::IO, mime::MIME"text/plain", m::LazyDiffVolume) = composite_show(io, mime, m) - -""" - diffvolume(m[, x]) - -Compute the differential volume (at a point `x`). If `J` is the Jacobian matrix, -possibly rectangular, then the differential volume is `sqrt(det(J'*J))`. - -If the map is square, then the differential volume is the absolute value of the -Jacobian determinant. -""" -function diffvolume(m) - if issquaremap(m) - if mapsize(m,1) > 1 - absmap(jacdet(m)) - else - jacdet(m) - end - else - LazyDiffVolume(m) - end -end - -function diffvolume(m, x) - if issquaremap(m) - if mapsize(m,1) > 1 - abs(jacdet(m, x)) - else - jacdet(m, x) - end - else - jac = jacobian(m, x) - sqrt(det(adjoint(jac)*jac)) - end -end - - - -const NumberLike = Union{Number,UniformScaling} - -""" - to_matrix(::Type{T}, A[, b]) - -Convert the `A` in the affine map `A*x` or `A*x+b` with domaintype `T` to a matrix. -""" -to_matrix(::Type{T}, A) where {T} = A -to_matrix(::Type{T}, A::AbstractMatrix) where {T} = A -to_matrix(::Type{T}, A::NumberLike) where {T<:Number} = A -to_matrix(::Type{<:StaticVector{N,T}}, A::Number) where {N,T} = A * one(SMatrix{N,N,T}) -to_matrix(::Type{<:StaticVector{N,T}}, A::UniformScaling) where {N,T} = A.λ * one(SMatrix{N,N,T}) -to_matrix(::Type{T}, A::Number) where {T<:AbstractVector} = A * I -to_matrix(::Type{T}, A::UniformScaling) where {T<:Number} = one(T) -to_matrix(::Type{T}, A::UniformScaling) where {T<:AbstractVector} = A - -to_matrix(::Type{T}, A, b) where {T} = A -to_matrix(::Type{T}, A::AbstractMatrix, b) where {T} = A -to_matrix(::Type{T}, A::Number, b::Number) where {T<:Number} = A -to_matrix(::Type{T}, A::UniformScaling{S}, b::Number) where {S,T<:Number} = - convert(promote_type(S,T,typeof(b)), A.λ) -to_matrix(::Type{<:StaticVector{N,T}}, A::NumberLike, b::StaticVector{N,T}) where {N,T} = A * one(SMatrix{N,N,T}) -to_matrix(::Type{T}, A::NumberLike, b::AbstractVector) where {S,T<:AbstractVector{S}} = - A * Array{S,2}(I, length(b), length(b)) - -""" - to_vector(::Type{T}, A[, b]) - -Convert the `b` in the affine map `A*x` or `A*x+b` with domaintype `T` to a vector. -""" -to_vector(::Type{T}, A) where {T} = zero(T) -to_vector(::Type{T}, A::StaticVector{M,S}) where {T,M,S} = zero(SVector{M,S}) -to_vector(::Type{T}, A::StaticMatrix{M,N,S}) where {T,M,N,S} = zero(SVector{M,S}) -to_vector(::Type{T}, A::AbstractArray) where {T} = zeros(eltype(T),size(A,1)) -to_vector(::Type{T}, A, b) where {T} = b - - -"Return a zero matrix of the same size as the map." -zeromatrix(m) = zeromatrix(m, domaintype(m), codomaintype(m)) -zeromatrix(m, ::Type{T}, ::Type{U}) where {T,U} = zeros(numtype(T),mapsize(m)) -zeromatrix(m, ::Type{T}, ::Type{U}) where {T<:Number,U<:Number} = zero(promote_type(T,U)) -zeromatrix(m, ::Type{T}, ::Type{U}) where {N,T<:StaticVector{N},U<:Number} = - transpose(zero(SVector{N,eltype(T)})) -zeromatrix(m, ::Type{T}, ::Type{U}) where {N,M,T<:StaticVector{N},U<:StaticVector{M}} = - zero(SMatrix{M,N,promote_type(eltype(T),eltype(U))}) -zeromatrix(m, ::Type{T}, ::Type{U}) where {T<:Number,M,U<:StaticVector{M}} = - zero(SVector{M,promote_type(T,eltype(U))}) -zeromatrix(m, ::Type{T}, ::Type{U}) where {T<:AbstractVector,U<:Number} = - transpose(zeros(promote_type(eltype(T),U), mapsize(m,1))) - - -"Return a zero vector of the same size as the codomain of the map." -zerovector(m) = zerovector(m, codomaintype(m)) -zerovector(m, ::Type{U}) where {U} = zero(U) -zerovector(m, ::Type{<:StaticVector{M,T}}) where {M,T} = zero(SVector{M,T}) -# If the output type is a vector, the map itself should store the size information. -zerovector(m, ::Type{<:AbstractVector{T}}) where {T} = zeros(T, mapsize(m,1)) - -"Return an identity matrix with the same size as the map." -identitymatrix(m) = identitymatrix(m, codomaintype(m)) -identitymatrix(m, ::Type{T}) where {T} = one(T) -identitymatrix(m, ::Type{<:StaticVector{N,T}}) where {N,T} = one(SMatrix{N,N,T}) -identitymatrix(m, ::Type{<:AbstractVector{T}}) where {T} = Diagonal{T}(ones(mapsize(m,1))) diff --git a/FunctionMaps/src/generic/lazy.jl b/FunctionMaps/src/generic/lazy.jl deleted file mode 100644 index a2870af..0000000 --- a/FunctionMaps/src/generic/lazy.jl +++ /dev/null @@ -1,53 +0,0 @@ - -""" -A lazy map has an action that is defined in terms of other maps. Those maps are -stored internally, and the action of the lazy map is computed on-the-fly and only -when invoked. -""" -abstract type LazyMap{T} <: Map{T} end - - -"A composite lazy map is defined in terms of several other maps." -abstract type CompositeLazyMap{T} <: LazyMap{T} end -"A simple lazy map derives from a single other map." -abstract type SimpleLazyMap{T} <: LazyMap{T} end - -supermap(m::SimpleLazyMap) = m.map -components(m::CompositeLazyMap) = m.maps - -Base.getindex(m::CompositeLazyMap, I::ComponentIndex...) = component(m, map(Indexing.to_index, I)...) - -isreal(m::SimpleLazyMap) = isreal(supermap(m)) -isreal(m::CompositeLazyMap) = all(map(isreal, components(m))) - -"A `DerivedMap` inherits all of its properties from another map, but has its own type." -abstract type DerivedMap{T} <: SimpleLazyMap{T} end - -applymap(m::DerivedMap, x) = supermap(m)(x) -appymap!(y, m::DerivedMap, x) = applymap!(y, supermap(m), x) - -jacobian(m::DerivedMap) = jacobian(supermap(m)) - - -"A `WrappedMap{T}` takes any object and turns it into a `Map{T}`." -struct WrappedMap{T,M} <: DerivedMap{T} - map :: M -end -WrappedMap{T}(map) where {T} = WrappedMap{T,typeof(map)}(map) -WrappedMap(map) = WrappedMap{Float64}(map) - -similarmap(m::WrappedMap, ::Type{T}) where {T} = WrappedMap{T}(m) - -convert(::Type{Map}, m::Map) = m -convert(::Type{Map}, m) = WrappedMap(m) -convert(::Type{Map{T}}, m) where {T} = WrappedMap{T}(m) -convert(::Type{Map{T}}, m::WrappedMap) where {T} = WrappedMap{T}(m.map) - -isequalmap(m1::WrappedMap, m2::Function) = m1.map == m2 -isequalmap(m1::Function, m2::WrappedMap) = m1 == m2.map - -map_hash(m::WrappedMap, h::UInt) = hashrec("WrappedMap", supermap(m), h) - -Display.displaystencil(m::WrappedMap{T}) where {T} = - ["WrappedMap{$T}(", supermap(m), ")"] -show(io::IO, mime::MIME"text/plain", m::WrappedMap) = composite_show(io, mime, m) diff --git a/FunctionMaps/src/generic/map.jl b/FunctionMaps/src/generic/map.jl deleted file mode 100644 index 82a782a..0000000 --- a/FunctionMaps/src/generic/map.jl +++ /dev/null @@ -1,176 +0,0 @@ - -"An `AbstractMap` represents a function `y=f(x)` of a single variable." -abstract type AbstractMap end - -"A `Map{T}` is a map of a single variable of type `T`." -abstract type Map{T} <: AbstractMap end - -Map(m) = convert(Map, m) -Map{T}(m) where {T} = convert(Map{T}, m) - -"A `TypedMap{T,U}` maps a variable of type `T` to a variable of type `U`." -abstract type TypedMap{T,U} <: Map{T} end - -const EuclideanMap{N,T} = Map{<:StaticVector{N,T}} -const VectorMap{T} = Map{Vector{T}} - -CompositeTypes.Display.displaysymbol(m::Map) = 'F' - -""" - domaintype(m) - -What is the expected type of a point in the domain of the function map `m`? -""" -domaintype(m) = domaintype(typeof(m)) -domaintype(::Type{M}) where {M} = Any -domaintype(::Type{<:Map{T}}) where {T} = T - -""" - codomaintype(m[, T]) - -What is the codomain type of the function map `m`, given that `T` is its domain type? -""" -codomaintype(m) = codomaintype(m, domaintype(m)) -codomaintype(m, ::Type{T}) where {T} = _codomaintype(typeof(m), T) -_codomaintype(::Type{M}, ::Type{T}) where {M,T} = Base.promote_op(applymap, M, T) -_codomaintype(::Type{M}, ::Type{Any}) where {M} = Any -_codomaintype(M::Type{<:TypedMap{T,U}}, ::Type{T}) where {T,U} = U -_codomaintype(M::Type{<:TypedMap{T,U}}, ::Type{Any}) where {T,U} = U - -prectype(::Type{<:Map{T}}) where T = prectype(T) -numtype(::Type{<:Map{T}}) where T = numtype(T) - -isreal(m::AbstractMap) = isrealtype(domaintype(m)) && isrealtype(codomaintype(m)) -isreal(::UniformScaling{T}) where {T} = isrealtype(T) -isreal(::Type{UniformScaling{T}}) where {T} = isrealtype(T) - -convert(::Type{AbstractMap}, m::AbstractMap) = m -convert(::Type{Map{T}}, m::Map{T}) where {T} = m -convert(::Type{Map{T}}, m::Map{S}) where {S,T} = similarmap(m, T) -convert(::Type{TypedMap{T,U}}, m::TypedMap{T,U}) where {T,U} = m -convert(::Type{TypedMap{T,U}}, m::TypedMap) where {T,U} = similarmap(m, T, U) - -convert_domaintype(::Type{T}, map::Map{T}) where {T} = map -convert_domaintype(::Type{U}, map::Map{T}) where {T,U} = convert(Map{U}, map) -convert_domaintype(::Type{Any}, map) = map -convert_domaintype(::Type{Any}, map::Map{T}) where T = map -convert_domaintype(::Type{T}, map) where T = _convert_domaintype(T, map, domaintype(map)) -_convert_domaintype(::Type{T}, map, ::Type{T}) where T = map -_convert_domaintype(::Type{U}, map, ::Type{T}) where {T,U} = Map{U}(map) - -convert_numtype(::Type{U}, map::Map{T}) where {T,U} = convert(Map{to_numtype(U,T)}, map) -convert_prectype(::Type{U}, map::Map{T}) where {T,U} = convert(Map{to_prectype(U,T)}, map) - -convert_codomaintype(::Type{U}, map) where U = - _convert_codomaintype(U, map, domaintype(map), codomaintype(map)) -# types match: we can return map -_convert_codomaintype(::Type{U}, map, ::Type{T}, ::Type{U}) where {T,U} = map -# ...or, outputs are vectors with the same eltype: that's good enough, further -# conversion might be unnecessary and expensive -_convert_codomaintype(::Type{<:AbstractVector{S}}, map, ::Type{T}, ::Type{<:AbstractVector{S}}) where {T,S} = map -# no match yet: we now try to convert the numtype -_convert_codomaintype(::Type{U}, map, ::Type{T}, ::Type{V}) where {T,U,V} = - _convert_codomaintype2(U, convert_numtype(numtype(U), map)) -_convert_codomaintype2(::Type{U}, map::Map) where U = - _convert_codomaintype2(U, map, domaintype(map), codomaintype(map)) -# better match this time around: we can return map -_convert_codomaintype2(::Type{U}, map, ::Type{T}, ::Type{U}) where {T,U} = map -_convert_codomaintype2(::Type{<:AbstractVector{S}}, map, ::Type{T}, ::Type{<:AbstractVector{S}}) where {T,S} = map -# still no match: we can't do this automatically -_convert_codomaintype2(::Type{U}, map, ::Type{T}, ::Type{V}) where {T,U,V} = - throw(ArgumentError("Don't know how to convert the codomain type of $(map) to $(U).")) - -convert_codomaintype(::Type{Any}, map) = map - - -# Users may call a map, concrete subtypes specialize the `applymap` function -(m::AbstractMap)(x) = applymap(m, x) - -# For Map{T}, we allow invocation with multiple arguments by conversion to T -(m::Map{T})(x) where {T} = promote_and_apply(m, x) -(m::Map{T})(x...) where {T} = promote_and_apply(m, convert(T, x)) - -"Promote map and point to compatible types." -promote_map_point_pair(m, x) = (m, x) -promote_map_point_pair(m::Map, x) = _promote_map_point_pair(m, x, promote_type(domaintype(m), typeof(x))) -# This is the line where we promote both the map and the point: -_promote_map_point_pair(m, x, ::Type{T}) where {T} = - convert(Map{T}, m), convert(T, x) -_promote_map_point_pair(m, x, ::Type{Any}) = m, x - -# Some exceptions: -# - types match, do nothing -promote_map_point_pair(m::Map{T}, x::T) where {T} = m, x -# - an Any map, do nothing -promote_map_point_pair(m::Map{Any}, x) = m, x -# - tuples: these are typically composite maps, promotion may happen later -promote_map_point_pair(m::Map{<:Tuple}, x::Tuple) = m, x -# - abstract vectors: promotion may be expensive -promote_map_point_pair(m::Map{<:AbstractVector}, x::AbstractVector) = m, x -# - SVector: promotion is likely cheap -promote_map_point_pair(m::EuclideanMap{N,T}, x::AbstractVector{S}) where {N,S,T} = - _promote_map_point_pair(m, x, SVector{N,promote_type(S,T)}) - -# For maps of type Map{T}, we call promote_map_point_pair and then applymap -promote_and_apply(m::Map, x) = applymap(promote_map_point_pair(m,x)...) - -applymap!(y, m, x) = y .= m(x) - -# Fallback for functions that are not of type AbstractMap -applymap(m, x) = m(x) -applymap(m::AbstractMap, x) = error("Please implement applymap for map $(m)") - - -"Can the maps be promoted to a common domain type without throwing an error?" -promotable_maps(maps...) = promotable_eltypes(map(domaintype, maps)...) - -"Promote the given maps to have a common domain type." -promote_maps() = () -promote_maps(m1) = m1 -function promote_maps(m1, m2) - T = promote_type(domaintype(m1), domaintype(m2)) - convert_domaintype(T, m1), convert_domaintype(T, m2) -end -function promote_maps(m1, m2, m3, maps...) - T = mapreduce(domaintype, promote_type, (m1,m2,m3,maps...)) - convert_domaintype.(T, (m1,m2,m3,maps...)) -end - -isvectorvalued_type(::Type{T}) where {T<:Number} = true -isvectorvalued_type(::Type{T}) where {T<:AbstractVector} = true -isvectorvalued_type(::Type{T}) where {T} = false - - -"Is the map a vector-valued function, i.e., a function from Rn to Rm?" -isvectorvalued(m) = - isvectorvalued_type(domaintype(m)) && isvectorvalued_type(codomaintype(m)) - -# mapsize should be defined for vector valued maps -# The size of a map equals the size of its jacobian -# The jacobian can be a number, a vector, an adjoint vector, or a matrix -mapsize(m, i) = _mapsize(m, i, mapsize(m)) -_mapsize(m, i, size::Tuple{Int,Int}) = i <= 2 ? size[i] : 1 -_mapsize(m, i, size::Tuple{Int}) = i <= 1 ? size[i] : 1 -_mapsize(m, i, size::Tuple{}) = 1 - -"Is the given map a square map?" -issquaremap(m) = isvectorvalued(m) && (mapsize(m,1) == mapsize(m,2)) - -isoverdetermined(m) = mapsize(m,1) >= mapsize(m,2) -isunderdetermined(m) = mapsize(m,1) <= mapsize(m,2) - -is_scalar_to_vector(m) = mapsize(m) isa Tuple{Int} -is_vector_to_scalar(m) = mapsize(m) isa Tuple{Int,Int} && codomaintype(m)<:Number -is_scalar_to_scalar(m) = mapsize(m) == () -is_vector_to_vector(m) = mapsize(m) isa Tuple{Int,Int} && !is_vector_to_scalar(m) - -# Display routines -map_stencil(m, x) = [Display.SymbolObject(m), '(', x, ')'] -map_stencil_broadcast(m, x) = [Display.SymbolObject(m), ".(", x, ')'] -map_stencil_broadcast(m::Function, x) = [repr(m), ".(", x, ')'] - -Display.object_parentheses(m::Map) = map_object_parentheses(m) -Display.stencil_parentheses(m::Map) = map_stencil_parentheses(m) - -map_object_parentheses(m) = false -map_stencil_parentheses(m) = false diff --git a/FunctionMaps/src/generic/product.jl b/FunctionMaps/src/generic/product.jl deleted file mode 100644 index 8cc6e03..0000000 --- a/FunctionMaps/src/generic/product.jl +++ /dev/null @@ -1,201 +0,0 @@ - -""" -A product map is diagonal and acts on each of the components of x separately: -`y = f(x)` becomes `y_i = f_i(x_i)`. -""" -abstract type ProductMap{T} <: CompositeLazyMap{T} end - -components(m::ProductMap) = m.maps -factors(d::ProductMap) = components(d) - -VcatMapElement = Union{Map{<:StaticVector},Map{<:Number}} - -ProductMap(maps::Tuple) = ProductMap(maps...) -ProductMap(maps::StaticVector) = ProductMap(maps...) -ProductMap(maps...) = TupleProductMap(maps...) -ProductMap(maps::VcatMapElement...) = VcatMap(maps...) -ProductMap(maps::AbstractVector) = VectorProductMap(maps) - -ProductMap{T}(maps...) where {T} = _TypedProductMap(T, maps...) -_TypedProductMap(::Type{T}, maps...) where {T<:Tuple} = TupleProductMap(maps...) -_TypedProductMap(::Type{SVector{N,T}}, maps...) where {N,T} = VcatMap{T}(maps...) -_TypedProductMap(::Type{T}, maps...) where {T<:AbstractVector} = VectorProductMap{T}(maps...) - -compatibleproductdims(d1::ProductMap, d2::ProductMap) = - mapsize(d1) == mapsize(d2) && - all(map(==, map(mapsize, components(d1)), map(mapsize, components(d2)))) - -isconstantmap(m::ProductMap) = mapreduce(isconstantmap, &, components(m)) -islinearmap(m::ProductMap) = mapreduce(islinearmap, &, components(m)) -isaffinemap(m::ProductMap) = mapreduce(isaffinemap, &, components(m)) - -affinematrix(m::ProductMap) = toexternalmatrix(m, map(affinematrix, components(m))) -affinevector(m::ProductMap) = toexternalpoint(m, map(affinevector, components(m))) -mapconstant(m::ProductMap) = toexternalpoint(m, map(mapconstant, components(m))) - -jacobian(m::ProductMap, x) = - toexternalmatrix(m, map(jacobian, components(m), tointernalpoint(m, x))) -function jacobian(m::ProductMap{T}) where {T} - if isaffinemap(m) - ConstantMap{T}(affinematrix(m)) - else - LazyJacobian(m) - end -end - -# diffvolume(m::ProductMap) = multiply_map(map(diffvolume, factors(m))...) - -similarmap(m::ProductMap, ::Type{T}) where {T} = ProductMap{T}(components(m)) - -tointernalpoint(m::ProductMap, x) = x -toexternalpoint(m::ProductMap, y) = y - -applymap(m::ProductMap, x) = - toexternalpoint(m, map(applymap, components(m), tointernalpoint(m, x))) - -productmap(m1, m2) = productmap1(m1, m2) -productmap1(m1, m2) = productmap2(m1, m2) -productmap2(m1, m2) = default_productmap(m1, m2) -default_productmap(m1, m2) = ProductMap(m1, m2) - -productmap(m1::ProductMap, m2::ProductMap) = - ProductMap(components(m1)..., components(m2)...) -productmap1(m1::ProductMap, m2) = ProductMap(components(m1)..., m2) -productmap2(m1, m2::ProductMap) = ProductMap(m1, components(m2)...) - -for op in (:inverse, :leftinverse, :rightinverse) - @eval $op(m::ProductMap) = ProductMap(map($op, components(m))) - @eval $op(m::ProductMap, x) = toexternalpoint(m, map($op, components(m), tointernalpoint(m, x))) -end - -function composedmap(m1::ProductMap, m2::ProductMap) - if compatibleproductdims(m1, m2) - ProductMap(map(composedmap, components(m1), components(m2))) - else - ComposedMap(m1,m2) - end -end - -mapsize(m::ProductMap) = (sum(t->mapsize(t,1), components(m)), sum(t->mapsize(t,2), components(m))) - -isequalmap(m1::ProductMap, m2::ProductMap) = all(map(isequalmap, components(m1), components(m2))) -map_hash(m::ProductMap, h::UInt) = hashrec("ProductMap", collect(components(m)), h) - -canonicalmap(m::ProductMap) = any(map(hascanonicalmap, factors(m))) ? - ProductMap(map(canonicalmap, factors(m))) : m -canonicalmap(::Equal, m::ProductMap) = any(map(hasequalmap, factors(m))) ? - ProductMap(map(equalmap, factors(m))) : m -canonicalmap(::Equivalent, m::ProductMap) = any(map(hasequivalentmap, factors(m))) ? - ProductMap(map(equivalentmap, factors(m))) : m - -Display.combinationsymbol(m::ProductMap) = Display.Symbol('⊗') -Display.displaystencil(m::ProductMap) = composite_displaystencil(m) -map_object_parentheses(m::ProductMap) = true -map_stencil_parentheses(m::ProductMap) = true -show(io::IO, mime::MIME"text/plain", m::ProductMap) = composite_show(io, mime, m) -show(io::IO, m::ProductMap) = composite_show_compact(io, m) - -""" -A `VcatMap` is a product map with domain and codomain vectors -concatenated (`vcat`) into a single vector. -""" -struct VcatMap{T,M,N,DIM1,DIM2,MAPS} <: ProductMap{SVector{N,T}} - maps :: MAPS -end - -VcatMap(maps::Union{Tuple,Vector}) = VcatMap(maps...) -VcatMap(maps...) = VcatMap{numtype(maps...)}(maps...) - -VcatMap{T}(maps::Union{Tuple,Vector}) where T = VcatMap{T}(maps...) -function VcatMap{T}(maps...) where T - M = sum(t->mapsize(t,1), maps) - N = sum(t->mapsize(t,2), maps) - VcatMap{T,M,N}(maps...) -end - -mapdim(map) = mapsize(map,2) - -VcatMap{T,M,N}(maps::Union{Tuple,Vector}) where {T,M,N} = VcatMap{T,M,N}(maps...) -function VcatMap{T,M,N}(maps...) where {T,M,N} - DIM1 = map(t->mapsize(t,1), maps) - DIM2 = map(t->mapsize(t,2), maps) - VcatMap{T,M,N,DIM1,DIM2}(convert_numtype.(Ref(T), maps)...) -end - -VcatMap{T,M,N,DIM1,DIM2}(maps...) where {T,M,N,DIM1,DIM2} = - VcatMap{T,M,N,DIM1,DIM2,typeof(maps)}(maps) - -mapsize(m::VcatMap{T,M,N}) where {T,M,N} = (M,N) - -tointernalpoint(m::VcatMap{T,M,N,DIM1,DIM2}, x) where {T,M,N,DIM1,DIM2} = - convert_fromcartesian(x, Val{DIM2}()) -toexternalpoint(m::VcatMap{T,M,N,DIM1,DIM2}, y) where {T,M,N,DIM1,DIM2} = - convert_tocartesian(y, Val{DIM1}()) - -size_as_matrix(A::AbstractArray) = size(A) -size_as_matrix(A::Number) = (1,1) - -# The Jacobian is block-diagonal -function toexternalmatrix(m::VcatMap{T,M,N}, matrices) where {T,M,N} - A = zeros(T, M, N) - k = 0 - l = 0 - for el in matrices - m,n = size_as_matrix(el) - A[k+1:k+m,l+1:l+n] .= el - k += m - l += n - end - SMatrix{M,N}(A) -end - - -""" -A `VectorProductMap` is a product map where all components are univariate maps, -with inputs and outputs collected into a `Vector`. -""" -struct VectorProductMap{T<:AbstractVector,M} <: ProductMap{T} - maps :: Vector{M} -end - -VectorProductMap(maps...) = VectorProductMap(maps) -VectorProductMap(maps) = VectorProductMap(collect(maps)) -function VectorProductMap(maps::Vector) - T = mapreduce(numtype, promote_type, maps) - VectorProductMap{Vector{T}}(maps) -end - -VectorProductMap{T}(maps...) where {T} = VectorProductMap{T}(maps) -VectorProductMap{T}(maps) where {T} = VectorProductMap{T}(collect(maps)) -function VectorProductMap{T}(maps::Vector) where {T} - Tmaps = convert.(Map{eltype(T)}, maps) - VectorProductMap{T,eltype(Tmaps)}(Tmaps) -end - -# the Jacobian is a diagonal matrix -toexternalmatrix(m::VectorProductMap, matrices) = Diagonal(matrices) - -mapsize(m::VectorProductMap) = (length(m.maps), length(m.maps)) - -""" -A `TupleProductMap` is a product map with all components collected in a tuple. -There is no vector-valued function associated with this map. -""" -struct TupleProductMap{T<:Tuple,MM} <: ProductMap{T} - maps :: MM -end - -TupleProductMap(maps::Vector) = TupleProductMap(maps...) -TupleProductMap(maps...) = TupleProductMap(maps) -function TupleProductMap(maps::Tuple) - T = Tuple{map(domaintype, maps)...} - TupleProductMap{T}(maps) -end - -TupleProductMap{T}(maps::Vector) where {T<:Tuple} = TupleProductMap{T}(maps...) -TupleProductMap{T}(maps...) where {T<:Tuple} = TupleProductMap{T}(maps) -function TupleProductMap{T}(maps::NTuple{N,<:AbstractMap}) where {N,T <: Tuple} - Tmaps = map((t,d) -> convert(Map{t},d), tuple(T.parameters...), maps) - TupleProductMap{T,typeof(Tmaps)}(Tmaps) -end -TupleProductMap{T}(maps) where {T <: Tuple} = TupleProductMap{T,typeof(maps)}(maps) diff --git a/FunctionMaps/src/types/affine.jl b/FunctionMaps/src/types/affine.jl deleted file mode 100644 index d19b7f6..0000000 --- a/FunctionMaps/src/types/affine.jl +++ /dev/null @@ -1,586 +0,0 @@ - -""" -An affine map has the general form `y = A*x + b`. - -We use `affinematrix(m)` and `affinevector(m)` to denote `A` and `b` respectively. Concrete -subtypes include linear maps of the form `y = A*x` and translations of the -form `y = x + b`. -""" -abstract type AbstractAffineMap{T} <: Map{T} end - -unsafe_matrix(m::AbstractAffineMap) = m.A -unsafe_vector(m::AbstractAffineMap) = m.b - -"Return the matrix `A` in the affine map `Ax+b`." -affinematrix(m::AbstractAffineMap) = to_matrix(domaintype(m), unsafe_matrix(m), unsafe_vector(m)) - -"Return the vector `b` in the affine map `Ax+b`." -affinevector(m::AbstractAffineMap) = to_vector(domaintype(m), unsafe_matrix(m), unsafe_vector(m)) - -applymap(m::AbstractAffineMap, x) = _affine_applymap(m, x, unsafe_matrix(m), unsafe_vector(m)) -_affine_applymap(m, x, A, b) = A*x + b - -applymap!(y, m::AbstractAffineMap, x) = _affine_applymap!(y, m, x, unsafe_matrix(m), unsafe_vector(m)) -function _affine_applymap!(y, m, x, A, b) - mul!(y, A, x) - y .+= b - y -end - -isreal(m::AbstractAffineMap) = _affine_isreal(m, unsafe_matrix(m), unsafe_vector(m)) -_affine_isreal(m, A, b) = isreal(A) && isreal(b) - -jacobian(m::AbstractAffineMap{T}) where {T} = ConstantMap{T}(affinematrix(m)) -jacobian(m::AbstractAffineMap, x) = affinematrix(m) - -jacdet(m::AbstractAffineMap, x) = _affine_jacdet(m, x, unsafe_matrix(m)) -_affine_jacdet(m, x, A) = det(A) -_affine_jacdet(m, x::Number, A::UniformScaling) = A.λ -_affine_jacdet(m, x::AbstractVector, A::Number) = A^length(x) -_affine_jacdet(m, x::AbstractVector, A::UniformScaling) = A.λ^length(x) - -function diffvolume(m::AbstractAffineMap{T}) where T - J = jacobian(m) - c = sqrt(det(affinevector(J)'*affinevector(J))) - ConstantMap{T}(c) -end - -islinearmap(m::AbstractMap) = false -islinearmap(m::AbstractAffineMap) = _affine_islinearmap(m, unsafe_vector(m)) -_affine_islinearmap(m, b) = all(b .== 0) - -isaffinemap(m) = false -isaffinemap(m::Map) = islinearmap(m) || isconstantmap(m) -isaffinemap(m::AbstractAffineMap) = true - -isequalmap(m1::AbstractAffineMap, m2::AbstractAffineMap) = - affinematrix(m1) == affinematrix(m2) && affinevector(m1) == affinevector(m2) - -isequalmap(m1::AbstractAffineMap, m2::IdentityMap) = - islinearmap(m1) && affinematrix(m1) == affinematrix(m2) -isequalmap(m1::IdentityMap, m2::AbstractAffineMap) = isequalmap(m2, m1) - -map_hash(m::AbstractAffineMap, h::UInt) = hashrec("AbstractAffineMap", affinematrix(m), affinevector(m), h) - -mapsize(m::AbstractAffineMap) = _affine_mapsize(m, domaintype(m), unsafe_matrix(m), unsafe_vector(m)) -_affine_mapsize(m, T, A::AbstractArray, b) = size(A) -_affine_mapsize(m, T, A::AbstractVector, b::AbstractVector) = (length(A),) -_affine_mapsize(m, T, A::Number, b::Number) = () -_affine_mapsize(m, T, A::Number, b::AbstractVector) = (length(b),length(b)) -_affine_mapsize(m, T, A::UniformScaling, b::Number) = () -_affine_mapsize(m, T, A::UniformScaling, b) = (length(b),length(b)) - - -Display.displaystencil(m::AbstractAffineMap) = vcat(["x -> "], map_stencil(m, 'x')) -show(io::IO, mime::MIME"text/plain", m::AbstractAffineMap) = composite_show(io, mime, m) - -map_stencil(m::AbstractAffineMap, x) = _affine_map_stencil(m, x, unsafe_matrix(m), unsafe_vector(m)) -_affine_map_stencil(m, x, A, b) = [A, " * ", x, " + ", b] -_affine_map_stencil(m, x, A, b::Real) = - b >= 0 ? [A, " * ", x, " + ", b] : [A, " * ", x, " - ", abs(b)] - -map_stencil_broadcast(m::AbstractAffineMap, x) = _affine_map_stencil_broadcast(m, x, unsafe_matrix(m), unsafe_vector(m)) -_affine_map_stencil_broadcast(m, x, A, b) = [A, " .* ", x, " .+ ", b] -_affine_map_stencil_broadcast(m, x, A::Number, b) = [A, " * ", x, " .+ ", b] - -map_object_parentheses(m::AbstractAffineMap) = true -map_stencil_parentheses(m::AbstractAffineMap) = true - -######################## -# Linear maps: y = A*x -######################## - -""" -The supertype of all linear maps `y = A*x`. -Concrete subtypes may differ in how `A` is represented. -""" -abstract type LinearMap{T} <: AbstractAffineMap{T} end - -mapsize(m::LinearMap) = _linearmap_size(m, domaintype(m), unsafe_matrix(m)) -_linearmap_size(m, T, A) = size(A) -_linearmap_size(m, ::Type{T}, A::Number) where {T<:Number} = () -_linearmap_size(m, ::Type{T}, A::AbstractVector) where {T<:Number} = (length(A),) -_linearmap_size(m, ::Type{T}, A::Number) where {N,T<:StaticVector{N}} = (N,N) - -affinematrix(m::LinearMap) = to_matrix(domaintype(m), unsafe_matrix(m)) -affinevector(m::LinearMap) = to_vector(domaintype(m), unsafe_matrix(m)) - -LinearMap(A::Number) = ScalarLinearMap(A) -LinearMap(A::SMatrix) = StaticLinearMap(A) -LinearMap(A::MMatrix) = StaticLinearMap(A) -LinearMap(A::Matrix) = VectorLinearMap(A) -LinearMap(A) = GenericLinearMap(A) - -LinearMap{T}(A::Number) where {T<:Number} = _LinearMap(A, T, promote_type(T,typeof(A))) -_LinearMap(A::Number, ::Type{T}, ::Type{T}) where {T<:Number} = ScalarLinearMap{T}(A) -_LinearMap(A::Number, ::Type{T}, ::Type{S}) where {T<:Number,S} = GenericLinearMap{T}(A) -LinearMap{T}(A::SMatrix{M,N}) where {M,N,S,T <: SVector{N,S}} = StaticLinearMap{S}(A) -LinearMap{T}(A::MMatrix{M,N}) where {M,N,S,T <: SVector{N,S}} = StaticLinearMap{S}(A) -LinearMap{T}(A::Matrix) where {S,T <: Vector{S}} = VectorLinearMap{S}(A) -LinearMap{T}(A) where {T} = GenericLinearMap{T}(A) - -# convenience functions -LinearMap(a::Number...) = LinearMap(promote(a...)) -LinearMap(a::NTuple{N,T}) where {N,T <: Number} = LinearMap{SVector{N,T}}(Diagonal(SVector{N,T}(a))) - -applymap(m::LinearMap, x) = _linear_applymap(m, x, unsafe_matrix(m)) -_linear_applymap(m, x, A) = A*x - -applymap!(y, m::LinearMap, x) = _linear_applymap!(y, m, x, unsafe_matrix(m)) -_linear_applymap!(y, m, x, A) = mul!(y, A, x) - -islinearmap(m::LinearMap) = true - -isreal(m::LinearMap) = _linear_isreal(m, unsafe_matrix(m)) -_linear_isreal(m, A) = isreal(A) - -isequalmap(m1::LinearMap, m2::LinearMap) = affinematrix(m1) == affinematrix(m2) - -# inverse should be called only on square maps, otherwise use -# leftinverse or rightinverse in order to use pinv instead of inv -inverse(m::LinearMap) = (@assert issquaremap(m); LinearMap(inv(m.A))) -inverse(m::LinearMap, x) = (@assert issquaremap(m); m.A \ x) - -function leftinverse(m::LinearMap) - @assert isoverdetermined(m) - LinearMap(matrix_pinv(m.A)) -end -function rightinverse(m::LinearMap) - @assert isunderdetermined(m) - LinearMap(matrix_pinv(m.A)) -end -function leftinverse(m::LinearMap, x) - @assert isoverdetermined(m) - m.A \ x -end -function rightinverse(m::LinearMap, x) - @assert isunderdetermined(m) - m.A \ x -end - -similarmap(m::LinearMap, ::Type{T}) where {T} = LinearMap{T}(m.A) - -convert(::Type{Map}, a::Number) = LinearMap(a) -convert(::Type{Map{T}}, a::Number) where {T} = LinearMap{T}(a) - -convert(::Type{LinearMap}, ::StaticIdentityMap{T}) where {T} = LinearMap{T}(1) -convert(::Type{LinearMap{T}}, ::StaticIdentityMap) where {T} = LinearMap{T}(1) -convert(::Type{AbstractAffineMap}, m::StaticIdentityMap) = convert(LinearMap, m) -convert(::Type{AbstractAffineMap{T}}, m::StaticIdentityMap) where {T} = convert(LinearMap, m) - - -map_stencil(m::LinearMap, x) = [unsafe_matrix(m), " * ", x] -map_stencil_broadcast(m::LinearMap, x) = _linear_map_stencil_broadcast(m, x, unsafe_matrix(m)) -_linear_map_stencil_broadcast(m, x, A) = [A, " .* ", x] -_linear_map_stencil_broadcast(m, x, A::Number) = [A, " * ", x] - - -"A `GenericLinearMap` is a linear map `y = A*x` for any type of `A`." -struct GenericLinearMap{T,AA} <: LinearMap{T} - A :: AA -end - -""" -What is the suggested domaintype for a generic linear map `A*x` with -the given argument 'A'? -""" -glm_domaintype(A) = Any -glm_domaintype(A::Number) = typeof(A) -glm_domaintype(A::AbstractMatrix{T}) where T = Vector{T} -glm_domaintype(A::StaticMatrix{M,N,T}) where {M,N,T} = SVector{N,T} -glm_domaintype(A::AbstractVector{T}) where {T} = T -glm_domaintype(A::Diagonal{T,<:StaticVector{N,T}}) where {N,T} = SVector{N,T} - -GenericLinearMap(A) = GenericLinearMap{glm_domaintype(A)}(A) - -# Allow any A -GenericLinearMap{T}(A) where {T} = GenericLinearMap{T,typeof(A)}(A) - -# Promote some eltypes if applicable -GenericLinearMap{T}(A::Number) where {T <: Number} = - _GenericLinearMap(A, T, promote_type(T,typeof(A))) -_GenericLinearMap(A::Number, ::Type{T}, ::Type{T}) where {T <: Number} = - GenericLinearMap{T,T}(A) -_GenericLinearMap(A::Number, ::Type{T}, ::Type{S}) where {S,T<:Number} = - GenericLinearMap{T,typeof(A)}(A) - -GenericLinearMap{T}(A::Number) where {T <: AbstractVector} = - _GenericLinearMap(A, T, promote_type(eltype(T),typeof(A))) -_GenericLinearMap(A::Number, ::Type{T}, ::Type{S}) where {S,T<:AbstractVector{S}} = - GenericLinearMap{T,S}(A) -_GenericLinearMap(A::Number, ::Type{T}, ::Type{S}) where {S,T<:AbstractVector} = - GenericLinearMap{T,typeof(A)}(A) - -GenericLinearMap{T}(A::AbstractMatrix{S}) where {S,T <: AbstractVector{S}} = - GenericLinearMap{T,typeof(A)}(A) -GenericLinearMap{T}(A::AbstractMatrix{U}) where {S,T <: AbstractVector{S},U} = - GenericLinearMap{T}(convert(AbstractMatrix{S}, A)) -GenericLinearMap{T}(A::AbstractVector{T}) where {T} = - GenericLinearMap{T,typeof(A)}(A) -GenericLinearMap{T}(A::AbstractVector{S}) where {S,T} = - GenericLinearMap{T}(convert(AbstractVector{T}, A)) -GenericLinearMap{T}(A::AbstractMatrix) where {T<:Number} = - throw(ArgumentError("Linear map with matrix A can not have scalar domaintype.")) - -# Preserve the action on vectors with a number type -inverse(m::GenericLinearMap{T,AA}) where {T<:AbstractVector,AA<:Number} = - LinearMap{T}(inv(m.A)) -leftinverse(m::GenericLinearMap{T,AA}) where {T<:AbstractVector,AA<:Number} = - LinearMap{T}(inv(m.A)) -rightinverse(m::GenericLinearMap{T,AA}) where {T<:AbstractVector,AA<:Number} = - LinearMap{T}(inv(m.A)) - -convert(::Type{Map}, A::UniformScaling) = GenericLinearMap{Vector{Any}}(A) -convert(::Type{Map{T}}, A::UniformScaling) where {T} = GenericLinearMap{T}(A) - -"A `ScalarLinearMap` is a linear map `y = A*x` for scalars." -struct ScalarLinearMap{T} <: LinearMap{T} - A :: T -end - -ScalarLinearMap(A::Number) = ScalarLinearMap{typeof(A)}(A) - -isreal(m::ScalarLinearMap{T}) where {T} = isrealtype(T) - -show(io::IO, m::ScalarLinearMap) = show_scalar_linear_map(io, m.A) -show_scalar_linear_map(io, A::Real) = print(io, "x -> $(A) * x") -show_scalar_linear_map(io, A::Complex) = print(io, "x -> ($(A)) * x") - - -"A `VectorLinearMap` is a linear map `y = A*x` using vectors and matrices." -struct VectorLinearMap{T} <: LinearMap{Vector{T}} - A :: Matrix{T} -end - -VectorLinearMap(A::AbstractMatrix{T}) where {T} = - VectorLinearMap{T}(A) - - -"A `StaticLinearMap` is a linear map `y = A*x` using static arrays." -struct StaticLinearMap{T,N,M,L} <: LinearMap{SVector{N,T}} - A :: SMatrix{M,N,T,L} -end - -StaticLinearMap(A::AbstractMatrix{T}) where {T} = - StaticLinearMap{T}(A) - -StaticLinearMap{T}(A::StaticMatrix{M,N,S}) where {M,N,T,S} = - StaticLinearMap{T}(convert(AbstractMatrix{T}, A)) -StaticLinearMap{T}(A::StaticMatrix{M,N,T}) where {M,N,T} = - StaticLinearMap{T,M,N}(A) -StaticLinearMap{T,N,M}(A::AbstractMatrix) where {T,N,M} = - StaticLinearMap{T,N,M,M*N}(A) - -convert(::Type{Map{SVector{N,T}}}, m::VectorLinearMap) where {N,T} = StaticLinearMap{T,N,N}(m.A) - - -# Implement the interface for abstract arrays, -# representing the linear map x->A*x -MapStyle(A::AbstractArray) = IsMap() - -Map(A::AbstractArray) = LinearMap(A) -Map{T}(A::AbstractArray) where T = LinearMap{T}(A) - -domaintype(A::AbstractArray) = domaintype(Map(A)) - -applymap(A::AbstractArray, x) = A*x -mapsize(A::AbstractArray) = size(A) - -islinearmap(A::AbstractArray) = true -isaffinemap(A::AbstractArray) = true -affinematrix(A::AbstractArray) = A -affinevector(A::AbstractArray) = zerovector(A) - -inverse(A::AbstractMatrix) = inv(A) -inverse(A::AbstractMatrix, x) = A \ x - -jacobian(A::AbstractMatrix) = ConstantMap{glm_domaintype(A)}(A) -jacobian(A::AbstractMatrix, x) = A - -canonicalmap(A::AbstractArray) = Map(A) -canonicalmap(::Equal, A::AbstractArray) = Map(A) - - - - -########################## -# Translations: y = x + b -########################## - - -"A `Translation` represents the map `y = x + b`." -abstract type Translation{T} <: AbstractAffineMap{T} end - -# unsafe_matrix(m::Translation) = I -affinematrix(m::Translation) = to_matrix(domaintype(m), LinearAlgebra.I, unsafe_vector(m)) -affinevector(m::Translation) = unsafe_vector(m) - -mapsize(m::Translation) = _translation_mapsize(m, domaintype(m), unsafe_vector(m)) -_translation_mapsize(m, ::Type{T}, b::Number) where {T} = () -_translation_mapsize(m, ::Type{T}, b) where {T} = (length(b),length(b)) - -map_stencil(m::Translation, x) = _translation_map_stencil(m, x, unsafe_vector(m)) -_translation_map_stencil(m, x, b) = [x, " + ", b] -_translation_map_stencil(m, x, b::Real) = - b >= 0 ? [x, " + ", b] : [x, " - ", abs(b)] -map_stencil_broadcast(m::Translation, x) = _translation_map_stencil_broadcast(m, x, unsafe_vector(m)) -_translation_map_stencil_broadcast(m, x, b) = [x, " .+ ", b] - -"Translation by a scalar value." -struct ScalarTranslation{T} <: Translation{T} - b :: T -end - -isreal(m::ScalarTranslation{T}) where {T} = isrealtype(T) - -show(io::IO, m::ScalarTranslation) = show_scalar_translation(io, m.b) -show_scalar_translation(io, b::Real) = print(io, "x -> x", b < 0 ? " - " : " + ", abs(b)) -show_scalar_translation(io, b) = print(io, "x -> x + ", b) - - -"Translation by a static vector." -struct StaticTranslation{T,N} <: Translation{SVector{N,T}} - b :: SVector{N,T} -end - -"Translation by a vector." -struct VectorTranslation{T} <: Translation{Vector{T}} - b :: Vector{T} -end - -"Translation by a generic vectorlike object." -struct GenericTranslation{T,B} <: Translation{T} - b :: B -end - -Translation(b::Number) = ScalarTranslation(b) -Translation(b::StaticVector) = StaticTranslation(b) -Translation(b::Vector) = VectorTranslation(b) -Translation(b) = GenericTranslation(b) - -Translation{T}(b::Number) where {T<:Number} = ScalarTranslation{T}(b) -Translation{T}(b::AbstractVector) where {N,S,T<:StaticVector{N,S}} = StaticTranslation{S,N}(b) -Translation{T}(b::Vector) where {S,T<:Vector{S}} = VectorTranslation{S}(b) -Translation{T}(b) where {T} = GenericTranslation{T}(b) - -jacdet(m::Translation{T}) where {T} = UnityMap{T,eltype(T)}() -jacdet(m::Translation{T}, x) where {T} = one(eltype(T)) - -isreal(m::Translation) = isreal(unsafe_vector(m)) - -isequalmap(m1::Translation, m2::Translation) = unsafe_vector(m1)==unsafe_vector(m2) - -similarmap(m::Translation, ::Type{T}) where {T} = Translation{T}(m.b) - -applymap(m::Translation, x) = _translation_applymap(m, x, unsafe_vector(m)) -_translation_applymap(m, x, b) = x + b -applymap!(y, m::Translation, x) = _translation_applymap!(y, m, x, unsafe_vector(m)) -_translation_applymap!(y, m, x, b) = y .= x .+ m.b - -inverse(m::Translation{T}) where {T} = Translation{T}(-m.b) -inverse(m::Translation, x) = x - m.b - -ScalarTranslation(b::Number) = ScalarTranslation{typeof(b)}(b) - -StaticTranslation(b::AbstractVector{T}) where {T} = StaticTranslation{T}(b) - -StaticTranslation{T}(b::StaticVector{N}) where {N,T} = - StaticTranslation{T,N}(b) - -VectorTranslation(b::AbstractVector{T}) where {T} = VectorTranslation{T}(b) - -GenericTranslation(b) = GenericTranslation{typeof(b)}(b) -GenericTranslation(b::AbstractVector{T}) where {T} = - GenericTranslation{Vector{T}}(b) -GenericTranslation(b::StaticVector{N,T}) where {N,T} = - GenericTranslation{SVector{N,T}}(b) - -GenericTranslation{T}(b) where {T} = GenericTranslation{T,typeof(b)}(b) -GenericTranslation{T}(b::Number) where {T<:Number} = - GenericTranslation{T,T}(b) - - -############################ -# Affine maps: y = A*x + b -############################ - - -""" -The supertype of all affine maps that store `A` and `b`. -Concrete subtypes differ in how `A` and `b` are represented. -""" -abstract type AffineMap{T} <: AbstractAffineMap{T} end - -AffineMap(A::Number, b::Number) = ScalarAffineMap(A, b) -AffineMap(A::StaticMatrix, b::StaticVector) = StaticAffineMap(A, b) -AffineMap(A::Matrix, b::Vector) = VectorAffineMap(A, b) -AffineMap(A::UniformScaling{Bool}, b::Number) = ScalarAffineMap(one(b), b) -AffineMap(A, b) = GenericAffineMap(A, b) - -AffineMap{T}(A::Number, b::Number) where {T<:Number} = ScalarAffineMap{T}(A, b) -AffineMap{T}(A::AbstractMatrix, b::AbstractVector) where {N,S,T<:SVector{N,S}} = StaticAffineMap{S,N}(A, b) -AffineMap{T}(A::Matrix, b::Vector) where {S,T<:Vector{S}} = VectorAffineMap{S}(A, b) -AffineMap{T}(A::UniformScaling{Bool}, b::Number) where {T} = ScalarAffineMap{T}(one(T), b) -AffineMap{T}(A, b) where {T} = GenericAffineMap{T}(A, b) - -similarmap(m::AffineMap, ::Type{T}) where {T} = AffineMap{T}(m.A, m.b) - -convert(::Type{AffineMap}, m) = (@assert isaffinemap(m); AffineMap(affinematrix(m), affinevector(m))) -convert(::Type{AffineMap{T}}, m) where {T} = (@assert isaffinemap(m); AffineMap{T}(affinematrix(m), affinevector(m))) -# avoid ambiguity errors with convert(::Type{T}, x::T) in Base: -convert(::Type{AffineMap}, m::AffineMap) = m -convert(::Type{AffineMap{T}}, m::AffineMap{T}) where T = m - -# If y = A*x+b, then x = inv(A)*(y-b) = inv(A)*y - inv(A)*b -inverse(m::AffineMap) = (@assert issquaremap(m); AffineMap(inv(m.A), -inv(m.A)*m.b)) -inverse(m::AffineMap, x) = (@assert issquaremap(m); m.A \ (x-m.b)) - -function leftinverse(m::AffineMap) - @assert isoverdetermined(m) - pA = matrix_pinv(m.A) - AffineMap(pA, -pA*m.b) -end -function rightinverse(m::AffineMap) - @assert isunderdetermined(m) - pA = matrix_pinv(m.A) - AffineMap(pA, -pA*m.b) -end -function leftinverse(m::AffineMap, x) - @assert isoverdetermined(m) - m.A \ (x-m.b) -end -function rightinverse(m::AffineMap, x) - @assert isunderdetermined(m) - m.A \ (x-m.b) -end - - -"An affine map for any combination of types of `A` and `b`." -struct GenericAffineMap{T,AA,B} <: AffineMap{T} - A :: AA - b :: B -end - -GenericAffineMap(A, b) = GenericAffineMap{typeof(b)}(A, b) -GenericAffineMap(A::AbstractVector{S}, b::AbstractVector{T}) where {S,T} = - GenericAffineMap{promote_type(S,T)}(A, b) -GenericAffineMap(A::AbstractArray{S}, b::AbstractVector{T}) where {S,T} = - GenericAffineMap{Vector{promote_type(S,T)}}(A, b) -GenericAffineMap(A::StaticMatrix{M,N,S}, b::StaticVector{M,T}) where {M,N,S,T} = - GenericAffineMap{SVector{N,promote_type(S,T)}}(A, b) -GenericAffineMap(A::StaticMatrix{M,N,S}, b::AbstractVector{T}) where {M,N,S,T} = - GenericAffineMap{SVector{N,promote_type(S,T)}}(A, b) -GenericAffineMap(A::S, b::AbstractVector{T}) where {S<:Number,T} = - GenericAffineMap{Vector{promote_type(S,T)}}(A, b) -GenericAffineMap(A::S, b::StaticVector{N,T}) where {S<:Number,N,T} = - GenericAffineMap{SVector{N,promote_type(S,T)}}(A, b) -GenericAffineMap(A::UniformScaling{Bool}, b) = - GenericAffineMap(UniformScaling{eltype(b)}(1), b) - - -# Fallback routine for generic A and b, special cases follow -GenericAffineMap{T}(A, b) where {T} = GenericAffineMap{T,typeof(A),typeof(b)}(A, b) - -GenericAffineMap{T}(A::AbstractVector{S}, b::AbstractVector{U}) where {T<:Number,S,U} = - GenericAffineMap{T}(convert(AbstractVector{T}, A), convert(AbstractVector{T}, b)) -GenericAffineMap{T}(A::AbstractVector{T}, b::AbstractVector{T}) where {T<:Number} = - GenericAffineMap{T,typeof(A),typeof(b)}(A, b) -GenericAffineMap{T}(A::Number, b) where {T} = GenericAffineMap{T,eltype(T),typeof(b)}(A, b) -GenericAffineMap{T}(A::Number, b::AbstractVector) where {N,S,T <: StaticVector{N,S}} = - GenericAffineMap{T,S,SVector{N,S}}(A, b) -# Promote element types of abstract arrays -GenericAffineMap{T}(A::AbstractArray, b::AbstractVector) where {S,T<:AbstractVector{S}} = - GenericAffineMap{T}(convert(AbstractArray{eltype(T)},A), convert(AbstractVector{eltype(T)}, b)) -GenericAffineMap{T}(A::AbstractArray{S}, b::AbstractVector{S}) where {S,T<:AbstractVector{S}} = - GenericAffineMap{T,typeof(A),typeof(b)}(A, b) -GenericAffineMap{T}(A::UniformScaling{Bool}, b::AbstractVector) where {S,T<:AbstractVector{S}} = - GenericAffineMap{T}(A*one(S), convert(AbstractVector{S}, b)) -GenericAffineMap{T}(A::UniformScaling{S}, b::AbstractVector{S}) where {S,T<:AbstractVector{S}} = - GenericAffineMap{T,typeof(A),typeof(b)}(A, b) - - -similarmap(m::GenericAffineMap, ::Type{T}) where {T} = AffineMap{T}(m.A, m.b) - -convert(::Type{GenericAffineMap{T}}, m::GenericAffineMap) where {T} = - GenericAffineMap{T}(m.A, m.b) - - - -"An affine map with scalar representation." -struct ScalarAffineMap{T} <: AffineMap{T} - A :: T - b :: T -end - -ScalarAffineMap(A, b) = ScalarAffineMap(promote(A, b)...) - -isreal(m::ScalarAffineMap{T}) where {T} = isrealtype(T) - -show(io::IO, m::ScalarAffineMap) = show_scalar_affine_map(io, m.A, m.b) -show_scalar_affine_map(io, A::Real, b::Real) = print(io, "x -> $(A) * x", b < 0 ? " - " : " + ", abs(b)) -show_scalar_affine_map(io, A::Complex, b::Complex) = print(io, "x -> ($(A)) * x + ", b) -show_scalar_affine_map(io, A, b) = print(io, "x -> ($(A)) * x + $(b)") - - -convert(::Type{ScalarAffineMap{T}}, m::ScalarAffineMap) where {T} = - ScalarAffineMap{T}(m.A, m.b) - -"An affine map with array and vector representation." -struct VectorAffineMap{T} <: AffineMap{Vector{T}} - A :: Matrix{T} - b :: Vector{T} -end - -VectorAffineMap(A::AbstractArray{T}, b::AbstractVector{T}) where {T} = - VectorAffineMap{T}(A, b) -function VectorAffineMap(A::AbstractArray{S}, b::AbstractVector{T}) where {S,T} - U = promote_type(S,T) - VectorAffineMap(convert(AbstractArray{U}, A), convert(AbstractVector{U}, b)) -end - -convert(::Type{VectorAffineMap{T}}, m::VectorAffineMap) where {T} = - VectorAffineMap{T}(m.A, m.b) - - - -"An affine map with representation using static arrays." -struct StaticAffineMap{T,N,M,L} <: AffineMap{SVector{N,T}} - A :: SMatrix{M,N,T,L} - b :: SVector{M,T} -end - -# Constructors: -# - first, we deduce T -StaticAffineMap(A::AbstractMatrix{T}, b::AbstractVector{T}) where {T} = - StaticAffineMap{T}(A, b) -function StaticAffineMap(A::AbstractMatrix{S}, b::AbstractVector{T}) where {S,T} - U = promote_type(S,T) - StaticAffineMap(convert(AbstractMatrix{U}, A), convert(AbstractVector{U}, b)) -end - -StaticAffineMap{T}(A::AbstractMatrix, b::AbstractVector) where {T} = - StaticAffineMap{T}(convert(AbstractMatrix{T}, A), convert(AbstractVector{T}, b)) - -# - then, we determine N and/or M, from the arguments -function StaticAffineMap{T}(A::AbstractMatrix{T}, b::StaticVector{M,T}) where {T,M} - @assert size(A) == (M,M) - StaticAffineMap{T,M,M}(A, b) -end -StaticAffineMap{T}(A::StaticMatrix{M,N,T}, b::AbstractVector) where {T,N,M} = - StaticAffineMap{T,N,M}(A, b) -StaticAffineMap{T}(A::StaticMatrix{M,N,T}, b::StaticVector{M,T}) where {T,N,M} = - StaticAffineMap{T,N,M}(A, b) -# line below catches ambiguity error -StaticAffineMap{T}(A::StaticMatrix{M1,N,T}, b::StaticVector{M2,T}) where {T,N,M1,M2} = - throw(ArgumentError("Non-matching dimensions")) -StaticAffineMap{T,N}(A::AbstractMatrix, b::AbstractVector) where {T,N} = - StaticAffineMap{T,N,N}(A, b) -StaticAffineMap{T,N}(A::StaticMatrix{M,N}, b::AbstractVector) where {T,N,M} = - StaticAffineMap{T,N,M}(A, b) - -# - finally invoke the constructor (and implicitly convert the data if necessary) -StaticAffineMap{T,N,M}(A::AbstractMatrix, b::AbstractVector) where {T,N,M} = - StaticAffineMap{T,N,M,M*N}(A, b) - -convert(::Type{Map{SVector{N,T}}}, m::VectorAffineMap) where {N,T} = - StaticAffineMap{T,N}(m.A, m.b) diff --git a/FunctionMaps/src/types/arithmetics.jl b/FunctionMaps/src/types/arithmetics.jl deleted file mode 100644 index d768b12..0000000 --- a/FunctionMaps/src/types/arithmetics.jl +++ /dev/null @@ -1,141 +0,0 @@ - -# Simplifications go here - -# Basic maps - -composedmap1(m1::IdentityMap, m2) = m2 -composedmap1(m1::IdentityMap{T}, m2::Map{T}) where {T} = m2 -composedmap1(m1::IdentityMap{T}, m2::Map{S}) where {S,T} = convert(Map{T}, m2) - -composedmap2(m1, m2::IdentityMap) = m1 -composedmap2(m1::Map{T}, m2::IdentityMap{T}) where {T} = m1 -composedmap2(m1::Map{S}, m2::IdentityMap{T}) where {S,T} = convert(Map{T}, m1) - -composedmap2(m1, m2::ConstantMap) = m2 -composedmap2(m1::Map{T}, m2::ConstantMap) where {T} = ConstantMap{T}(mapconstant(m2)) -composedmap1(m1::ConstantMap{T}, m2) where {T} = ConstantMap{T}(m2(mapconstant(m1))) - -composedmap2(m1, m2::ZeroMap) = m2 -composedmap2(m1::Map{T}, m2::ZeroMap{S,U}) where {S,T,U} = ZeroMap{T,U}() -composedmap1(m1::ConstantMap{T}, m2::ZeroMap{S,U}) where {S,T,U} = ZeroMap{T,U}() - -multiply_map1(m1::ZeroMap, m2) = m1 -multiply_map2(m1, m2::ZeroMap) = m2 -multiply_map2(m1, m2::ConstantMap) = _constant_multiply_map(m1, m2) -_constant_multiply_map(m1::ConstantMap{T}, m2::ConstantMap{S}) where {S,T} = - ConstantMap{promote_type(S,T)}(mapconstant(m1)*mapconstant(m2)) -_constant_multiply_map(m1, m2::ConstantMap) = default_multiply_map(m1, m2) - -multiply_map(m1::Function, m2::Function) = t -> m1(t)*m2(t) - -sum_map1(m1::ZeroMap, m2) = m2 -sum_map2(m1, m2::ZeroMap) = m1 -sum_map2(m1, m2::ConstantMap) = _constant_sum_map(m1, m2) -_constant_sum_map(m1::ConstantMap{T}, m2::ConstantMap{S}) where {S,T} = - ConstantMap{promote_type(S,T)}(mapconstant(m1)+mapconstant(m2)) -_constant_sum_map(m1, m2::ConstantMap) = default_sum_map(m1, m2) - -## Affine maps - -composedmap(m1::AbstractAffineMap, m2::AbstractAffineMap) = affine_composition(m1, m2) - -""" -Compute the affine map that represents map2 after map1, that is: -`y = a2*(a1*x+b1)+b2 = a2*a1*x + a2*b1 + b2`. -""" -affine_composition(map1::AbstractAffineMap, map2::AbstractAffineMap) = - AffineMap(affinematrix(map2) * affinematrix(map1), - affinematrix(map2)*affinevector(map1) + affinevector(map2)) - -affine_composition(map1::AffineMap, map2::AffineMap) = - AffineMap(unsafe_matrix(map2) * unsafe_matrix(map1), - unsafe_matrix(map2)*unsafe_vector(map1) + unsafe_vector(map2)) - -affine_composition(map1::LinearMap, map2::LinearMap) = - LinearMap(unsafe_matrix(map2) * unsafe_matrix(map1)) - -affine_composition(map1::LinearMap, map2::AffineMap) = - AffineMap(unsafe_matrix(map2) * unsafe_matrix(map1), unsafe_vector(map2)) - -affine_composition(map1::AffineMap, map2::LinearMap) = - AffineMap(unsafe_matrix(map2) * unsafe_matrix(map1), unsafe_matrix(map2)*unsafe_vector(map1)) - -affine_composition(map1::Translation, map2::Translation) = - Translation(unsafe_vector(map2) + unsafe_vector(map1)) - -affine_composition(map1::Translation, map2::LinearMap) = - AffineMap(unsafe_matrix(map2), unsafe_matrix(map2)*unsafe_vector(map1)) - -affine_composition(map1::LinearMap, map2::Translation) = - AffineMap(unsafe_matrix(map1), unsafe_vector(map2)) - -# The sum of two affine maps is again an affine map -sum_map(map1::AbstractAffineMap, map2::AbstractAffineMap) = - AffineMap(affinematrix(map1)+affinematrix(map2), affinevector(map1)+affinevector(map2)) - - -isequalmap(m1::ProductMap, m2::IdentityMap) = all(map(isidentitymap, components(m1))) -isequalmap(m1::IdentityMap, m2::ProductMap) = isequalmap(m2, m1) - - -""" -Map the interval `[a,b]` to the interval `[c,d]`. - -This function deals with infinite intervals, and the type of the -map returned may depend on the value (finiteness) of the given endpoints. -""" -interval_map(a, b, c, d) = interval_map(promote(a,b,c,d)...) - -function interval_map(a::T, b::T, c::T, d::T) where {T} - FT = float(T) - if isfinite(a) && isfinite(b) && isfinite(c) && isfinite(d) - bounded_interval_map(a, b, c, d) - elseif isfinite(a) && !isfinite(b) && isfinite(c) && !isfinite(d) - # (a,Inf) to (c,Inf) - AffineMap(one(FT), c-a) - elseif isfinite(a) && !isfinite(b) && !isfinite(c) && isfinite(d) - # (a,Inf) to (Inf,d) - AffineMap(-one(FT), d+a) - elseif !isfinite(a) && isfinite(b) && isfinite(c) && !isfinite(d) - # (Inf,b) to (c,Inf) - AffineMap(-one(FT), c+b) - elseif !isfinite(a) && isfinite(b) && !isfinite(c) && isfinite(d) - # (Inf,b) to (Inf,d) - AffineMap(one(FT), d-b) - elseif !isfinite(a) && !isfinite(b) && !isfinite(c) && !isfinite(d) - if (a < 0) && (b > 0) && (c < 0) && (d > 0) - # (-Inf,Inf) to (-Inf,Inf) - StaticIdentityMap{FT}() - elseif (a < 0) && (b > 0) && (c > 0) && (d < 0) - # (-Inf,Inf) to (Inf,-Inf) - LinearMap(-one(FT)) - elseif (a > 0) && (b < 0) && (c < 0) && (d > 0) - # (Inf,-Inf) to (-Inf,Inf) - LinearMap(-one(FT)) - elseif (a > 0) && (b < 0) && (c > 0) && (d < 0) - # (Inf,-Inf) to (Inf,-Inf) - StaticIdentityMap{FT}() - elseif (a > 0) && (b > 0) && (c > 0) && (d > 0) - # (Inf,Inf) to (Inf,Inf) - StaticIdentityMap{FT}() - elseif (a < 0) && (b < 0) && (c < 0) && (d < 0) - # (-Inf,-Inf) to (-Inf,-Inf) - StaticIdentityMap{FT}() - elseif (a < 0) && (b < 0) && (c > 0) && (d > 0) - # (-Inf,-Inf) to (Inf,Inf) - LinearMap(-one(FT)) - elseif (a > 0) && (b > 0) && (c < 0) && (d < 0) - # (Inf,Inf) to (-Inf,-Inf) - LinearMap(-one(FT)) - else - throw(ArgumentError("Requested affine map is unbounded")) - end - else - throw(ArgumentError("Requested affine map is unbounded")) - end -end - -"Like interval_map, but guaranteed to return a scalar affine map." -bounded_interval_map(a, b, c, d) = bounded_interval_map(promote(a,b,c,d)...) -bounded_interval_map(a::T, b::T, c::T, d::T) where {T} = - AffineMap((d-c)/(b-a), c - a*(d-c)/(b-a)) diff --git a/FunctionMaps/src/types/basic.jl b/FunctionMaps/src/types/basic.jl deleted file mode 100644 index c9728b3..0000000 --- a/FunctionMaps/src/types/basic.jl +++ /dev/null @@ -1,160 +0,0 @@ - -"Supertype of identity maps." -abstract type IdentityMap{T} <: Map{T} end - -IdentityMap(n::Int) = DynamicIdentityMap(n) -IdentityMap() = StaticIdentityMap() -IdentityMap(::Val{N}) where {N} = StaticIdentityMap(Val(N)) - -IdentityMap{T}(n::Int) where {T} = DynamicIdentityMap{T}(n) -IdentityMap{T}(n::Int) where {T<:StaticTypes} = StaticIdentityMap{T}() -IdentityMap{T}(::Val{N}) where {N,T} = StaticIdentityMap{T}(Val(N)) -IdentityMap{T}() where {T} = StaticIdentityMap{T}() - -applymap(map::IdentityMap, x) = x -applymap!(y, map::IdentityMap, x) = y .= x - -inverse(m::IdentityMap) = m -inverse(m::IdentityMap, x) = x - -islinearmap(::IdentityMap) = true -isreal(::IdentityMap{T}) where {T} = isrealtype(T) - -isidentitymap(::IdentityMap) = true -isidentitymap(m::Map{T}) where {T} = m == StaticIdentityMap{T}() - -mapsize(m::IdentityMap{T}) where {T<:Number} = () -mapsize(m::IdentityMap{T}) where {T} = (euclideandimension(T),euclideandimension(T)) - -affinematrix(m::IdentityMap) = identitymatrix(m) -affinevector(m::IdentityMap) = zerovector(m) - -jacobian(m::IdentityMap) = ConstantMap(affinematrix(m)) -jacobian(m::IdentityMap, x) = affinematrix(m) - -jacdet(m::IdentityMap, x) = 1 - -determinantmap(m::IdentityMap{T}) where {T} = UnityMap{T,prectype(T)}() - -mapcompose(m1::IdentityMap) = m1 -mapcompose(m1::IdentityMap, maps...) = mapcompose(maps...) -mapcompose2(m1, m2::IdentityMap, maps...) = mapcompose(m1, maps...) - -show(io::IO, m::IdentityMap{T}) where {T} = print(io, "x -> x") -map_object_parentheses(m::IdentityMap) = true - -"The identity map for variables of type `T`." -struct StaticIdentityMap{T} <: IdentityMap{T} -end - -StaticIdentityMap() = StaticIdentityMap{Float64}() -StaticIdentityMap(::Val{N}) where {N} = StaticIdentityMap{SVector{N,Float64}}() - -StaticIdentityMap{T}(n::Int) where {T} = - (@assert n == euclideandimension(T); StaticIdentityMap{T}()) -StaticIdentityMap{T}(::Val{N}) where {N,T} = - (@assert N == euclideandimension(T); StaticIdentityMap{T}()) - -similarmap(m::StaticIdentityMap, ::Type{T}) where {T<:StaticTypes} = StaticIdentityMap{T}() -similarmap(m::StaticIdentityMap, ::Type{T}) where {T} = - DynamicIdentityMap{T}(euclideandimension(T)) - -convert(::Type{StaticIdentityMap{T}}, ::StaticIdentityMap) where {T} = StaticIdentityMap{T}() - -isequalmap(m1::StaticIdentityMap, m2::StaticIdentityMap) = true -map_hash(m::StaticIdentityMap, h::UInt) = hash("StaticIdentityMap", h) - -"Identity map with dynamic size determined by a dimension field." -struct DynamicIdentityMap{T} <: IdentityMap{T} - dimension :: Int -end - -const EuclideanIdentityMap{N,T} = StaticIdentityMap{SVector{N,T}} -const VectorIdentityMap{T} = DynamicIdentityMap{Vector{T}} - -DynamicIdentityMap(dimension::Int) = VectorIdentityMap(dimension) -VectorIdentityMap(dimension::Int) = VectorIdentityMap{Float64}(dimension) - -mapsize(m::DynamicIdentityMap) = (m.dimension, m.dimension) - -similarmap(m::DynamicIdentityMap, ::Type{T}) where {T} = - DynamicIdentityMap{T}(m.dimension) -similarmap(m::DynamicIdentityMap, ::Type{T}) where {T<:StaticTypes} = - StaticIdentityMap{T}() - -isequalmap(m1::DynamicIdentityMap, m2::DynamicIdentityMap) = m1.dimension == m2.dimension -map_hash(m::DynamicIdentityMap, h::UInt) = hashrec("DynamicIdentityMap", m.dimension, h) - - -"The supertype of constant maps from `T` to `U`." -abstract type ConstantMap{T,U} <: TypedMap{T,U} end - -applymap(m::ConstantMap, x) = mapconstant(m) - -isconstantmap(m::AbstractMap) = false -isconstantmap(m::ConstantMap) = true - -isreal(m::ConstantMap{T,U}) where {T,U} = - isrealtype(T) && isrealtype(U) && isreal(mapconstant(m)) - -mapsize(m::ConstantMap) = _constant_mapsize(m, mapconstant(m)) -_constant_mapsize(m::ConstantMap{T,U}, c) where {T<:Number,U<:Number} = () -_constant_mapsize(m::ConstantMap{T,U}, c) where {T<:Number,U} = (length(c),) -_constant_mapsize(m::ConstantMap{T,U}, c) where {T,U<:Number} = (1,euclideandimension(T)) -_constant_mapsize(m::ConstantMap{T,U}, c) where {T,U} = (length(c), euclideandimension(T)) - -affinematrix(m::ConstantMap) = zeromatrix(m) -affinevector(m::ConstantMap) = mapconstant(m) - -jacobian(m::ConstantMap{T}) where {T} = ConstantMap{T}(affinematrix(m)) -jacobian(m::ConstantMap, x) = affinematrix(m) - -jacdet(::ConstantMap, x) = 0 - -determinantmap(m::ConstantMap{T}) where {T} = ConstantMap{T}(det(mapconstant(m))) -absmap(m::ConstantMap{T}) where {T} = ConstantMap{T}(abs(mapconstant(m))) - -diffvolume(m::ConstantMap{T,U}) where {T,U} = ZeroMap{T,U}() - -isequalmap(m1::ConstantMap, m2::ConstantMap) = mapconstant(m1)==mapconstant(m2) -map_hash(m::ConstantMap, h::UInt) = hashrec("ConstantMap", mapconstant(m), h) - -similarmap(m::ConstantMap, ::Type{T}) where {T} = ConstantMap{T}(mapconstant(m)) -similarmap(m::ConstantMap, ::Type{T}, ::Type{U}) where {T,U} = ConstantMap{T,U}(m.c) - -ConstantMap() = ConstantMap{Float64}() -ConstantMap(c) = FixedConstantMap(c) -ConstantMap{T}() where {T} = UnityMap{T}() -ConstantMap{T}(c) where {T} = FixedConstantMap{T}(c) -ConstantMap{T,U}() where {T,U} = UnityMap{T,U}() -ConstantMap{T,U}(c) where {T,U} = FixedConstantMap{T,U}(c) - -show(io::IO, m::ConstantMap{T}) where {T} = print(io, "x -> $(mapconstant(m))") -map_object_parentheses(m::ConstantMap) = true - - -"The zero map `f(x) = 0`." -struct ZeroMap{T,U} <: ConstantMap{T,U} -end -ZeroMap{T}() where {T} = ZeroMap{T,T}() -mapconstant(m::ZeroMap{T,U}) where {T,U} = zero(U) -similarmap(m::ZeroMap{S,U}, ::Type{T}) where {T,S,U} = ZeroMap{T,U}() -similarmap(m::ZeroMap, ::Type{T}, ::Type{U}) where {T,U} = ZeroMap{T,U}() - - -"The unity map `f(x) = 1`." -struct UnityMap{T,U} <: ConstantMap{T,U} -end -UnityMap{T}() where {T} = UnityMap{T,real(numtype(T))}() -mapconstant(m::UnityMap{T,U}) where {T,U} = one(U) -similarmap(m::UnityMap{S,U}, ::Type{T}) where {T,S,U} = UnityMap{T,U}() -similarmap(m::UnityMap, ::Type{T}, ::Type{U}) where {T,U} = UnityMap{T,U}() - - -"The constant map `f(x) = c`." -struct FixedConstantMap{T,U} <: ConstantMap{T,U} - c :: U -end -FixedConstantMap{T}(c::U) where {T,U} = FixedConstantMap{T,U}(c) -FixedConstantMap(c::T) where {T} = FixedConstantMap{T}(c) -mapconstant(m::FixedConstantMap) = m.c diff --git a/FunctionMaps/src/types/coordinates.jl b/FunctionMaps/src/types/coordinates.jl deleted file mode 100644 index daf9b76..0000000 --- a/FunctionMaps/src/types/coordinates.jl +++ /dev/null @@ -1,138 +0,0 @@ - -# This file contains a few specific maps, included mainly to test -# functionality on non-trivial maps - -""" -A Cartesion to Polar map. First dimension is interpreted as radial distance, -second as an angle. -The unit circle is mapped to the square `[-1,1]x[-1,1]`. -""" -struct CartToPolarMap{T} <: Map{SVector{2,T}} -end - -CartToPolarMap() = CartToPolarMap{Float64}() - -mapsize(m::CartToPolarMap) = (2,2) - -applymap(map::CartToPolarMap{T}, x) where {T} = - SVector{2,T}(sqrt(x[1]^2+x[2]^2)*2-1, atan(x[2],x[1])/pi) - -function jacobian(m::CartToPolarMap{T}, x) where {T} - d = sqrt(x[1]^2+x[2]^2) - r = 1/(1+(x[2]/x[1])^2) - SMatrix{2,2,T}(2*x[1]/d, 1/pi*r*(-x[2]/x[1]^2), 2*x[2]/d, 1/pi*r*1/x[1]) -end - -inverse(m::CartToPolarMap{T}) where {T} = PolarToCartMap{T}() -inverse(m::CartToPolarMap, x) = inverse(m)(x) - -isreal(m::CartToPolarMap) = true - -convert(::Type{Map{SVector{2,T}}}, ::CartToPolarMap) where {T} = CartToPolarMap{T}() - -isequalmap(m1::CartToPolarMap, m2::CartToPolarMap) = true - - -""" -A Polar to Cartesian map. The angle is mapped to the second dimension, -radius to the first. -The square `[-1,1]x[-1,1]` is mapped to the unit circle. -""" -struct PolarToCartMap{T} <: Map{SVector{2,T}} -end - -PolarToCartMap() = PolarToCartMap{Float64}() - -mapsize(m::PolarToCartMap) = (2,2) - -applymap(map::PolarToCartMap{T}, x) where {T} = SVector{2,T}((x[1]+1)/2*cos(pi*x[2]), (x[1]+1)/2*sin(pi*x[2])) - -jacobian(m::PolarToCartMap{T}, x) where {T} = - SMatrix{2,2,T}(cos(pi*x[2])/2, sin(pi*x[2])/2, -pi*(x[1]+1)/2*sin(pi*x[2]), pi*(x[1]+1)/2*cos(pi*x[2])) - -inverse(m::PolarToCartMap{T}) where {T} = CartToPolarMap{T}() -inverse(m::PolarToCartMap, x) = inverse(m)(x) - -isreal(m::PolarToCartMap) = true - -convert(::Type{Map{SVector{2,T}}}, ::PolarToCartMap) where {T} = PolarToCartMap{T}() - -isequalmap(m1::PolarToCartMap, m2::PolarToCartMap) = true - - -""" -The map `[cos(2πt), sin(2πt)]` from `[0,1]` to the unit circle in `ℝ^2`. -""" -struct UnitCircleMap{T} <: Map{T} end - -UnitCircleMap() = UnitCircleMap{Float64}() - -mapsize(m::UnitCircleMap) = (2,) - -applymap(m::UnitCircleMap{T}, t) where {T} = SVector(cos(2*T(pi)*t), sin(2*T(pi)*t)) - -function jacobian(m::UnitCircleMap{T}, t) where {T} - a = 2*T(pi) - SVector(-a*sin(a*t), a*cos(a*t)) -end - -# we know that the differential volume is 2*pi -diffvolume(m::UnitCircleMap{T}) where {T} = ConstantMap{T}(2*T(pi)) - -""" -`AngleMap` is a left inverse of `UnitCircleMap`. A 2D vector `x` is projected onto -the intersection point with the unit circle of the line connecting `x` to the -origin. The angle of this point, scaled to the interval `[0,1)`, is the result. -""" -struct AngleMap{T} <: Map{SVector{2,T}} -end - -AngleMap() = AngleMap{Float64}() - -function applymap(m::AngleMap{T}, x) where {T} - twopi = 2*convert(T, pi) - θ = atan(x[2],x[1]) - if θ < 0 - # atan2 returns an angle in (-π,π], convert to [0,2π) using periodicity. - θ += twopi - end - # And divide by 2π to scale to [0,1) - θ / twopi -end - -mapsize(m::AngleMap) = (1,2) - -function jacobian(m::AngleMap{T}, t) where {T} - x = t[1]; y = t[2] - twopi = 2*convert(T, pi) - v = SVector(one(T) / (1+(y/x)^2) * (-y/x^2) * 1/twopi, - one(T) / (1+(y/x)^2) * one(T)/x * 1/twopi) - transpose(v) -end - - -leftinverse(m::UnitCircleMap{T}) where {T} = AngleMap{T}() -leftinverse(m::UnitCircleMap, x) = leftinverse(m)(x) -rightinverse(m::AngleMap{T}) where {T} = UnitCircleMap{T}() -rightinverse(m::AngleMap, x) = rightinverse(m)(x) - - -""" -The map `r*[cos(2πt), sin(2πt)]` from `[0,1]^2` to the unit disk in `ℝ^2`. -""" -struct UnitDiskMap{T} <: Map{SVector{2,T}} end - -UnitDiskMap() = UnitDiskMap{Float64}() - -mapsize(m::UnitDiskMap) = (2,2) - -applymap(m::UnitDiskMap{T}, x) where {T} = - SVector(x[1]*cos(2*T(pi)*x[2]), x[1]*sin(2*T(pi)*x[2])) - -function jacobian(m::UnitDiskMap{T}, x) where {T} - a = 2*T(pi) - SMatrix{2,2}(cos(a*x[2]), sin(a*x[2]), -a*x[1]*sin(a*x[2]), a*x[1]*cos(a*x[2])) -end - -# we know that the differential volume is 2*r*pi -diffvolume(m::UnitDiskMap, x) = 2*x[1]*pi diff --git a/FunctionMaps/src/util/common.jl b/FunctionMaps/src/util/common.jl deleted file mode 100644 index 0c962d6..0000000 --- a/FunctionMaps/src/util/common.jl +++ /dev/null @@ -1,195 +0,0 @@ - -isrealtype(::Type{<:Real}) = true -isrealtype(::Type{<:Complex}) = false -isrealtype(::Type{T}) where {T} = isrealtype(eltype(T)) - -const StaticTypes = Union{Number,<:StaticVector{N} where N,<:NTuple{N,Any} where N} - -"Apply the `hash` function recursively to the given arguments." -hashrec(x) = hash(x) -hashrec(x, args...) = hash(x, hashrec(args...)) - -# Workaround for #88, manually compute the hash of an array using all its elements -hashrec() = zero(UInt) -function hashrec(A::AbstractArray, args...) - h = hash(size(A)) - for x in A - h = hash(x, h) - end - hash(h, hashrec(args...)) -end - -"What is the euclidean dimension of the given type (if applicable)?" -euclideandimension(::Type{T}) where {T <: Number} = 1 -euclideandimension(::Type{T}) where {N,T <: StaticVector{N}} = N -euclideandimension(::Type{T}) where {N,T <: NTuple{N,Any}} = N -euclideandimension(::Type{T}) where {T} = - throw(ArgumentError("Don't know the euclidean dimension of $(T).")) - -################# -# Element type -################# - -""" - convert_eltype(T, x) - -Convert `x` such that its `eltype` equals `T`. -""" -convert_eltype(::Type{T}, d::AbstractArray) where {T} = convert(AbstractArray{T}, d) -convert_eltype(::Type{T}, d::AbstractRange) where {T} = map(T, d) -convert_eltype(::Type{T}, d::Set) where {T} = convert(Set{T}, d) -convert_eltype(::Type{T}, d::Number) where {T} = convert(T, d) - -"Are the given element types promotable to a concrete supertype?" -promotable_eltypes(types...) = isconcretetype(promote_type(types...)) -promotable_eltypes(::Type{S}, ::Type{T}) where {S<:AbstractVector,T<:AbstractVector} = - promotable_eltypes(eltype(S), eltype(T)) - -################# -# Precision type -################# - -""" - prectype(x[, ...]) - -The floating point precision type associated with the argument(s). -""" -prectype(x) = prectype(typeof(x)) -prectype(::Type{T}) where T = Any # fallback -prectype(::Type{T}) where {T<:AbstractFloat} = T -prectype(::Type{T}) where {T<:Number} = prectype(float(T)) -prectype(::Type{<:AbstractArray{T}}) where T = prectype(T) -prectype(::Type{<:Complex{T}}) where {T} = prectype(T) -prectype(::Type{NTuple{N,T}}) where {N,T} = prectype(T) -prectype(::Type{Tuple{A}}) where {A} = prectype(A) -prectype(::Type{Tuple{A,B}}) where {A,B} = prectype(A,B) -prectype(::Type{Tuple{A,B,C}}) where {A,B,C} = prectype(A,B,C) -@generated function prectype(T::Type{<:Tuple}) - quote $(promote_type(map(prectype, T.parameters[1].parameters)...)) end -end - -prectype(a...) = promote_type(map(prectype, a)...) - -""" - convert_prectype(T, x) - -Convert `x` such that its `prectype` equals `T`. -""" -convert_prectype(::Type{T}, x) where {T} = - convert_eltype(to_prectype(T, eltype(x)), x) - -""" - to_prectype(U, T) - -Return the type to which `T` can be converted, such that the `prectype` becomes `U`. -""" -to_prectype(::Type{U}, ::Type{T}) where {T,U} = throw(ArgumentError("Don't know how to convert the prectype of $(T) to $(U).")) -to_prectype(::Type{U}, ::Type{T}) where {T <: Real,U <: Real} = U -to_prectype(::Type{U}, ::Type{Complex{T}}) where {T <: Real,U <: Real} = Complex{U} -to_prectype(::Type{U}, ::Type{Vector{T}}) where {T,U} = Vector{to_prectype(U,T)} -to_prectype(::Type{U}, ::Type{<:StaticVector{N,T}}) where {N,T,U} = SVector{N,to_prectype(U,T)} -to_prectype(::Type{U}, ::Type{MVector{N,T}}) where {N,T,U} = MVector{N,to_prectype(U,T)} -to_prectype(::Type{U}, ::Type{<:StaticMatrix{M,N,T}}) where {M,N,T,U} = SMatrix{M,N,to_prectype(U,T)} -to_prectype(::Type{U}, ::Type{MMatrix{M,N,T}}) where {M,N,T,U} = MMatrix{M,N,to_prectype(U,T)} - -""" - promote_prectype(a, b[, ...]) - -Promote the precision types of the arguments to a joined supertype. -""" -promote_prectype(a) = a -promote_prectype(a, b) = _promote_prectype(prectype(a,b), a, b) -promote_prectype(a, b, c...) = _promote_prectype(prectype(a,b,c...), a, b, c...) -_promote_prectype(U, a) = convert_prectype(U, a) -_promote_prectype(U, a, b) = convert_prectype(U, a), convert_prectype(U, b) -_promote_prectype(U, a, b, c...) = - (convert_prectype(U, a), convert_prectype(U, b), _promote_prectype(U, c...)...) - -################# -# Numeric type -################# - -"The numeric element type of `x` in a Euclidean space." -numtype(x) = numtype(typeof(x)) -numtype(::Type{T}) where T = Any -numtype(::Type{T}) where {T<:Number} = T -numtype(::Type{<:AbstractArray{T}}) where T = T -numtype(::Type{NTuple{N,T}}) where {N,T} = T -numtype(::Type{Tuple{A,B}}) where {A,B} = promote_type(numtype(A), numtype(B)) -numtype(::Type{Tuple{A,B,C}}) where {A,B,C} = promote_type(numtype(A), numtype(B), numtype(C)) -numtype(::Type{Tuple{A,B,C,D}}) where {A,B,C,D} = promote_type(numtype(A), numtype(B), numtype(C), numtype(D)) -@generated function numtype(T::Type{<:Tuple}) - quote $(promote_type(T.parameters[1].parameters...)) end -end - -numtype(a...) = promote_type(map(numtype, a)...) - -""" - convert_numtype(T, x) - -Convert `x` such that its `numtype` equals `T`. -""" -convert_numtype(::Type{T}, x) where {T} = - convert_eltype(to_numtype(T, eltype(x)), x) - -""" - to_numtype(U, T) - -Return the type to which `T` can be converted, such that the `numtype` becomes `U`. -""" -to_numtype(::Type{U}, ::Type{T}) where {T,U} = throw(ArgumentError("Don't know how to convert the numtype of $(T) to $(U).")) -to_numtype(::Type{U}, ::Type{T}) where {T <: Number,U <: Number} = U -to_numtype(::Type{U}, ::Type{Vector{T}}) where {T,U} = Vector{to_numtype(U,T)} -to_numtype(::Type{U}, ::Type{<:StaticVector{N,T}}) where {N,T,U} = SVector{N,to_numtype(U,T)} -to_numtype(::Type{U}, ::Type{MVector{N,T}}) where {N,T,U} = MVector{N,to_numtype(U,T)} -to_numtype(::Type{U}, ::Type{<:StaticMatrix{M,N,T}}) where {M,N,T,U} = SMatrix{M,N,to_numtype(U,T)} -to_numtype(::Type{U}, ::Type{MMatrix{M,N,T}}) where {M,N,T,U} = MMatrix{M,N,to_numtype(U,T)} - -""" - promote_numtype(a, b[, ...]) - -Promote the numeric types of the arguments to a joined supertype. -""" -promote_numtype(a) = a -promote_numtype(a, b) = _promote_numtype(numtype(a,b), a, b) -promote_numtype(a, b, c...) = _promote_numtype(numtype(a,b,c...), a, b, c...) -_promote_numtype(U, a) = convert_numtype(U, a) -_promote_numtype(U, a, b) = convert_numtype(U, a), convert_numtype(U, b) -_promote_numtype(U, a, b, c...) = - (convert_numtype(U, a), convert_numtype(U, b), _promote_numtype(U, c...)...) - -## Conversion from nested vectors to flat vectors and back - -""" -Convert a vector from a cartesian format to a nested tuple according to the -given dimensions. - -For example: -`convert_fromcartesian([1,2,3,4,5], Val{(2,2,1)}()) -> ([1,2],[3,4],5)` -""" -@generated function convert_fromcartesian(x::AbstractVector, ::Val{DIM}) where {DIM} - dimsum = [0; cumsum([d for d in DIM])] - E = Expr(:tuple, [ (dimsum[i+1]-dimsum[i] > 1 ? Expr(:call, :SVector, [:(x[$j]) for j = dimsum[i]+1:dimsum[i+1]]...) : :(x[$(dimsum[i+1])])) for i in 1:length(DIM)]...) - return quote $(E) end -end - -"The inverse function of `convert_fromcartesian`." -@generated function convert_tocartesian(x, ::Val{DIM}) where {DIM} - dimsum = [0; cumsum([d for d in DIM])] - E = vcat([[:(x[$i][$j]) for j in 1:DIM[i]] for i in 1:length(DIM)]...) - quote SVector($(E...)) end -end - -# we use matrix_pinv rather than pinv to preserve static matrices -matrix_pinv(A) = pinv(A) -# matrix_pinv(A::SMatrix{M,N}) where {M,N} = SMatrix{N,M}(pinv(A)) -matrix_pinv(A::StaticVector{N}) where {N} = convert(Transpose{Float64, SVector{N,Float64}}, pinv(A)) - -## Composite objects - -"Factors of a product-like composite object (equivalent to `components(d)`)." -function factors end -"The number of factors of a product-like composite object." -nfactors(d) = length(factors(d)) -"Factor `I...` of a product-like composite object." -factor(d, I...) = getindex(factors(d), I...) diff --git a/FunctionMaps/test/aqua.jl b/FunctionMaps/test/aqua.jl deleted file mode 100644 index 4f3dad2..0000000 --- a/FunctionMaps/test/aqua.jl +++ /dev/null @@ -1,14 +0,0 @@ -using Aqua - -# Aqua.test_all(FunctionMaps) - -Aqua.test_ambiguities(FunctionMaps) -# The unbound tests annoyingly flag several valid uses of NTuple{N,T} -# Aqua.test_unbound_args(FunctionMaps) -Aqua.test_undefined_exports(FunctionMaps) -Aqua.test_project_extras(FunctionMaps) -Aqua.test_stale_deps(FunctionMaps) -Aqua.test_deps_compat(FunctionMaps) -# Piracies to be removed in next breaking release -# Aqua.test_piracies(FunctionMaps) -Aqua.test_persistent_tasks(FunctionMaps) diff --git a/FunctionMaps/test/runtests.jl b/FunctionMaps/test/runtests.jl deleted file mode 100644 index 5253792..0000000 --- a/FunctionMaps/test/runtests.jl +++ /dev/null @@ -1,20 +0,0 @@ -module FunctionMapsTests - -using Test, LinearAlgebra, StaticArrays -using CompositeTypes, CompositeTypes.Indexing - -include("using_fmaps.jl") - -include("aqua.jl") - -include("test_common.jl") -include("test_interface.jl") - -include("test_generic.jl") -include("test_basic.jl") -include("test_affine.jl") -include("test_product.jl") -include("test_maps.jl") -include("test_arithmetics.jl") - -end diff --git a/FunctionMaps/test/runtests2.jl b/FunctionMaps/test/runtests2.jl deleted file mode 100644 index d9850d3..0000000 --- a/FunctionMaps/test/runtests2.jl +++ /dev/null @@ -1,18 +0,0 @@ -module FunctionMapsTests - -using Test, LinearAlgebra, StaticArrays -using CompositeTypes, CompositeTypes.Indexing - -include("using_fmaps2.jl") - -include("test_common.jl") -include("test_interface.jl") - -include("test_generic.jl") -include("test_basic.jl") -include("test_affine.jl") -include("test_product.jl") -include("test_maps.jl") -include("test_arithmetics.jl") - -end diff --git a/FunctionMaps/test/test_affine.jl b/FunctionMaps/test/test_affine.jl deleted file mode 100644 index 10a9d4d..0000000 --- a/FunctionMaps/test/test_affine.jl +++ /dev/null @@ -1,237 +0,0 @@ -function test_affine_maps(T) - A = rand(T,2,2) - @test FunctionMaps.to_matrix(Vector{T}, A) == A - @test FunctionMaps.to_matrix(T, 2) == 2 - @test FunctionMaps.to_matrix(SVector{2,T}, 2) == SMatrix{2,2}(2,0,0,2) - @test FunctionMaps.to_matrix(SVector{2,T}, LinearAlgebra.I) == SMatrix{2,2}(1,0,0,1) - @test FunctionMaps.to_matrix(Vector{T}, 2) == UniformScaling(2) - @test FunctionMaps.to_matrix(Vector{T}, LinearAlgebra.I) == LinearAlgebra.I - # test fallback with nonsensical call - @test FunctionMaps.to_matrix(Tuple{Int}, 2) == 2 - - @test FunctionMaps.to_matrix(T, A, 2) == A - @test FunctionMaps.to_matrix(T, 2, 3) == 2 - @test FunctionMaps.to_matrix(T, UniformScaling(2), 3) == 2 - @test FunctionMaps.to_matrix(T, LinearAlgebra.I, zero(T)) isa T - @test FunctionMaps.to_matrix(SVector{2,T}, 2, SVector(1,1)) == SMatrix{2,2}(2,0,0,2) - @test FunctionMaps.to_matrix(Vector{T}, 2, [1,2]) == [2 0 ; 0 2] - - @test FunctionMaps.to_vector(T, 2) == 0 - @test FunctionMaps.to_vector(SVector{2,T}, 2) == SVector(0,0) - @test FunctionMaps.to_vector(Vector{T}, A) == [0,0] - @test FunctionMaps.to_vector(T, 2, 3) == 3 - - if T != BigFloat # BigFloat's make pinv fail for StaticArrays - @test FunctionMaps.matrix_pinv(SMatrix{2,2}(rand(T),rand(T),rand(T),rand(T))) isa SMatrix{2,2} - @test FunctionMaps.matrix_pinv(SVector(rand(T),rand(T))) isa Transpose{T,SVector{2,T}} - end - - test_linearmap(T) - test_translation(T) - test_affinemap(T) - test_abstractarraymap(T) -end - - - -function test_linearmap(T) - m1 = LinearMap(2one(T)) - @test m1 isa LinearMap{T} - @test domaintype(m1) == T - @test islinearmap(m1) - @test isaffinemap(m1) - @test m1(3) == 6 - @test m1(3one(T)) == 6 - @test m1(3one(widen(T))) isa widen(T) - @test m1(SVector(1,2)) == SVector(2,4) - @test m1([1.0,2.0]) == [2.0,4.0] - mw1 = convert(Map{widen(T)}, m1) - @test mw1 isa Map{widen(T)} - @test mw1.A isa widen(T) - @test jacobian(m1) isa ConstantMap{T} - @test jacobian(m1, 1) == 2 - @test LinearMap(one(T)) == StaticIdentityMap{T}() - - m2 = LinearMap(2) - @test domaintype(m2) == Int - @test m2 === ScalarLinearMap(2) - @test m2(one(T)) isa T - @test m2 == convert(Map, 2) - @test jacobian(m2, 1) == 2 - @test jacobian(m2) isa ConstantMap{Int} - @test jacobian(m2, 1) == 2 - @test LinearMap(1) == StaticIdentityMap{T}() - - m3 = LinearMap(SMatrix{2,2}(one(T), 2one(T), 3one(T), 4one(T))) - @test m3 isa LinearMap{SVector{2,T}} - @test m3 === StaticLinearMap(m3.A) - @test m3(SVector(1,2)) == SVector(7, 10) - @test m3(SVector{2,T}(1,2)) == SVector{2,T}(7, 10) - @test m3 ∘ m3 isa LinearMap - @test LinearMap(SA[1 0; 0 1]) == StaticIdentityMap{domaintype(m3)}() - - A = rand(T,2,2) - m4 = LinearMap(A) - @test m4 isa LinearMap{Vector{T}} - @test m4([1,2]) == A * [1,2] - y = zeros(T,2) - @test (FunctionMaps.applymap!(y, m4, [1,2]); y == A * [1,2]) - @test jacobian(m4, [1,2]) == A - - m5 = LinearMap{Vector{T}}(UniformScaling(2*one(T))) - @test m5 isa LinearMap{Vector{T}} - @test m5([1,2]) == 2 * [1,2] - @test jacobian(m5, [1,2]) == UniformScaling(2) - - m6 = LinearMap{SVector{2,T}}(2) - @test m6 isa GenericLinearMap{SVector{2,T},T} - @test m6(SA[one(T),one(T)]) == [2,2] - @test leftinverse(m6) isa LinearMap{SVector{2,T}} - @test leftinverse(m6)([2,2]) == [1,1] - - # Test construction and conversion - @test LinearMap{T}(1) isa ScalarLinearMap{T} - @test LinearMap{Vector{T}}(rand(Int,5,5)) isa VectorLinearMap{T} - @test LinearMap{SVector{2,T}}(SMatrix{2,2}(1, 2, 3, 4)) isa StaticLinearMap{T,2,2,4} - - @test GenericLinearMap(one(T)) isa GenericLinearMap{T,T} - @test GenericLinearMap(SMatrix{2,2,T}(1,1,1,1)) isa GenericLinearMap{SVector{2,T},SMatrix{2,2,T,4}} - - @test VectorLinearMap(rand(T,3,3)) isa VectorLinearMap{T} - @test StaticLinearMap(SMatrix{2,2,T}(1,2,3,4)) isa StaticLinearMap{T,2,2,4} - - @test convert(Map{SVector{2,T}}, LinearMap(rand(T,2,2))) isa StaticLinearMap{T,2,2,4} - @test convert(Map{T}, 1) isa LinearMap{T} -end - - - -function test_translation(T) - v = randvec(T,3) - m = Translation(v) - @test !islinearmap(m) - @test isaffinemap(m) - @test inverse(inverse(m)) == m - @test jacobian(m) isa ConstantMap - @test affinevector(jacobian(m)) == [1 0 0; 0 1 0; 0 0 1] - @test jacdet(m) isa UnityMap{SVector{3,T},T} - - @test mapsize(Translation(one(T))) == () - - # Test construction and conversion - @test Translation(one(T)) isa ScalarTranslation{T} - @test Translation{T}(1) isa ScalarTranslation{T} - @test Translation(SVector{2,T}(1,2)) isa StaticTranslation{T,2} - @test Translation{SVector{2,T}}(rand(T,2)) isa StaticTranslation{T,2} - @test Translation(rand(Int,2)) isa VectorTranslation{Int} - @test Translation{Vector{T}}(rand(Int,2)) isa VectorTranslation{T} - @test Translation(1:10) isa GenericTranslation{Vector{Int}} - @test Translation{Vector{T}}(1:10) isa GenericTranslation{Vector{T}} - - @test StaticTranslation(MVector(one(T),2)) isa StaticTranslation{T,2} - @test StaticTranslation{T}(SVector(1,2)) isa StaticTranslation{T,2} - - @test VectorTranslation(1:10) isa VectorTranslation{Int} - @test VectorTranslation{T}(1:10) isa VectorTranslation{T} - - @test GenericTranslation(1:10) isa GenericTranslation{Vector{Int}} - @test GenericTranslation{Vector{T}}(1:10) isa GenericTranslation{Vector{T}} - - @test GenericTranslation(one(T)) isa GenericTranslation{T,T} - @test GenericTranslation{T}(1) isa GenericTranslation{T,T} - @test_throws InexactError GenericTranslation{Int}(one(T)/2) - - @test convert(Map{SVector{2,T}}, Translation(rand(T,2))) isa StaticTranslation{T,2} -end - - -function test_affinemap(T) - m1 = AffineMap(T(2), T(3)) - @test m1 isa ScalarAffineMap{T} - @test !islinearmap(m1) - @test isaffinemap(m1) - @test m1(2) == 7 - - @test convert(ScalarAffineMap{BigFloat}, m1) == m1 - - @test m1 ∘ m1 isa AffineMap - @test (m1 ∘ m1)(2) == 2*(2*2+3)+3 - - m2 = AffineMap(T(2), SVector{2,T}(1,2)) - @test m2 isa GenericAffineMap{SVector{2,T}} - @test m2(SVector(1,2)) == SVector(3,6) - @test mapsize(m2) == (2,2) - - m3 = AffineMap(UniformScaling(2*one(T)), [one(T),2*one(T)]) - @test m3 isa AffineMap{Vector{T}} - @test mapsize(m3) == (2,2) - @test m3([1,2]) == 2 * [1,2] + [1,2] - y = zeros(T,2) - @test (FunctionMaps.applymap!(y, m3, [1,2]); y == m3([1,2])) - @test jacobian(m3, [1,2]) == [2 0; 0 2] - @test jacdet(m3, [1,2]) == 4 - - @test AffineMap(LinearAlgebra.I, 2.0) isa ScalarAffineMap - m4 = GenericAffineMap(LinearAlgebra.I, 2.0) - @test m4.A isa UniformScaling{Float64} - @test jacdet(m4, 3.0) == 1.0 - - @test mapsize(AffineMap(rand(3),rand(3))) == (3,) - @test mapsize(AffineMap(LinearAlgebra.I,2)) == () - - # Test construction and conversion - @test AffineMap(1, 2*one(T)) isa ScalarAffineMap{T} - @test AffineMap{BigFloat}(1, 2.0) isa ScalarAffineMap{BigFloat} - @test AffineMap(rand(T,2),rand(T,2)) isa GenericAffineMap{T} - @test AffineMap(rand(T,2,2),rand(T,2)) isa VectorAffineMap{T} - @test AffineMap{Vector{T}}(rand(Int,2,2),rand(Int,2)) isa VectorAffineMap{T} - @test AffineMap{SVector{2,T}}(rand(T,2,2),rand(T,2)) isa StaticAffineMap{T,2,2,4} - @test AffineMap{SVector{2,T}}(rand(Int,2,2),rand(T,2)) isa StaticAffineMap{T,2,2,4} - @test AffineMap{SVector{2,T}}(rand(Int,2,2),rand(Int,2)) isa StaticAffineMap{T,2,2,4} - @test AffineMap{SVector{2,T}}(1,rand(Int,2)) isa GenericAffineMap{SVector{2,T},T,SVector{2,T}} - @test AffineMap{SVector{2,T}}(SMatrix{3,2}(3,2,1,4,5,6),SVector(1,2,3)) isa StaticAffineMap{T,2,3,6} - - @test GenericAffineMap(SMatrix{2,2,T}(1,2,3,4), [1,2]) isa GenericAffineMap{SVector{2,T},SMatrix{2,2,T,4}} - @test GenericAffineMap(SMatrix{1,1,Int}(1),SVector{1,T}(1)) isa GenericAffineMap{SVector{1,T}} - @test GenericAffineMap(rand(T,2,2),rand(Int,2)) isa GenericAffineMap{Vector{T},Matrix{T},Vector{T}} - @test GenericAffineMap{Vector{T}}(1, [1,2]) isa GenericAffineMap{Vector{T},T,Vector{Int}} - @test GenericAffineMap{Vector{T}}(rand(Int,2,2), rand(Int,2)) isa GenericAffineMap{Vector{T},Matrix{T},Vector{T}} - @test convert(GenericAffineMap{Vector{T}}, GenericAffineMap(rand(Int,2,2), rand(Int,2))) isa GenericAffineMap{Vector{T}} - @test GenericAffineMap(rand(T,2),rand(T,2)) isa GenericAffineMap{T} - @test GenericAffineMap{T}(rand(Int,2),rand(Int,2)) isa GenericAffineMap{T} - - @test VectorAffineMap(rand(T,2,2),rand(T,2)) isa VectorAffineMap{T} - - @test StaticAffineMap(SMatrix{1,1,T}(1),SVector{1,T}(1)) isa StaticAffineMap{T} - @test StaticAffineMap(SMatrix{1,1,T}(1),SVector{1,Int}(1)) isa StaticAffineMap{T} - - @test StaticAffineMap(rand(T,2,2),SVector{2,T}(1,1)) isa StaticAffineMap{T} - @test_throws AssertionError StaticAffineMap(rand(T,3,2),SVector{2,T}(1,1)) - @test StaticAffineMap(SMatrix{2,2,T}(1,2,3,4),[1,1]) isa StaticAffineMap{T,2} - @test StaticAffineMap(SMatrix{3,2,T}(1,2,3,4,5,6),SVector{3,T}(1,2,3)) isa StaticAffineMap{T,2,3,6} - - @test convert(Map{SVector{2,T}}, AffineMap(rand(2,2),rand(2))) isa StaticAffineMap{T,2,2,4} -end - -function test_abstractarraymap(T) - A = rand(T, 3, 3) - @test FunctionMaps.MapStyle(A) isa FunctionMaps.IsMap - @test FunctionMaps.checkmap(A) == A - @test Map(A) isa LinearMap - @test Map(A) isa VectorLinearMap - @test Map(SA[1 2; 3 4]) isa StaticLinearMap - @test Map(Diagonal(rand(3))) isa GenericLinearMap - @test mapsize(A) == size(A) - @test applymap(A, 1:3) == A*(1:3) - @test islinearmap(A) - @test isaffinemap(A) - @test affinematrix(A) == A - @test affinevector(A) == [0,0,0] - @test jacobian(A) isa ConstantMap - @test jacobian(A, 1:3) == A - @test inverse(A) ≈ inv(A) - @test inverse(A, 1:3) ≈ A \ (1:3) - @test canonicalmap(A) isa LinearMap - @test FunctionMaps.equalmap(A) isa LinearMap - @test isequalmap(A, LinearMap(A)) -end diff --git a/FunctionMaps/test/test_arithmetics.jl b/FunctionMaps/test/test_arithmetics.jl deleted file mode 100644 index d1503ac..0000000 --- a/FunctionMaps/test/test_arithmetics.jl +++ /dev/null @@ -1,22 +0,0 @@ -@testset "map interval" begin - @test interval_map(1.0, Inf, 2.0, Inf) == AffineMap(1.0, 1.0) - @test interval_map(1.0, Inf, Inf, 2.0) == AffineMap(-1.0, 3.0) - @test interval_map(Inf, 1.0, 2.0, Inf) == AffineMap(-1.0, 3.0) - @test interval_map(Inf, 1.0, Inf, 2.0) == AffineMap(1.0, 1.0) - @test interval_map(-Inf, Inf, -Inf, Inf) == IdentityMap() - @test interval_map(-Inf, Inf, Inf, -Inf) == LinearMap(-1) - @test interval_map(Inf, -Inf, -Inf, Inf) == LinearMap(-1) - @test interval_map(Inf, -Inf, Inf, -Inf) == IdentityMap() - @test interval_map(Inf, Inf, Inf, Inf) == IdentityMap() - @test interval_map(-Inf, -Inf, -Inf, -Inf) == IdentityMap() - @test interval_map(Inf, Inf, -Inf, -Inf) == LinearMap(-1) - @test interval_map(-Inf, -Inf, Inf, Inf) == LinearMap(-1) - @test_throws ArgumentError interval_map(-Inf, Inf, -Inf, -Inf) - @test_throws ArgumentError interval_map(-1, 1, -1, Inf) - @test_throws ArgumentError interval_map(-1, 1, -1, Inf) - - a = 1.0; b = 2.0; c = 3.0; d = 5; - @test interval_map(a, b, c, d)(a) == c - @test interval_map(a, b, c, d)(b) == d - @test interval_map(a, b, c, d) == bounded_interval_map(a, b, c, d) -end diff --git a/FunctionMaps/test/test_basic.jl b/FunctionMaps/test/test_basic.jl deleted file mode 100644 index a24c35f..0000000 --- a/FunctionMaps/test/test_basic.jl +++ /dev/null @@ -1,114 +0,0 @@ -function test_isomorphisms(T) - m1 = FunctionMaps.VectorToNumber{T}() - @test m1(SA[1.0]) == 1.0 - m1b = FunctionMaps.NumberToVector{T}() - @test inverse(m1) == m1b - @test inverse(m1, 0.4) == m1b(0.4) - @test inverse(m1b) == m1 - @test inverse(m1b, SA[0.4]) == m1(SA[0.4]) - @test m1b(1.0) == SA[1.0] - @test mapsize(m1) == (1,1) - @test mapsize(m1b) == (1,) - @test mapsize(m1 ∘ inverse(m1)) == () - @test mapsize(inverse(m1) ∘ m1) == (1,1) - @test jacobian(m1 ∘ inverse(m1), 1) == one(T) - @test jacobian(inverse(m1) ∘ m1, [1]) == ones(T,1,1) - @test jacobian(m1) isa ConstantMap - @test jacobian(m1b) isa ConstantMap - - m2 = FunctionMaps.VectorToComplex{T}() - @test m2(SA[one(T), one(T)]) == 1 + im - m2b = FunctionMaps.ComplexToVector{T}() - @test inverse(m2) == m2b - @test inverse(m2, zero(T)+0im) == m2b(zero(T)+0im) - @test inverse(m2b) == m2 - @test inverse(m2b, SA[one(T),one(T)]) == one(T)+one(T)*im - @test m2b(one(T)+one(T)*im) == SA[one(T),one(T)] - @test mapsize(m2) == (1,2) - @test mapsize(m2b) == (2,) - @test mapsize(m2 ∘ m2b) == () - @test mapsize(m2b ∘ m2) == (2,2) - - m3 = FunctionMaps.VectorToTuple{2,T}() - m3b = FunctionMaps.TupleToVector{2,T}() - @test m3(SA[one(T), one(T)]) == (one(T),one(T)) - @test inverse(m3b, SA[one(T), one(T)]) == (one(T),one(T)) - @test inverse(m3) == m3b - @test inverse(m3b) == m3 - @test m3b( (one(T),one(T)) ) == SA[one(T),one(T)] - @test inverse(m3, (one(T),one(T)) ) == SA[one(T),one(T)] - - m4 = FunctionMaps.NestedToFlat{3,T,Tuple{Tuple{T,T},T},(2,1)}() - m4b = FunctionMaps.FlatToNested{3,T,Tuple{Tuple{T,T},T},(2,1)}() - x4 = ([T(1),T(2)],T(3)) - y4 = T[1,2,3] - @test m4(x4) == y4 - @test m4b(y4) == x4 - @test inverse(m4) == m4b - @test inverse(m4b) == m4 - @test inverse(m4, y4) == x4 - @test inverse(m4b, x4) == y4 -end - -function test_scaling_maps(T) - test_generic_map(LinearMap(T(2))) - test_generic_map(LinearMap(T(2), T(3))) - test_generic_map(LinearMap(T(2), T(3), T(4))) - test_generic_map(LinearMap(T(2), T(3), T(4), T(5))) -end - -function test_identity_map(T) - i1 = StaticIdentityMap{T}() - i2 = StaticIdentityMap{SVector{2,T}}() - test_generic_map(i1) - test_generic_map(i2) - @test i1 == i2 - @test hash(i1) == hash(i2) - @test islinearmap(i1) - @test isaffinemap(i1) - @test convert(StaticIdentityMap{SVector{2,T}}, i1) === i2 - @test jacobian(i1) isa ConstantMap - @test jacobian(i1, 1) == 1 - @test jacdet(i1, 1) == 1 - m1 = convert(FunctionMaps.AbstractAffineMap{T}, i1) - @test m1 isa LinearMap{T} - @test jacdet(m1, 1) == 1 - @test convert(FunctionMaps.AbstractAffineMap, i1) isa LinearMap{T} - m2 = convert(FunctionMaps.LinearMap{T}, i1) - @test m2 isa LinearMap{T} - @test jacdet(m2, 1) == 1 - - i3 = VectorIdentityMap{T}(10) - test_generic_map(i3) - r = rand(T, 10) - @test i3(r) ≈ r - @test hash(i3) == hash(VectorIdentityMap{Int}(10)) - - @test IdentityMap() ∘ LinearMap(2) == LinearMap(2.0) -end - -function test_basic_maps(T) - @test UnityMap{SVector{2,Float64}}() == UnityMap{SVector{2,Float64},Float64}() - @test hash(ConstantMap(2)) == hash(ConstantMap(2.0)) - @test FunctionMaps.absmap(ConstantMap(-2)) == ConstantMap(2) -end - -function test_wrapped_maps(T) - m1 = WrappedMap{T}(cos) - m2 = WrappedMap{T}(sin) - @test m1(one(T)) ≈ cos(one(T)) - @test m2(one(T)) ≈ sin(one(T)) - @test isequalmap(m1, cos) - @test isequalmap(sin, m2) - m3 = m1 ∘ m2 - @test m3(one(T)) ≈ cos(sin(one(T))) - - @test WrappedMap(cos) isa WrappedMap{Float64} - @test convert(Map, cos) isa WrappedMap - @test convert(Map{BigFloat}, m1) isa WrappedMap{BigFloat} - - @test convert(Map{T}, cos) isa FunctionMaps.WrappedMap{T,typeof(cos)} - - @test convert(Map, LinearAlgebra.I) isa GenericLinearMap{Vector{Any}} - @test convert(Map{Vector{T}}, LinearAlgebra.I) isa GenericLinearMap{Vector{T}} -end diff --git a/FunctionMaps/test/test_common.jl b/FunctionMaps/test/test_common.jl deleted file mode 100644 index b573fee..0000000 --- a/FunctionMaps/test/test_common.jl +++ /dev/null @@ -1,97 +0,0 @@ -function test_dimension() - @test euclideandimension(Int) == 1 - @test euclideandimension(Float64) == 1 - @test euclideandimension(ComplexF64) == 1 - @test euclideandimension(SVector{2,Float64}) == 2 - @test euclideandimension(MVector{2,Float64}) == 2 - @test euclideandimension(Tuple{Int,Int}) == 2 - @test euclideandimension(Tuple{Int,Float64}) == 2 - @test_throws ArgumentError euclideandimension(Vector{Float64}) -end - -function test_prectype() - @test prectype(1.0) == Float64 - @test prectype(big(1.0)) == BigFloat - @test prectype(1) == typeof(float(1)) - @test prectype(SVector(1,2)) == typeof(float(1)) - @test prectype(SVector(1,big(2))) == typeof(float(big(2))) - @test prectype(1.0+2.0im) == Float64 - @test prectype([1.0+2.0im, 3.0]) == Float64 - @test prectype(NTuple{2,Int}) == Float64 - @test prectype(typeof((1.0,))) == Float64 - @test prectype(typeof((1.0, 2.0))) == Float64 - @test prectype(typeof((1.0, 2.0, 3.0))) == Float64 - @test prectype(typeof((1.0, big(2.0), 3.0+im))) == BigFloat - @test prectype(NTuple{4,Int}) == Float64 - @test @inferred(prectype(1, 2.0)) == Float64 - @test @inferred(prectype(typeof((1, 2.0, 3, 40+im)))) == Float64 - - @test convert_prectype(Float64, 2) == 2 - @test convert_prectype(Float64, 2) isa Float64 - @test convert_prectype(BigFloat, 1.0+im) == 1+im - @test convert_prectype(BigFloat, 1.0+im) isa Complex{BigFloat} - @test convert_prectype(Float64, SA[1,2]) == SA[1.0,2.0] - @test convert_prectype(Float64, SA[1,2]) isa SVector{2,Float64} - @test convert_prectype(Float64, MVector(1,2)) == MVector(1.0,2.0) - @test convert_prectype(Float64, MVector(1,2)) isa MVector - @test convert_prectype(Float64, SA[1 2; 3 4]) == SA[1.0 2.0; 3.0 4.0] - @test convert_prectype(Float64, SA[1 2; 3 4]) isa SMatrix{2,2,Float64} - @test convert_prectype(Float64, MMatrix{2,2}(1, 2, 3, 4)) == SA[1.0 3.0; 2.0 4.0] - @test convert_prectype(Float64, MMatrix{2,2}(1, 2, 3, 4)) isa MMatrix{2,2,Float64} - @test convert_prectype(BigFloat, SA[1,2+im]) isa SVector{2,Complex{BigFloat}} - @test_throws ArgumentError convert_prectype(BigFloat, "a") - - @test promote_prectype(2) == 2 - @test promote_prectype(2, 3.0) isa Tuple{Float64,Float64} - @test promote_prectype(2, 3.0+im, big(4)) isa Tuple{BigFloat,Complex{BigFloat},BigFloat} -end - - -function test_numtype() - @test numtype(1.0) == Float64 - @test numtype(big(1.0)) == BigFloat - @test numtype(1) == Int - @test numtype([1,2,3], [4,5,6]) == Int - @test numtype(SVector(1,2)) == Int - @test numtype(SVector(1,big(2))) == BigInt - @test numtype(Array{Float64,2}) == Float64 - @test numtype(1.0+2.0im) == Complex{Float64} - @test numtype([1.0+2.0im, 3.0]) == Complex{Float64} - @test numtype(NTuple{2,Complex{Int}}) == Complex{Int} - @test numtype(typeof((1.0,))) == Float64 - @test numtype(typeof((1.0, 2.0))) == Float64 - @test numtype((1.0, 2.0, 3.0)) == Float64 - @test numtype((1.0, 2.0, 3.0, 4.0)) == Float64 - @test numtype(1.0, big(2.0), 3.0+im) == Complex{BigFloat} - @test numtype(typeof((1.0, big(2.0), 3.0+im))) == Complex{BigFloat} - @test @inferred(numtype(1, 2.0)) == Float64 - @test @inferred(numtype(typeof((1, 2.0, 3, 40+im)))) == Complex{Float64} - - @test convert_numtype(Float64, 2) == 2 - @test convert_numtype(Float64, 2) isa Float64 - @test convert_numtype(Float64, SA[1,2]) == SA[1.0,2.0] - @test convert_numtype(Float64, SA[1,2]) isa SVector{2,Float64} - @test convert_numtype(Float64, MVector(1,2)) == MVector(1.0,2.0) - @test convert_numtype(Float64, MVector(1,2)) isa MVector - @test convert_numtype(Float64, SA[1 2; 3 4]) == SA[1.0 2.0; 3.0 4.0] - @test convert_numtype(Float64, SA[1 2; 3 4]) isa SMatrix{2,2,Float64} - @test convert_numtype(Float64, MMatrix{2,2}(1, 2, 3, 4)) == SA[1.0 3.0; 2.0 4.0] - @test convert_numtype(Float64, MMatrix{2,2}(1, 2, 3, 4)) isa MMatrix{2,2,Float64} - @test_throws ArgumentError convert_numtype(BigFloat, "a") - - @test promote_numtype(2) == 2 - @test promote_numtype(2, 3.0) isa Tuple{Float64,Float64} - @test promote_numtype(2, 3.0+im, big(4)) isa Tuple{Complex{BigFloat},Complex{BigFloat},Complex{BigFloat}} -end - -@testset "common functionality" begin - @testset "dimension" begin - test_dimension() - end - @testset "prectype" begin - test_prectype() - end - @testset "numtype" begin - test_numtype() - end -end diff --git a/FunctionMaps/test/test_generic.jl b/FunctionMaps/test/test_generic.jl deleted file mode 100644 index 2a62783..0000000 --- a/FunctionMaps/test/test_generic.jl +++ /dev/null @@ -1,145 +0,0 @@ -function generic_map_tests(T) - for map in maps_to_test(T) - @test prectype(map) == T - test_generic_map(map) - end -end - -# generic map test suite -function test_generic_map(m) - @test convert(Map{domaintype(m)}, m) == m - - x = suitable_point_to_map(m) - @test applymap(m, x) == m(x) - y = m(x) - @test isreal(m) == (isreal(x) && isreal(y)) - - S = domaintype(m) - U = codomaintype(m) - @test x isa S - @test y isa U - - if isaffinemap(m) && !isconstantmap(m) && !(prectype(m) == BigFloat) - test_generic_inverse(m) - else - try - # The map may not support an inverse, let's try - test_generic_inverse(m) - catch e - # make sure it failed because the inverse is not defined - @test e isa MethodError || prectype(m) == BigFloat - end - end - - if isaffinemap(m) - M = affinematrix(m) - @test size(M) == mapsize(m) - test_generic_jacobian(m) - else - try # jacobian may not be implemented - test_generic_jacobian(m) - catch e - end - end - - if domaintype(m) == Float64 - @test convert(Map{BigFloat}, m) isa Map{BigFloat} - @test convert(Map{BigFloat}, m) == m - end - if prectype(m) == Float64 - U = widertype(domaintype(m)) - @test convert(Map{U}, m) isa Map{U} - @test convert(Map{U}, m) == m - end - - if isaffinemap(m) - A = affinematrix(m) - b = affinevector(m) - @test m(x) == A*x+b - end - - @test hash(m) == map_hash(m) -end - -function test_generic_inverse(m) - M = mapsize(m,1) - N = mapsize(m,2) - x = suitable_point_to_map(m) - y = m(x) - - if M==N - inverse(m, y) # trigger exception outside test if not implemented - - minv = inverse(m) - @test minv(y) ≈ x - @test inverse(m)(y) ≈ x - @test inverse(m, y) ≈ x - @test m\y ≈ x - @test LazyInverse(m)(y) ≈ inverse(m, y) - end - if M >= N && numtype(m)!=BigFloat - leftinverse(m, y) # trigger exception outside test if not implemented - - mli = leftinverse(m) - @test mli(y) ≈ x - @test leftinverse(m, y) ≈ x - end - if M <= N && numtype(m)!=BigFloat - rightinverse(m, y) # trigger exception outside test if not implemented - - mri = rightinverse(m) - @test m(mri(y)) ≈ y - @test m(rightinverse(m, y)) ≈ y - end -end - -function test_generic_jacobian(m) - x = suitable_point_to_map(m) - # Trigger exception outside of test if jacobian(m, x) is not implemented. - # test_generic_jacobian(m) should be called within a try/catch block - jacobian(m, x) - - δ = sqrt(eps(prectype(m))) - x2 = x .+ δ - # we intentionally test jacobian(m, x) before testing jacobian(m) - if !(m isa ProductMap) - @test norm(m(x2) - m(x) + jacobian(m, x)*(x-x2)) < 100δ - end - jac = jacobian(m) - @test jac(x) ≈ jacobian(m, x) - j = jacobian(m, x) - @test size(j) == mapsize(m) - if j isa AbstractArray{Float64} - y = similar(j) - jacobian!(y, m, x) - @test y ≈ jac(x) - end - if issquarematrix(jac(x)) - @test jacdet(m, x) ≈ det(jacobian(m, x)) - end - if mapsize(m) == () - d = diffvolume(m, x) - s = sqrt(det(transpose(jacobian(m,x))*jacobian(m,x))) - @test d ≈ s || d ≈ -s - else - @test abs(diffvolume(m, x) - sqrt(det(transpose(jacobian(m,x))*jacobian(m,x)))) < 1e-7 - end - if isaffinemap(m) - @test diffvolume(m) isa ConstantMap - end -end - -function test_generic_functionality() - m = LinearMap(2) - @test convert_domaintype(Float64, m) isa ScalarLinearMap{Float64} - @test convert_domaintype(Complex{Float64}, m) isa ScalarLinearMap{Complex{Float64}} - @test convert_domaintype(SVector{2,BigFloat}, m) isa GenericLinearMap{SVector{2, BigFloat}, BigFloat} - @test convert_domaintype(Any, m) === m - @test convert_codomaintype(Float64, m) isa ScalarLinearMap{Float64} - @test convert_codomaintype(Complex{Float64}, m) isa ScalarLinearMap{Complex{Float64}} - @test_throws ArgumentError convert_codomaintype(SVector{2,BigFloat}, m) - - m2 = LinearMap{SVector{2,Float64}}(2) - @test convert_domaintype(SVector{2,BigFloat}, m2) isa GenericLinearMap{SVector{2, BigFloat}, BigFloat} - @test convert_codomaintype(SVector{2,BigFloat}, m2) isa GenericLinearMap{SVector{2, BigFloat}, BigFloat} -end diff --git a/FunctionMaps/test/test_interface.jl b/FunctionMaps/test/test_interface.jl deleted file mode 100644 index 6920ee2..0000000 --- a/FunctionMaps/test/test_interface.jl +++ /dev/null @@ -1,16 +0,0 @@ - -struct MySimpleMap end -FunctionMaps.MapStyle(::Type{MySimpleMap}) = IsMap() - -@testset "map interface" begin - @test MapStyle(2.0) isa NotMap - @test MapStyle(2.0) == MapStyle(typeof(2.0)) - @test_throws ArgumentError checkmap(2.0) - m = LinearMap(2) - @test MapStyle(m) isa IsMap - @test functionmap(m) === m - @test checkmap(m) === m - m2 = MySimpleMap() - @test MapStyle(m2) isa IsMap - @test checkmap(m2) == m2 -end diff --git a/FunctionMaps/test/test_maps.jl b/FunctionMaps/test/test_maps.jl deleted file mode 100644 index 264be85..0000000 --- a/FunctionMaps/test/test_maps.jl +++ /dev/null @@ -1,191 +0,0 @@ - -maps_to_test(T) = [ - StaticIdentityMap{T}(), - VectorIdentityMap{T}(10), - ConstantMap{T}(one(T)), - ConstantMap{T}(SVector{2,T}(1,2)), - ZeroMap{T}(), - UnityMap{T}(), - AffineMap(T(1.2), T(2.4)), # scalar map - AffineMap(-T(1.2), T(2.4)), # scalar map with negative A - AffineMap(randvec(T, 2, 2), randvec(T, 2)), # static map - AffineMap(randvec(T, 3, 2), randvec(T, 3)), # static map, rectangular - AffineMap(rand(T, 2, 2), rand(T, 2)), # vector map - AffineMap(rand(T, 3, 2), rand(T, 3)), # vector map, rectangular - AffineMap(LinearAlgebra.I, one(T)/2), # use UniformScaling object as A - AffineMap(LinearAlgebra.I, randvec(T,2)), # use UniformScaling object as A - GenericAffineMap(randvec(T, 2, 2), randvec(T, 2)), - GenericAffineMap(T(1.2), randvec(T, 2)), - GenericAffineMap(randvec(T, 3, 2), randvec(T, 3)), - Translation(randvec(T, 3)), - LinearMap(randvec(T, 2, 2)), - LinearMap(randvec(T, 2)), - LinearMap(randvec(T, 2, 2)) ∘ AffineMap(T(1.2), randvec(T, 2)), - AffineMap(5.0, 2.0) ∘ VectorToComplex{T}() ∘ UnitCircleMap{T}(), - LinearMap(SMatrix{2,2}(1,2,3,T(4))) ∘ CartToPolarMap{T}() ∘ LinearMap(SMatrix{2,2}(1,2,3,T(4))), - # Interval{Any}(0.0, 1.0) -] - -randvec(T,n) = SVector{n,T}(rand(n)) -randvec(T,m,n) = SMatrix{m,n,T}(rand(m,n)) - -suitable_point_to_map(m::Map) = suitable_point_to_map(m, domaintype(m)) - -suitable_point_to_map(m::Map, ::Type{SVector{N,T}}) where {N,T} = SVector{N,T}(rand(N)) -suitable_point_to_map(m::Map, ::Type{T}) where {T<:Number} = rand(T) -suitable_point_to_map(m::Map, ::Type{<:AbstractVector{T}}) where {T} = rand(T, mapsize(m,2)) - -suitable_point_to_map(m::ProductMap) = - map(suitable_point_to_map, components(m)) -suitable_point_to_map(m::VcatMap{T,M,N}) where {T,M,N} = - SVector{N,T}(rand(T,N)) - -suitable_point_to_map(::CartToPolarMap{T}) where {T} = randvec(T,2) -suitable_point_to_map(::PolarToCartMap{T}) where {T} = randvec(T,2) - -widertype(T) = widen(T) -widertype(::Type{SVector{N,T}}) where {N,T} = SVector{N,widen(T)} -widertype(::Type{Vector{T}}) where {T} = Vector{widen(T)} -widertype(::Type{Tuple{A}}) where {A} = Tuple{widen(A)} -widertype(::Type{Tuple{A,B}}) where {A,B} = Tuple{widen(A),widen(B)} -widertype(::Type{Tuple{A,B,C}}) where {A,B,C} = Tuple{widen(A),widen(B),widen(C)} -widertype(::Type{NTuple{N,T}}) where {N,T} = NTuple{N,widen(T)} - -issquarematrix(A) = false -issquarematrix(A::AbstractArray) = size(A,1)==size(A,2) - - -function test_maps() - @testset "generic functionality" begin - test_generic_functionality() - end - - @testset "generic map tests" begin - generic_map_tests(Float64) - generic_map_tests(BigFloat) - end - - # Test special maps - @testset "identity map" begin - test_identity_map(Float64) - test_identity_map(BigFloat) - end - @testset "basic maps" begin - test_basic_maps(Float64) - test_basic_maps(BigFloat) - end - @testset "affine maps" begin - test_affine_maps(Float64) - test_affine_maps(BigFloat) - end - @testset "composite maps" begin - test_composite_map(Float64) - test_composite_map(BigFloat) - end - @testset "product maps" begin - test_product_map(Float64) - test_product_map(BigFloat) - end - @testset "wrapped maps" begin - test_wrapped_maps(Float64) - test_wrapped_maps(BigFloat) - end - @testset "scaling maps" begin - test_scaling_maps(Float64) - test_scaling_maps(BigFloat) - end - @testset "isomorphisms" begin - test_isomorphisms(Float64) - test_isomorphisms(BigFloat) - end - @testset "Mixed maps" begin - test_mixed_maps() - end -end - -function test_composite_map(T) - a = T(0) - b = T(1) - c = T(2) - d = T(3) - ma = StaticIdentityMap{T}() - mb = interval_map(a, b, c, d) - - r = suitable_point_to_map(ma) - m1 = ma∘mb - test_generic_map(m1) - @test m1(r) ≈ ma(mb(r)) - m2 = m1∘mb - test_generic_map(m2) - @test m2(r) ≈ m1(mb(r)) - m3 = mb∘m2 - test_generic_map(m3) - @test m3(r) ≈ mb(m2(r)) - m = m2∘m3 - test_generic_map(m) - @test m(r) ≈ m2(m3(r)) - - m5 = ComposedMap(LinearMap(rand(T,2,2)), AffineMap(rand(T,2,2),rand(T,2))) - test_generic_map(m5) - @test jacobian(m5) isa ConstantMap - @test m5[Component(1)] isa LinearMap - @test m5[Component(2)] isa AffineMap - @test ComposedMap(m5[Component(1:2)]...) == m5 - @test_throws BoundsError m5[Component(3)] - - m6 = multiply_map(ma,ma) - @test m6(one(T)/2) == one(T)/4 - @test jacobian(m6) isa SumMap - @test jacobian(m6)(one(T)) == 2 - @test jacobian(m6, one(T)) == 2 - - m7 = sum_map(ma,ma) - @test m7(one(T)) == 2 - @test jacobian(m7) isa ConstantMap - @test jacobian(m7, one(T)) == 2 - - @test mapsize(ComposedMap(LinearMap(2),LinearMap(rand(T,2)),LinearMap(rand(T,2,2)))) == (2,) - @test mapsize(ComposedMap(LinearMap(rand(T,2)),LinearMap(rand(T,2)'),LinearMap(rand(T,2)))) == (2,) - @test mapsize(ComposedMap(LinearMap(rand(T,2,2)),LinearMap(rand(T,2)'),LinearMap(2))) == (1,2) - @test mapsize(ComposedMap(LinearMap(one(T)),LinearMap(one(T)))) == () - - @test composedmap() == () - @test composedmap(ma) == ma - @test composedmap(ma,ma) == ma - @test composedmap(ma,ma,ma) == ma - - @test composite_jacobian(ma) == jacobian(ma) - - @test multiply_map() == () - @test multiply_map(ma) == ma - @test multiply_map(ma,ma)(2*one(T)) == 4 - @test multiply_map(ma,ma,ma)(2*one(T)) == 8 - - @test sum_jacobian() == () - @test sum_jacobian(ma) == jacobian(ma) -end - -function test_mixed_maps() - m1 = composedmap(cos, sin) - @test domaintype(m1) == Any - @test m1(0.4) == sin(cos(0.4)) - - m2 = composedmap(AffineMap(2.0, 3.0), cos) - @test domaintype(m2) == Float64 - @test m2(0.4) ≈ cos(2*0.4+3) - @inferred m2(0.4) - @test repr(m2) == "cos ∘ (x -> 2.0 * x + 3.0)" - - m3 = composedmap(sin, AffineMap(2.0, 3.0), cos) - @test m3(0.5) ≈ cos(2*sin(0.5)+3) - @test repr(m3) == "cos ∘ (x -> 2.0 * x + 3.0) ∘ sin" - @test domaintype(m3) == Any - @inferred m3(0.5) - - m4 = productmap(sin, cos) - @test m4 isa TupleProductMap - @test domaintype(m4) == Tuple{Any,Any} - @test m4(0.3,0.5) == (sin(0.3), cos(0.5)) -end - -test_maps() diff --git a/FunctionMaps/test/test_product.jl b/FunctionMaps/test/test_product.jl deleted file mode 100644 index db053c0..0000000 --- a/FunctionMaps/test/test_product.jl +++ /dev/null @@ -1,48 +0,0 @@ -function test_product_map(T) - ma = StaticIdentityMap{T}() - mb = interval_map(T(0), T(1), T(2), T(3)) - - r1 = suitable_point_to_map(ma) - r2 = suitable_point_to_map(ma) - r3 = suitable_point_to_map(ma) - r4 = suitable_point_to_map(ma) - r5 = suitable_point_to_map(ma) - - m1 = productmap(ma,mb) - test_generic_map(m1) - @test m1(SVector(r1,r2)) ≈ SVector(ma(r1),mb(r2)) - m2 = productmap(m1,mb) - test_generic_map(m2) - @test m2(SVector(r1,r2,r3)) ≈ SVector(ma(r1),mb(r2),mb(r3)) - m3 = productmap(mb,m2) - test_generic_map(m3) - @test m3(SVector(r1,r2,r3,r4)) ≈ SVector(mb(r1),ma(r2),mb(r3),mb(r4)) - m4 = productmap(m1,m2) - test_generic_map(m4) - @test m4(SVector(r1,r2,r3,r4,r5)) ≈ SVector(m1(SVector(r1,r2))...,m2(SVector(r3,r4,r5))...) - - @test ProductMap(ma, mb) == ProductMap([ma,mb]) - @test hash(ProductMap(ma, mb)) == hash(ProductMap([ma,mb])) - - m5 = productmap(AffineMap(SMatrix{2,2,T}(1.0,2,3,4), SVector{2,T}(1,3)), LinearMap{T}(2.0)) - @test domaintype(component(m5,1)) == SVector{2,T} - @test domaintype(component(m5,2)) == T - test_generic_map(m5) - x = SVector{3,T}(rand(T,3)) - @test m5(x) ≈ SVector(component(m5,1)(SVector(x[1],x[2]))...,component(m5,2)(x[3])) - - m6 = ProductMap([ma,mb]) - @test m6 isa VectorProductMap - @test mapsize(m6) == (2,2) - @test convert(Map{SVector{2,T}}, m6) isa VcatMap - test_generic_map(m6) - - @test components(productmap(LinearMap(1), LinearMap(one(T)))) isa Tuple{<:Map{T},<:Map{T}} - - # test a non-rectangular product map - m_rect = AffineMap(SA[one(T) 2; 3 4; 5 6], SA[one(T),one(T),zero(T)]) - pm = productmap(m_rect, m_rect) - @test pm isa VcatMap{T,6,4,(3,3),(2,2)} - @test jacobian(pm, SA[one(T),0,0,0]) isa SMatrix{6,4,T} - @test diffvolume(pm, SA[one(T),0,0,0]) == 24 -end diff --git a/FunctionMaps/test/using_fmaps.jl b/FunctionMaps/test/using_fmaps.jl deleted file mode 100644 index 1058146..0000000 --- a/FunctionMaps/test/using_fmaps.jl +++ /dev/null @@ -1,53 +0,0 @@ -using FunctionMaps - -# for test_common.jl -using FunctionMaps: - convert_numtype, - convert_prectype, - promote_numtype, - promote_prectype, - euclideandimension - -# for test_interface.jl -using FunctionMaps: - MapStyle, IsMap, NotMap, - functionmap, checkmap - -# for test_generic.jl -using FunctionMaps: - convert_domaintype, convert_codomaintype, - map_hash, - LazyInverse, jacobian! - -# for test_maps.jl -using FunctionMaps: ScalarAffineMap, - VectorAffineMap, - StaticAffineMap, - GenericAffineMap, - ScalarLinearMap, - VectorLinearMap, - StaticLinearMap, - GenericLinearMap, - ScalarTranslation, - VectorTranslation, - StaticTranslation, - GenericTranslation, - ProductMap, - TupleProductMap, VcatMap, VectorProductMap, - WrappedMap, - interval_map, multiply_map, - SumMap, sum_map, - composedmap, ComposedMap, - composite_jacobian, sum_jacobian, - CartToPolarMap, PolarToCartMap, - UnitCircleMap, AngleMap, UnitDiskMap, - VectorToComplex - -# for test_affine.jl -using FunctionMaps: - to_matrix, to_vector, - matrix_pinv - -# for test_arithmetics.jl -using FunctionMaps: - interval_map, bounded_interval_map diff --git a/FunctionMaps/test/using_fmaps2.jl b/FunctionMaps/test/using_fmaps2.jl deleted file mode 100644 index 45ea364..0000000 --- a/FunctionMaps/test/using_fmaps2.jl +++ /dev/null @@ -1,53 +0,0 @@ -using DomainSets, DomainSets.FunctionMaps - -# for test_common.jl -using DomainSets.FunctionMaps: - convert_numtype, - convert_prectype, - promote_numtype, - promote_prectype, - euclideandimension - -# for test_interface.jl -using DomainSets.FunctionMaps: - MapStyle, IsMap, NotMap, - functionmap, checkmap - -# for test_generic.jl -using DomainSets.FunctionMaps: - convert_domaintype, convert_codomaintype, - map_hash, - LazyInverse, jacobian! - -# for test_maps.jl -using DomainSets.FunctionMaps: ScalarAffineMap, - VectorAffineMap, - StaticAffineMap, - GenericAffineMap, - ScalarLinearMap, - VectorLinearMap, - StaticLinearMap, - GenericLinearMap, - ScalarTranslation, - VectorTranslation, - StaticTranslation, - GenericTranslation, - ProductMap, - TupleProductMap, VcatMap, VectorProductMap, - WrappedMap, - interval_map, multiply_map, - SumMap, sum_map, - composedmap, ComposedMap, - composite_jacobian, sum_jacobian, - CartToPolarMap, PolarToCartMap, - UnitCircleMap, AngleMap, UnitDiskMap, - VectorToComplex - -# for test_affine.jl -using DomainSets.FunctionMaps: - to_matrix, to_vector, - matrix_pinv - -# for test_arithmetics.jl -using DomainSets.FunctionMaps: - interval_map, bounded_interval_map diff --git a/Project.toml b/Project.toml index 13b2638..6f93241 100644 --- a/Project.toml +++ b/Project.toml @@ -1,9 +1,10 @@ name = "DomainSets" uuid = "5b8099bc-c8ec-5219-889f-1d9e522a28bf" -version = "0.7.15" +version = "0.8.0" [deps] CompositeTypes = "b152e2b5-7a66-4b01-a709-34e65c35f657" +FunctionMaps = "a85aefff-f8ca-4649-a888-c8e5398bc76c" IntervalSets = "8197267c-284f-5f27-9208-e0e47529a953" LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e" Random = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c" @@ -18,6 +19,7 @@ DomainSetsMakieExt = "Makie" [compat] Aqua = "0.8" CompositeTypes = "0.1.4" +FunctionMaps = "0.1" IntervalSets = "0.7.4" LinearAlgebra = "1" Makie = "0.19, 0.20, 0.21, 0.22" diff --git a/src/DomainSets.jl b/src/DomainSets.jl index 557ddaf..2a35170 100644 --- a/src/DomainSets.jl +++ b/src/DomainSets.jl @@ -5,8 +5,7 @@ using IntervalSets using LinearAlgebra using Random using StaticArrays - -include("../FunctionMaps/src/FunctionMaps.jl") +using FunctionMaps ################################ ## Exhaustive list of imports @@ -17,14 +16,25 @@ import Base: # Set operations setdiff, in, isempty, issubset, intersect, union, &, \, # Arrays - eltype, hash, isreal, + eltype, hash, # Types, promotions and conversions convert, promote, - # for maps (both deprecated and to be removed) - size, inv, # Display show +import FunctionMaps: + convert_eltype, + convert_prectype, + convert_numtype, + prectype, + numtype, + factors, + isrealtype, + tointernalpoint, + toexternalpoint, + compatibleproductdims + + import CompositeTypes: component, components import IntervalSets: (..), endpoints, Domain, AbstractInterval, TypedEndpointsInterval, @@ -33,8 +43,27 @@ import IntervalSets: (..), endpoints, Domain, AbstractInterval, TypedEndpointsIn infimum, supremum export .. -import LinearAlgebra: cross, ×, pinv - +using FunctionMaps: + CanonicalType, + Equal, + AbstractAffineMap, + GenericLinearMap, + GenericAffineMap, + interval_map, + UnitCircleMap, + AngleMap, + UnitDiskMap, + StaticTypes, + hashrec, + euclideandimension, + convert_eltype, + promotable_eltypes, + promote_prectype, + promote_numtype, + convert_fromcartesian, + convert_tocartesian, + nfactors, + factor ################################ ## Exhaustive list of exports @@ -48,32 +77,6 @@ export prectype, numtype, convert_prectype, promote_prectype, iscomposite, component, components, ncomponents -## Maps - -include("maps_imports.jl") - -# from maps/map.jl -export AbstractMap, Map, MapRef, TypedMap, - applymap, isequalmap, - domaintype, codomaintype, - inverse, leftinverse, rightinverse, - mapsize, jacobian, jacdet, diffvolume -# from maps/composite.jl -export ComposedMap, composedmap -# from maps/product.jl -export ProductMap, productmap -# from maps/basic.jl -export IdentityMap, - StaticIdentityMap, VectorIdentityMap, - ZeroMap, UnityMap, ConstantMap, - isconstant, constant, - isconstantmap, mapconstant -# from maps/affine.jl -export AffineMap, Translation, LinearMap, - affinematrix, affinevector, - islinear, isaffine, - islinearmap, isaffinemap - ## Generic domains # from generic/domain.jl @@ -83,6 +86,7 @@ export Domain, EuclideanDomain, VectorDomain, isopenset, isclosedset, iscompact, boundary, ∂, interior, closure, + isrealdomain, choice # from generic/geometry.jl @@ -106,6 +110,7 @@ export superdomain # from generic/productdomain.jl export ProductDomain, productdomain, center, + cartesianproduct, VcatDomain, VectorProductDomain, TupleProductDomain, factors @@ -173,10 +178,6 @@ export LevelSet, ZeroSet, # from domain/indicator.jl export IndicatorFunction -## Applications -# from applications/rotation.jl -export rotation_map, - CartToPolarMap, PolarToCartMap include("util/common.jl") diff --git a/src/deprecated.jl b/src/deprecated.jl index 329a205..1d6bf85 100644 --- a/src/deprecated.jl +++ b/src/deprecated.jl @@ -1,38 +1,57 @@ -@deprecate convert_prectype(d::Domain, ::Type{T}) where {T} convert_prectype(T, d) -@deprecate convert_numtype(d::Domain, ::Type{T}) where {T} convert_numtype(T, d) - -@deprecate point_in_domain(d) choice(d) - -@deprecate broadcast_in(A, d::Domain) vectorized_in(A, d) -@deprecate broadcast_approx_in(A, d::Domain) vectorized_approx_in(A, d) -@deprecate broadcast_approx_in(A, d::Domain, tol) vectorized_approx_in(A, d, tol) - -@deprecate UnionDomain(domain::Domain) UnionDomain((domain,)) -@deprecate UnionDomain{T}(domain::Domain) where {T} UnionDomain{T}((domain,)) -@deprecate IntersectDomain(domain::Domain) IntersectDomain((domain,)) -@deprecate IntersectDomain{T}(domain::Domain) where T IntersectDomain((domain,)) - -@deprecate issubset(d1::AnyDomain, d2) issubset(d1, DomainRef(d2)) -@deprecate union(d1::AnyDomain, d2, domains...) union(d1, DomainRef(d2), domains...) -@deprecate union(d1, d2::AnyDomain, domains...) union(DomainRef(d1), d2, domains...) -@deprecate intersect(d1, d2::AnyDomain) intersect(DomainRef(d1), d2) -@deprecate intersect(d1::AnyDomain, d2) intersect(d1, DomainRef(d2)) -@deprecate \(d1::AnyDomain, d2) d1 \ DomainRef(d2) -@deprecate \(d1, d2::AnyDomain) DomainRef(d1) \ d2 -@deprecate setdiff(d1::AnyDomain, d2) setdiff(d1, DomainRef(d2)) -@deprecate setdiff(d1, d2::AnyDomain) setdiff(DomainRef(d1), d2) -@deprecate (&)(d1::AnyDomain, d2) d1 & DomainRef(d2) -@deprecate (&)(d1, d2::AnyDomain) DomainRef(d1) & d2 - -@deprecate isequal1(d1,d2) isequaldomain1(d1,d2) -@deprecate isequal2(d1,d2) isequaldomain2(d1,d2) - -# AsDomain was defined in versions 0.7 and 0.7.1, and then replaced by DomainRef -@deprecate AsDomain(d) DomainRef(d) - -@deprecate inv(m::AbstractMap) inverse(m) - -@deprecate size(m::AbstractMap) mapsize(m) -@deprecate size(m::AbstractMap, i) mapsize(m, i) -@deprecate issquare(m::AbstractMap) issquaremap(m) +# New deprecations in 0.8: +import Base: isreal +@deprecate isreal(d::Domain) isrealdomain(d) + +LinearAlgebra.cross(d1::Domain, domains...) = linearalgebra_x_becomes_domainsets_x(d1, domains...) +@deprecate linearalgebra_x_becomes_domainsets_x(d1::Domain, domains...) cartesianproduct(d1, domains...) + +#################### +# Change in exports +#################### + +# Quite a few functions were exported in DomainSets 0.7.15 but are no longer +# exported in DomainSets 0.8. Of those, some have moved to FunctionMaps, some +# are just internal to DomainSets. + +# Any missing exports not in the two lists below have been deprecated and removed. + +# A) Functions that now live in FunctionMaps (exported or unexported) +# "AffineMap" +# "ConstantMap" +# "IdentityMap" +# "LinearMap" +# "Map" +# "MapRef" +# "ProductMap" +# "Translation" +# "UnityMap" +# "ZeroMap" +# "affinematrix" +# "affinevector" +# "applymap" +# "codomaintype" +# "composedmap" +# "diffvolume" +# "domaintype" +# "inverse" +# "isaffinemap" +# "isconstantmap" +# "isequalmap" +# "islinearmap" +# "jacdet" +# "jacobian" +# "leftinverse" +# "mapconstant" +# "mapsize" +# "productmap" +# "rightinverse" + +# B) Functions that are no longer exported but remain available in DomainSets +# "CartToPolarMap" +# "ComposedMap" +# "PolarToCartMap" +# "StaticIdentityMap" +# "TypedMap" +# "VectorIdentityMap" +# "rotation_map" diff --git a/src/domains/complex.jl b/src/domains/complex.jl index a9de5d4..42daec6 100644 --- a/src/domains/complex.jl +++ b/src/domains/complex.jl @@ -33,14 +33,16 @@ show(io::IO, d::ComplexUnitDisk{T,:open}) where {T} = print(io, "ComplexUnitDisk canonicaldomain(::Parameterization, d::ComplexUnitCircle{T}) where {T} = UnitInterval{T}() mapfrom_canonical(::Parameterization, d::ComplexUnitCircle{T}) where {T} = - VectorToComplex{T}() ∘ UnitCircleMap{T}() + FunctionMaps.VectorToComplex{T}() ∘ UnitCircleMap{T}() canonicaldomain(::Isomorphic, d::ComplexUnitCircle{T}) where {T} = UnitCircle{T}() -mapfrom_canonical(::Isomorphic, d::ComplexUnitCircle{T}) where {T} = VectorToComplex{T}() +mapfrom_canonical(::Isomorphic, d::ComplexUnitCircle{T}) where {T} = + FunctionMaps.VectorToComplex{T}() canonicaldomain(::Parameterization, d::ComplexUnitDisk{T}) where {T} = UnitSquare{T}() mapfrom_canonical(::Parameterization, d::ComplexUnitDisk{T}) where {T} = - VectorToComplex{T}() ∘ UnitDiskMap{T}() + FunctionMaps.VectorToComplex{T}() ∘ UnitDiskMap{T}() canonicaldomain(::Isomorphic, d::ComplexUnitDisk{T}) where {T} = UnitCircle{T}() -mapfrom_canonical(::Isomorphic, d::ComplexUnitDisk{T}) where {T} = VectorToComplex{T}() +mapfrom_canonical(::Isomorphic, d::ComplexUnitDisk{T}) where {T} = + FunctionMaps.VectorToComplex{T}() diff --git a/src/domains/interval.jl b/src/domains/interval.jl index dd4d6da..625b803 100644 --- a/src/domains/interval.jl +++ b/src/domains/interval.jl @@ -149,7 +149,7 @@ identitymap(d::TypedEndpointsInterval{L,R,T}) where {L,R,T<:Integer} = IdentityM canonicaldomain(d::TypedEndpointsInterval{:closed,:closed,T}) where {T} = ChebyshevInterval{T}() mapfrom_canonical(d::TypedEndpointsInterval{:closed,:closed}) = - bounded_interval_map(-1, 1, endpoints(d)...) + FunctionMaps.bounded_interval_map(-1, 1, endpoints(d)...) canonicaldomain(d::FixedInterval) = d canonicaldomain(ctype::Equal, d::FixedInterval) = d diff --git a/src/domains/trivial.jl b/src/domains/trivial.jl index 5afc18f..7bfbe19 100644 --- a/src/domains/trivial.jl +++ b/src/domains/trivial.jl @@ -51,9 +51,9 @@ isequaldomain2(d1, d2::EmptySpace) = isempty(d1) domainhash(d::EmptySpace, h::UInt) = hash("EmptySpace", h) mapto1(d1::EmptySpace{T}, d2) where {T} = - isempty(d2) ? StaticIdentityMap{T}() : throw(ArgumentError("Can't map empty to non-empty domain")) + isempty(d2) ? IdentityMap{T}() : throw(ArgumentError("Can't map empty to non-empty domain")) mapto2(d1, d2::EmptySpace{T}) where {T} = - isempty(d1) ? StaticIdentityMap{T}() : throw(ArgumentError("Can't map non-empty to empty domain")) + isempty(d1) ? IdentityMap{T}() : throw(ArgumentError("Can't map non-empty to empty domain")) """ diff --git a/src/generic/broadcast.jl b/src/generic/broadcast.jl index 0bdc6f1..06b6077 100644 --- a/src/generic/broadcast.jl +++ b/src/generic/broadcast.jl @@ -41,7 +41,7 @@ broadcasted(::DomainSetStyle, ::typeof(/), d::AnyDomain, a::Number) = broadcasted(::DomainSetStyle, ::typeof(\), a::Number, d::AnyDomain) = mapped_domain(LinearMap(a), domain(d)) -broadcasted(::DomainSetStyle, m::AbstractMap, d::AnyDomain) = +broadcasted(::DomainSetStyle, m::Map, d::AnyDomain) = map_domain(m, domain(d)) broadcasted(::DomainSetStyle, fun::Function, d::AnyDomain) = diff --git a/src/generic/canonical.jl b/src/generic/canonical.jl index b71e9d3..6efd4a8 100644 --- a/src/generic/canonical.jl +++ b/src/generic/canonical.jl @@ -141,12 +141,12 @@ mapto_canonical(::Isomorphic, d) = leftinverse(mapfrom_canonical(Isomorphic(), d canonicaldomain(::Isomorphic, d::Domain{<:StaticVector{1,T}}) where {T} = convert(Domain{T}, d) mapfrom_canonical(::Isomorphic, d::Domain{<:StaticVector{1,T}}) where {T} = - NumberToVector{T}() + FunctionMaps.NumberToVector{T}() canonicaldomain(::Isomorphic, d::Domain{NTuple{N,T}}) where {N,T} = convert(Domain{SVector{N,T}}, d) mapfrom_canonical(::Isomorphic, d::Domain{NTuple{N,T}}) where {N,T} = - VectorToTuple{N,T}() + FunctionMaps.VectorToTuple{N,T}() diff --git a/src/generic/domain.jl b/src/generic/domain.jl index 851a2f6..b195b61 100644 --- a/src/generic/domain.jl +++ b/src/generic/domain.jl @@ -141,7 +141,7 @@ approx_indomain(x, d, tol) = in(x, d) isapprox(d1::Domain, d2::Domain; kwds...) = d1 == d2 -isreal(d::Domain) = isrealtype(eltype(d)) +isrealdomain(d::Domain) = isrealtype(eltype(d)) "Return a point from the given domain." function choice(d) end diff --git a/src/generic/mapped.jl b/src/generic/mapped.jl index ee75d19..6904069 100644 --- a/src/generic/mapped.jl +++ b/src/generic/mapped.jl @@ -37,7 +37,7 @@ corners(d::AbstractMappedDomain) = [forward_map(d, x) for x in corners(superdoma show(io::IO, mime::MIME"text/plain", d::AbstractMappedDomain) = composite_show(io, mime, d) Display.displaystencil(d::AbstractMappedDomain) = - map_stencil_broadcast(forward_map(d), superdomain(d)) + FunctionMaps.map_stencil_broadcast(forward_map(d), superdomain(d)) Display.object_parentheses(d::AbstractMappedDomain) = FunctionMaps.map_object_parentheses(forward_map(d)) Display.stencil_parentheses(d::AbstractMappedDomain) = @@ -56,11 +56,13 @@ _promote_map_domain_pair(map, domain, ::Type{T}, ::Type{T}) where T = map, domai _promote_map_domain_pair(map, domain, ::Type{S}, ::Type{T}) where {S,T} = _promote_map_domain_pair(map, domain, S, T, promote_type(S,T)) _promote_map_domain_pair(map, domain, ::Type{S}, ::Type{T}, ::Type{U}) where {S,T,U} = - convert_codomaintype(U, map), convert_eltype(U, domain) + FunctionMaps.convert_codomaintype(U, map), convert_eltype(U, domain) _promote_map_domain_pair(map, domain, ::Type{Any}, ::Type{T}, ::Type{Any}) where T = map, domain _promote_map_domain_pair(map, domain, ::Type{S}, ::Type{Any}, ::Type{Any}) where S = map, domain +_promote_map_domain_pair(map, domain, ::Type{Any}, ::Type{Any}, ::Type{Any}) = + map, domain # ensure dimensions of static arrays match _promote_map_domain_pair(map, domain, S::Type{<:StaticVector{N}}, T::Type{<:StaticVector{N}}) where N = _promote_map_domain_pair(map, domain, S, T, promote_type(S,T)) @@ -84,7 +86,7 @@ _MappedDomain(invmap, domain) = _MappedDomainT(mapped_domain_eltype(invmap,domain), invmap, domain) MappedDomain{T}(invmap, domain) where T = - _MappedDomainT(T, convert_domaintype(T, invmap), convert_eltype(T, domain)) + _MappedDomainT(T, FunctionMaps.convert_domaintype(T, invmap), convert_eltype(T, domain)) _MappedDomainT(::Type{T}, invmap, domain) where T = MappedDomain{T,typeof(invmap),typeof(domain)}(invmap, domain) diff --git a/src/generic/productdomain.jl b/src/generic/productdomain.jl index 4f6840d..c52f725 100644 --- a/src/generic/productdomain.jl +++ b/src/generic/productdomain.jl @@ -104,9 +104,24 @@ productdomain(d1::ProductDomain, d2::ProductDomain) = productdomain1(d1::ProductDomain, d2) = ProductDomain(factors(d1)..., d2) productdomain2(d1, d2::ProductDomain) = ProductDomain(d1, factors(d2)...) -# Only override cross for variables of type Domain, it may have a different -# meaning for other variables (like the vector cross product) -cross(d::Domain...) = productdomain(d...) +""" + cartesianproduct(domains...) + +Return the cartesian product of the given domains. + +The function may simplify the result in some cases. Use the `ProductDomain` +constructor to obtain an explicit representation as a product. + +# Examples +```julia +julia> using DomainSets: × + +julia> (0..1) × (2..3) +(0 .. 1) × (2 .. 3) +```` +""" +cartesianproduct(d...) = productdomain(d...) +× = cartesianproduct Base.:^(d::Domain, n::Int) = productdomain(ntuple(i->d, n)...) @@ -136,12 +151,12 @@ end # multiplication with a scalar number function map_domain(linmap::GenericLinearMap{<:StaticVector{N,S},A}, domain::ProductDomain{<:StaticVector{N,T}}) where {N,S,T,A<:Number} - c = unsafe_matrix(linmap) + c = FunctionMaps.unsafe_matrix(linmap) ProductDomain{SVector{N,promote_type(S,T)}}(map(d -> c .* d, components(domain))) end function map_domain(transmap::Translation{<:StaticVector{N,S}}, domain::ProductDomain{<:StaticVector{N,T}}) where {N,S,T} - vec = unsafe_vector(transmap) + vec = FunctionMaps.unsafe_vector(transmap) ProductDomain{SVector{N,promote_type(S,T)}}( map( (d,v) -> d .+ v, components(domain), tointernalpoint(domain, vec))) end @@ -242,4 +257,4 @@ end TupleProductDomain{T}(domains::Tuple) where {T} = TupleProductDomain{T,typeof(domains)}(domains) -matching_product_map(d::TupleProductDomain, maps) = TupleProductMap(maps) +matching_product_map(d::TupleProductDomain, maps) = FunctionMaps.TupleProductMap(maps) diff --git a/src/maps_imports.jl b/src/maps_imports.jl deleted file mode 100644 index c3d5715..0000000 --- a/src/maps_imports.jl +++ /dev/null @@ -1,105 +0,0 @@ - -# In this file we list all functions defined in FunctionMaps that -# were defined in the namespace of DomainSets - -using .FunctionMaps: - MapRef, - # canonical.jl - CanonicalType, CanonicalExtensionType, - canonicalextensiontype, - Equal, Equivalent, - # composite.jl - ComposedMap, similarmap, codomaintype, applymap, - applymap_rec, mapsize, jacobian, backpropagate, - inverse, leftinverse, rightinverse, - leftinverse_rec, rightinverse_rec, - composedmap, - isequalmap, map_hash, - MulMap, multiply_map, multiply_map1, multiply_map2, - SumMap, sum_map, sum_map1, sum_map2, - composite_jacobian, - mul_jacobian, sum_jacobian, - # inverse.jl - LazyInverse, implements_inverse, - # isomorhpism - Isomorphism, VectorToNumber, NumberToVector, - VectorToComplex, ComplexToVector, - VectorToTuple, TupleToVector, - NestedToFlat, FlatToNested, - # jacobian.jl - LazyJacobian, DeterminantMap, determinantmap, jacdet, - AbsMap, absmap, - LazyDiffVolume, diffvolume, - NumberLike, - to_matrix, to_vector, zeromatrix, zerovector, identitymatrix, - # lazy.jl - LazyMap, CompositeLazyMap, SimpleLazyMap, - supermap, - DerivedMap, WrappedMap, - # map.jl - AbstractMap, Map, TypedMap, EuclideanMap, VectorMap, - domaintype, codomaintype, prectype, numtype, - convert_domaintype, convert_codomaintype, convert_prectype, convert_numtype, - promote_map_point_pair, - promote_and_apply, promote_maps, - isvectorvalued_type, isvectorvalued, - issquaremap, isoverdetermined, isunderdetermined, - isreal, - is_scalar_to_vector, is_scalar_to_scalar, - is_vector_to_scalar, is_vector_to_vector, - isequalmap1, isequalmap2, default_isequalmap, - map_stencil, map_stencil_broadcast, - # product.jl - ProductMap, VcatMapElement, - compatibleproductdims, - tointernalpoint, toexternalpoint, - productmap, productmap1, productmap2, - VcatMap, mapdim, size_as_matrix, toexternalmatrix, - VectorProductMap, TupleProductMap, - # affine.jl - AbstractAffineMap, - unsafe_matrix, unsafe_vector, matrix, vector, - affinematrix, affinevector, - applymap!, - islinear, isaffine, - islinearmap, isaffinemap, - LinearMap, GenericLinearMap, ScalarLinearMap, - VectorLinearMap, StaticLinearMap, - Translation, ScalarTranslation, StaticTranslation, - VectorTranslation, GenericTranslation, - AffineMap, GenericAffineMap, ScalarAffineMap, - VectorAffineMap, StaticAffineMap, - # arithmetics.jl - affine_composition, - interval_map, bounded_interval_map, - # basic.jl - IdentityMap, isidentity, isidentitymap, - StaticIdentityMap, DynamicIdentityMap, - EuclideanIdentityMap, VectorIdentityMap, - ConstantMap, isconstantmap, mapconstant, - isconstant, constant, - ZeroMap, UnityMap, FixedConstantMap, - # coordinates.jl - CartToPolarMap, PolarToCartMap, - UnitCircleMap, AngleMap, UnitDiskMap, - # common.jl - isrealtype, - StaticTypes, - hashrec, - euclideandimension, - convert_eltype, promotable_eltypes, - prectype, convert_prectype, to_prectype, promote_prectype, - numtype, convert_numtype, to_numtype, promote_numtype, - convert_fromcartesian, convert_tocartesian, - matrix_pinv, - factors, nfactors, factor - -import .FunctionMaps: - convert_eltype, - convert_prectype, - convert_numtype, - prectype, - numtype, - factors, - isreal, - tointernalpoint, toexternalpoint, compatibleproductdims diff --git a/test/aqua.jl b/test/aqua.jl index 690eba0..ab1733b 100644 --- a/test/aqua.jl +++ b/test/aqua.jl @@ -4,7 +4,7 @@ using Aqua # A few ambiguities remain in deprecated function calls. Tests to be activated # in next breaking release when those deprecations are removed: -# Aqua.test_ambiguities(DomainSets) +Aqua.test_ambiguities(DomainSets) # The unbound tests annoyingly flag several valid uses of NTuple{N,T} # Aqua.test_unbound_args(DomainSets) diff --git a/test/runtests.jl b/test/runtests.jl index 4828fde..adef6c3 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -1,25 +1,18 @@ using Test, LinearAlgebra, StaticArrays, Random, StableRNGs -using DomainSets +using DomainSets, FunctionMaps using CompositeTypes.Indexing +using DomainSets: × + const io = IOBuffer() const textmime = MIME"text/plain"() println("#############################") println("# Aqua automated tests") println("#############################") -include("aqua.jl") - -# First run the test suite of submodule FunctionMaps.jl -println("#############################") -println("# Tests of FunctionMaps.jl") -println("#############################") -include("../FunctionMaps/test/runtests2.jl") +# include("aqua.jl") -println("#############################") -println("# Tests of DomainSets.jl") -println("#############################") include("test_common.jl") include("test_generic_domain.jl") include("test_specific_domains.jl") diff --git a/test/test_applications.jl b/test/test_applications.jl index fb7f849..4cfa26d 100644 --- a/test/test_applications.jl +++ b/test/test_applications.jl @@ -1,5 +1,23 @@ -using DomainSets: isreal +# these functions are copied from the FunctionMaps.jl tests: +# --- start of copy +randvec(T,n) = SVector{n,T}(rand(n)) +randvec(T,m,n) = SMatrix{m,n,T}(rand(m,n)) +suitable_point_to_map(m::Map) = suitable_point_to_map(m, domaintype(m)) +suitable_point_to_map(m::Map, ::Type{SVector{N,T}}) where {N,T} = SVector{N,T}(rand(N)) +suitable_point_to_map(m::Map, ::Type{T}) where {T<:Number} = rand(T) +suitable_point_to_map(m::Map, ::Type{<:AbstractVector{T}}) where {T} = rand(T, FunctionMaps.mapsize(m,2)) +suitable_point_to_map(m::FunctionMaps.ProductMap) = + map(suitable_point_to_map, components(m)) +suitable_point_to_map(m::FunctionMaps.VcatMap{T,M,N}) where {T,M,N} = + SVector{N,T}(rand(T,N)) +# --- end of copy + +suitable_point_to_map(::FunctionMaps.CartToPolarMap{T}) where {T} = randvec(T,2) +suitable_point_to_map(::FunctionMaps.PolarToCartMap{T}) where {T} = randvec(T,2) + +using DomainSets: rotation_map +using FunctionMaps: CartToPolarMap, PolarToCartMap function test_rotation_map(T) ϕ = T(pi)/4 @@ -23,28 +41,28 @@ function test_rotation_map(T) phi = T(rand()) psi = T(rand()) m2 = rotation_map(theta) - FunctionMapsTests.test_generic_map(m2) + # FunctionMapsTests.test_generic_map(m2) m3 = rotation_map(phi, theta, psi) - FunctionMapsTests.test_generic_map(m3) + # FunctionMapsTests.test_generic_map(m3) - r = FunctionMapsTests.suitable_point_to_map(m2) + r = suitable_point_to_map(m2) @test norm(m2(r))≈norm(r) - r = FunctionMapsTests.suitable_point_to_map(m3) + r = suitable_point_to_map(m3) @test norm(m3(r))≈norm(r) @test islinearmap(m3) end function test_cart_polar_map(T) m1 = CartToPolarMap{T}() - FunctionMapsTests.test_generic_map(m1) + # FunctionMapsTests.test_generic_map(m1) @test !islinearmap(m1) - @test isreal(m1) + @test isrealmap(m1) m2 = PolarToCartMap{T}() - FunctionMapsTests.test_generic_map(m2) + # FunctionMapsTests.test_generic_map(m2) @test !islinearmap(m2) - @test isreal(m2) + @test isrealmap(m2) @test inverse(m1) == m2 @test inverse(m2) == m1 @@ -141,5 +159,6 @@ end m1 = LinearMap(SMatrix{2,2}(1,2,3,4.0)) m2 = CartToPolarMap() cmap = m1 ∘ m2 ∘ m1 - FunctionMapsTests.test_generic_map(cmap) + @test cmap isa Map + # FunctionMapsTests.test_generic_map(cmap) end diff --git a/test/test_canonical.jl b/test/test_canonical.jl index e4aac1c..b5fad01 100644 --- a/test/test_canonical.jl +++ b/test/test_canonical.jl @@ -40,8 +40,8 @@ end d2 = UnitBall(Val(1)) @test canonicaldomain(iso, d2) == UnitBall{Float64}() - @test mapfrom_canonical(iso, d2) isa DomainSets.NumberToVector - @test mapto_canonical(iso, d2) isa DomainSets.VectorToNumber + @test mapfrom_canonical(iso, d2) isa FunctionMaps.NumberToVector + @test mapto_canonical(iso, d2) isa FunctionMaps.VectorToNumber @test_throws ArgumentError mapto(UnitCircle(), UnitDisk()) @test mapto(UnitInterval(), UnitInterval()) isa IdentityMap @@ -73,13 +73,6 @@ end @test_throws ArgumentError mapfrom_canonical(Interval{:closed,:open}(-Inf,2.0)) @test_throws ArgumentError mapfrom_canonical(Interval{:open,:closed}(2.0,Inf)) - @test DomainSets.interval_map(-Inf,Inf,-Inf,Inf) isa StaticIdentityMap - @test DomainSets.interval_map(Inf,-Inf,Inf,-Inf) isa StaticIdentityMap - @test DomainSets.interval_map(-Inf,Inf,Inf,-Inf) == LinearMap(-1) - @test DomainSets.interval_map(Inf,-Inf,-Inf,Inf) == LinearMap(-1) - @test DomainSets.interval_map(Inf,Inf,Inf,Inf) isa StaticIdentityMap - @test DomainSets.interval_map(-Inf,-Inf,-Inf,-Inf) isa StaticIdentityMap - @testset "canonical types" begin struct MyCanonicalType <: DomainSets.CanonicalType end diff --git a/test/test_domain_ball.jl b/test/test_domain_ball.jl index 4a1f352..83ef63a 100644 --- a/test/test_domain_ball.jl +++ b/test/test_domain_ball.jl @@ -294,7 +294,7 @@ function test_spheres() @test hasparameterization(UnitSphere(2)) p = parameterization(C) - @test mapsize(p) == (2,) + @test FunctionMaps.mapsize(p) == (2,) x = applymap(p, 1/2) @test jacobian(p, 0.4) ≈ SA[-2pi*sin(2pi*0.4), 2pi*cos(2pi*0.4)] @test diffvolume(p, 0.4) ≈ 2*pi @@ -306,7 +306,7 @@ function test_spheres() @test applymap(q, -x) ≈ 1 @test rightinverse(q) == p @test rightinverse(q)(0) ≈ rightinverse(q, 0) - @test jacobian(q) isa DomainSets.LazyJacobian + @test jacobian(q) isa FunctionMaps.LazyJacobian @test jacobian(q, x) isa LinearAlgebra.Transpose{Float64,SVector{2,Float64}} @test boundingbox(C) == ProductDomain(ChebyshevInterval(), ChebyshevInterval()) @@ -332,7 +332,7 @@ function test_spheres() @test affinematrix(mapfrom_canonical(C)) == [2 0; 0 2] @test affinevector(mapfrom_canonical(C)) == [1; 1] @test parameterdomain(C) == UnitInterval() - @test mapfrom_parameterdomain(C) isa ComposedMap + @test mapfrom_parameterdomain(C) isa FunctionMaps.ComposedMap @test mapfrom_parameterdomain(C)(0.5) ≈ [-1; 1] @test mapfrom_parameterdomain(C, 0.5) ≈ [-1; 1] @test mapto_parameterdomain(C)([-1; 1]) ≈ 0.5 diff --git a/test/test_domain_product.jl b/test/test_domain_product.jl index 7896356..784deda 100644 --- a/test/test_domain_product.jl +++ b/test/test_domain_product.jl @@ -47,7 +47,7 @@ function test_product_domains() @test SA[-1.1,0.3] ∉ d2 d3 = VcatDomain(1.05 * UnitDisk(), -1.0 .. 1.0) - @inferred(cross(1.05 * UnitDisk(), -1.0 .. 1.0)) === d3 + @inferred(cartesianproduct(1.05 * UnitDisk(), -1.0 .. 1.0)) === d3 @test d3 isa VcatDomain @test eltype(d3) == SVector{3,Float64} @test SA[0.5,0.5,0.8] ∈ d3 @@ -375,14 +375,14 @@ function test_product_domains() @testset "product mapto" begin d1 = ProductDomain(1.0..2.0, 1.0..2.0) d2 = ProductDomain(2.0..4.0, 2.0..4.0) - @test mapto(d1, d2) isa DomainSets.VcatMap - @test jacobian(mapto(d1,d2)) isa ConstantMap - @test jacdet(mapto(d1,d2)) == DomainSets.ConstantMap{SVector{2,Float64}}(4.0) + @test mapto(d1, d2) isa FunctionMaps.VcatMap + @test jacobian(mapto(d1,d2)) isa FunctionMaps.ConstantMap + @test jacdet(mapto(d1,d2)) == FunctionMaps.ConstantMap{SVector{2,Float64}}(4.0) d1v = ProductDomain([1.0..2.0, 1.0..2.0]) d2v = ProductDomain([2.0..4.0, 2.0..4.0]) - @test mapto(d1v, d2v) isa DomainSets.VectorProductMap - @test jacobian(mapto(d1v,d2v)) isa DomainSets.ConstantMap - @test jacdet(mapto(d1v,d2v)) == DomainSets.ConstantMap{Vector{Float64}}(4.0) + @test mapto(d1v, d2v) isa FunctionMaps.VectorProductMap + @test jacobian(mapto(d1v,d2v)) isa FunctionMaps.ConstantMap + @test jacdet(mapto(d1v,d2v)) == FunctionMaps.ConstantMap{Vector{Float64}}(4.0) end @testset "VcatDomain == bug" begin diff --git a/test/test_generic_domain.jl b/test/test_generic_domain.jl index 42498d6..b1a3a21 100644 --- a/test/test_generic_domain.jl +++ b/test/test_generic_domain.jl @@ -10,8 +10,8 @@ widen_eltype(::Type{Vector{T}}) where {T<:Number} = Vector{widen(T)} # These tests check whether the given domain correctly implements the # interface of a domain. function test_generic_domain(d) - @test isreal(d) == isrealtype(domaineltype(d)) - @test isreal(d) == isrealtype(domain_numtype(d)) + @test isrealdomain(d) == isrealtype(domaineltype(d)) + @test isrealdomain(d) == isrealtype(domain_numtype(d)) if d isa Domain @test convert(Domain{eltype(d)}, d) == d