mirror of
https://github.com/kalmarek/PropertyT.jl.git
synced 2024-11-27 01:10:28 +01:00
Merge branch 'AutF4' of git.wmi.amu.edu.pl:kalmar/PropertyT.jl into AutF4
# Conflicts: # src/PropertyT.jl # src/sdps.jl
This commit is contained in:
commit
883f6dcd18
@ -4,6 +4,7 @@ using JLD
|
|||||||
using GroupRings
|
using GroupRings
|
||||||
using Memento
|
using Memento
|
||||||
|
|
||||||
|
using Groups
|
||||||
import Nemo: Group, GroupElem
|
import Nemo: Group, GroupElem
|
||||||
|
|
||||||
const logger = Memento.config("info", fmt="{msg}")
|
const logger = Memento.config("info", fmt="{msg}")
|
||||||
@ -45,23 +46,21 @@ function ΔandSDPconstraints(name::String, G::Group)
|
|||||||
return Δ, sdp_constraints
|
return Δ, sdp_constraints
|
||||||
end
|
end
|
||||||
|
|
||||||
function ΔandSDPconstraints{T<:GroupElem}(name::String, S::Vector{T}, radius::Int)
|
function ΔandSDPconstraints{T<:GroupElem}(name::String, S::Vector{T}, Id::T; radius::Int=2)
|
||||||
S, Id = generating_set()
|
|
||||||
info(logger, "Computing pm, Δ, sdp_constraints...")
|
info(logger, "Computing pm, Δ, sdp_constraints...")
|
||||||
t = @timed Δ, sdp_constraints = ΔandSDPconstraints(S, radius)
|
Δ, sdp_constraints = ΔandSDPconstraints(S, Id, radius=radius)
|
||||||
info(logger, timed_msg(t))
|
|
||||||
pm_fname, Δ_fname = pmΔfilenames(name)
|
pm_fname, Δ_fname = pmΔfilenames(name)
|
||||||
save(pm_fname, "pm", parent(Δ).pm)
|
save(pm_fname, "pm", parent(Δ).pm)
|
||||||
save(Δ_fname, "Δ", Δ.coeffs)
|
save(Δ_fname, "Δ", Δ.coeffs)
|
||||||
|
return Δ, sdp_constraints
|
||||||
end
|
end
|
||||||
|
|
||||||
function ΔandSDPconstraints{T<:GroupElem}(S::Vector{T}, r::Int=2)
|
function ΔandSDPconstraints{T<:GroupElem}(S::Vector{T}, Id::T; radius::Int=2)
|
||||||
Id = parent(S[1])()
|
B, sizes = Groups.generate_balls(S, Id, radius=2*radius)
|
||||||
B, sizes = Groups.generate_balls(S, Id, radius=2*r)
|
|
||||||
info(logger, "Generated balls of sizes $sizes")
|
info(logger, "Generated balls of sizes $sizes")
|
||||||
|
|
||||||
info(logger, "Creating product matrix...")
|
info(logger, "Creating product matrix...")
|
||||||
t = @timed pm = GroupRings.create_pm(B, GroupRings.reverse_dict(B), sizes[r]; twisted=true)
|
t = @timed pm = GroupRings.create_pm(B, GroupRings.reverse_dict(B), sizes[radius]; twisted=true)
|
||||||
info(logger, timed_msg(t))
|
info(logger, timed_msg(t))
|
||||||
|
|
||||||
info(logger, "Creating sdp_constratints...")
|
info(logger, "Creating sdp_constratints...")
|
||||||
@ -70,7 +69,7 @@ function ΔandSDPconstraints{T<:GroupElem}(S::Vector{T}, r::Int=2)
|
|||||||
|
|
||||||
RG = GroupRing(parent(Id), B, pm)
|
RG = GroupRing(parent(Id), B, pm)
|
||||||
|
|
||||||
Δ = splaplacian(RG, S, B[1:sizes[r]], sizes[2*r])
|
Δ = splaplacian(RG, S, Id, sizes[2*radius])
|
||||||
return Δ, sdp_constraints
|
return Δ, sdp_constraints
|
||||||
end
|
end
|
||||||
|
|
||||||
@ -143,6 +142,7 @@ end
|
|||||||
Kazhdan_from_sgap(λ,N) = sqrt(2*λ/N)
|
Kazhdan_from_sgap(λ,N) = sqrt(2*λ/N)
|
||||||
|
|
||||||
function setup_logging(name::String)
|
function setup_logging(name::String)
|
||||||
|
isdir(name) || mkdir(name)
|
||||||
|
|
||||||
Memento.add_handler(logger,
|
Memento.add_handler(logger,
|
||||||
Memento.DefaultHandler(joinpath(name,"full_$(string((now()))).log"),
|
Memento.DefaultHandler(joinpath(name,"full_$(string((now()))).log"),
|
||||||
@ -155,20 +155,16 @@ function setup_logging(name::String)
|
|||||||
end
|
end
|
||||||
|
|
||||||
|
|
||||||
function check_property_T(name::String, S, solver, upper_bound, tol, radius)
|
function check_property_T(name::String, S, Id, solver, upper_bound, tol, radius)
|
||||||
|
|
||||||
if !isdir(name)
|
isdir(name) || mkdir(name)
|
||||||
mkdir(name)
|
|
||||||
end
|
|
||||||
|
|
||||||
setup_logging(name)
|
|
||||||
|
|
||||||
if all(isfile.(pmΔfilenames(name)))
|
if all(isfile.(pmΔfilenames(name)))
|
||||||
# cached
|
# cached
|
||||||
Δ, sdp_constraints = ΔandSDPconstraints(name, parent(S[1]))
|
Δ, sdp_constraints = ΔandSDPconstraints(name, parent(S[1]))
|
||||||
else
|
else
|
||||||
# compute
|
# compute
|
||||||
Δ, sdp_constraints = ΔandSDPconstraints(name, S, radius)
|
Δ, sdp_constraints = ΔandSDPconstraints(name, S, Id, radius=radius)
|
||||||
end
|
end
|
||||||
|
|
||||||
info(logger, "|S| = $(length(S))")
|
info(logger, "|S| = $(length(S))")
|
||||||
@ -202,7 +198,7 @@ function check_property_T(name::String, S, solver, upper_bound, tol, radius)
|
|||||||
end
|
end
|
||||||
if sgap > 0
|
if sgap > 0
|
||||||
info(logger, "λ ≥ $(Float64(trunc(sgap,12)))")
|
info(logger, "λ ≥ $(Float64(trunc(sgap,12)))")
|
||||||
Kazhdan_κ = Kazhdan_from_sgap(sgap, S)
|
Kazhdan_κ = Kazhdan_from_sgap(sgap, length(S))
|
||||||
Kazhdan_κ = Float64(trunc(Kazhdan_κ, 12))
|
Kazhdan_κ = Float64(trunc(Kazhdan_κ, 12))
|
||||||
info(logger, "κ($name, S) ≥ $Kazhdan_κ: Group HAS property (T)!")
|
info(logger, "κ($name, S) ≥ $Kazhdan_κ: Group HAS property (T)!")
|
||||||
return true
|
return true
|
||||||
|
@ -81,7 +81,7 @@ function distance_to_cone{T<:Rational, S<:Interval}(λ::T, sqrt_matrix::Array{S,
|
|||||||
SOS = compute_SOS(sqrt_matrix, Δ)
|
SOS = compute_SOS(sqrt_matrix, Δ)
|
||||||
info(logger, "ɛ(∑ξᵢ*ξᵢ) ∈ $(GroupRings.augmentation(SOS))")
|
info(logger, "ɛ(∑ξᵢ*ξᵢ) ∈ $(GroupRings.augmentation(SOS))")
|
||||||
λ_int = @interval(λ)
|
λ_int = @interval(λ)
|
||||||
Δ_int = GroupRingElem([@interval(c) for c in Δ.coeffs], parent(Δ).pm)
|
Δ_int = GroupRingElem([@interval(c) for c in Δ.coeffs], parent(Δ))
|
||||||
SOS_diff = EOI(Δ_int, λ_int) - SOS
|
SOS_diff = EOI(Δ_int, λ_int) - SOS
|
||||||
eoi_SOS_L1_dist = norm(SOS_diff,1)
|
eoi_SOS_L1_dist = norm(SOS_diff,1)
|
||||||
|
|
||||||
@ -91,7 +91,7 @@ function distance_to_cone{T<:Rational, S<:Interval}(λ::T, sqrt_matrix::Array{S,
|
|||||||
info(logger, "ɛ(Δ² - λΔ - ∑ξᵢ*ξᵢ) ∈ $(ɛ_dist)")
|
info(logger, "ɛ(Δ² - λΔ - ∑ξᵢ*ξᵢ) ∈ $(ɛ_dist)")
|
||||||
info(logger, "‖Δ² - λΔ - ∑ξᵢ*ξᵢ‖₁ ∈ $(eoi_SOS_L1_dist)")
|
info(logger, "‖Δ² - λΔ - ∑ξᵢ*ξᵢ‖₁ ∈ $(eoi_SOS_L1_dist)")
|
||||||
|
|
||||||
distance_to_cone = λ - 2^(len-1)*eoi_SOS_L₁_dist
|
distance_to_cone = λ - 2^(len-1)*eoi_SOS_L1_dist
|
||||||
return distance_to_cone
|
return distance_to_cone
|
||||||
end
|
end
|
||||||
|
|
||||||
@ -99,14 +99,14 @@ function distance_to_cone{T<:AbstractFloat}(λ::T, sqrt_matrix::Array{T,2}, Δ::
|
|||||||
SOS = compute_SOS(sqrt_matrix, Δ)
|
SOS = compute_SOS(sqrt_matrix, Δ)
|
||||||
|
|
||||||
SOS_diff = EOI(Δ, λ) - SOS
|
SOS_diff = EOI(Δ, λ) - SOS
|
||||||
eoi_SOS_L₁_dist = norm(SOS_diff,1)
|
eoi_SOS_L1_dist = norm(SOS_diff,1)
|
||||||
|
|
||||||
info(logger, "λ = $λ")
|
info(logger, "λ = $λ")
|
||||||
ɛ_dist = GroupRings.augmentation(SOS_diff)
|
ɛ_dist = GroupRings.augmentation(SOS_diff)
|
||||||
info(logger, "ɛ(Δ² - λΔ - ∑ξᵢ*ξᵢ) ≈ $(@sprintf("%.10f", ɛ_dist))")
|
info(logger, "ɛ(Δ² - λΔ - ∑ξᵢ*ξᵢ) ≈ $(@sprintf("%.10f", ɛ_dist))")
|
||||||
info(logger, "‖Δ² - λΔ - ∑ξᵢ*ξᵢ‖₁ ≈ $(@sprintf("%.10f", eoi_SOS_L₁_dist))")
|
info(logger, "‖Δ² - λΔ - ∑ξᵢ*ξᵢ‖₁ ≈ $(@sprintf("%.10f", eoi_SOS_L1_dist))")
|
||||||
|
|
||||||
distance_to_cone = λ - 2^(len-1)*eoi_SOS_L₁_dist
|
distance_to_cone = λ - 2^(len-1)*eoi_SOS_L1_dist
|
||||||
return distance_to_cone
|
return distance_to_cone
|
||||||
end
|
end
|
||||||
|
|
||||||
@ -115,7 +115,7 @@ function check_distance_to_positive_cone(Δ::GroupRingElem, λ, P;
|
|||||||
|
|
||||||
isapprox(eigvals(P), abs(eigvals(P)), atol=tol) ||
|
isapprox(eigvals(P), abs(eigvals(P)), atol=tol) ||
|
||||||
warn("The solution matrix doesn't seem to be positive definite!")
|
warn("The solution matrix doesn't seem to be positive definite!")
|
||||||
@assert P == Symmetric(P)
|
# @assert P == Symmetric(P)
|
||||||
Q = real(sqrtm(P))
|
Q = real(sqrtm(P))
|
||||||
|
|
||||||
info(logger, "------------------------------------------------------------")
|
info(logger, "------------------------------------------------------------")
|
||||||
@ -130,6 +130,7 @@ function check_distance_to_positive_cone(Δ::GroupRingElem, λ, P;
|
|||||||
return fp_distance
|
return fp_distance
|
||||||
end
|
end
|
||||||
|
|
||||||
|
info(logger, "")
|
||||||
info(logger, "Projecting columns of rationalized Q to the augmentation ideal...")
|
info(logger, "Projecting columns of rationalized Q to the augmentation ideal...")
|
||||||
δ = eps(λ)
|
δ = eps(λ)
|
||||||
Q_ℚ = ℚ(Q, δ)
|
Q_ℚ = ℚ(Q, δ)
|
||||||
@ -140,20 +141,20 @@ function check_distance_to_positive_cone(Δ::GroupRingElem, λ, P;
|
|||||||
|
|
||||||
info(logger, "Checking in interval arithmetic")
|
info(logger, "Checking in interval arithmetic")
|
||||||
Q_ℚω_int = Float64.(Q_ℚω) ± δ
|
Q_ℚω_int = Float64.(Q_ℚω) ± δ
|
||||||
t = @timed Interval_dist_to_Σ² = distance_to_cone(λ_ℚ, Q_ℚω_int, Δ_ℚ, len=len)
|
t = @timed Interval_dist_to_ΣSq = distance_to_cone(λ_ℚ, Q_ℚω_int, Δ_ℚ, len=len)
|
||||||
info(logger, timed_msg(t))
|
info(logger, timed_msg(t))
|
||||||
info(logger, "The Augmentation-projected actual distance (to positive cone) ∈ $(Interval_dist_to_Σ²)")
|
info(logger, "The Augmentation-projected actual distance (to positive cone) ∈ $(Interval_dist_to_ΣSq)")
|
||||||
info(logger, "------------------------------------------------------------")
|
info(logger, "------------------------------------------------------------")
|
||||||
|
|
||||||
if Interval_dist_to_Σ².lo ≤ 0 || !rational
|
if Interval_dist_to_ΣSq.lo ≤ 0 || !rational
|
||||||
return Interval_dist_to_Σ²
|
return Interval_dist_to_ΣSq
|
||||||
else
|
else
|
||||||
info(logger, "Checking Projected SOS decomposition in exact rational arithmetic...")
|
info(logger, "Checking Projected SOS decomposition in exact rational arithmetic...")
|
||||||
t = @timed ℚ_dist_to_Σ² = distance_to_cone(λ_ℚ, Q_ℚω, Δ_ℚ, len=len)
|
t = @timed ℚ_dist_to_ΣSq = distance_to_cone(λ_ℚ, Q_ℚω, Δ_ℚ, len=len)
|
||||||
info(logger, timed_msg(t))
|
info(logger, timed_msg(t))
|
||||||
@assert isa(ℚ_dist_to_Σ², Rational)
|
@assert isa(ℚ_dist_to_ΣSq, Rational)
|
||||||
info(logger, "Augmentation-projected rational distance (to positive cone) ≥ $(Float64(trunc(ℚ_dist_to_Σ²,8)))")
|
info(logger, "Augmentation-projected rational distance (to positive cone) ≥ $(Float64(trunc(ℚ_dist_to_ΣSq,8)))")
|
||||||
info(logger, "------------------------------------------------------------")
|
info(logger, "------------------------------------------------------------")
|
||||||
return ℚ_dist_to_Σ²
|
return ℚ_dist_to_ΣSq
|
||||||
end
|
end
|
||||||
end
|
end
|
||||||
|
@ -13,9 +13,9 @@ function constraints_from_pm(pm, total_length=maximum(pm))
|
|||||||
return constraints
|
return constraints
|
||||||
end
|
end
|
||||||
|
|
||||||
function splaplacian(RG::GroupRing, S, basis, n=length(basis), T::Type=Int)
|
function splaplacian(RG::GroupRing, S, Id=RG.group(), n=length(basis),T::Type=Int)
|
||||||
result = RG(spzeros(T, n))
|
result = RG(spzeros(T, n))
|
||||||
result[RG.group()] = T(length(S))
|
result[Id] = T(length(S))
|
||||||
for s in S
|
for s in S
|
||||||
result[s] -= one(T)
|
result[s] -= one(T)
|
||||||
end
|
end
|
||||||
@ -27,7 +27,7 @@ function create_SDP_problem(Δ::GroupRingElem, matrix_constraints; upper_bound=I
|
|||||||
Δ² = Δ*Δ
|
Δ² = Δ*Δ
|
||||||
@assert length(Δ.coeffs) == length(matrix_constraints)
|
@assert length(Δ.coeffs) == length(matrix_constraints)
|
||||||
m = JuMP.Model();
|
m = JuMP.Model();
|
||||||
JuMP.@variable(m, P[1:N, 1:N], SDP)
|
JuMP.@variable(m, P[1:N, 1:N])
|
||||||
JuMP.@SDconstraint(m, P >= 0)
|
JuMP.@SDconstraint(m, P >= 0)
|
||||||
JuMP.@constraint(m, sum(P[i] for i in eachindex(P)) == 0)
|
JuMP.@constraint(m, sum(P[i] for i in eachindex(P)) == 0)
|
||||||
|
|
||||||
|
Loading…
Reference in New Issue
Block a user