forked from StudioInfinity/dyna3
Use module names as filenames
You're right: this naming convention seems to be standard for Julia modules now.
This commit is contained in:
Aaron Fenyes
parent
4e02ee16fc
commit
65d23fb667
2 changed files with 1 additions and 1 deletions
207
engine-proto/Engine.jl
Normal file
207
engine-proto/Engine.jl
Normal file
|
|
@ -0,0 +1,207 @@
|
|||
include("HittingSet.jl")
|
||||
|
||||
module Engine
|
||||
|
||||
export Construction, mprod
|
||||
|
||||
import Subscripts
|
||||
using LinearAlgebra
|
||||
using AbstractAlgebra
|
||||
using Groebner
|
||||
using ..HittingSet
|
||||
|
||||
# --- commutative algebra ---
|
||||
|
||||
# as of version 0.36.6, AbstractAlgebra only supports ideals in multivariate
|
||||
# polynomial rings when the coefficients are integers. we use Groebner to extend
|
||||
# support to rationals and to finite fields of prime order
|
||||
Generic.reduce_gens(I::Generic.Ideal{U}) where {T <: FieldElement, U <: MPolyRingElem{T}} =
|
||||
Generic.Ideal{U}(base_ring(I), groebner(gens(I)))
|
||||
|
||||
function codimension(I::Generic.Ideal{U}, maxdepth = Inf) where {T <: RingElement, U <: MPolyRingElem{T}}
|
||||
leading = [exponent_vector(f, 1) for f in gens(I)]
|
||||
targets = [Set(findall(.!iszero.(exp_vec))) for exp_vec in leading]
|
||||
length(HittingSet.solve(HittingSetProblem(targets), maxdepth))
|
||||
end
|
||||
|
||||
dimension(I::Generic.Ideal{U}, maxdepth = Inf) where {T <: RingElement, U <: MPolyRingElem{T}} =
|
||||
length(gens(base_ring(I))) - codimension(I, maxdepth)
|
||||
|
||||
# --- primitve elements ---
|
||||
|
||||
abstract type Element{T} end
|
||||
|
||||
mutable struct Point{T} <: Element{T}
|
||||
coords::Vector{MPolyRingElem{T}}
|
||||
vec::Union{Vector{MPolyRingElem{T}}, Nothing}
|
||||
rel::Nothing
|
||||
|
||||
## [to do] constructor argument never needed?
|
||||
Point{T}(
|
||||
coords::Vector{MPolyRingElem{T}} = MPolyRingElem{T}[],
|
||||
vec::Union{Vector{MPolyRingElem{T}}, Nothing} = nothing
|
||||
) where T = new(coords, vec, nothing)
|
||||
end
|
||||
|
||||
function buildvec!(pt::Point)
|
||||
coordring = parent(pt.coords[1])
|
||||
pt.vec = [one(coordring), dot(pt.coords, pt.coords), pt.coords...]
|
||||
end
|
||||
|
||||
mutable struct Sphere{T} <: Element{T}
|
||||
coords::Vector{MPolyRingElem{T}}
|
||||
vec::Union{Vector{MPolyRingElem{T}}, Nothing}
|
||||
rel::Union{MPolyRingElem{T}, Nothing}
|
||||
|
||||
## [to do] constructor argument never needed?
|
||||
Sphere{T}(
|
||||
coords::Vector{MPolyRingElem{T}} = MPolyRingElem{T}[],
|
||||
vec::Union{Vector{MPolyRingElem{T}}, Nothing} = nothing,
|
||||
rel::Union{MPolyRingElem{T}, Nothing} = nothing
|
||||
) where T = new(coords, vec, rel)
|
||||
end
|
||||
|
||||
function buildvec!(sph::Sphere)
|
||||
coordring = parent(sph.coords[1])
|
||||
sph.vec = sph.coords
|
||||
sph.rel = mprod(sph.coords, sph.coords) + one(coordring)
|
||||
end
|
||||
|
||||
const coordnames = IdDict{Symbol, Vector{Union{Symbol, Nothing}}}(
|
||||
nameof(Point) => [nothing, nothing, :xₚ, :yₚ, :zₚ],
|
||||
nameof(Sphere) => [:rₛ, :sₛ, :xₛ, :yₛ, :zₛ]
|
||||
)
|
||||
|
||||
coordname(elem::Element, index) = coordnames[nameof(typeof(elem))][index]
|
||||
|
||||
function pushcoordname!(coordnamelist, indexed_elem::Tuple{Any, Element}, coordindex)
|
||||
elemindex, elem = indexed_elem
|
||||
name = coordname(elem, coordindex)
|
||||
if !isnothing(name)
|
||||
subscript = Subscripts.sub(string(elemindex))
|
||||
push!(coordnamelist, Symbol(name, subscript))
|
||||
end
|
||||
end
|
||||
|
||||
function takecoord!(coordlist, indexed_elem::Tuple{Any, Element}, coordindex)
|
||||
elem = indexed_elem[2]
|
||||
if !isnothing(coordname(elem, coordindex))
|
||||
push!(elem.coords, popfirst!(coordlist))
|
||||
end
|
||||
end
|
||||
|
||||
# --- primitive relations ---
|
||||
|
||||
abstract type Relation{T} end
|
||||
|
||||
mprod(v, w) = (v[1]*w[2] + w[1]*v[2]) / 2 - dot(v[3:end], w[3:end])
|
||||
|
||||
# elements: point, sphere
|
||||
struct LiesOn{T} <: Relation{T}
|
||||
elements::Vector{Element{T}}
|
||||
|
||||
LiesOn{T}(pt::Point{T}, sph::Sphere{T}) where T = new{T}([pt, sph])
|
||||
end
|
||||
|
||||
equation(rel::LiesOn) = mprod(rel.elements[1].vec, rel.elements[2].vec)
|
||||
|
||||
# elements: sphere, sphere
|
||||
struct AlignsWithBy{T} <: Relation{T}
|
||||
elements::Vector{Element{T}}
|
||||
cos_angle::T
|
||||
|
||||
AlignsWithBy{T}(sph1::Sphere{T}, sph2::Sphere{T}, cos_angle::T) where T = new{T}([sph1, sph2], cos_angle)
|
||||
end
|
||||
|
||||
equation(rel::AlignsWithBy) = mprod(rel.elements[1].vec, rel.elements[2].vec) - rel.cos_angle
|
||||
|
||||
# --- constructions ---
|
||||
|
||||
mutable struct Construction{T}
|
||||
elements::Set{Element{T}}
|
||||
relations::Set{Relation{T}}
|
||||
|
||||
function Construction{T}(; elements = Set{Element{T}}(), relations = Set{Relation{T}}()) where T
|
||||
allelements = union(elements, (rel.elements for rel in relations)...)
|
||||
new{T}(allelements, relations)
|
||||
end
|
||||
end
|
||||
|
||||
function Base.push!(ctx::Construction{T}, elem::Element{T}) where T
|
||||
push!(ctx.elements, elem)
|
||||
end
|
||||
|
||||
function Base.push!(ctx::Construction{T}, rel::Relation{T}) where T
|
||||
push!(ctx.relations, rel)
|
||||
union!(ctx.elements, rel.elements)
|
||||
end
|
||||
|
||||
function realize(ctx::Construction{T}) where T
|
||||
# collect coordinate names
|
||||
coordnamelist = Symbol[]
|
||||
elemenum = enumerate(ctx.elements)
|
||||
for coordindex in 1:5
|
||||
for indexed_elem in elemenum
|
||||
pushcoordname!(coordnamelist, indexed_elem, coordindex)
|
||||
end
|
||||
end
|
||||
|
||||
# construct coordinate ring
|
||||
coordring, coordqueue = polynomial_ring(parent_type(T)(), coordnamelist, ordering = :degrevlex)
|
||||
|
||||
# retrieve coordinates
|
||||
for (_, elem) in elemenum
|
||||
empty!(elem.coords)
|
||||
end
|
||||
for coordindex in 1:5
|
||||
for indexed_elem in elemenum
|
||||
takecoord!(coordqueue, indexed_elem, coordindex)
|
||||
end
|
||||
end
|
||||
|
||||
# construct coordinate vectors
|
||||
for (_, elem) in elemenum
|
||||
buildvec!(elem)
|
||||
end
|
||||
|
||||
# turn relations into equations
|
||||
eqns = vcat(
|
||||
equation.(ctx.relations),
|
||||
[elem.rel for elem in ctx.elements if !isnothing(elem.rel)]
|
||||
)
|
||||
Generic.Ideal(coordring, eqns)
|
||||
end
|
||||
|
||||
end
|
||||
|
||||
# ~~~ sandbox setup ~~~
|
||||
|
||||
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 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 = 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 freeom")
|
||||
Loading…
Add table
Add a link
Reference in a new issue