Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

some methods for computing orth. discriminants #2748

Merged
merged 3 commits into from
Sep 8, 2023
Merged
Show file tree
Hide file tree
Changes from 2 commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
1 change: 1 addition & 0 deletions experimental/OrthogonalDiscriminants/docs/doc.main
Original file line number Diff line number Diff line change
Expand Up @@ -2,6 +2,7 @@
"Orthogonal discriminants" => [
"introduction.md",
"access.md",
"compute.md",
"misc.md",
],
]
16 changes: 16 additions & 0 deletions experimental/OrthogonalDiscriminants/docs/src/compute.md
Original file line number Diff line number Diff line change
@@ -0,0 +1,16 @@
```@meta
CurrentModule = Oscar.OrthogonalDiscriminants
DocTestSetup = quote
using Oscar
end
```

# Criteria for computing orthogonal discriminants

## Character-theoretical criteria

```@docs
od_from_order
od_from_eigenvalues
od_for_specht_module
```
Original file line number Diff line number Diff line change
Expand Up @@ -17,8 +17,10 @@ import Oscar.Partition
import Oscar.partition

# The following code can be loaded at compile time.
include("utils.jl")
include("data.jl")
include("gram_det.jl")
include("theoretical.jl")
include("exports.jl")

end # module
Expand Down
21 changes: 19 additions & 2 deletions experimental/OrthogonalDiscriminants/src/data.jl
Original file line number Diff line number Diff line change
Expand Up @@ -89,11 +89,28 @@ function orthogonal_discriminants(tbl::Oscar.GAPGroupCharacterTable)
end


function comment_matches(str)
error("dummy function")
# Return the character described by `d`.
function character_of_entry(d::Dict)
@req haskey(d, :groupname) "the dictionary has no :groupname"
tbl = character_table(d[:groupname])
@req haskey(d, :characteristic) "the dictionary has no :characteristic"
p = d[:characteristic]
if p != 0
tbl = mod(tbl, p)
end
@req haskey(d, :charpos) "the dictionary has no :charpos"
return tbl[d[:charpos]]
end


# Return `true` if `str` occurs as an entry in `d[:comment]`
function comment_matches(d::Dict, str::String)
@req haskey(d, :comment) "the dictionary has no :comment"
return str in d[:comment]
end


# Compare two character fields, by comparing their embeddings.
function is_equal_field(emb1, emb2)
dom1 = domain(emb1)
dom2 = domain(emb2)
Expand Down
174 changes: 174 additions & 0 deletions experimental/OrthogonalDiscriminants/src/theoretical.jl
Original file line number Diff line number Diff line change
@@ -0,0 +1,174 @@

# character-theoretical methods


@doc raw"""
od_from_order(chi::GAPGroupClassFunction)

Return `(flag, val)` where `flag` is `true` if the order of the group
of `chi` divides only one of the orders of the two orthogonal groups
[`omega_group`](@ref)`(+/-1, d, q)`, where `d` is the degree of `chi`
and `q` is the order of the field of definition of `chi`.

In this case, `val` is `"O+"` or `"O-"`.

```jldoctest
julia> t = character_table("L3(2)");

julia> Oscar.OrthogonalDiscriminants.od_from_order(mod(t, 3)[4])
(true, "O-")

julia> Oscar.OrthogonalDiscriminants.od_from_order(mod(t, 2)[4])
(false, "")
```
"""
function od_from_order(chi::GAPGroupClassFunction)
characteristic(chi) == 0 && return (false, "")
d = numerator(degree(chi))
q = order_field_of_definition(chi)
tbl = ordinary_table(chi.table)
ord = order(ZZRingElem, tbl)

# Compute the order of the subgroup that shall embed into
# the perfect group `omega_group(epsilon, d, q)`.
n = sum(class_lengths(tbl)[class_positions_of_solvable_residuum(tbl)])
flag1, flag2 = order_omega_mod_N(d, q, n)
if flag1
return flag2 ? (false, "") : (true, "O+")
else
return flag2 ? (true, "O-") : (false, "")
end
end


@doc raw"""
od_from_eigenvalues(chi::GAPGroupClassFunction)

Return `(flag, val)` where `flag` is `true` if there is a conjugacy class
on which representing matrices for `chi` have no eigenvalue $\pm 1$.
In this case, if `chi` is orthogonally stable (this is not checked here)
then `val` is a string that describes the orthogonal discriminant of `chi`.

If `flag` is `false` then `val` is equal to `""`.

This criterion works only if the characteristic of `chi` is not $2$,
`(false, "")` is returned if the characteristic is $2$.

# Examples
```jldoctest
julia> t = character_table("A5");

julia> Oscar.OrthogonalDiscriminants.od_from_eigenvalues(t[4])
(true, "5")

julia> Oscar.OrthogonalDiscriminants.od_from_eigenvalues(mod(t, 3)[4])
(true, "O-")

julia> Oscar.OrthogonalDiscriminants.od_from_eigenvalues(mod(t, 2)[4])
(false, "")
```
"""
function od_from_eigenvalues(chi::GAPGroupClassFunction)
p = characteristic(chi)
p == 2 && return (false, "")

tbl = chi.table
ord = orders_class_representatives(tbl)
for i in 2:length(chi)
n = ord[i]
ev = multiplicities_eigenvalues(chi, i)
if ev[end] != 0 || (iseven(n) && ev[divexact(n, 2)] != 0)
continue
end

F, z = cyclotomic_field(n)
od = prod(x -> x[1]^x[2], [(z^i-z^-i, ev[i]) for i in 1:n])
if mod(degree(chi), 4) == 2
od = -od
end

K, _ = abelian_closure(QQ)
if p == 0
# Coerce `od` into the character field of `chi`.
F, emb = character_field(chi)
od = preimage(emb, K(od))

# Reduce the representative `od` mod obvious squares
# in the character field.
od = reduce_mod_squares(F, od)

# Embed this value into the alg. closure.
od = emb(od)

# Turn the value into a string (using Atlas notation).
str = atlas_description(od)
else
# Decide if the reduction mod `p` is a square in the char. field.
str = is_square(reduce(K(od), character_field(chi)[1])) ? "O+" : "O-"
end

return true, str
end

return false, ""
end

@doc raw"""
od_for_specht_module(chi::GAPGroupClassFunction)

Return `(flag, val)` where `flag` is `true` if `chi` is an ordinary
irreducible character of a symmetric group or of an alternating group
such that `chi` extends to the corresponding symmetric group.
In this case, if `chi` is orthogonally stable (this is not checked here)
then `val` is a string that describes the orthogonal discriminant of `chi`;
the discriminant is computed using the Jantzen-Schaper formula,
via [`gram_determinant_specht_module`](@ref).

`(false, "")` is returned in all cases where this criterion is not applicable.

# Examples
```jldoctest
julia> t = character_table("A5");

julia> Oscar.OrthogonalDiscriminants.od_for_specht_module(t[4])
(true, "5")

julia> Oscar.OrthogonalDiscriminants.od_for_specht_module(mod(t, 3)[4])
(false, "")
```
"""
function od_for_specht_module(chi::GAPGroupClassFunction)
characteristic(chi) == 0 || return (false, "")

# Find out to which alternating or symmetric group `chi` belongs.
tbl = chi.table
name = identifier(tbl)
startswith(name, "A") || return (false, "")
pos = findfirst('.', name)
if pos == nothing
n = parse(Int, name[2:end])
else
n = parse(Int, name[2:(pos-1)])
end
n == nothing && return (false, "")
name == "A$n" || name == "A$n.2" || name == "A6.2_1" || return (false, "")

chipos = findfirst(isequal(chi), tbl)
chipos == nothing && return (false, "")
para = character_parameters(tbl)[chipos]
isa(para, Vector{Int}) || return (false, "")

# Now we know that `chi` belongs to Sym(n) or extends to Sym(n)
gramdet = gram_determinant_specht_module(partition(para))
res = 1
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Is this type stable?

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Better ZZ(1).
(However, for the return value of the function, this does not matter because we call string(res).)

for pair in gramdet
if is_odd(pair[2])
res = res * pair[1]
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
res = res * pair[1]
res *= pair[1]

end
end
if mod(degree(ZZRingElem, chi), 4) == 2
res = - res
end

return true, string(res)
end
103 changes: 103 additions & 0 deletions experimental/OrthogonalDiscriminants/src/utils.jl
Original file line number Diff line number Diff line change
@@ -0,0 +1,103 @@
@doc raw"""
order_omega_mod_N(d::IntegerUnion, q::IntegerUnion, N::IntegerUnion) -> Pair{Bool, Bool}

Return `(flag_plus, flag_minus)` where `flag_plus` and `flag_minus`
are `true` or `false`, depending on whether `N` divides the order
of the orthogonal groups $\Omega^+(d, q)$ and $\Omega^-(d, q)$.

# Examples
```jldoctest
julia> Oscar.OrthogonalDiscriminants.order_omega_mod_N(4, 2, 60)
(false, true)

julia> Oscar.OrthogonalDiscriminants.order_omega_mod_N(4, 5, 60)
(true, true)
```
"""
function order_omega_mod_N(d::IntegerUnion, q::IntegerUnion, N::IntegerUnion)
@req is_even(d) "d must be even"
m = div(d, 2)
exp = 0
while mod(N, q) == 0
exp = exp + 1
fingolfin marked this conversation as resolved.
Show resolved Hide resolved
N = div(N, q)
end
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Actually, I think it might be possible to replace this by

Suggested change
exp = 0
while mod(N, q) == 0
exp = exp + 1
N = div(N, q)
end
exp, N = remove(N, q)

facts = collect(factor(q))
p = facts[1][1]
if mod(N, p) == 0
exp = exp + 1
fingolfin marked this conversation as resolved.
Show resolved Hide resolved
while mod(N, p) == 0
N = div(N, p)
end
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
while mod(N, p) == 0
N = div(N, p)
end
_, N = remove(N, p)

end
if m*(m-1) < exp
# A group of order `N` does not embed in any candidate.
return (false, false)
end

q2 = ZZ(q)^2
q2i = ZZ(1)
for i in 1:(m-1)
q2i = q2 * q2i
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
q2i = q2 * q2i
q2i *= q2

if i == 1 && is_odd(q)
g = gcd(N, div(q2i-1, 2))
else
g = gcd(N, q2i-1)
end
N = div(N, g)
if N == 1
# A group of order N may embed in both candidates.
return (true, true)
end
end

# embeds in + type?, embeds in - type?
return (mod(q^m-1, N) == 0, mod(q^m+1, N) == 0)
end


@doc raw"""
reduce_mod_squares(F::AnticNumberField, val::nf_elem)

Return an element of `F` that is equal to `val` modulo squares in `F`.

If `val` describes an integer then the result corresponds to the
squarefree part of this integer.
Otherwise the coefficients of the result have a squarefree g.c.d.

# Examples
```jldoctest
julia> F, z = cyclotomic_field(4);

julia> Oscar.OrthogonalDiscriminants.reduce_mod_squares(F, 4*z^0)
1

julia> Oscar.OrthogonalDiscriminants.reduce_mod_squares(F, -8*z^0)
-2
```
"""
function reduce_mod_squares(F::AnticNumberField, val::nf_elem)
@req parent(val) === F "val must have parent F"
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Why is F even an argument when the only place using it is this check?

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Yes.
I still have to get used to the fact that the elements carry the relevant context.

is_zero(val) && return val
d = denominator(val)
if ! isone(d)
val = val * d^2
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
val = val * d^2
val *= d^2

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Do you want to say that this syntax is preferable?
If yes then should this be stated in the "Deveoper Style Guide"? (There are quite some places in the code where one could switch to this syntax.)

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

From my POV either is fine, this is a matter of taste. Pick whichever you prefer :-)

end
if is_integer(val)
val = ZZ(val)
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

It may not matter, but just FYI, this breaks type stability (and I think possibly even of the whole function?) for val if it is not already a ZZRingElem.

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I do not understand this comment. In the end, an element of F is returned.

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The code inside this function is unstable: Julia cannot determine a concrete type for val anymore: val was initially an nf_elem, but now is a ZZRingElem -- so the best it can do is to say that val is a Union{nf_elem, ZZRingElem} in this function. You can see this in @code_warntype.

This kind of issue can percolate: Julia may or may not be able to deduce which sign method is called, and thus may or may not be able to predict the type of sgn, and then good, and finally for the return value.

In this case, it seems to be OK as the union splitting code in the Julia compiler can deal with it. But in general I would recommend to avoid re-using a variable for a value of a different type (instead call it val2 or whatever)

sgn = sign(val)
good = filter(x -> is_odd(x[2]), collect(factor(val)))
return F(prod([x[1] for x in good], init = sgn))
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The following might avoid the second allocation. Dunno if it makes things better or not, though...

Suggested change
good = filter(x -> is_odd(x[2]), collect(factor(val)))
return F(prod([x[1] for x in good], init = sgn))
good = [x[1] for x in collect(factor(val)) if is_odd(x[2])]
return F(prod(good, init = sgn))

end
# Just get rid of the square part of the gcd of the coefficients.
c = map(numerator, coefficients(val))
s = 1
for (p, e) in collect(factor(gcd(c)))
if iseven(e)
s = s * p^e
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
s = s * p^e
s *= p^e

elseif e > 1
s = s * p^(e-1)
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
s = s * p^(e-1)
s *= p^(e-1)

end
end
return F(c // s)
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I am confused, isn't c a vector at this point?

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Yes, we want the vector where each entry is divided by s.
Perhaps it is better to divide val by s.

end
2 changes: 2 additions & 0 deletions experimental/OrthogonalDiscriminants/test/runtests.jl
Original file line number Diff line number Diff line change
Expand Up @@ -2,3 +2,5 @@ using Oscar
using Test

include("gram_det.jl")
include("utils.jl")
include("theoretical.jl")
Loading
Loading