Merge pull request #389 from oscar-system/PolyhedralGeometry/scalar-t…

…ypes

src/Polymake.jl

@@ -161,6 +161,7 @@ include("convert.jl")
 
 include("integers.jl")
 include("rationals.jl")
+include("quadraticextension.jl")
 include("sets.jl")
 include("std/lists.jl")
 include("std/pairs.jl")

src/arrays.jl

@@ -1,5 +1,5 @@
 const Array_suppT = Union{Int64, CxxWrap.CxxLong,
-                        Integer, Rational,
+                        Integer, Rational, QuadraticExtension{Rational},
                         String, CxxWrap.StdString,
                         StdPair{CxxWrap.CxxLong,CxxWrap.CxxLong},
                         StdList{StdPair{CxxWrap.CxxLong,CxxWrap.CxxLong}},

src/convert.jl

@@ -96,14 +96,16 @@ NodeMap{Dir, T}(g::Graph{Dir}) where Dir<:DirType where T<:Set{<:Union{Int64, Cx
 
 
 convert_to_pm_type(::Type{HomologyGroup{T}}) where T<:Integer = HomologyGroup{T}
+convert_to_pm_type(::Type{<:QuadraticExtension{T}}) where T<:Rational = QuadraticExtension{Rational}
 # convert_to_pm_type(::Type{<:Union{AbstractSet, Set}}) = Set
 
 # specific converts for container types we wrap:
 convert_to_pm_type(::Type{<:Set{<:Base.Integer}}) = Set{Int64}
 convert_to_pm_type(::Type{<:Base.AbstractSet{<:Base.Integer}}) = Set{Int64}
 
 for (pmT, jlT) in [(Integer, Base.Integer),
-                   (Rational, Union{Base.Rational, Rational})]
+                   (Rational, Union{Base.Rational, Rational}),
+                   (QuadraticExtension{Polymake.Rational}, QuadraticExtension{Polymake.Rational})]
     @eval begin
         convert_to_pm_type(::Type{<:AbstractMatrix{T}}) where T<:$jlT = Matrix{$pmT}
         convert_to_pm_type(::Type{<:AbstractVector{T}}) where T<:$jlT = Vector{$pmT}

src/meta.jl

@@ -1,7 +1,7 @@
 module Meta
 import JSON
 import Polymake: appname_module_dict, module_appname_dict, shell_context_help
-import Polymake: Rational, PolymakeType, PropertyValue, OptionSet
+import Polymake: Rational, PolymakeType, PropertyValue, OptionSet, QuadraticExtension
 
 struct UnparsablePolymakeFunction <: Exception
     msg::String
@@ -33,6 +33,7 @@ translate_type_to_pm_string(::Type{<:Base.Rational}) = "Rational"
 translate_type_to_pm_string(::Type{<:Base.Integer}) = "Integer"
 translate_type_to_pm_string(::typeof(min)) = "Min"
 translate_type_to_pm_string(::typeof(max)) = "Max"
+translate_type_to_pm_string(::Type{<:QuadraticExtension{T}}) where T = string("QuadraticExtension<", translate_type_to_pm_string(T), ">")
 
 translate_type_to_pm_string(T) = throw(DomainError(T, "$T has been passed as a type parameter but no translation to a C++ template was defined. You may define such translation by appropriately extending
     `Polymake.Meta.translate_type_to_pm_string`."))

src/quadraticextension.jl

@@ -0,0 +1,46 @@
+const qe_suppT = Union{Polymake.Rational}
+
+QuadraticExtension{T}(a::Number, b::Number, r::Number) where T<:qe_suppT =
+    QuadraticExtension{T}(convert(T, a), convert(T, b), convert(T, r))
+
+QuadraticExtension{T}(a::Number) where T<:qe_suppT = QuadraticExtension{T}(a, 0, 0)
+
+QuadraticExtension(x...) = QuadraticExtension{Rational}(x...)
+
+Base.zero(::Type{<:QuadraticExtension{T}}) where T<:qe_suppT = QuadraticExtension{T}(0)
+Base.zero(::QuadraticExtension{T}) where T<:qe_suppT = QuadraticExtension{T}(0)
+Base.one(::Type{<:QuadraticExtension{T}}) where T<:qe_suppT = QuadraticExtension{T}(1)
+Base.one(::QuadraticExtension{T}) where T<:qe_suppT = QuadraticExtension{T}(1)
+
+generating_field_elements(qe::QuadraticExtension{T}) where T<:qe_suppT = (a = _a(qe), b = _b(qe),  r =_r(qe))
+
+# we might need to split this up
+# if we should wrap `QuadraticExtension` with another scalar than `Rational`
+function Base.promote_rule(::Type{<:QuadraticExtension{T}},
+    ::Type{<:Union{T, Base.Integer, Base.Rational{<:Base.Integer}}}) where T<:qe_suppT
+    return QuadraticExtension{T}
+end
+
+import Base: <, //, <=
+# defining for `Real` to avoid disambiguities
+Base.:<(x::Real, y::QuadraticExtension{T}) where T<:qe_suppT = convert(QuadraticExtension{T}, x) < y
+Base.:<(x::QuadraticExtension{T}, y::Real) where T<:qe_suppT = x < convert(QuadraticExtension{T}, y)
+Base.://(x::Real, y::QuadraticExtension{T}) where T<:qe_suppT = convert(QuadraticExtension{T}, x) // y
+Base.://(x::QuadraticExtension{T}, y::Real) where T<:qe_suppT = x // convert(QuadraticExtension{T}, y)
+
+Base.:<=(x::QuadraticExtension{T}, y::QuadraticExtension{T}) where T<:qe_suppT = x < y || x == y
+Base.:/(x::QuadraticExtension{T}, y::QuadraticExtension{T}) where T<:qe_suppT = x // y
+
+# no-copy convert
+convert(::Type{<:QuadraticExtension{T}}, qe::QuadraticExtension{T}) where T<:qe_suppT = qe
+
+function convert(to::Type{<:Number}, qe::QuadraticExtension)
+    !iszero(_b(qe)) && !iszero(_r(qe)) && throw(DomainError("Given QuadraticExtension not trivial."))
+    return convert(to, _a(qe))
+end
+
+# compatibility with Float64
+Float64(x::QuadraticExtension{T}) where T<:qe_suppT = Float64(_a(x)) + Float64(_b(x)) * sqrt(Float64(_r(x)))
+Base.promote_rule(::Type{<:QuadraticExtension{Rational}}, ::Type{<:AbstractFloat}) = Float64
+
+convert(to::Type{<:AbstractFloat}, qe::QuadraticExtension) = convert(to, Float64(qe))

src/setup_types.jl

@@ -15,10 +15,11 @@ const SmallObject = Union{
     Map,
     Graph,
     HomologyGroup,
+    QuadraticExtension,
     EdgeMap,
     NodeMap,
 }
-const VecOrMat_eltypes = Union{Int64, Integer, Rational, Float64, CxxWrap.CxxLong}
+const VecOrMat_eltypes = Union{Int64, Integer, Rational, Float64, QuadraticExtension{Rational}, CxxWrap.CxxLong}
 
 const TypeConversionFunctions = Dict(
     Symbol("Int") => to_int,

test/arrays.jl

@@ -3,6 +3,7 @@
         (Int64, 2),
         (Polymake.Integer, Polymake.Integer(2)),
         (Polymake.Rational, Polymake.Rational(2 //3)),
+        (Polymake.QuadraticExtension{Polymake.Rational}, Polymake.QuadraticExtension{Polymake.Rational}(1, 2, 3)),
         (String, "a"),
         (Polymake.Set{Polymake.to_cxx_type(Int)}, Polymake.Set{Int64}([1,2,1])),
         (Polymake.Matrix{Polymake.Integer}, Polymake.Matrix{Polymake.Integer}([1 0; 2 1])),

test/convert.jl

@@ -21,7 +21,7 @@
 
     @testset "convert_to_pm_type(PolymakeType)" begin
         for T in [Polymake.Integer, Polymake.Rational, Polymake.Array, Polymake.IncidenceMatrix,
-            Polymake.Matrix, Polymake.Set{Int64}, Polymake.SparseMatrix, Polymake.TropicalNumber, Polymake.Vector]
+            Polymake.Matrix, Polymake.Set{Int64}, Polymake.SparseMatrix, Polymake.TropicalNumber, Polymake.Vector, Polymake.QuadraticExtension{Polymake.Rational}]
             @test T == Polymake.convert_to_pm_type(T)
         end
     end

test/matrices.jl

@@ -2,7 +2,7 @@
     IntTypes = [Int32, Int64, UInt64, BigInt]
     FloatTypes = [Float32, Float64, BigFloat]
 
-    for T in [Int64, Polymake.Integer, Polymake.Rational, Float64]
+    for T in [Int64, Polymake.Integer, Polymake.Rational, Float64, Polymake.QuadraticExtension{Polymake.Rational}]
         @test Polymake.Matrix{T} <: AbstractMatrix
         @test Polymake.Matrix{T}(undef, 3,4) isa AbstractMatrix
         @test Polymake.Matrix{T}(undef, 3,4) isa Polymake.Matrix
@@ -26,7 +26,7 @@
             @test Polymake.Matrix(T.(jl_m)) isa
                 Polymake.Matrix{T<:Union{Int32,Int64} ? Polymake.to_cxx_type(Int64) : Polymake.Integer}
 
-            for ElType in [Polymake.Integer, Polymake.Rational, Float64]
+            for ElType in [Polymake.Integer, Polymake.Rational, Float64, Polymake.QuadraticExtension{Polymake.Rational}]
                 for m in [jl_m, jl_m//T(1), jl_m/T(1)]
                     @test Polymake.Matrix{ElType}(m) isa Polymake.Matrix{ElType}
                     @test convert(Polymake.Matrix{ElType}, m) isa Polymake.Matrix{ElType}
@@ -150,6 +150,34 @@
                 @test string(V) == "pm::Matrix<pm::Rational>\n5/3 2 3\n10/3 5 6\n"
             end
         end
+
+        @testset "Polymake.Matrix{Polymake.QuadraticExtension{Polymake.Rational}}" begin
+            V = Polymake.Matrix{Polymake.QuadraticExtension{Polymake.Rational}}(jl_m)
+            # linear indexing:
+            @test V[1] == 1//1
+            @test V[2] == 4//1
+
+            @test eltype(V) == Polymake.QuadraticExtension{Polymake.Rational}
+
+            @test_throws BoundsError V[0, 1]
+            @test_throws BoundsError V[2, 5]
+            @test_throws BoundsError V[3, 1]
+
+            @test length(V) == 6
+            @test size(V) == (2,3)
+
+            for T in [IntTypes; Polymake.Integer]
+                V = Polymake.Matrix{Polymake.QuadraticExtension{Polymake.Rational}}(jl_m) # local copy
+                @test setindex!(V, T(5)//T(3), 1, 1) isa Polymake.Matrix{Polymake.QuadraticExtension{Polymake.Rational}}
+                @test V[T(1), 1] isa Polymake.QuadraticExtension{Polymake.Rational}
+                @test V[1, T(1)] == 5//3
+                # testing the return value of brackets operator
+                @test V[2] = T(10)//T(3) isa typeof(T(10)//T(3))
+                V[2] = T(10)//T(3)
+                @test V[2] == 10//3
+                @test string(V) == "pm::Matrix<pm::QuadraticExtension<pm::Rational> >\n5/3 2 3\n10/3 5 6\n"
+            end
+        end
 
         @testset "Polymake.Matrix{Float64}" begin
             V = Polymake.Matrix{Float64}(jl_m)
@@ -187,6 +215,7 @@
                 V = Polymake.Matrix{Polymake.Integer}(undef, 2, 3)
                 W = Polymake.Matrix{Polymake.Rational}(undef, 2, 3)
                 U = Polymake.Matrix{Float64}(undef, 2, 3)
+                Y = Polymake.Matrix{Polymake.QuadraticExtension{Polymake.Rational}}(undef, 2, 3)
 
                 @test (X .= T.(jl_m)) isa Polymake.Matrix{Polymake.to_cxx_type(Int64)}
                 @test (X .= T.(jl_m).//1) isa Polymake.Matrix{Polymake.to_cxx_type(Int64)}
@@ -199,8 +228,11 @@
 
                 @test (U .= T.(jl_m)) isa Polymake.Matrix{Float64}
                 @test (U .= T.(jl_m).//1) isa Polymake.Matrix{Float64}
+
+                @test (Y .= T.(jl_m)) isa Polymake.Matrix{Polymake.QuadraticExtension{Polymake.Rational}}
+                @test (Y .= T.(jl_m).//1) isa Polymake.Matrix{Polymake.QuadraticExtension{Polymake.Rational}}
 
-                @test X == U == V == W
+                @test X == U == V == W == Y
 
                 # TODO:
                 # @test (V .== jl_m) isa BitPolymake.Array
@@ -232,6 +264,9 @@
         W = Polymake.Matrix{Polymake.Rational}(jl_w)
         jl_u = jl_m/4
         U = Polymake.Matrix{Float64}(jl_u)
+        sr2 = Polymake.QuadraticExtension{Polymake.Rational}(0, 1, 2)
+        jl_y = sr2 * jl_m
+        Y = Polymake.Matrix{Polymake.QuadraticExtension{Polymake.Rational}}(jl_y)
 
         @test -X isa Polymake.Polymake.MatrixAllocated{Polymake.to_cxx_type(Int64)}
         @test -X == -jl_m
@@ -244,6 +279,9 @@
 
         @test -U isa Polymake.Polymake.MatrixAllocated{Float64}
         @test -U == -jl_u
+
+        @test -Y isa Polymake.Matrix{Polymake.QuadraticExtension{Polymake.Rational}}
+        @test -Y == -jl_y
 
         int_scalar_types = [IntTypes; Polymake.Integer]
         rational_scalar_types = [[Base.Rational{T} for T in IntTypes]; Polymake.Rational]
@@ -252,7 +290,7 @@
         @test Int32(2)X isa Polymake.Matrix{Polymake.to_cxx_type(Int64)}
 
         for T in int_scalar_types
-            for (mat, ElType) in [(V, Polymake.Integer), (W, Polymake.Rational), (U, Float64)]
+            for (mat, ElType) in [(V, Polymake.Integer), (W, Polymake.Rational), (U, Float64), (Y, Polymake.QuadraticExtension{Polymake.Rational})]
                 op = *
                 @test op(T(2), mat) isa Polymake.Matrix{ElType}
                 @test op(mat, T(2)) isa Polymake.Matrix{ElType}
@@ -277,6 +315,13 @@
                     @test broadcast(op, mat, T(2)) isa Polymake.Matrix{ElType}
                 end
             end
+            let (op, ElType) = (//, Polymake.QuadraticExtension{Polymake.Rational})
+
+                @test op(Y, T(2)) isa Polymake.Matrix{ElType}
+                @test broadcast(op, T(2), Y) isa Polymake.Matrix{ElType}
+                @test broadcast(op, Y, T(2)) isa Polymake.Matrix{ElType}
+
+            end
             let (op, ElType) = (/, Float64)
                 mat = U
                 @test op(mat, T(2)) isa Polymake.Matrix{ElType}
@@ -286,7 +331,7 @@
         end
 
         for T in rational_scalar_types
-            for (mat, ElType) in [(V, Polymake.Rational), (W, Polymake.Rational), (U, Float64)]
+            for (mat, ElType) in [(V, Polymake.Rational), (W, Polymake.Rational), (U, Float64), (Y, Polymake.QuadraticExtension{Polymake.Rational})]
 
                 op = *
                 @test op(T(2), mat) isa Polymake.Matrix{ElType}
@@ -340,12 +385,37 @@
                 @test broadcast(op, mat, T(2)) isa Polymake.Matrix{ElType}
             end
         end
+
+        let T = Polymake.QuadraticExtension{Polymake.Rational}, mat = Y, ElType = Polymake.QuadraticExtension{Polymake.Rational}
+            op = *
+            @test op(T(2), mat) isa Polymake.Matrix{ElType}
+            @test op(mat, T(2)) isa Polymake.Matrix{ElType}
+            @test broadcast(op, T(2), mat) isa Polymake.Matrix{ElType}
+            @test broadcast(op, mat, T(2)) isa Polymake.Matrix{ElType}
+
+            op = +
+            @test op(mat, T.(jl_m)) isa Polymake.Matrix{ElType}
+            @test op(T.(jl_m), mat) isa Polymake.Matrix{ElType}
+
+            @test broadcast(op, mat, T.(jl_m)) isa Polymake.Matrix{ElType}
+            @test broadcast(op, T.(jl_m), mat) isa Polymake.Matrix{ElType}
+
+            @test broadcast(op, T(2), mat) isa Polymake.Matrix{ElType}
+            @test broadcast(op, mat, T(2)) isa Polymake.Matrix{ElType}
+
+            op = //
+            # @test op(T(2), mat) isa Polymake.Matrix{ElType}
+            @test op(mat, T(2)) isa Polymake.Matrix{ElType}
+            @test broadcast(op, T(2), mat) isa Polymake.Matrix{ElType}
+            @test broadcast(op, mat, T(2)) isa Polymake.Matrix{ElType}
+        end
 
-        for T in [int_scalar_types; rational_scalar_types; FloatTypes]
+        for T in [int_scalar_types; rational_scalar_types; FloatTypes; Polymake.QuadraticExtension{Polymake.Rational}]
             @test T(2)*X == X*T(2) == T(2) .* X == X .* T(2) == 2jl_m
             @test T(2)*V == V*T(2) == T(2) .* V == V .* T(2) == 2jl_m
             @test T(2)*W == W*T(2) == T(2) .* W == W .* T(2) == 2jl_w
             @test T(2)*U == U*T(2) == T(2) .* U == U .* T(2) == 2jl_u
+            @test T(2)*Y == Y*T(2) == T(2) .* Y == Y .* T(2) == 2jl_y
 
             @test X + T.(jl_m) == T.(jl_m) + X == X .+ T.(jl_m) == T.(jl_m) .+ X == 2jl_m
 
@@ -354,6 +424,8 @@
             @test W + T.(4jl_w) == T.(4jl_w) + W == W .+ T.(4jl_w) == T.(4jl_w) .+ W == 5jl_w
 
             @test U + T.(4jl_u) == T.(4jl_u) + U == U .+ T.(4jl_u) == T.(4jl_u) .+ U == 5jl_u
+
+            @test Y + T.(2 * jl_m) == T.(2 * jl_m) + Y == Y .+ T.(2 * jl_m) == T.(2 * jl_m) .+ Y == (1 + sr2) * jl_y
         end
     end
 end

test/perlobj.jl

@@ -136,6 +136,10 @@
 
         c = polytope.cube(3, 1//4, -1//4)
         @test c.VERTICES[1,2] == -1//4
+
+        i = Polymake.polytope.icosahedron()
+        @test i.VERTICES[1, :] == [1, 0, Polymake.QuadraticExtension{Polymake.Rational}(1//4, 1//4, 5), 1//2]
+        @test i.VOLUME == Polymake.QuadraticExtension{Polymake.Rational}(5//4, 5//12, 5)
     end
 
     @testset "attachments" begin

test/polynomial.jl

@@ -4,7 +4,7 @@ using SparseArrays
     FloatTypes = [Float32, Float64, BigFloat]
     RationalTypes = [Rational{I} for I in IntTypes]
 
-    for C in [Int64, Polymake.Integer, Polymake.Rational, Float64]
+    for C in [Int64, Polymake.Integer, Polymake.Rational, Float64, Polymake.QuadraticExtension{Polymake.Rational}]
         @test Polymake.Polynomial{C,Int64} <: Any
         @test Polymake.Polynomial{C,Int64}([1, 2],[3 4 5; 6 7 8]) isa Any
         @test Polymake.Polynomial{C,Int64}([1, 2],[3 4 5; 6 7 8]) isa Polymake.Polynomial
@@ -22,7 +22,7 @@ using SparseArrays
     end
 
     @testset "Low-level operations" begin
-        for (C,s) in [(Int64, "long"), (Polymake.Integer, "pm::Integer"), (Polymake.Rational, "pm::Rational"), (Float64, "double")]
+        for (C,s) in [(Int64, "long"), (Polymake.Integer, "pm::Integer"), (Polymake.Rational, "pm::Rational"), (Float64, "double"), (Polymake.QuadraticExtension{Polymake.Rational}, "pm::QuadraticExtension<pm::Rational>")]
             p = Polymake.Polynomial(C.(jl_v),jl_m)
             @test Polymake.nvars(p) == size(jl_m)[2]
             @test Polymake.nvars(p) isa Int
@@ -48,17 +48,17 @@ using SparseArrays
     end
 
     @testset "Equality" begin
-        for C1 in [Int64, Polymake.Integer, Polymake.Rational, Float64], C2 in [Int64, Polymake.Integer, Polymake.Rational, Float64]
+        for C1 in [Int64, Polymake.Integer, Polymake.Rational, Float64, Polymake.QuadraticExtension{Polymake.Rational}], C2 in [Int64, Polymake.Integer, Polymake.Rational, Float64, Polymake.QuadraticExtension{Polymake.Rational}]
             @test Polymake.Polynomial{C1}(jl_v,jl_m) == Polymake.Polynomial{C2}(jl_v,jl_m)
         end
     end
 
     @testset "Arithmetic" begin
         jl_v2 = [5, 6]
         jl_m2 = [3 4 5; 6 7 8]
-        for C1 in [Int64, Polymake.Integer, Polymake.Rational, Float64]
+        for C1 in [Int64, Polymake.Integer, Polymake.Rational, Float64, Polymake.QuadraticExtension{Polymake.Rational}]
             p = Polymake.Polynomial{C1}(jl_v,jl_m)
-            for C2 in [Int64, Polymake.Integer, Polymake.Rational, Float64]
+            for C2 in [Int64, Polymake.Integer, Polymake.Rational, Float64, Polymake.QuadraticExtension{Polymake.Rational}]
                 q = Polymake.Polynomial{C2}(jl_v2,jl_m2)
                 @test p + q isa Polymake.Polynomial{Polymake.to_cxx_type(promote_type(C1,C2))}
                 @test p + q == Polymake.Polynomial([6, 2, 6],[3 4 5; 6 7 0; 6 7 8])