diff --git a/engine-proto/Engine.jl b/engine-proto/Engine.jl index 41d3ed7..ac9ed35 100644 --- a/engine-proto/Engine.jl +++ b/engine-proto/Engine.jl @@ -2,13 +2,12 @@ include("HittingSet.jl") module Engine -export Construction, mprod, codimension, dimension +export Construction, mprod import Subscripts using LinearAlgebra using AbstractAlgebra using Groebner -using HomotopyContinuation: Variable, Expression, System using ..HittingSet # --- commutative algebra --- @@ -28,34 +27,6 @@ end dimension(I::Generic.Ideal{U}, maxdepth = Inf) where {T <: RingElement, U <: MPolyRingElem{T}} = length(gens(base_ring(I))) - codimension(I, maxdepth) -# hat tip Sascha Timme -# https://github.com/JuliaHomotopyContinuation/HomotopyContinuation.jl/issues/520#issuecomment-1317681521 -function Base.convert(::Type{Expression}, f::MPolyRingElem) - variables = Variable.(symbols(parent(f))) - f_data = zip(coefficients(f), exponent_vectors(f)) - sum(cf * prod(variables .^ exp_vec) for (cf, exp_vec) in f_data) -end - -# create a ModelKit.System from an ideal in a multivariate polynomial ring. the -# variable ordering is taken from the polynomial ring -function System(I::Generic.Ideal) - eqns = Expression.(gens(I)) - variables = Variable.(symbols(base_ring(I))) - System(eqns, variables = variables) -end - -## [to do] not needed right now -# create a ModelKit.System from a list of elements of a multivariate polynomial -# ring. the variable ordering is taken from the polynomial ring -##function System(eqns::AbstractVector{MPolyRingElem}) -## if isempty(eqns) -## return System([]) -## else -## variables = Variable.(symbols(parent(f))) -## return System(Expression.(eqns), variables = variables) -## end -##end - # --- primitve elements --- abstract type Element{T} end @@ -218,75 +189,39 @@ function realize(ctx::Construction{T}) where T append!(eqns, [sum(sph.coords[k] for sph in ctx.spheres) for k in 3:4]) end - (Generic.Ideal(coordring, eqns), eqns) + Generic.Ideal(coordring, eqns) end end # ~~~ sandbox setup ~~~ -using AbstractAlgebra -using HomotopyContinuation - CoeffType = Rational{Int64} a = Engine.Point{CoeffType}() s = Engine.Sphere{CoeffType}() a_on_s = Engine.LiesOn{CoeffType}(a, s) ctx = Engine.Construction{CoeffType}(elements = Set([a]), relations= Set([a_on_s])) -##ideal_a_s = Engine.realize(ctx) -##println("A point on a sphere: ", Engine.dimension(ideal_a_s), " degrees of freedom") +ideal_a_s = Engine.realize(ctx) +println("A point on a sphere: ", Engine.dimension(ideal_a_s), " degrees of freeom") b = Engine.Point{CoeffType}() b_on_s = Engine.LiesOn{CoeffType}(b, s) Engine.push!(ctx, b) Engine.push!(ctx, s) Engine.push!(ctx, b_on_s) -ideal_ab_s, eqns_ab_s = Engine.realize(ctx) -println("Two points on a sphere: ", Engine.dimension(ideal_ab_s), " degrees of freedom") +ideal_ab_s = Engine.realize(ctx) +println("Two points on a sphere: ", Engine.dimension(ideal_ab_s), " degrees of freeom") -##spheres = [Engine.Sphere{CoeffType}() for _ in 1:3] -##tangencies = [ -## Engine.AlignsWithBy{CoeffType}( -## spheres[n], -## spheres[mod1(n+1, length(spheres))], -## CoeffType(-1//1) -## ) -## for n in 1:3 -##] -##ctx_tan_sph = Engine.Construction{CoeffType}(elements = Set(spheres), relations = Set(tangencies)) -##ideal_tan_sph = Engine.realize(ctx_tan_sph) -##println("Three mutually tangent spheres: ", Engine.dimension(ideal_tan_sph), " degrees of freedom") - -# --- test rational cut --- - -cut = [ - sum(vcat(a.coords, (s.coords - [0, 0, 0, 0, 1]))) - sum(vcat([2, 1, 1] .* a.coords, [1, 2, 1, 1, 1] .* (s.coords - [0, 0, 0, 0, 1]))) - sum(vcat([1, 2, 0] .* a.coords, [1, 1, 0, 1, 2] .* (s.coords - [0, 0, 0, 0, 1]))) +spheres = [Engine.Sphere{CoeffType}() for _ in 1:3] +tangencies = [ + Engine.AlignsWithBy{CoeffType}( + spheres[n], + spheres[mod1(n+1, length(spheres))], + CoeffType(-1//1) + ) + for n in 1:3 ] -cut_ideal_ab_s = Generic.Ideal(base_ring(ideal_ab_s), [gens(ideal_ab_s); cut]) -cut_dim = Engine.dimension(cut_ideal_ab_s) -println("Two points on a sphere, after cut: ", cut_dim, " degrees of freedom") -if cut_dim == 0 - vbls = Variable.(symbols(base_ring(ideal_ab_s))) - cut_system = System([eqns_ab_s; cut], variables = vbls) - cut_result = HomotopyContinuation.solve(cut_system) - println("non-singular solutions:") - for soln in solutions(cut_result) - display(soln) - end - println("singular solutions:") - for sing in singular(cut_result) - display(sing.solution) - end - - # test corresponding witness set - cut_matrix = [1 1 1 1 0 1 1 0 1 1 0; 1 2 1 2 0 1 1 0 1 1 0; 1 1 0 1 0 1 2 0 2 0 0] - cut_subspace = LinearSubspace(cut_matrix, [1, 1, 2]) - witness = witness_set(System(eqns_ab_s, variables = vbls), cut_subspace) - println("witness solutions:") - for wtns in solutions(witness) - display(wtns) - end -end \ No newline at end of file +ctx_tan_sph = Engine.Construction{CoeffType}(elements = Set(spheres), relations = Set(tangencies)) +ideal_tan_sph = Engine.realize(ctx_tan_sph) +println("Three mutually tangent spheres: ", Engine.dimension(ideal_tan_sph), " degrees of freeom") \ No newline at end of file