From 6a3f7c5868bceb02fddd39d09efaf260a7454e69 Mon Sep 17 00:00:00 2001 From: Claus Fieker Date: Wed, 6 Nov 2024 09:02:53 +0100 Subject: [PATCH 1/9] add test --- test/NumField/NfAbs/NfAbs.jl | 14 ++++++++++++++ 1 file changed, 14 insertions(+) diff --git a/test/NumField/NfAbs/NfAbs.jl b/test/NumField/NfAbs/NfAbs.jl index 8b3460f127..1e32a3c5dc 100644 --- a/test/NumField/NfAbs/NfAbs.jl +++ b/test/NumField/NfAbs/NfAbs.jl @@ -85,3 +85,17 @@ end @test !is_normal(K) @test length(automorphism_list(K)) == 2 end + + +@testset "NumField/NfAbs/MultDep" begin + k, a = wildanger_field(5,13); + zk = lll(maximal_order(k)) + class_group(zk) + h = zk.__attrs[:ClassGrpCtx] + r = vcat(h.R_gen, h.R_rel); + r = [x for x = r if isa(x, AbsSimpleNumFieldElem)] + q = Hecke.syzygies(r) + @test all(isone, evaluate(FacElem(r, q[i, :])) for i=1:nrows(q)) +end + + From 4016d4038374a06c2250c56261a4f33f074e8ae1 Mon Sep 17 00:00:00 2001 From: Claus Fieker Date: Wed, 6 Nov 2024 08:59:14 +0100 Subject: [PATCH 2/9] include --- src/NumField/NfAbs.jl | 1 + 1 file changed, 1 insertion(+) diff --git a/src/NumField/NfAbs.jl b/src/NumField/NfAbs.jl index 4aa2f43f60..93f456e9b3 100644 --- a/src/NumField/NfAbs.jl +++ b/src/NumField/NfAbs.jl @@ -21,3 +21,4 @@ include("NfAbs/MPolyFactor.jl") include("NfAbs/MPolyAbsFact.jl") include("NfAbs/ConjugatesNS.jl") include("NfAbs/Cyclotomic.jl") +include("NfAbs/MultDep.jl") From 6ce5d74dd31a72d487e6d7bfea61f23ff7a1eade Mon Sep 17 00:00:00 2001 From: Claus Fieker Date: Wed, 6 Nov 2024 08:58:43 +0100 Subject: [PATCH 3/9] move into 'proper' position --- src/NumField/NfAbs/MultDep.jl | 497 ++++++++++++++++++++++++++++++++++ 1 file changed, 497 insertions(+) create mode 100644 src/NumField/NfAbs/MultDep.jl diff --git a/src/NumField/NfAbs/MultDep.jl b/src/NumField/NfAbs/MultDep.jl new file mode 100644 index 0000000000..538ae3e0e6 --- /dev/null +++ b/src/NumField/NfAbs/MultDep.jl @@ -0,0 +1,497 @@ +module MultDep + +using Hecke +import Base.* + +""" +Given A[i] elements in K, find matrices I and U s.th. +A^I is a basis for + /Units < K^*/Units +actually a sub-group of the S-unit group for a coprime base +generated from the elements in A. + +A^U is a generating set for the relations: + < Units +and every u s.th. A^u is a unit is in the span of U + +The coprime base is returned as well, the columns of I correspond to the +ordeing of the base. +""" +function syzygies_sunits_mod_units(A::Vector{AbsSimpleNumFieldElem}; use_ge::Bool = false, max_ord::Union{Nothing, AbsSimpleNumFieldOrder}=nothing) + k = parent(A[1]) + @assert all(i->parent(i) === k, A) + if !use_ge + if max_ord === nothing + zk = maximal_order(k) + else + zk = max_ord + end + cp = coprime_base(A, zk) + else + cp = coprime_base(A) + end + sort!(cp, lt = (a,b) -> norm(a) > norm(b)) + M = sparse_matrix(FlintZZ) + for a = A + T = Tuple{Int, ZZRingElem}[] + for i = 1:length(cp) + p = cp[i] + v = valuation(a, p) + if v != 0 + push!(T, (i, ZZRingElem(v))) + end +# isone(I) && break + end + push!(M, sparse_row(FlintZZ, T)) + end + h, t = Hecke.hnf_with_transform(matrix(M)) + h = h[1:rank(h), :] + return h, t[nrows(h)+1:end, :], cp + # THINK! do we want or not... + # - M is naturally sparse, hence it makes sense + # - for this application we need all units, hence the complete + # kernel - which is huge and dense... + # - cp seems to not be used for anything afterwards. + # It will be when we actually create the group, in the DiscLog + h, t = Hecke.hnf_kannan_bachem(M, Val(true), truncate = true) + return h, t, cp +end + +#= +function units(h::SMat, t, b::Vector{AbsSimpleNumFieldElem}) + u = FacElem{AbsSimpleNumFieldElem, AbsSimpleNumField}[] + for i=nrows(h)+1:length(b) + k = [ZZRingElem(0) for i=b] + k[i] = 1 + for i in length(t):-1:1 + Hecke.apply_right!(k, t[i]) + end + push!(u, FacElem(b, k)) + end + return u +end +=# + +#= For Ge's PhD + - smallest common overorder (non-coprime case) + - ideals in non-maximal orders with checking if the operation + worked - if not automatically redoing with larger order +=# +function *(O1::AbsNumFieldOrder, O2::AbsNumFieldOrder) + k = number_field(O1) + @assert k === number_field(O2) + b1 = basis(O1, k) + n = length(b1) + b2 = basis(O2, k) + p = [x*y for (x,y) in Base.Iterators.ProductIterator((b1, b2))] + d = reduce(lcm, [denominator(x) for x = p]) + M = zero_matrix(FlintZZ, n*n, n) + z = ZZRingElem() + for i = 1:n*n + a = p[i]*d + Hecke.elem_to_mat_row!(M, i, z, a) + end + h = hnf(M) + b = AbsSimpleNumFieldElem[] + for i=1:n + push!(b, Hecke.elem_from_mat_row(k, h, i, d)) + end + return Hecke.Order(k, b) +end + +mutable struct GeIdeal + a::AbsNumFieldOrderIdeal + function GeIdeal(a::AbsNumFieldOrderIdeal) + o =order(a) + if o.is_maximal == 1 + return new(a) + end + #keep track of known maximality and use this to speed things up + o = ring_of_multipliers(a) + if o === order(a) + return new(a) + else + return new(extend(a, o)) + end + end + function GeIdeal(a::AbsSimpleNumFieldElem) + o = Hecke.any_order(parent(a)) + d = denominator(a, o) + return new(o(a*d)*o), new(o(d)*o) + end +end + +import Hecke.gcd, Hecke.isone, Hecke.*, Hecke.gcd_into!, Hecke.copy, Hecke.divexact, + Hecke.is_unit, Hecke.coprime_base, Hecke.valuation + +function make_compatible!(a::GeIdeal, b::GeIdeal) + o1 = order(a.a) + o2 = order(b.a) + if o1 === o2 + return + end + o = o1*o2 + a.a = extend(a.a, o) + b.a = extend(b.a, o) +end + +#is_unit, divexact, gcd, isone, *, gcd_into!, copy +isone(a::GeIdeal) = isone(a.a) +is_unit(a::GeIdeal) = isone(a) +copy(a::GeIdeal) = GeIdeal(a.a) + +function gcd(a::GeIdeal, b::GeIdeal) + make_compatible!(a, b) + return GeIdeal(a.a + b.a) +end + +function gcd_into!(a::GeIdeal, b::GeIdeal, c::GeIdeal) + make_compatible!(b, c) + a.a = b.a + c.a + return a +end + +function *(a::GeIdeal, b::GeIdeal) + make_compatible!(a, b) + return GeIdeal(a.a * b.a) +end + +function divexact(a::GeIdeal, b::GeIdeal; check::Bool=true) + make_compatible!(a, b) + return GeIdeal(divexact(a.a, b.a; check=check)) +end + +Hecke.norm(a::GeIdeal) = norm(a.a) +#XXX: not sure if those 2 should get "exported", Ge does not seem to be +# overly fast.... +#TODO: do an integer coprime base first, then refine. Possibly also do +# a p-maximality for all p from the integer coprime base. +#TODO: find examples where Ge is useful +function coprime_base(A::Vector{AbsSimpleNumFieldElem}) + c = Vector{GeIdeal}() + for a = A + n,d = GeIdeal(a) + isone(n) || push!(c, n) + isone(d) || push!(c, d) + end + return Hecke.AbstractAlgebra.coprime_base_steel(c) +end + +function coprime_base(A::Vector{AbsSimpleNumFieldElem}, O::AbsNumFieldOrder) + c = Vector{AbsNumFieldOrderIdeal{AbsSimpleNumField, AbsSimpleNumFieldElem}}() + for a = A + n,d = integral_split(a*O) + isone(n) || push!(c, n) + isone(d) || push!(c, d) + end + return coprime_base(c) +end + + +function valuation(a::AbsSimpleNumFieldElem, p::GeIdeal) + return valuation(a, p.a) +end + +""" +reduce the rows of v modulo the lattice spanned by the rows of M. +M should be LLL reduced. +""" +function size_reduce(M::ZZMatrix, v::ZZMatrix) + s = gram_schmidt_orthogonalisation(QQ.(transpose(M))) + d = diagonal(transpose(s)*s) + for i=1:nrows(v) + for j=nrows(M):-1:1 + x = ZZ(round((v[i:i, :]* s[:, j:j]/d[j])[1,1])) + v[i:i, :] -= x*M[j:j, :] + end + end + return v +end + +""" +A a vector of units in fac-elem form. Find matrices U and V s.th. +A^U is a basis for /Tor +and +A^V is a generating system for the relations of A in Units/Tor +""" +function syzygies_units_mod_tor(A::Vector{FacElem{AbsSimpleNumFieldElem, AbsSimpleNumField}}) + p = next_prime(100) + K = base_ring(parent(A[1])) + m = maximum(degree, keys(factor(GF(p), K.pol).fac)) + while m > 4 + p = next_prime(p) + m = maximum(degree, keys(factor(GF(p), K.pol).fac)) + end + #experimentally, the runtime is dominated by log + u = FacElem{AbsSimpleNumFieldElem, AbsSimpleNumField}[] + prec = 10 + + r1, r2 = signature(K) + r = r1+r2 -1 + n = degree(K) + C = Hecke.qAdicConj(K, p) + la = transpose(matrix(conjugates_log(A[1], C, prec, all = false, flat = true))) + lu = zero_matrix(base_ring(la), 0, n) + uu = [] + indep = Int[] + dep = Int[] + for ia = 1:length(A) + a = A[ia] + while true + @vtime :qAdic 1 la = transpose(matrix(conjugates_log(a, C, prec, all = false, flat = true))) + if iszero(la) + @vtime :qAdic 1 @hassert :qAdic 1 verify_gamma([a], [ZZRingElem(1)], ZZRingElem(p)^prec) + @vprint :qAdic 1 println("torsion found") + gamma = vcat([ZZ(0) for i=1:r+length(uu)], [ZZ(1)]) + push!(uu, (a, gamma)) + push!(dep, ia) + break + end + lv = vcat(lu, la) + #check_precision and change + if false && any(x->precision(x) < prec, lv) + println("loss of precision - not sure what to do") + for i=1:rows(lv) + for j = cols(lv) #seems to not do anything + lv[i, j] = setprecision(lv[i, j], min_p) + @assert precision(lv[i,j]) == min_p + end + end + end + @vtime :qAdic 1 k = kernel(lv, side = :left) + @assert nrows(k) < 2 + if nrows(k) == 0 + @vprint :qAdic 1 "new independent unit found\n" + push!(u, a) + push!(indep, ia) + lu = vcat(lu, la) + @assert length(u) <= r + else # length == 1 extend the module + @vprint :qAdic 1 "dependent unit found, looking for relation\n" + s = QQFieldElem[] + for x in k[1, :] + @vtime :qAdic 1 y = lift_reco(FlintQQ, x, reco = true) + if y === nothing + prec *= 2 + @vprint :qAdic 1 "increase prec to ", prec + lu = matrix([conjugates_log(x, C, prec, all = false, flat = true) for x = u])' + break + end + push!(s, y) + end + if length(s) < ncols(k) + continue + end + d = reduce(lcm, map(denominator, s)) + gamma = ZZRingElem[FlintZZ(x*d)::ZZRingElem for x = s] + @assert reduce(gcd, gamma) == 1 # should be a primitive relation + if !verify_gamma(push!(copy(u), a), gamma, ZZRingElem(p)^prec) + prec *= 2 + @vprint :qAdic 1 "increase prec to ", prec + lu = matrix([conjugates_log(x, C, prec, all = false, flat = true) for x = u])' + continue + end + @assert length(gamma) == length(u)+1 + gamma = vcat(gamma[1:length(u)], [0 for i=length(u)+1:r+length(uu)], [gamma[end]]) + push!(uu, (a, gamma)) + push!(dep, ia) + end + break + end + end + #= + let u_1, .., u_n be units and + has rank s and + r_i in Z^n be such that + prod u_i^r_i = 1 (OK, sum of the logs is zero) + rank = s as well #wrong: + #Input: u^2, u^3, u^7 in rank r = 2 + #U = [-3 0 2 0; -7 0 0 2] + the r_i form a Q-generating basis for the relations. + They are indep as they are growing in length (new columns per new element, non + are torsion, so the final non-zero entries are moving right.) + Essentially, the gamma of above are the r_i + Let [H|0] = [r_i|i]*T be the hnf with trafo, so T in Gl(n, Z) + Then + = <[u_i|i] T> + [r_i|i] * [u_i|i]^t = 0 (by construction) + [r_i|i] T inv(T) [u[i] | i] = 0 + Careful! H may not have full rank, hence the shape is different + [H | 0] [v_i | i] = 0 + so, since H is triangular(!!) v_1, ... v_n-s = 0 + and = + + for the case of n=s+1 this is mostly the "normal" construction. + Note: as a side, the relations do not have to be primitive. + If they are, (and n=s+1), then H = 1 + We deal with that by saturation + =# + + for i=1:length(uu)-1 + append!(uu[i][2], zeros(FlintZZ, length(uu[end][2])-length(uu[i][2]))) + end + if length(uu) == 0 #all torsion + return [], A + else + U = matrix(FlintZZ, length(uu), length(uu[end][2]), reduce(vcat, [x[2] for x = uu])) + U = hcat(U[:, 1:length(u)], U[:, r+1:ncols(U)]) + end + U = saturate(U) + _, U = hnf_with_transform(transpose(U)) + k = base_ring(A[1]) + + U = inv(U) + V = sub(U, 1:nrows(U), 1:ncols(U)-length(u)) #the torsion part + U = sub(U, 1:nrows(U), ncols(U)-length(u)+1:ncols(U)) + #U can be reduced modulo V... + V = lll(transpose(V)) + U = size_reduce(V, transpose(U)) + U = lll(U) + + #so basis is A[indep] cat A[dep] ^U + #rels: A[tor], .. * V + nt = zero_matrix(ZZ, length(A), length(dep) + length(indep)) + for i=1:length(indep) + nt[i, indep[i]] = 1 + end + for i=1:length(dep) + nt[i+length(indep), dep[i]] = 1 + end + rel = nt*transpose(V) + return nt*transpose(U), rel +end + + +function verify_gamma(a::Vector{FacElem{AbsSimpleNumFieldElem, AbsSimpleNumField}}, g::Vector{ZZRingElem}, v::ZZRingElem) + #knowing that sum g[i] log(a[i]) == 0 mod v, prove that prod a[i]^g[i] is + #torsion + #= I claim N(1-a) > v^n for n the field degree: + Let K be one of the p-adic fields involved, set b = a^g + then log(K(b)) = 0 (v = p^l) by assumption + so val(log(K(b))) >= l, but + val(X) = 1/deg(K) val(norm(X)) for p-adics + This is true for all completions involved, and sum degrees is n + =# + + t = prod([a[i]^g[i] for i=1:length(a)]) + # t is either 1 or 1-t is large, norm(1-t) is div. by p^ln + #in this case T2(1-t) is large, by the arguments above: T2>= (np^l)^2=:B + # and, see the bottom, \|Log()\|_2^2 >= 1/4 arcosh((B-2)/2)^2 + B = ArbField(nbits(v)*2)(v)^2 + B = 1/2 *acosh((B-2)/2)^2 + p = Hecke.upper_bound(ZZRingElem, log(B)/log(parent(B)(2))) + @vprint :qAdic 1 "using", p, nbits(v)*2 + b = conjugates_arb_log(t, max(-Int(div(p, 2)), 2)) +# @show B , sum(x*x for x = b), is_torsion_unit(t)[1] + @hassert :qAdic 1 (B > sum(x*x for x = b)) == is_torsion_unit(t)[1] + return B > sum(x*x for x = b) +end + +""" +A padic number a is internally written as + p^v*u +for a unit u given in some precision. +This returns u, v and the precision +""" +function canonical_split(a::PadicFieldElem) + u = ZZ() + ccall((:fmpz_set, Hecke.libflint), Nothing, (Ref{ZZRingElem}, Ref{Int}), u, a.u) + return u, a.v, a.N +end + +#TODO: different name ... +function lift_reco(::QQField, a::PadicFieldElem; reco::Bool = false) + if reco + u, v, N = canonical_split(a) + R = parent(a) + fl, c, d = rational_reconstruction(u, prime(R, N-v)) + !fl && return nothing + + x = FlintQQ(c, d) + if v < 0 + return x//prime(R, -v) + else + return x*prime(R, v) + end + else + return lift(FlintQQ, a) + end +end + +""" +Given torsion units A[i] in factored form, find a generator +for the group as well as a basis for the relations +""" +function syzygies_tor(A::Vector{FacElem{AbsSimpleNumFieldElem, AbsSimpleNumField}}) + K = base_ring(parent(A[1])) + bnd = Int(Hecke._torsion_group_order_divisor(K)) + #bnd is a multiple of the torsion in the field + #so any residue field where bnd divides q-1 can be used to find the + #relations, Tor is a subgroup of k^* in this case. + O = any_order(K) + #XXX: bad, but prime ideals in non-max orders seem a bit dodgy + # in particular extend_easy gives garbage + #= + k, a = wildanger_field(5,13); + O = any_order(k); + P = prime_ideals_over(O, 71)[1] + F, mF = residue_field(O, P) + mF(O(gen(k))) # => 25 + mmF = Hecke.extend_easy(mF, k); + mmF(gen(k)) # => 70 + degree(P) # => 0 + P.prim_elem # => 21 but should be a root... + =# + + O = maximal_order(K) + for p = PrimesSet(bnd, -1, bnd, 1) + if discriminant(O) % p == 0 + continue + end + P = prime_ideals_over(O, p)[1] + if degree(P) > 1 + continue + end + F, mF = Hecke.ResidueFieldSmallDegree1(O, P) + mF = Hecke.extend_easy(mF, K) + a = try + map(mF, A) + catch e + if isa(e, Hecke.BadPrime) + #if some element in the base of the fac elems is not a unit mod P + #use the next prime... + continue + end + rethrow(e) + end + U, mU = unit_group(F) + b = map(pseudo_inv(mU), a) + B = free_abelian_group(length(A)) + h = hom(B, U, b) + k, mk = kernel(h) + i, mi = image(h) + @assert ngens(i) == 1 + return preimage(h, mi(i[1])).coeff, vcat([mk(x).coeff for x = gens(k)]...) + end +end + +""" +Given non-zero elements A[i] in K, find a basis for the relations, returned + as a matrix. +""" +function syzygies(A::Vector{AbsSimpleNumFieldElem}; use_ge::Bool = false, max_ord::Union{Nothing, AbsSimpleNumFieldOrder} = nothing) + _, t, _ = syzygies_sunits_mod_units(A; use_ge, max_ord) + u = [FacElem(A, t[i, :]) for i = 1:nrows(t)] + _, tt = syzygies_units_mod_tor(u) + u = Hecke._transform(u, tt) + _, ttt = syzygies_tor(u) + return ttt*transpose(tt)*t +end + +export syzygies + +end + +using .MultDep + + From d8f5126be735b50d2c6026e40af2eaa7a622e118 Mon Sep 17 00:00:00 2001 From: Claus Fieker Date: Wed, 6 Nov 2024 08:55:03 +0100 Subject: [PATCH 4/9] fix index problem, returned too many elementary divisors Try elementary_divisors(sparse_matrix(identity_matrix(ZZ, 2))) without will give 3 EDs, the last one being 0 --- src/Sparse/UpperTriangular.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/Sparse/UpperTriangular.jl b/src/Sparse/UpperTriangular.jl index 830ff289d4..889c65d955 100644 --- a/src/Sparse/UpperTriangular.jl +++ b/src/Sparse/UpperTriangular.jl @@ -92,7 +92,7 @@ function elementary_divisors(A::SMat{ZZRingElem}) e[i] *= p^v[i] end end - for i = length(e):min(nrows(A), ncols(A)) + for i = length(e)+1:min(nrows(A), ncols(A)) push!(e, ZZRingElem(0)) end return e From d6de3b84b6bc18179a62ae82c1a91bdebb9ead0e Mon Sep 17 00:00:00 2001 From: Claus Fieker Date: Wed, 6 Nov 2024 08:49:54 +0100 Subject: [PATCH 5/9] update MultDep - it seems to be working "fine" now. Needs further work on the Ge stuff and some of the internal interface are bad. Example k, a = wildanger_field(5,13); zk = lll(maximal_order(k)) class_group(zk) h = zk.__attrs[:ClassGrpCtx] r = vcat(h.R_gen, h.R_rel); r = [x for x = r if isa(x, AbsSimpleNumFieldElem)] q = Main.MultDep.syzygies(r) Then the rows of q are the relations --- examples/MultDep.jl | 282 +++++++++++++++++++++++++++++++------------- 1 file changed, 203 insertions(+), 79 deletions(-) diff --git a/examples/MultDep.jl b/examples/MultDep.jl index dfe934b2bf..79006f738d 100644 --- a/examples/MultDep.jl +++ b/examples/MultDep.jl @@ -3,11 +3,29 @@ module MultDep using Hecke import Base.* -function multiplicative_group_mod_units_fac_elem(A::Vector{AbsSimpleNumFieldElem}; use_max_ord::Bool = false) +""" +Given A[i] elements in K, find matrices I and U s.th. +A^I is a basis for + /Units < K^*/Units +actually a sub-group of the S-unit group for a coprime base +generated from the elements in A. + +A^U is a generating set for the relations: + < Units +and every u s.th. A^u is a unit is in the span of U + +The coprime base is returned as well, the columns of I correspond to the +ordeing of the base. +""" +function syzygies_sunits_mod_units(A::Vector{AbsSimpleNumFieldElem}; use_ge::Bool = false, max_ord::Union{Nothing, AbsSimpleNumFieldOrder}=nothing) k = parent(A[1]) @assert all(i->parent(i) === k, A) - if use_max_ord - zk = maximal_order(k) + if !use_ge + if max_ord === nothing + zk = maximal_order(k) + else + zk = max_ord + end cp = coprime_base(A, zk) else cp = coprime_base(A) @@ -26,10 +44,20 @@ function multiplicative_group_mod_units_fac_elem(A::Vector{AbsSimpleNumFieldElem end push!(M, sparse_row(FlintZZ, T)) end + h, t = Hecke.hnf_with_transform(matrix(M)) + h = h[1:rank(h), :] + return h, t[nrows(h)+1:end, :], cp + # THINK! do we want or not... + # - M is naturally sparse, hence it makes sense + # - for this application we need all units, hence the complete + # kernel - which is huge and dense... + # - cp seems to not be used for anything afterwards. + # It will be when we actually create the group, in the DiscLog h, t = Hecke.hnf_kannan_bachem(M, Val(true), truncate = true) return h, t, cp end +#= function units(h::SMat, t, b::Vector{AbsSimpleNumFieldElem}) u = FacElem{AbsSimpleNumFieldElem, AbsSimpleNumField}[] for i=nrows(h)+1:length(b) @@ -42,51 +70,16 @@ function units(h::SMat, t, b::Vector{AbsSimpleNumFieldElem}) end return u end +=# -function unit_group_mod_torsion_fac_elem(O::AbsNumFieldOrder, u::Vector{FacElem{AbsSimpleNumFieldElem, AbsSimpleNumField}}) - U = Hecke._unit_group_init(O) - s = signature(O) - r = s[1] + s[2] - 1 - for y = u - is_tors, p = is_torsion_unit(y, false, U.tors_prec) - U.tors_prec = max(p, U.tors_prec) - if is_tors - continue - end - Hecke._add_unit(U, y) - if length(U.units) >= r - break - end - end - if length(U.units) < r - # maybe use pAdic stuff here... - error("not complete yet") - end - - U.full_rank = true - - U.units = Hecke.reduce(U.units, U.tors_prec) - - for y = u - is_tors, p = is_torsion_unit(y, false, U.tors_prec) - U.tors_prec = max(p, U.tors_prec) - if is_tors - continue - end - x = Hecke.reduce_mod_units([y], U)[1] - is_tors, p = is_torsion_unit(x, false, U.tors_prec) - U.tors_prec = max(p, U.tors_prec) - if is_tors - continue - end - Hecke._add_dependent_unit(U, x) - end - return U -end - +#= For Ge's PhD + - smallest common overorder (non-coprime case) + - ideals in non-maximal orders with checking if the operation + worked - if not automatically redoing with larger order +=# function *(O1::AbsNumFieldOrder, O2::AbsNumFieldOrder) - k = nf(O1) - @assert k === nf(O2) + k = number_field(O1) + @assert k === number_field(O2) b1 = basis(O1, k) n = length(b1) b2 = basis(O2, k) @@ -169,7 +162,11 @@ function divexact(a::GeIdeal, b::GeIdeal; check::Bool=true) end Hecke.norm(a::GeIdeal) = norm(a.a) - +#XXX: not sure if those 2 should get "exported", Ge does not seem to be +# overly fast.... +#TODO: do an integer coprime base first, then refine. Possibly also do +# a p-maximality for all p from the integer coprime base. +#TODO: find examples where Ge is useful function coprime_base(A::Vector{AbsSimpleNumFieldElem}) c = Vector{GeIdeal}() for a = A @@ -177,7 +174,7 @@ function coprime_base(A::Vector{AbsSimpleNumFieldElem}) isone(n) || push!(c, n) isone(d) || push!(c, d) end - return coprime_base(c) + return Hecke.AbstractAlgebra.coprime_base_steel(c) end function coprime_base(A::Vector{AbsSimpleNumFieldElem}, O::AbsNumFieldOrder) @@ -195,8 +192,29 @@ function valuation(a::AbsSimpleNumFieldElem, p::GeIdeal) return valuation(a, p.a) end +""" +reduce the rows of v modulo the lattice spanned by the rows of M. +M should be LLL reduced. +""" +function size_reduce(M::ZZMatrix, v::ZZMatrix) + s = gram_schmidt_orthogonalisation(QQ.(transpose(M))) + d = diagonal(transpose(s)*s) + for i=1:nrows(v) + for j=nrows(M):-1:1 + x = ZZ(round((v[i:i, :]* s[:, j:j]/d[j])[1,1])) + v[i:i, :] -= x*M[j:j, :] + end + end + return v +end -function mult_syzygies_units(A::Vector{FacElem{AbsSimpleNumFieldElem, AbsSimpleNumField}}) +""" +A a vector of units in fac-elem form. Find matrices U and V s.th. +A^U is a basis for /Tor +and +A^V is a generating system for the relations of A in Units/Tor +""" +function syzygies_units_mod_tor(A::Vector{FacElem{AbsSimpleNumFieldElem, AbsSimpleNumField}}) p = next_prime(100) K = base_ring(parent(A[1])) m = maximum(degree, keys(factor(GF(p), K.pol).fac)) @@ -211,16 +229,22 @@ function mult_syzygies_units(A::Vector{FacElem{AbsSimpleNumFieldElem, AbsSimpleN r1, r2 = signature(K) r = r1+r2 -1 n = degree(K) - C = qAdicConj(K, p) - la = matrix(conjugates_log(A[1], C, prec, all = false, flat = true))' + C = Hecke.qAdicConj(K, p) + la = transpose(matrix(conjugates_log(A[1], C, prec, all = false, flat = true))) lu = zero_matrix(base_ring(la), 0, n) uu = [] - for a = A + indep = Int[] + dep = Int[] + for ia = 1:length(A) + a = A[ia] while true - @vtime :qAdic 1 la = matrix(conjugates_log(a, C, prec, all = false, flat = true))' + @vtime :qAdic 1 la = transpose(matrix(conjugates_log(a, C, prec, all = false, flat = true))) if iszero(la) @vtime :qAdic 1 @hassert :qAdic 1 verify_gamma([a], [ZZRingElem(1)], ZZRingElem(p)^prec) @vprint :qAdic 1 println("torsion found") + gamma = vcat([ZZ(0) for i=1:r+length(uu)], [ZZ(1)]) + push!(uu, (a, gamma)) + push!(dep, ia) break end lv = vcat(lu, la) @@ -237,11 +261,13 @@ function mult_syzygies_units(A::Vector{FacElem{AbsSimpleNumFieldElem, AbsSimpleN @vtime :qAdic 1 k = kernel(lv, side = :left) @assert nrows(k) < 2 if nrows(k) == 0 - println("new ") + @vprint :qAdic 1 "new independent unit found\n" push!(u, a) + push!(indep, ia) lu = vcat(lu, la) @assert length(u) <= r else # length == 1 extend the module + @vprint :qAdic 1 "dependent unit found, looking for relation\n" s = QQFieldElem[] for x in k[1, :] @vtime :qAdic 1 y = lift_reco(FlintQQ, x, reco = true) @@ -259,7 +285,7 @@ function mult_syzygies_units(A::Vector{FacElem{AbsSimpleNumFieldElem, AbsSimpleN d = reduce(lcm, map(denominator, s)) gamma = ZZRingElem[FlintZZ(x*d)::ZZRingElem for x = s] @assert reduce(gcd, gamma) == 1 # should be a primitive relation - @time if !verify_gamma(push!(copy(u), a), gamma, ZZRingElem(p)^prec) + if !verify_gamma(push!(copy(u), a), gamma, ZZRingElem(p)^prec) prec *= 2 @vprint :qAdic 1 "increase prec to ", prec lu = matrix([conjugates_log(x, C, prec, all = false, flat = true) for x = u])' @@ -268,6 +294,7 @@ function mult_syzygies_units(A::Vector{FacElem{AbsSimpleNumFieldElem, AbsSimpleN @assert length(gamma) == length(u)+1 gamma = vcat(gamma[1:length(u)], [0 for i=length(u)+1:r+length(uu)], [gamma[end]]) push!(uu, (a, gamma)) + push!(dep, ia) end break end @@ -277,14 +304,19 @@ function mult_syzygies_units(A::Vector{FacElem{AbsSimpleNumFieldElem, AbsSimpleN has rank s and r_i in Z^n be such that prod u_i^r_i = 1 (OK, sum of the logs is zero) - rank = s as well - so the r_i form a Q-basis for the relations. + rank = s as well #wrong: + #Input: u^2, u^3, u^7 in rank r = 2 + #U = [-3 0 2 0; -7 0 0 2] + the r_i form a Q-generating basis for the relations. + They are indep as they are growing in length (new columns per new element, non + are torsion, so the final non-zero entries are moving right.) Essentially, the gamma of above are the r_i Let [H|0] = [r_i|i]*T be the hnf with trafo, so T in Gl(n, Z) Then = <[u_i|i] T> [r_i|i] * [u_i|i]^t = 0 (by construction) [r_i|i] T inv(T) [u[i] | i] = 0 + Careful! H may not have full rank, hence the shape is different [H | 0] [v_i | i] = 0 so, since H is triangular(!!) v_1, ... v_n-s = 0 and = @@ -292,35 +324,44 @@ function mult_syzygies_units(A::Vector{FacElem{AbsSimpleNumFieldElem, AbsSimpleN for the case of n=s+1 this is mostly the "normal" construction. Note: as a side, the relations do not have to be primitive. If they are, (and n=s+1), then H = 1 + We deal with that by saturation =# for i=1:length(uu)-1 append!(uu[i][2], zeros(FlintZZ, length(uu[end][2])-length(uu[i][2]))) end - if length(uu) == 0 - U = matrix(FlintZZ, length(uu), length(uu[end][2]), reduce(vcat, [x[2] for x = uu])) + if length(uu) == 0 #all torsion + return [], A else U = matrix(FlintZZ, length(uu), length(uu[end][2]), reduce(vcat, [x[2] for x = uu])) + U = hcat(U[:, 1:length(u)], U[:, r+1:ncols(U)]) end - _, U = hnf_with_transform(U') - if false - U = inv(U) - V = sub(U, 1:rows(U), 1:cols(U)-length(u)) - U = sub(U, 1:rows(U), cols(U)-length(u)+1:cols(U)) - #U can be reduced modulo V... - Z = zero_matrix(FlintZZ, cols(V), cols(U)) - I = identity_matrix(FlintZZ, cols(U)) * p^(2*prec) - k = base_ring(A[1]) - A = [ Z V'; I U'] - l = lll(A) - U = sub(l, cols(V)+1:rows(l), cols(U)+1:cols(l)) - U = lll(U) - else - U = lll(U') + U = saturate(U) + _, U = hnf_with_transform(transpose(U)) + k = base_ring(A[1]) + + U = inv(U) + V = sub(U, 1:nrows(U), 1:ncols(U)-length(u)) #the torsion part + U = sub(U, 1:nrows(U), ncols(U)-length(u)+1:ncols(U)) + #U can be reduced modulo V... + V = lll(transpose(V)) + U = size_reduce(V, transpose(U)) + U = lll(U) + + #so basis is A[indep] cat A[dep] ^U + #rels: A[tor], .. * V + nt = zero_matrix(ZZ, length(A), length(dep) + length(indep)) + for i=1:length(indep) + nt[i, indep[i]] = 1 end - return Hecke._transform(vcat(u, FacElem{AbsSimpleNumFieldElem,AbsSimpleNumField}[FacElem(k(1)) for i=length(u)+1:r], [x[1] for x = uu]), U') + for i=1:length(dep) + nt[i+length(indep), dep[i]] = 1 + end + rel = nt*transpose(V) + return nt*transpose(U), rel end + function verify_gamma(a::Vector{FacElem{AbsSimpleNumFieldElem, AbsSimpleNumField}}, g::Vector{ZZRingElem}, v::ZZRingElem) #knowing that sum g[i] log(a[i]) == 0 mod v, prove that prod a[i]^g[i] is #torsion @@ -338,7 +379,7 @@ function verify_gamma(a::Vector{FacElem{AbsSimpleNumFieldElem, AbsSimpleNumField # and, see the bottom, \|Log()\|_2^2 >= 1/4 arcosh((B-2)/2)^2 B = ArbField(nbits(v)*2)(v)^2 B = 1/2 *acosh((B-2)/2)^2 - p = Hecke.upper_bound(log(B)/log(parent(B)(2)), ZZRingElem) + p = Hecke.upper_bound(ZZRingElem, log(B)/log(parent(B)(2))) @vprint :qAdic 1 "using", p, nbits(v)*2 b = conjugates_arb_log(t, max(-Int(div(p, 2)), 2)) # @show B , sum(x*x for x = b), is_torsion_unit(t)[1] @@ -346,9 +387,22 @@ function verify_gamma(a::Vector{FacElem{AbsSimpleNumFieldElem, AbsSimpleNumField return B > sum(x*x for x = b) end +""" +A padic number a is internally written as + p^v*u +for a unit u given in some precision. +This returns u, v and the precision +""" +function canonical_split(a::PadicFieldElem) + u = ZZ() + ccall((:fmpz_set, Hecke.libflint), Nothing, (Ref{ZZRingElem}, Ref{Int}), u, a.u) + return u, a.v, a.N +end + +#TODO: different name ... function lift_reco(::QQField, a::PadicFieldElem; reco::Bool = false) if reco - u, v, N = getUnit(a) + u, v, N = canonical_split(a) R = parent(a) fl, c, d = rational_reconstruction(u, prime(R, N-v)) !fl && return nothing @@ -364,6 +418,76 @@ function lift_reco(::QQField, a::PadicFieldElem; reco::Bool = false) end end -Hecke.number_of_rows(A::Matrix{T}) where {T} = size(A)[1] -Hecke.number_of_columns(A::Matrix{T}) where {T} = size(A)[2] +""" +Given torsion units A[i] in factored form, find a generator +for the group as well as a basis for the relations +""" +function syzygies_tor(A::Vector{FacElem{AbsSimpleNumFieldElem, AbsSimpleNumField}}) + K = base_ring(parent(A[1])) + bnd = Int(Hecke._torsion_group_order_divisor(K)) + #bnd is a multiple of the torsion in the field + #so any residue field where bnd divides q-1 can be used to find the + #relations, Tor is a subgroup of k^* in this case. + O = any_order(K) + #XXX: bad, but prime ideals in non-max orders seem a bit dodgy + # in particular extend_easy gives garbage + #= + k, a = wildanger_field(5,13); + O = any_order(k); + P = prime_ideals_over(O, 71)[1] + F, mF = residue_field(O, P) + mF(O(gen(k))) # => 25 + mmF = Hecke.extend_easy(mF, k); + mmF(gen(k)) # => 70 + degree(P) # => 0 + P.prim_elem # => 21 but should be a root... + =# + O = maximal_order(K) + for p = PrimesSet(bnd, -1, bnd, 1) + if discriminant(O) % p == 0 + continue + end + P = prime_ideals_over(O, p)[1] + if degree(P) > 1 + continue + end + F, mF = Hecke.ResidueFieldSmallDegree1(O, P) + mF = Hecke.extend_easy(mF, K) + a = try + map(mF, A) + catch e + if isa(e, Hecke.BadPrime) + #if some element in the base of the fac elems is not a unit mod P + #use the next prime... + continue + end + rethrow(e) + end + U, mU = unit_group(F) + b = map(pseudo_inv(mU), a) + B = free_abelian_group(length(A)) + h = hom(B, U, b) + k, mk = kernel(h) + i, mi = image(h) + @assert ngens(i) == 1 + return preimage(h, mi(i[1])).coeff, vcat([mk(x).coeff for x = gens(k)]...) + end +end + +""" +Given non-zero elements A[i] in K, find a basis for the relations, returned + as a matrix. +""" +function syzygies(A::Vector{AbsSimpleNumFieldElem}; use_ge::Bool = false, max_ord::Union{Nothing, AbsSimpleNumFieldOrder} = nothing) + _, t, _ = syzygies_sunits_mod_units(A; use_ge, max_ord) + u = [FacElem(A, t[i, :]) for i = 1:nrows(t)] + _, tt = syzygies_units_mod_tor(u) + u = Hecke._transform(u, tt) + _, ttt = syzygies_tor(u) + return ttt*transpose(tt)*t +end + +export syzygies + +end From cb8f62df02bf134b954eeb47aa09c164396eb8de Mon Sep 17 00:00:00 2001 From: Claus Fieker Date: Wed, 6 Nov 2024 09:08:19 +0100 Subject: [PATCH 6/9] delete --- examples/MultDep.jl | 493 -------------------------------------------- 1 file changed, 493 deletions(-) delete mode 100644 examples/MultDep.jl diff --git a/examples/MultDep.jl b/examples/MultDep.jl deleted file mode 100644 index 79006f738d..0000000000 --- a/examples/MultDep.jl +++ /dev/null @@ -1,493 +0,0 @@ -module MultDep - -using Hecke -import Base.* - -""" -Given A[i] elements in K, find matrices I and U s.th. -A^I is a basis for - /Units < K^*/Units -actually a sub-group of the S-unit group for a coprime base -generated from the elements in A. - -A^U is a generating set for the relations: - < Units -and every u s.th. A^u is a unit is in the span of U - -The coprime base is returned as well, the columns of I correspond to the -ordeing of the base. -""" -function syzygies_sunits_mod_units(A::Vector{AbsSimpleNumFieldElem}; use_ge::Bool = false, max_ord::Union{Nothing, AbsSimpleNumFieldOrder}=nothing) - k = parent(A[1]) - @assert all(i->parent(i) === k, A) - if !use_ge - if max_ord === nothing - zk = maximal_order(k) - else - zk = max_ord - end - cp = coprime_base(A, zk) - else - cp = coprime_base(A) - end - sort!(cp, lt = (a,b) -> norm(a) > norm(b)) - M = sparse_matrix(FlintZZ) - for a = A - T = Tuple{Int, ZZRingElem}[] - for i = 1:length(cp) - p = cp[i] - v = valuation(a, p) - if v != 0 - push!(T, (i, ZZRingElem(v))) - end -# isone(I) && break - end - push!(M, sparse_row(FlintZZ, T)) - end - h, t = Hecke.hnf_with_transform(matrix(M)) - h = h[1:rank(h), :] - return h, t[nrows(h)+1:end, :], cp - # THINK! do we want or not... - # - M is naturally sparse, hence it makes sense - # - for this application we need all units, hence the complete - # kernel - which is huge and dense... - # - cp seems to not be used for anything afterwards. - # It will be when we actually create the group, in the DiscLog - h, t = Hecke.hnf_kannan_bachem(M, Val(true), truncate = true) - return h, t, cp -end - -#= -function units(h::SMat, t, b::Vector{AbsSimpleNumFieldElem}) - u = FacElem{AbsSimpleNumFieldElem, AbsSimpleNumField}[] - for i=nrows(h)+1:length(b) - k = [ZZRingElem(0) for i=b] - k[i] = 1 - for i in length(t):-1:1 - Hecke.apply_right!(k, t[i]) - end - push!(u, FacElem(b, k)) - end - return u -end -=# - -#= For Ge's PhD - - smallest common overorder (non-coprime case) - - ideals in non-maximal orders with checking if the operation - worked - if not automatically redoing with larger order -=# -function *(O1::AbsNumFieldOrder, O2::AbsNumFieldOrder) - k = number_field(O1) - @assert k === number_field(O2) - b1 = basis(O1, k) - n = length(b1) - b2 = basis(O2, k) - p = [x*y for (x,y) in Base.Iterators.ProductIterator((b1, b2))] - d = reduce(lcm, [denominator(x) for x = p]) - M = zero_matrix(FlintZZ, n*n, n) - z = ZZRingElem() - for i = 1:n*n - a = p[i]*d - Hecke.elem_to_mat_row!(M, i, z, a) - end - h = hnf(M) - b = AbsSimpleNumFieldElem[] - for i=1:n - push!(b, Hecke.elem_from_mat_row(k, h, i, d)) - end - return Hecke.Order(k, b) -end - -mutable struct GeIdeal - a::AbsNumFieldOrderIdeal - function GeIdeal(a::AbsNumFieldOrderIdeal) - o =order(a) - if o.is_maximal == 1 - return new(a) - end - #keep track of known maximality and use this to speed things up - o = ring_of_multipliers(a) - if o === order(a) - return new(a) - else - return new(extend(a, o)) - end - end - function GeIdeal(a::AbsSimpleNumFieldElem) - o = Hecke.any_order(parent(a)) - d = denominator(a, o) - return new(o(a*d)*o), new(o(d)*o) - end -end - -import Hecke.gcd, Hecke.isone, Hecke.*, Hecke.gcd_into!, Hecke.copy, Hecke.divexact, - Hecke.is_unit, Hecke.coprime_base, Hecke.valuation - -function make_compatible!(a::GeIdeal, b::GeIdeal) - o1 = order(a.a) - o2 = order(b.a) - if o1 === o2 - return - end - o = o1*o2 - a.a = extend(a.a, o) - b.a = extend(b.a, o) -end - -#is_unit, divexact, gcd, isone, *, gcd_into!, copy -isone(a::GeIdeal) = isone(a.a) -is_unit(a::GeIdeal) = isone(a) -copy(a::GeIdeal) = GeIdeal(a.a) - -function gcd(a::GeIdeal, b::GeIdeal) - make_compatible!(a, b) - return GeIdeal(a.a + b.a) -end - -function gcd_into!(a::GeIdeal, b::GeIdeal, c::GeIdeal) - make_compatible!(b, c) - a.a = b.a + c.a - return a -end - -function *(a::GeIdeal, b::GeIdeal) - make_compatible!(a, b) - return GeIdeal(a.a * b.a) -end - -function divexact(a::GeIdeal, b::GeIdeal; check::Bool=true) - make_compatible!(a, b) - return GeIdeal(divexact(a.a, b.a; check=check)) -end - -Hecke.norm(a::GeIdeal) = norm(a.a) -#XXX: not sure if those 2 should get "exported", Ge does not seem to be -# overly fast.... -#TODO: do an integer coprime base first, then refine. Possibly also do -# a p-maximality for all p from the integer coprime base. -#TODO: find examples where Ge is useful -function coprime_base(A::Vector{AbsSimpleNumFieldElem}) - c = Vector{GeIdeal}() - for a = A - n,d = GeIdeal(a) - isone(n) || push!(c, n) - isone(d) || push!(c, d) - end - return Hecke.AbstractAlgebra.coprime_base_steel(c) -end - -function coprime_base(A::Vector{AbsSimpleNumFieldElem}, O::AbsNumFieldOrder) - c = Vector{AbsNumFieldOrderIdeal{AbsSimpleNumField, AbsSimpleNumFieldElem}}() - for a = A - n,d = integral_split(a*O) - isone(n) || push!(c, n) - isone(d) || push!(c, d) - end - return coprime_base(c) -end - - -function valuation(a::AbsSimpleNumFieldElem, p::GeIdeal) - return valuation(a, p.a) -end - -""" -reduce the rows of v modulo the lattice spanned by the rows of M. -M should be LLL reduced. -""" -function size_reduce(M::ZZMatrix, v::ZZMatrix) - s = gram_schmidt_orthogonalisation(QQ.(transpose(M))) - d = diagonal(transpose(s)*s) - for i=1:nrows(v) - for j=nrows(M):-1:1 - x = ZZ(round((v[i:i, :]* s[:, j:j]/d[j])[1,1])) - v[i:i, :] -= x*M[j:j, :] - end - end - return v -end - -""" -A a vector of units in fac-elem form. Find matrices U and V s.th. -A^U is a basis for /Tor -and -A^V is a generating system for the relations of A in Units/Tor -""" -function syzygies_units_mod_tor(A::Vector{FacElem{AbsSimpleNumFieldElem, AbsSimpleNumField}}) - p = next_prime(100) - K = base_ring(parent(A[1])) - m = maximum(degree, keys(factor(GF(p), K.pol).fac)) - while m > 4 - p = next_prime(p) - m = maximum(degree, keys(factor(GF(p), K.pol).fac)) - end - #experimentally, the runtime is dominated by log - u = FacElem{AbsSimpleNumFieldElem, AbsSimpleNumField}[] - prec = 10 - - r1, r2 = signature(K) - r = r1+r2 -1 - n = degree(K) - C = Hecke.qAdicConj(K, p) - la = transpose(matrix(conjugates_log(A[1], C, prec, all = false, flat = true))) - lu = zero_matrix(base_ring(la), 0, n) - uu = [] - indep = Int[] - dep = Int[] - for ia = 1:length(A) - a = A[ia] - while true - @vtime :qAdic 1 la = transpose(matrix(conjugates_log(a, C, prec, all = false, flat = true))) - if iszero(la) - @vtime :qAdic 1 @hassert :qAdic 1 verify_gamma([a], [ZZRingElem(1)], ZZRingElem(p)^prec) - @vprint :qAdic 1 println("torsion found") - gamma = vcat([ZZ(0) for i=1:r+length(uu)], [ZZ(1)]) - push!(uu, (a, gamma)) - push!(dep, ia) - break - end - lv = vcat(lu, la) - #check_precision and change - if false && any(x->precision(x) < prec, lv) - println("loss of precision - not sure what to do") - for i=1:rows(lv) - for j = cols(lv) #seems to not do anything - lv[i, j] = setprecision(lv[i, j], min_p) - @assert precision(lv[i,j]) == min_p - end - end - end - @vtime :qAdic 1 k = kernel(lv, side = :left) - @assert nrows(k) < 2 - if nrows(k) == 0 - @vprint :qAdic 1 "new independent unit found\n" - push!(u, a) - push!(indep, ia) - lu = vcat(lu, la) - @assert length(u) <= r - else # length == 1 extend the module - @vprint :qAdic 1 "dependent unit found, looking for relation\n" - s = QQFieldElem[] - for x in k[1, :] - @vtime :qAdic 1 y = lift_reco(FlintQQ, x, reco = true) - if y === nothing - prec *= 2 - @vprint :qAdic 1 "increase prec to ", prec - lu = matrix([conjugates_log(x, C, prec, all = false, flat = true) for x = u])' - break - end - push!(s, y) - end - if length(s) < ncols(k) - continue - end - d = reduce(lcm, map(denominator, s)) - gamma = ZZRingElem[FlintZZ(x*d)::ZZRingElem for x = s] - @assert reduce(gcd, gamma) == 1 # should be a primitive relation - if !verify_gamma(push!(copy(u), a), gamma, ZZRingElem(p)^prec) - prec *= 2 - @vprint :qAdic 1 "increase prec to ", prec - lu = matrix([conjugates_log(x, C, prec, all = false, flat = true) for x = u])' - continue - end - @assert length(gamma) == length(u)+1 - gamma = vcat(gamma[1:length(u)], [0 for i=length(u)+1:r+length(uu)], [gamma[end]]) - push!(uu, (a, gamma)) - push!(dep, ia) - end - break - end - end - #= - let u_1, .., u_n be units and - has rank s and - r_i in Z^n be such that - prod u_i^r_i = 1 (OK, sum of the logs is zero) - rank = s as well #wrong: - #Input: u^2, u^3, u^7 in rank r = 2 - #U = [-3 0 2 0; -7 0 0 2] - the r_i form a Q-generating basis for the relations. - They are indep as they are growing in length (new columns per new element, non - are torsion, so the final non-zero entries are moving right.) - Essentially, the gamma of above are the r_i - Let [H|0] = [r_i|i]*T be the hnf with trafo, so T in Gl(n, Z) - Then - = <[u_i|i] T> - [r_i|i] * [u_i|i]^t = 0 (by construction) - [r_i|i] T inv(T) [u[i] | i] = 0 - Careful! H may not have full rank, hence the shape is different - [H | 0] [v_i | i] = 0 - so, since H is triangular(!!) v_1, ... v_n-s = 0 - and = - - for the case of n=s+1 this is mostly the "normal" construction. - Note: as a side, the relations do not have to be primitive. - If they are, (and n=s+1), then H = 1 - We deal with that by saturation - =# - - for i=1:length(uu)-1 - append!(uu[i][2], zeros(FlintZZ, length(uu[end][2])-length(uu[i][2]))) - end - if length(uu) == 0 #all torsion - return [], A - else - U = matrix(FlintZZ, length(uu), length(uu[end][2]), reduce(vcat, [x[2] for x = uu])) - U = hcat(U[:, 1:length(u)], U[:, r+1:ncols(U)]) - end - U = saturate(U) - _, U = hnf_with_transform(transpose(U)) - k = base_ring(A[1]) - - U = inv(U) - V = sub(U, 1:nrows(U), 1:ncols(U)-length(u)) #the torsion part - U = sub(U, 1:nrows(U), ncols(U)-length(u)+1:ncols(U)) - #U can be reduced modulo V... - V = lll(transpose(V)) - U = size_reduce(V, transpose(U)) - U = lll(U) - - #so basis is A[indep] cat A[dep] ^U - #rels: A[tor], .. * V - nt = zero_matrix(ZZ, length(A), length(dep) + length(indep)) - for i=1:length(indep) - nt[i, indep[i]] = 1 - end - for i=1:length(dep) - nt[i+length(indep), dep[i]] = 1 - end - rel = nt*transpose(V) - return nt*transpose(U), rel -end - - -function verify_gamma(a::Vector{FacElem{AbsSimpleNumFieldElem, AbsSimpleNumField}}, g::Vector{ZZRingElem}, v::ZZRingElem) - #knowing that sum g[i] log(a[i]) == 0 mod v, prove that prod a[i]^g[i] is - #torsion - #= I claim N(1-a) > v^n for n the field degree: - Let K be one of the p-adic fields involved, set b = a^g - then log(K(b)) = 0 (v = p^l) by assumption - so val(log(K(b))) >= l, but - val(X) = 1/deg(K) val(norm(X)) for p-adics - This is true for all completions involved, and sum degrees is n - =# - - t = prod([a[i]^g[i] for i=1:length(a)]) - # t is either 1 or 1-t is large, norm(1-t) is div. by p^ln - #in this case T2(1-t) is large, by the arguments above: T2>= (np^l)^2=:B - # and, see the bottom, \|Log()\|_2^2 >= 1/4 arcosh((B-2)/2)^2 - B = ArbField(nbits(v)*2)(v)^2 - B = 1/2 *acosh((B-2)/2)^2 - p = Hecke.upper_bound(ZZRingElem, log(B)/log(parent(B)(2))) - @vprint :qAdic 1 "using", p, nbits(v)*2 - b = conjugates_arb_log(t, max(-Int(div(p, 2)), 2)) -# @show B , sum(x*x for x = b), is_torsion_unit(t)[1] - @hassert :qAdic 1 (B > sum(x*x for x = b)) == is_torsion_unit(t)[1] - return B > sum(x*x for x = b) -end - -""" -A padic number a is internally written as - p^v*u -for a unit u given in some precision. -This returns u, v and the precision -""" -function canonical_split(a::PadicFieldElem) - u = ZZ() - ccall((:fmpz_set, Hecke.libflint), Nothing, (Ref{ZZRingElem}, Ref{Int}), u, a.u) - return u, a.v, a.N -end - -#TODO: different name ... -function lift_reco(::QQField, a::PadicFieldElem; reco::Bool = false) - if reco - u, v, N = canonical_split(a) - R = parent(a) - fl, c, d = rational_reconstruction(u, prime(R, N-v)) - !fl && return nothing - - x = FlintQQ(c, d) - if v < 0 - return x//prime(R, -v) - else - return x*prime(R, v) - end - else - return lift(FlintQQ, a) - end -end - -""" -Given torsion units A[i] in factored form, find a generator -for the group as well as a basis for the relations -""" -function syzygies_tor(A::Vector{FacElem{AbsSimpleNumFieldElem, AbsSimpleNumField}}) - K = base_ring(parent(A[1])) - bnd = Int(Hecke._torsion_group_order_divisor(K)) - #bnd is a multiple of the torsion in the field - #so any residue field where bnd divides q-1 can be used to find the - #relations, Tor is a subgroup of k^* in this case. - O = any_order(K) - #XXX: bad, but prime ideals in non-max orders seem a bit dodgy - # in particular extend_easy gives garbage - #= - k, a = wildanger_field(5,13); - O = any_order(k); - P = prime_ideals_over(O, 71)[1] - F, mF = residue_field(O, P) - mF(O(gen(k))) # => 25 - mmF = Hecke.extend_easy(mF, k); - mmF(gen(k)) # => 70 - degree(P) # => 0 - P.prim_elem # => 21 but should be a root... - =# - - O = maximal_order(K) - for p = PrimesSet(bnd, -1, bnd, 1) - if discriminant(O) % p == 0 - continue - end - P = prime_ideals_over(O, p)[1] - if degree(P) > 1 - continue - end - F, mF = Hecke.ResidueFieldSmallDegree1(O, P) - mF = Hecke.extend_easy(mF, K) - a = try - map(mF, A) - catch e - if isa(e, Hecke.BadPrime) - #if some element in the base of the fac elems is not a unit mod P - #use the next prime... - continue - end - rethrow(e) - end - U, mU = unit_group(F) - b = map(pseudo_inv(mU), a) - B = free_abelian_group(length(A)) - h = hom(B, U, b) - k, mk = kernel(h) - i, mi = image(h) - @assert ngens(i) == 1 - return preimage(h, mi(i[1])).coeff, vcat([mk(x).coeff for x = gens(k)]...) - end -end - -""" -Given non-zero elements A[i] in K, find a basis for the relations, returned - as a matrix. -""" -function syzygies(A::Vector{AbsSimpleNumFieldElem}; use_ge::Bool = false, max_ord::Union{Nothing, AbsSimpleNumFieldOrder} = nothing) - _, t, _ = syzygies_sunits_mod_units(A; use_ge, max_ord) - u = [FacElem(A, t[i, :]) for i = 1:nrows(t)] - _, tt = syzygies_units_mod_tor(u) - u = Hecke._transform(u, tt) - _, ttt = syzygies_tor(u) - return ttt*transpose(tt)*t -end - -export syzygies - -end From 52c9e9dd3fe9a81aaa0b253c8412b7ce1e4661fa Mon Sep 17 00:00:00 2001 From: Claus Fieker Date: Wed, 6 Nov 2024 10:12:08 +0100 Subject: [PATCH 7/9] reduce exponents --- src/Map/NfOrd.jl | 1 + 1 file changed, 1 insertion(+) diff --git a/src/Map/NfOrd.jl b/src/Map/NfOrd.jl index 5d01396592..021f35fbc1 100644 --- a/src/Map/NfOrd.jl +++ b/src/Map/NfOrd.jl @@ -1087,6 +1087,7 @@ function image(mF::NfToGFMor_easy, a::FacElem{AbsSimpleNumFieldElem, AbsSimpleNu p = mF.defining_pol q = one(Fq) t = mF.t + quo = gcd(quo, order(Fq)-1) for (k, v) = a.fac vv = v if quo != 0 From 5edb8446f1078d9cd2036e78a7c0122b9284d487 Mon Sep 17 00:00:00 2001 From: Claus Fieker Date: Thu, 7 Nov 2024 08:50:03 +0100 Subject: [PATCH 8/9] I hate lin. alg. --- src/NumField/NfAbs/MultDep.jl | 10 +++++++--- 1 file changed, 7 insertions(+), 3 deletions(-) diff --git a/src/NumField/NfAbs/MultDep.jl b/src/NumField/NfAbs/MultDep.jl index 538ae3e0e6..ea8af64577 100644 --- a/src/NumField/NfAbs/MultDep.jl +++ b/src/NumField/NfAbs/MultDep.jl @@ -336,8 +336,11 @@ function syzygies_units_mod_tor(A::Vector{FacElem{AbsSimpleNumFieldElem, AbsSimp U = matrix(FlintZZ, length(uu), length(uu[end][2]), reduce(vcat, [x[2] for x = uu])) U = hcat(U[:, 1:length(u)], U[:, r+1:ncols(U)]) end + U = saturate(U) + _, U = hnf_with_transform(transpose(U)) + k = base_ring(A[1]) U = inv(U) @@ -350,13 +353,14 @@ function syzygies_units_mod_tor(A::Vector{FacElem{AbsSimpleNumFieldElem, AbsSimp #so basis is A[indep] cat A[dep] ^U #rels: A[tor], .. * V - nt = zero_matrix(ZZ, length(A), length(dep) + length(indep)) + nt = zero_matrix(ZZ, length(A), length(A)) for i=1:length(indep) - nt[i, indep[i]] = 1 + nt[indep[i], i] = 1 end for i=1:length(dep) - nt[i+length(indep), dep[i]] = 1 + nt[dep[i], i+length(indep)] = 1 end + @assert matrix([collect(1:length(A))]) * nt == matrix([vcat(indep, dep)]) rel = nt*transpose(V) return nt*transpose(U), rel end From a59230a505a4e9c11c926a11587ebe6cf9438a45 Mon Sep 17 00:00:00 2001 From: Claus Fieker Date: Fri, 8 Nov 2024 09:02:39 +0100 Subject: [PATCH 9/9] deal with 0 properly in lift --- src/NumField/NfAbs/MultDep.jl | 7 +++++-- 1 file changed, 5 insertions(+), 2 deletions(-) diff --git a/src/NumField/NfAbs/MultDep.jl b/src/NumField/NfAbs/MultDep.jl index ea8af64577..2c30422420 100644 --- a/src/NumField/NfAbs/MultDep.jl +++ b/src/NumField/NfAbs/MultDep.jl @@ -274,7 +274,7 @@ function syzygies_units_mod_tor(A::Vector{FacElem{AbsSimpleNumFieldElem, AbsSimp if y === nothing prec *= 2 @vprint :qAdic 1 "increase prec to ", prec - lu = matrix([conjugates_log(x, C, prec, all = false, flat = true) for x = u])' + lu = transpose(matrix([conjugates_log(x, C, prec, all = false, flat = true) for x = u])) break end push!(s, y) @@ -405,6 +405,9 @@ end #TODO: different name ... function lift_reco(::QQField, a::PadicFieldElem; reco::Bool = false) + if iszero(a) + return QQ(0) + end if reco u, v, N = canonical_split(a) R = parent(a) @@ -418,7 +421,7 @@ function lift_reco(::QQField, a::PadicFieldElem; reco::Bool = false) return x*prime(R, v) end else - return lift(FlintQQ, a) + return lift(QQ, a) end end