def __init__(self, field = None):
"""
Create a linear feedback shift cryptosystem.
INPUT: A string monoid over a binary alphabet.
OUTPUT:
EXAMPLES::
sage: E = LFSRCryptosystem(FiniteField(2))
sage: E
LFSR cryptosystem over Finite Field of size 2
TESTS::
sage: E = LFSRCryptosystem(FiniteField(2))
sage: E == loads(dumps(E))
True
TODO: Implement LFSR cryptosystem for arbitrary rings. The current
implementation is limited to the finite field of 2 elements only
because of the dependence on binary strings.
"""
if field is None:
field = FiniteField(2)
if field.cardinality() != 2:
raise NotImplementedError("Not yet implemented.")
S = BinaryStrings()
P = PolynomialRing(FiniteField(2),'x')
SymmetricKeyCryptosystem.__init__(self, S, S, None)
self._field = field
def check_consistency(self, n):
"""
Check that the pseudo-Conway polynomials of degree dividing
`n` in this lattice satisfy the required compatibility
conditions.
EXAMPLES::
sage: from sage.rings.finite_rings.conway_polynomials import PseudoConwayLattice
sage: PCL = PseudoConwayLattice(2, use_database=False)
sage: PCL.check_consistency(6)
sage: PCL.check_consistency(60) # long time
"""
p = self.p
K = FiniteField(p**n, modulus = self.polynomial(n), names='a')
a = K.gen()
for m in n.divisors():
assert (a**((p**n-1)//(p**m-1))).minimal_polynomial() == self.polynomial(m)
def __init__(self, field = None):
"""
Create a shrinking generator cryptosystem.
INPUT: A string monoid over a binary alphabet.
OUTPUT:
EXAMPLES::
sage: E = ShrinkingGeneratorCryptosystem()
sage: E
Shrinking generator cryptosystem over Finite Field of size 2
"""
if field is None:
field = FiniteField(2)
if field.cardinality() != 2:
raise NotImplementedError("Not yet implemented.")
S = BinaryStrings()
SymmetricKeyCryptosystem.__init__(self, S, S, None)
self._field = field
def BIBD_5q_5_for_q_prime_power(q):
r"""
Return a `(5q,5,1)`-BIBD with `q\equiv 1\pmod 4` a prime power.
See Theorem 24 [ClaytonSmith]_.
INPUT:
- ``q`` (integer) -- a prime power such that `q\equiv 1\pmod 4`.
EXAMPLES::
sage: from sage.combinat.designs.bibd import BIBD_5q_5_for_q_prime_power
sage: for q in [25, 45, 65, 85, 125, 145, 185, 205, 305, 405, 605]: # long time
....: _ = BIBD_5q_5_for_q_prime_power(q/5) # long time
"""
from sage.rings.finite_rings.finite_field_constructor import FiniteField
if q%4 != 1 or not is_prime_power(q):
raise ValueError("q is not a prime power or q%4!=1.")
d = (q-1)/4
B = []
F = FiniteField(q,'x')
a = F.primitive_element()
L = {b:i for i,b in enumerate(F)}
for b in L:
B.append([i*q + L[b] for i in range(5)])
for i in range(5):
for j in range(d):
B.append([ i*q + L[b ],
((i+1)%5)*q + L[ a**j+b ],
((i+1)%5)*q + L[-a**j+b ],
((i+4)%5)*q + L[ a**(j+d)+b],
((i+4)%5)*q + L[-a**(j+d)+b],
])
return B
开发者ID:TaraFife,项目名称:sage,代码行数:38,代码来源:bibd.py
示例5: kirkman_triple_system
def kirkman_triple_system(v,existence=False):
r"""
Return a Kirkman Triple System on `v` points.
A Kirkman Triple System `KTS(v)` is a resolvable Steiner Triple System. It
exists if and only if `v\equiv 3\pmod{6}`.
INPUT:
- `n` (integer)
- ``existence`` (boolean; ``False`` by default) -- whether to build the
`KTS(n)` or only answer whether it exists.
.. SEEALSO::
:meth:`IncidenceStructure.is_resolvable`
EXAMPLES:
A solution to Kirkmman's original problem::
sage: kts = designs.kirkman_triple_system(15)
sage: classes = kts.is_resolvable(1)[1]
sage: names = '0123456789abcde'
sage: to_name = lambda (r,s,t): ' '+names[r]+names[s]+names[t]+' '
sage: rows = [' '.join(('Day {}'.format(i) for i in range(1,8)))]
sage: rows.extend(' '.join(map(to_name,row)) for row in zip(*classes))
sage: print '\n'.join(rows)
Day 1 Day 2 Day 3 Day 4 Day 5 Day 6 Day 7
07e 18e 29e 3ae 4be 5ce 6de
139 24a 35b 46c 05d 167 028
26b 03c 14d 257 368 049 15a
458 569 06a 01b 12c 23d 347
acd 7bd 78c 89d 79a 8ab 9bc
TESTS::
sage: for i in range(3,300,6):
....: _ = designs.kirkman_triple_system(i)
"""
if v%6 != 3:
if existence:
return False
raise ValueError("There is no KTS({}) as v!=3 mod(6)".format(v))
if existence:
return False
elif v == 3:
return BalancedIncompleteBlockDesign(3,[[0,1,2]],k=3,lambd=1)
elif v == 9:
classes = [[[0, 1, 5], [2, 6, 7], [3, 4, 8]],
[[1, 6, 8], [3, 5, 7], [0, 2, 4]],
[[1, 4, 7], [0, 3, 6], [2, 5, 8]],
[[4, 5, 6], [0, 7, 8], [1, 2, 3]]]
KTS = BalancedIncompleteBlockDesign(v,[tr for cl in classes for tr in cl],k=3,lambd=1,copy=False)
KTS._classes = classes
return KTS
# Construction 1.1 from [Stinson91] (originally Theorem 6 from [RCW71])
#
# For all prime powers q=1 mod 6, there exists a KTS(2q+1)
elif ((v-1)//2)%6 == 1 and is_prime_power((v-1)//2):
from sage.rings.finite_rings.finite_field_constructor import FiniteField as GF
q = (v-1)//2
K = GF(q,'x')
a = K.primitive_element()
t = (q-1)/6
# m is the solution of a^m=(a^t+1)/2
from sage.groups.generic import discrete_log
m = discrete_log((a**t+1)/2, a)
assert 2*a**m == a**t+1
# First parallel class
first_class = [[(0,1),(0,2),'inf']]
b0 = K.one(); b1 = a**t; b2 = a**m
first_class.extend([(b0*a**i,1),(b1*a**i,1),(b2*a**i,2)]
for i in range(t)+range(2*t,3*t)+range(4*t,5*t))
b0 = a**(m+t); b1=a**(m+3*t); b2=a**(m+5*t)
first_class.extend([[(b0*a**i,2),(b1*a**i,2),(b2*a**i,2)]
for i in range(t)])
# Action of K on the points
action = lambda v,x : (v+x[0],x[1]) if len(x) == 2 else x
# relabel to integer
relabel = {(p,x): i+(x-1)*q
for i,p in enumerate(K)
for x in [1,2]}
relabel['inf'] = 2*q
classes = [[[relabel[action(p,x)] for x in tr] for tr in first_class]
for p in K]
KTS = BalancedIncompleteBlockDesign(v,[tr for cl in classes for tr in cl],k=3,lambd=1,copy=False)
KTS._classes = classes
#.........这里部分代码省略.........
def polynomial(self, n):
r"""
Return the pseudo-Conway polynomial of degree `n` in this
lattice.
INPUT:
- ``n`` -- positive integer
OUTPUT:
- a pseudo-Conway polynomial of degree `n` for the prime `p`.
ALGORITHM:
Uses an algorithm described in [HL99]_, modified to find
pseudo-Conway polynomials rather than Conway polynomials. The
major difference is that we stop as soon as we find a
primitive polynomial.
REFERENCE:
.. [HL99] \L. Heath and N. Loehr (1999). New algorithms for
generating Conway polynomials over finite fields.
Proceedings of the tenth annual ACM-SIAM symposium on
discrete algorithms, pp. 429-437.
EXAMPLES::
sage: from sage.rings.finite_rings.conway_polynomials import PseudoConwayLattice
sage: PCL = PseudoConwayLattice(2, use_database=False)
sage: PCL.polynomial(3)
x^3 + x + 1
sage: PCL.polynomial(4)
x^4 + x^3 + 1
sage: PCL.polynomial(60)
x^60 + x^59 + x^58 + x^55 + x^54 + x^53 + x^52 + x^51 + x^48 + x^46 + x^45 + x^42 + x^41 + x^39 + x^38 + x^37 + x^35 + x^32 + x^31 + x^30 + x^28 + x^24 + x^22 + x^21 + x^18 + x^17 + x^16 + x^15 + x^14 + x^10 + x^8 + x^7 + x^5 + x^3 + x^2 + x + 1
"""
if n in self.nodes:
return self.nodes[n]
p = self.p
n = Integer(n)
if n == 1:
f = self.ring.gen() - FiniteField(p).multiplicative_generator()
self.nodes[1] = f
return f
# Work in an arbitrary field K of order p**n.
K = FiniteField(p**n, names='a')
# TODO: something like the following
# gcds = [n.gcd(d) for d in self.nodes.keys()]
# xi = { m: (...) for m in gcds }
xi = {q: self.polynomial(n//q).any_root(K, -n//q, assume_squarefree=True)
for q in n.prime_divisors()}
# The following is needed to ensure that in the concrete instantiation
# of the "new" extension all previous choices are compatible.
_frobenius_shift(K, xi)
# Construct a compatible element having order the lcm of orders
q, x = xi.popitem()
v = p**(n//q) - 1
for q, xitem in six.iteritems(xi):
w = p**(n//q) - 1
g, alpha, beta = v.xgcd(w)
x = x**beta * xitem**alpha
v = v.lcm(w)
r = p**n - 1
# Get the missing part of the order to be primitive
g = r // v
# Iterate through g-th roots of x until a primitive one is found
z = x.nth_root(g)
root = K.multiplicative_generator()**v
while z.multiplicative_order() != r:
z *= root
# The following should work but tries to create a huge list
# whose length overflows Python's ints for large parameters
#Z = x.nth_root(g, all=True)
#for z in Z:
# if z.multiplicative_order() == r:
# break
f = z.minimal_polynomial()
self.nodes[n] = f
return f
def T2starGeneralizedQuadrangleGraph(q, dual=False, hyperoval=None, field=None, check_hyperoval=True):
r"""
Return the collinearity graph of the generalized quadrangle `T_2^*(q)`, or of its dual
Let `q=2^k` and `\Theta=PG(3,q)`. `T_2^*(q)` is a generalized quadrangle [GQwiki]_
of order `(q-1,q+1)`, see 3.1.3 in [PT09]_. Fix a plane `\Pi \subset \Theta` and a
`hyperoval <http://en.wikipedia.org/wiki/Oval_(projective_plane)#Even_q>`__
`O \subset \Pi`. The points of `T_2^*(q):=T_2^*(O)` are the points of `\Theta`
outside `\Pi`, and the lines are the lines of `\Theta` outside `\Pi`
that meet `\Pi` in a point of `O`.
INPUT:
- ``q`` -- a power of two
- ``dual`` -- if ``False`` (default), return the graph of `T_2^*(O)`.
Otherwise return the graph of the dual `T_2^*(O)`.
- ``hyperoval`` -- a hyperoval (i.e. a complete 2-arc; a set of points in the plane
meeting every line in 0 or 2 points) in the plane of points with 0th coordinate
0 in `PG(3,q)` over the field ``field``. Each point of ``hyperoval`` must be a length 4
vector over ``field`` with 1st non-0 coordinate equal to 1. By default, ``hyperoval`` and
``field`` are not specified, and constructed on the fly. In particular, ``hyperoval``
we build is the classical one, i.e. a conic with the point of intersection of its
tangent lines.
- ``field`` -- an instance of a finite field of order `q`, must be provided
if ``hyperoval`` is provided.
- ``check_hyperoval`` -- (default: ``True``) if ``True``,
check ``hyperoval`` for correctness.
EXAMPLES:
using the built-in construction::
sage: g=graphs.T2starGeneralizedQuadrangleGraph(4); g
T2*(O,4); GQ(3, 5): Graph on 64 vertices
sage: g.is_strongly_regular(parameters=True)
(64, 18, 2, 6)
sage: g=graphs.T2starGeneralizedQuadrangleGraph(4,dual=True); g
T2*(O,4)*; GQ(5, 3): Graph on 96 vertices
sage: g.is_strongly_regular(parameters=True)
(96, 20, 4, 4)
supplying your own hyperoval::
sage: F=GF(4,'b')
sage: O=[vector(F,(0,0,0,1)),vector(F,(0,0,1,0))]+map(lambda x: vector(F, (0,1,x^2,x)),F)
sage: g=graphs.T2starGeneralizedQuadrangleGraph(4, hyperoval=O, field=F); g
T2*(O,4); GQ(3, 5): Graph on 64 vertices
sage: g.is_strongly_regular(parameters=True)
(64, 18, 2, 6)
TESTS::
sage: F=GF(4,'b') # repeating a point...
sage: O=[vector(F,(0,1,0,0)),vector(F,(0,0,1,0))]+map(lambda x: vector(F, (0,1,x^2,x)),F)
sage: graphs.T2starGeneralizedQuadrangleGraph(4, hyperoval=O, field=F)
Traceback (most recent call last):
...
RuntimeError: incorrect hyperoval size
sage: O=[vector(F,(0,1,1,0)),vector(F,(0,0,1,0))]+map(lambda x: vector(F, (0,1,x^2,x)),F)
sage: graphs.T2starGeneralizedQuadrangleGraph(4, hyperoval=O, field=F)
Traceback (most recent call last):
...
RuntimeError: incorrect hyperoval
"""
from sage.combinat.designs.incidence_structures import IncidenceStructure
from sage.combinat.designs.block_design import ProjectiveGeometryDesign as PG
from sage.modules.free_module_element import free_module_element as vector
p, k = is_prime_power(q,get_data=True)
if k==0 or p!=2:
raise ValueError('q must be a power of 2')
if field is None:
F = FiniteField(q, 'a')
else:
F = field
Theta = PG(3, 1, F, point_coordinates=1)
Pi = set(filter(lambda x: x[0]==F.zero(), Theta.ground_set()))
if hyperoval is None:
O = filter(lambda x: x[1]+x[2]*x[3]==0 or (x[1]==1 and x[2]==0 and x[3]==0), Pi)
O = set(O)
else:
map(lambda x: x.set_immutable(), hyperoval)
O = set(hyperoval)
if check_hyperoval:
if len(O) != q+2:
raise RuntimeError("incorrect hyperoval size")
for L in Theta.blocks():
if set(L).issubset(Pi):
if not len(O.intersection(L)) in [0,2]:
raise RuntimeError("incorrect hyperoval")
L = map(lambda z: filter(lambda y: not y in O, z),
filter(lambda x: len(O.intersection(x)) == 1, Theta.blocks()))
if dual:
G = IncidenceStructure(L).intersection_graph()
#.........这里部分代码省略.........
def hadamard_matrix_paleyI(n, normalize=True):
"""
Implements the Paley type I construction.
The Paley type I case corresponds to the case `p \cong 3 \mod{4}` for a
prime `p` (see [Hora]_).
INPUT:
- ``n`` -- the matrix size
- ``normalize`` (boolean) -- whether to normalize the result.
EXAMPLES:
We note that this method by default returns a normalised Hadamard matrix ::
sage: from sage.combinat.matrices.hadamard_matrix import hadamard_matrix_paleyI
sage: hadamard_matrix_paleyI(4)
[ 1 1 1 1]
[ 1 -1 1 -1]
[ 1 -1 -1 1]
[ 1 1 -1 -1]
Otherwise, it returns a skew Hadamard matrix `H`, i.e. `H=S+I`, with
`S=-S^\top` ::
sage: M=hadamard_matrix_paleyI(4, normalize=False); M
[ 1 1 1 1]
[-1 1 1 -1]
[-1 -1 1 1]
[-1 1 -1 1]
sage: S=M-identity_matrix(4); -S==S.T
True
TESTS::
sage: from sage.combinat.matrices.hadamard_matrix import is_hadamard_matrix
sage: test_cases = [x+1 for x in range(100) if is_prime_power(x) and x%4==3]
sage: all(is_hadamard_matrix(hadamard_matrix_paleyI(n),normalized=True,verbose=True)
....: for n in test_cases)
True
sage: all(is_hadamard_matrix(hadamard_matrix_paleyI(n,normalize=False),verbose=True)
....: for n in test_cases)
True
"""
p = n - 1
if not(is_prime_power(p) and (p % 4 == 3)):
raise ValueError("The order %s is not covered by the Paley type I construction." % n)
from sage.rings.finite_rings.finite_field_constructor import FiniteField
K = FiniteField(p,'x')
K_list = list(K)
K_list.insert(0,K.zero())
H = matrix(ZZ, [[(1 if (x-y).is_square() else -1)
for x in K_list]
for y in K_list])
for i in range(n):
H[i,0] = -1
H[0,i] = 1
if normalize:
for i in range(n):
H[i,i] = -1
H = normalise_hadamard(H)
return H
def GDD_4_2(q, existence=False, check=True):
r"""
Return a `(2q,\{4\},\{2\})`-GDD for `q` a prime power with `q\equiv 1\pmod{6}`.
This method implements Lemma VII.5.17 from [BJL99] (p.495).
INPUT:
- ``q`` (integer)
- ``existence`` (boolean) -- instead of building the design, return:
- ``True`` -- meaning that Sage knows how to build the design
- ``Unknown`` -- meaning that Sage does not know how to build the
design, but that the design may exist (see :mod:`sage.misc.unknown`).
- ``False`` -- meaning that the design does not exist.
- ``check`` -- (boolean) Whether to check that output is correct before
returning it. As this is expected to be useless (but we are cautious
guys), you may want to disable it whenever you want speed. Set to ``True``
by default.
EXAMPLE::
sage: from sage.combinat.designs.group_divisible_designs import GDD_4_2
sage: GDD_4_2(7,existence=True)
True
sage: GDD_4_2(7)
Group Divisible Design on 14 points of type 2^7
sage: GDD_4_2(8,existence=True)
Unknown
sage: GDD_4_2(8)
Traceback (most recent call last):
...
NotImplementedError
"""
if q <= 1 or q % 6 != 1 or not is_prime_power(q):
if existence:
return Unknown
raise NotImplementedError
if existence:
return True
from sage.rings.finite_rings.finite_field_constructor import FiniteField as GF
G = GF(q, "x")
w = G.primitive_element()
e = w ** ((q - 1) // 3)
# A first parallel class is defined. G acts on it, which yields all others.
first_class = [[(0, 0), (1, w ** i), (1, e * w ** i), (1, e * e * w ** i)] for i in range((q - 1) // 6)]
label = {p: i for i, p in enumerate(G)}
classes = [[[2 * label[x[1] + g] + (x[0] + j) % 2 for x in S] for S in first_class] for g in G for j in range(2)]
return GroupDivisibleDesign(
2 * q,
groups=[[i, i + 1] for i in range(0, 2 * q, 2)],
blocks=sum(classes, []),
K=[4],
G=[2],
check=check,
copy=False,
)
def DesarguesianProjectivePlaneDesign(n, point_coordinates=True, check=True):
r"""
Return the Desarguesian projective plane of order ``n`` as a 2-design.
The Desarguesian projective plane of order `n` can also be defined as the
projective plane over a field of order `n`. For more information, have a
look at :wikipedia:`Projective_plane`.
INPUT:
- ``n`` -- an integer which must be a power of a prime number
- ``point_coordinates`` (boolean) -- whether to label the points with their
homogeneous coordinates (default) or with integers.
- ``check`` -- (boolean) Whether to check that output is correct before
returning it. As this is expected to be useless (but we are cautious
guys), you may want to disable it whenever you want speed. Set to
``True`` by default.
.. SEEALSO::
:func:`ProjectiveGeometryDesign`
EXAMPLES::
sage: designs.DesarguesianProjectivePlaneDesign(2)
(7,3,1)-Balanced Incomplete Block Design
sage: designs.DesarguesianProjectivePlaneDesign(3)
(13,4,1)-Balanced Incomplete Block Design
sage: designs.DesarguesianProjectivePlaneDesign(4)
(21,5,1)-Balanced Incomplete Block Design
sage: designs.DesarguesianProjectivePlaneDesign(5)
(31,6,1)-Balanced Incomplete Block Design
sage: designs.DesarguesianProjectivePlaneDesign(6)
Traceback (most recent call last):
...
ValueError: the order of a finite field must be a prime power
"""
K = FiniteField(n, 'a')
n2 = n**2
relabel = {x:i for i,x in enumerate(K)}
Kiter = relabel # it is much faster to iterate throug a dict than through
# the finite field K
# we decompose the (equivalence class) of points [x:y:z] of the projective
# plane into an affine plane, an affine line and a point. At the same time,
# we relabel the points with the integers from 0 to n^2 + n as follows:
# - the affine plane is the set of points [x:y:1] (i.e. the third coordinate
# is non-zero) and gets relabeled from 0 to n^2-1
affine_plane = lambda x,y: relabel[x] + n * relabel[y]
# - the affine line is the set of points [x:1:0] (i.e. the third coordinate is
# zero but not the second one) and gets relabeld from n^2 to n^2 + n - 1
line_infinity = lambda x: n2 + relabel[x]
# - the point is [1:0:0] and gets relabeld n^2 + n
point_infinity = n2 + n
blcks = []
# the n^2 lines of the form "x = sy + az"
for s in Kiter:
for a in Kiter:
# points in the affine plane
blcks.append([affine_plane(s*y+a, y) for y in Kiter])
# point at infinity
blcks[-1].append(line_infinity(s))
# the n horizontals of the form "y = az"
for a in Kiter:
# points in the affine plane
blcks.append([affine_plane(x,a) for x in Kiter])
# point at infinity
blcks[-1].append(point_infinity)
# the line at infinity "z = 0"
blcks.append(range(n2,n2+n+1))
if check:
from .designs_pyx import is_projective_plane
if not is_projective_plane(blcks):
raise RuntimeError('There is a problem in the function DesarguesianProjectivePlane')
from .bibd import BalancedIncompleteBlockDesign
B = BalancedIncompleteBlockDesign(n2+n+1, blcks, check=check)
if point_coordinates:
zero = K.zero()
one = K.one()
d = {affine_plane(x,y): (x,y,one)
for x in Kiter
for y in Kiter}
d.update({line_infinity(x): (x,one,zero)
for x in Kiter})
d[n2+n]=(one,zero,zero)
B.relabel(d)
return B
#.........这里部分代码省略.........
returning it. As this is expected to be useless (but we are cautious
guys), you may want to disable it whenever you want speed. Set to
``True`` by default.
EXAMPLES::
sage: H = designs.HughesPlane(9)
sage: H
(91,10,1)-Balanced Incomplete Block Design
We prove in the following computations that the Desarguesian plane ``H`` is
not Desarguesian. Let us consider the two triangles `(0,1,10)` and `(57, 70,
59)`. We show that the intersection points `D_{0,1} \cap D_{57,70}`,
`D_{1,10} \cap D_{70,59}` and `D_{10,0} \cap D_{59,57}` are on the same line
while `D_{0,70}`, `D_{1,59}` and `D_{10,57}` are not concurrent::
sage: blocks = H.blocks()
sage: line = lambda p,q: next(b for b in blocks if p in b and q in b)
sage: b_0_1 = line(0, 1)
sage: b_1_10 = line(1, 10)
sage: b_10_0 = line(10, 0)
sage: b_57_70 = line(57, 70)
sage: b_70_59 = line(70, 59)
sage: b_59_57 = line(59, 57)
sage: set(b_0_1).intersection(b_57_70)
{2}
sage: set(b_1_10).intersection(b_70_59)
{73}
sage: set(b_10_0).intersection(b_59_57)
{60}
sage: line(2, 73) == line(73, 60)
True
sage: b_0_57 = line(0, 57)
sage: b_1_70 = line(1, 70)
sage: b_10_59 = line(10, 59)
sage: p = set(b_0_57).intersection(b_1_70)
sage: q = set(b_1_70).intersection(b_10_59)
sage: p == q
False
TESTS:
Some wrong input::
sage: designs.HughesPlane(5)
Traceback (most recent call last):
...
EmptySetError: No Hughes plane of non-square order exists.
sage: designs.HughesPlane(16)
Traceback (most recent call last):
...
EmptySetError: No Hughes plane of even order exists.
Check that it works for non-prime `q`::
sage: designs.HughesPlane(3**4) # not tested - 10 secs
(6643,82,1)-Balanced Incomplete Block Design
"""
if not q2.is_square():
raise EmptySetError("No Hughes plane of non-square order exists.")
if q2%2 == 0:
raise EmptySetError("No Hughes plane of even order exists.")
q = q2.sqrt()
K = FiniteField(q2, prefix='x')
F = FiniteField(q, prefix='y')
A = q3_minus_one_matrix(F)
A = A.change_ring(K)
m = K.list()
V = VectorSpace(K, 3)
zero = K.zero()
one = K.one()
points = [(x, y, one) for x in m for y in m] + \
[(x, one, zero) for x in m] + \
[(one, zero, zero)]
relabel = {tuple(p):i for i,p in enumerate(points)}
blcks = []
for a in m:
if a not in F or a == 1:
# build L(a)
aa = ~a
l = []
l.append(V((-a, one, zero)))
for x in m:
y = - aa * (x+one)
if not y.is_square():
y *= aa**(q-1)
l.append(V((x, y, one)))
# compute the orbit of L(a)
blcks.append([relabel[normalize_hughes_plane_point(p,q)] for p in l])
for i in range(q2 + q):
l = [A*j for j in l]
blcks.append([relabel[normalize_hughes_plane_point(p,q)] for p in l])
from .bibd import BalancedIncompleteBlockDesign
return BalancedIncompleteBlockDesign(q2**2+q2+1, blcks, check=check)
def BIBD_from_arc_in_desarguesian_projective_plane(n,k,existence=False):
r"""
Returns a `(n,k,1)`-BIBD from a maximal arc in a projective plane.
This function implements a construction from Denniston [Denniston69]_, who
describes a maximal :meth:`arc
<sage.combinat.designs.bibd.BalancedIncompleteBlockDesign.arc>` in a
:func:`Desarguesian Projective Plane
<sage.combinat.designs.block_design.DesarguesianProjectivePlaneDesign>` of
order `2^k`. From two powers of two `n,q` with `n<q`, it produces a
`((n-1)(q+1)+1,n,1)`-BIBD.
INPUT:
- ``n,k`` (integers) -- must be powers of two (among other restrictions).
- ``existence`` (boolean) -- whether to return the BIBD obtained through
this construction (default), or to merely indicate with a boolean return
value whether this method *can* build the requested BIBD.
EXAMPLES:
A `(232,8,1)`-BIBD::
sage: from sage.combinat.designs.bibd import BIBD_from_arc_in_desarguesian_projective_plane
sage: from sage.combinat.designs.bibd import BalancedIncompleteBlockDesign
sage: D = BIBD_from_arc_in_desarguesian_projective_plane(232,8)
sage: BalancedIncompleteBlockDesign(232,D)
(232,8,1)-Balanced Incomplete Block Design
A `(120,8,1)`-BIBD::
sage: D = BIBD_from_arc_in_desarguesian_projective_plane(120,8)
sage: BalancedIncompleteBlockDesign(120,D)
(120,8,1)-Balanced Incomplete Block Design
Other parameters::
sage: all(BIBD_from_arc_in_desarguesian_projective_plane(n,k,existence=True)
....: for n,k in
....: [(120, 8), (232, 8), (456, 8), (904, 8), (496, 16),
....: (976, 16), (1936, 16), (2016, 32), (4000, 32), (8128, 64)])
True
Of course, not all can be built this way::
sage: BIBD_from_arc_in_desarguesian_projective_plane(7,3,existence=True)
False
sage: BIBD_from_arc_in_desarguesian_projective_plane(7,3)
Traceback (most recent call last):
...
ValueError: This function cannot produce a (7,3,1)-BIBD
REFERENCE:
.. [Denniston69] R. H. F. Denniston,
Some maximal arcs in finite projective planes.
Journal of Combinatorial Theory 6, no. 3 (1969): 317-319.
http://dx.doi.org/10.1016/S0021-9800(69)80095-5
"""
q = (n-1)//(k-1)-1
if (k % 2 or
q % 2 or
q <= k or
n != (k-1)*(q+1)+1 or
not is_prime_power(k) or
not is_prime_power(q)):
if existence:
return False
raise ValueError("This function cannot produce a ({},{},1)-BIBD".format(n,k))
if existence:
return True
n = k
# From now on, the code assumes the notations of [Denniston69] for n,q, so
# that the BIBD returned by the method will have the requested parameters.
from sage.rings.finite_rings.finite_field_constructor import FiniteField as GF
from sage.libs.gap.libgap import libgap
from sage.matrix.constructor import Matrix
K = GF(q,'a')
one = K.one()
# An irreducible quadratic form over K[X,Y]
GO = libgap.GeneralOrthogonalGroup(-1,2,q)
M = libgap.InvariantQuadraticForm(GO)['matrix']
M = Matrix(M)
M = M.change_ring(K)
Q = lambda xx,yy : M[0,0]*xx**2+(M[0,1]+M[1,0])*xx*yy+M[1,1]*yy**2
# Here, the additive subgroup H (of order n) of K mentioned in
# [Denniston69] is the set of all elements of K of degree < log_n
# (seeing elements of K as polynomials in 'a')
K_iter = list(K) # faster iterations
log_n = is_prime_power(n,get_data=True)[1]
#.........这里部分代码省略.........
开发者ID:TaraFife,项目名称:sage,代码行数:101,代码来源:bibd.py
示例14: v_4_1_rbibd
def v_4_1_rbibd(v,existence=False):
r"""
Return a `(v,4,1)`-RBIBD.
INPUT:
- `n` (integer)
- ``existence`` (boolean; ``False`` by default) -- whether to build the
design or only answer whether it exists.
.. SEEALSO::
- :meth:`IncidenceStructure.is_resolvable`
- :func:`resolvable_balanced_incomplete_block_design`
.. NOTE::
A resolvable `(v,4,1)`-BIBD exists whenever `1\equiv 4\pmod(12)`. This
function, however, only implements a construction of `(v,4,1)`-BIBD such
that `v=3q+1\equiv 1\pmod{3}` where `q` is a prime power (see VII.7.5.a
from [BJL99]_).
EXAMPLE::
sage: rBIBD = designs.resolvable_balanced_incomplete_block_design(28,4)
sage: rBIBD.is_resolvable()
True
sage: rBIBD.is_t_design(return_parameters=True)
(True, (2, 28, 4, 1))
TESTS::
sage: for q in prime_powers(2,30):
....: if (3*q+1)%12 == 4:
....: _ = designs.resolvable_balanced_incomplete_block_design(3*q+1,4) # indirect doctest
"""
# Volume 1, VII.7.5.a from [BJL99]_
if v%3 != 1 or not is_prime_power((v-1)//3):
if existence:
return Unknown
raise NotImplementedError("I don't know how to build a ({},{},1)-RBIBD!".format(v,4))
from sage.rings.finite_rings.finite_field_constructor import FiniteField as GF
q = (v-1)//3
nn = (q-1)//4
G = GF(q,'x')
w = G.primitive_element()
e = w**(nn)
assert e**2 == -1
first_class = [[(w**i,j),(-w**i,j),(e*w**i,j+1),(-e*w**i,j+1)]
for i in range(nn) for j in range(3)]
first_class.append([(0,0),(0,1),(0,2),'inf'])
label = {p:i for i,p in enumerate(G)}
classes = [[[v-1 if x=='inf' else (x[1]%3)*q+label[x[0]+g] for x in S]
for S in first_class]
for g in G]
BIBD = BalancedIncompleteBlockDesign(v,
blocks = sum(classes,[]),
k=4,
check=True,
copy=False)
BIBD._classes = classes
assert BIBD.is_resolvable()
return BIBD
def DuadicCodeOddPair(F,S1,S2):
"""
Constructs the "odd pair" of duadic codes associated to the
"splitting" S1, S2 of n.
.. warning::
Maybe the splitting should be associated to a sum of
q-cyclotomic cosets mod n, where q is a *prime*.
EXAMPLES::
sage: from sage.coding.code_constructions import _is_a_splitting
sage: n = 11; q = 3
sage: C = Zmod(n).cyclotomic_cosets(q); C
[[0], [1, 3, 4, 5, 9], [2, 6, 7, 8, 10]]
sage: S1 = C[1]
sage: S2 = C[2]
sage: _is_a_splitting(S1,S2,11)
True
sage: codes.DuadicCodeOddPair(GF(q),S1,S2)
([11, 6] Cyclic Code over GF(3),
[11, 6] Cyclic Code over GF(3))
This is consistent with Theorem 6.1.3 in [HP2003]_.
"""
from .cyclic_code import CyclicCode
n = len(S1) + len(S2) + 1
if not _is_a_splitting(S1,S2,n):
raise TypeError("%s, %s must be a splitting of %s."%(S1,S2,n))
q = F.order()
k = Mod(q,n).multiplicative_order()
FF = GF(q**k,"z")
z = FF.gen()
zeta = z**((q**k-1)/n)
P1 = PolynomialRing(FF,"x")
x = P1.gen()
g1 = prod([x-zeta**i for i in S1+[0]])
g2 = prod([x-zeta**i for i in S2+[0]])
j = sum([x**i/n for i in range(n)])
P2 = PolynomialRing(F,"x")
x = P2.gen()
coeffs1 = [_lift2smallest_field(c)[0] for c in (g1+j).coefficients(sparse=False)]
coeffs2 = [_lift2smallest_field(c)[0] for c in (g2+j).coefficients(sparse=False)]
gg1 = P2(coeffs1)
gg2 = P2(coeffs2)
gg1 = gcd(gg1, x**n - 1)
gg2 = gcd(gg2, x**n - 1)
C1 = CyclicCode(length = n, generator_pol = gg1)
C2 = CyclicCode(length = n, generator_pol = gg2)
return C1,C2
请发表评论