Difference families¶
This module gathers everything related to difference families. One can build a
difference family (or check that it can be built) with difference_family()
:
sage: G,F = designs.difference_family(13,4,1)
It defines the following functions:
is_difference_family() |
Check if the input is a (k , l )-difference family. |
difference_family() |
Return a (k , l )-difference family on an Abelian group of size v . |
radical_difference_family() |
Return a radical difference family. |
radical_difference_set() |
Return a radical difference set. |
singer_difference_set() |
Return a difference set associated to hyperplanes in a projective space. |
df_q_6_1() |
Return a difference family with parameter ![]() |
one_radical_difference_family() |
Return a radical difference family using an exhaustive search. |
twin_prime_powers_difference_set() |
Return a twin prime powers difference family. |
REFERENCES:
[BJL99-1] | T. Beth, D. Jungnickel, H. Lenz “Design theory Vol. I.” Second edition. Encyclopedia of Mathematics and its Applications, 69. Cambridge University Press, (1999). |
[BLJ99-2] | T. Beth, D. Jungnickel, H. Lenz “Design theory Vol. II.” Second edition. Encyclopedia of Mathematics and its Applications, 78. Cambridge University Press, (1999). |
[Bo39] | R. C. Bose, “On the construction of balanced incomplete block designs”, Ann. Eugenics, vol. 9, (1939), 353–399. |
[Bu95] | (1, 2) M. Buratti “On simple radical difference families”, J. of Combinatorial Designs, vol. 3, no. 2 (1995) |
[Wi72] | (1, 2, 3) R. M. Wilson “Cyclotomy and difference families in elementary Abelian groups”, J. of Num. Th., 4 (1972), pp. 17-47. |
Functions¶
-
sage.combinat.designs.difference_family.
block_stabilizer
(G, B)¶ Compute the left stabilizer of the block
B
under the action ofG
.This function return the list of all
such that
(as a set).
INPUT:
G
– a group (additive or multiplicative).B
– a subset ofG
.
EXAMPLES:
sage: from sage.combinat.designs.difference_family import block_stabilizer sage: Z8 = Zmod(8) sage: block_stabilizer(Z8, [Z8(0),Z8(2),Z8(4),Z8(6)]) [0, 2, 4, 6] sage: block_stabilizer(Z8, [Z8(0),Z8(2)]) [0] sage: C = cartesian_product([Zmod(4),Zmod(3)]) sage: block_stabilizer(C, [C((0,0)),C((2,0)),C((0,1)),C((2,1))]) [(0, 0), (2, 0)] sage: b = map(Zmod(45),[1, 3, 7, 10, 22, 25, 30, 35, 37, 38, 44]) sage: block_stabilizer(Zmod(45),b) [0]
-
sage.combinat.designs.difference_family.
df_q_6_1
(K, existence=False, check=True)¶ Return a
-difference family over the finite field
.
The construction uses Theorem 11 of [Wi72].
EXAMPLES:
sage: from sage.combinat.designs.difference_family import is_difference_family, df_q_6_1 sage: prime_powers = [v for v in xrange(31,500,30) if is_prime_power(v)] sage: parameters = [v for v in prime_powers if df_q_6_1(GF(v,'a'), existence=True)] sage: print parameters [31, 151, 181, 211, 241, 271, 331, 361, 421] sage: for v in parameters: ....: K = GF(v, 'a') ....: df = df_q_6_1(K, check=True) ....: assert is_difference_family(K, df, v, 6, 1)
-
sage.combinat.designs.difference_family.
difference_family
(v, k, l=1, existence=False, explain_construction=False, check=True)¶ Return a (
k
,l
)-difference family on an Abelian group of cardinalityv
.Let
be a finite Abelian group. For a given subset
of
, we define
to be the multi-set of differences
. A
-difference family is a collection of
-subsets of
,
such that the union of the difference sets
for
, seen as a multi-set, contains each element of
exactly
-times.
When there is only one block, i.e.
, then a
-difference family is also called a difference set.
See also Wikipedia article Difference_set.
If there is no such difference family, an
EmptySetError
is raised and if there is no construction at the momentNotImplementedError
is raised.INPUT:
v,k,l
– parameters of the difference family. Ifl
is not provided it is assumed to be1
.existence
– ifTrue
, then return eitherTrue
if Sage knows how to build such design,Unknown
if it does not andFalse
if it knows that the design does not exist.explain_construction
– instead of returning a difference family, returns a string that explains the construction used.check
– boolean (default:True
). IfTrue
then the result of the computation is checked before being returned. This should not be needed but ensures that the output is correct.
OUTPUT:
A pair
(G,D)
made of a groupand a difference family
on that group. Or, if
existence
isTrue
a troolean or ifexplain_construction
isTrue
a string.EXAMPLES:
sage: G,D = designs.difference_family(73,4) sage: G Finite Field of size 73 sage: D [[0, 1, 5, 18], [0, 3, 15, 54], [0, 9, 45, 16], [0, 27, 62, 48], [0, 8, 40, 71], [0, 24, 47, 67]] sage: print designs.difference_family(73, 4, explain_construction=True) The database contains a (73,4)-evenly distributed set sage: G,D = designs.difference_family(15,7,3) sage: G The cartesian product of (Finite Field of size 3, Finite Field of size 5) sage: D [[(1, 1), (1, 4), (2, 2), (2, 3), (0, 0), (1, 0), (2, 0)]] sage: print designs.difference_family(15,7,3,explain_construction=True) Twin prime powers difference family sage: print designs.difference_family(91,10,1,explain_construction=True) Singer difference set
For
we look at the set of small prime powers for which a construction is available:
sage: def prime_power_mod(r,m): ....: k = m+r ....: while True: ....: if is_prime_power(k): ....: yield k ....: k += m sage: from itertools import islice sage: l6 = {True:[], False: [], Unknown: []} sage: for q in islice(prime_power_mod(1,30), 60): ....: l6[designs.difference_family(q,6,existence=True)].append(q) sage: l6[True] [31, 121, 151, 181, 211, ..., 3061, 3121, 3181] sage: l6[Unknown] [61] sage: l6[False] [] sage: l7 = {True: [], False: [], Unknown: []} sage: for q in islice(prime_power_mod(1,42), 60): ....: l7[designs.difference_family(q,7,existence=True)].append(q) sage: l7[True] [169, 337, 379, 421, 463, 547, 631, 673, 757, 841, 883, 967, ..., 4621, 4957, 5167] sage: l7[Unknown] [43, 127, 211, 2017, 2143, 2269, 2311, 2437, 2521, 2647, ..., 4999, 5041, 5209] sage: l7[False] []
List available constructions:
sage: for v in xrange(2,100): ....: constructions = [] ....: for k in xrange(2,10): ....: for l in xrange(1,10): ....: if designs.difference_family(v,k,l,existence=True): ....: constructions.append((k,l)) ....: _ = designs.difference_family(v,k,l) ....: if constructions: ....: print "%2d: %s"%(v, ', '.join('(%d,%d)'%(k,l) for k,l in constructions)) 3: (2,1) 4: (3,2) 5: (2,1), (4,3) 6: (5,4) 7: (2,1), (3,1), (3,2), (4,2), (6,5) 8: (7,6) 9: (2,1), (4,3), (8,7) 10: (9,8) 11: (2,1), (4,6), (5,2), (5,4), (6,3) 13: (2,1), (3,1), (3,2), (4,1), (4,3), (5,5), (6,5) 15: (3,1), (4,6), (5,6), (7,3) 16: (3,2), (5,4), (6,2) 17: (2,1), (4,3), (5,5), (8,7) 19: (2,1), (3,1), (3,2), (4,2), (6,5), (9,4), (9,8) 21: (3,1), (4,3), (5,1), (6,3), (6,5) 22: (4,2), (6,5), (7,4), (8,8) 23: (2,1) 25: (2,1), (3,1), (3,2), (4,1), (4,3), (6,5), (7,7), (8,7) 27: (2,1), (3,1) 28: (3,2), (6,5) 29: (2,1), (4,3), (7,3), (7,6), (8,4), (8,6) 31: (2,1), (3,1), (3,2), (4,2), (5,2), (5,4), (6,1), (6,5) 33: (3,1), (5,5), (6,5) 34: (4,2) 35: (5,2) 37: (2,1), (3,1), (3,2), (4,1), (4,3), (6,5), (9,2), (9,8) 39: (3,1), (6,5) 40: (3,2), (4,1) 41: (2,1), (4,3), (5,1), (5,4), (6,3), (8,7) 43: (2,1), (3,1), (3,2), (4,2), (6,5), (7,2), (7,3), (7,6), (8,4) 45: (3,1), (5,1) 46: (4,2), (6,2) 47: (2,1) 49: (2,1), (3,1), (3,2), (4,1), (4,3), (6,5), (8,7), (9,3) 51: (3,1), (5,2), (6,3) 52: (4,1) 53: (2,1), (4,3) 55: (3,1), (9,4) 57: (3,1), (7,3), (8,1) 59: (2,1) 61: (2,1), (3,1), (3,2), (4,1), (4,3), (5,1), (5,4), (6,2), (6,3), (6,5) 63: (3,1) 64: (3,2), (4,1), (7,2), (7,6), (9,8) 65: (5,1) 67: (2,1), (3,1), (3,2), (6,5) 69: (3,1) 71: (2,1), (5,2), (5,4), (7,3), (7,6), (8,4) 73: (2,1), (3,1), (3,2), (4,1), (4,3), (6,5), (8,7), (9,1), (9,8) 75: (3,1), (5,2) 76: (4,1) 79: (2,1), (3,1), (3,2), (6,5) 81: (2,1), (3,1), (4,3), (5,1), (5,4), (8,7) 83: (2,1) 85: (4,1), (7,2), (7,3), (8,2) 89: (2,1), (4,3), (8,7) 91: (6,1) 97: (2,1), (3,1), (3,2), (4,1), (4,3), (6,5), (8,7), (9,3)
TESTS:
Check more of the Wilson constructions from [Wi72]:
sage: Q5 = [241, 281,421,601,641, 661, 701, 821,881] sage: Q9 = [73, 1153, 1873, 2017] sage: Q15 = [76231] sage: Q4 = [13, 73, 97, 109, 181, 229, 241, 277, 337, 409, 421, 457] sage: Q8 = [1009, 3137, 3697] sage: for Q,k in [(Q4,4),(Q5,5),(Q8,8),(Q9,9),(Q15,15)]: ....: for q in Q: ....: assert designs.difference_family(q,k,1,existence=True) is True ....: _ = designs.difference_family(q,k,1)
Check Singer difference sets:
sage: sgp = lambda q,d: ((q**(d+1)-1)//(q-1), (q**d-1)//(q-1), (q**(d-1)-1)//(q-1)) sage: for q in range(2,10): ....: if is_prime_power(q): ....: for d in [2,3,4]: ....: v,k,l = sgp(q,d) ....: assert designs.difference_family(v,k,l,existence=True) is True ....: _ = designs.difference_family(v,k,l)
Check twin primes difference sets:
sage: for p in [3,5,7,9,11]: ....: v = p*(p+2); k = (v-1)/2; lmbda = (k-1)/2 ....: G,D = designs.difference_family(v,k,lmbda)
Check the database:
sage: from sage.combinat.designs.database import DF,EDS sage: for v,k,l in DF: ....: assert designs.difference_family(v,k,l,existence=True) is True ....: df = designs.difference_family(v,k,l,check=True) sage: for k in EDS: ....: for v in EDS[k]: ....: assert designs.difference_family(v,k,1,existence=True) is True ....: df = designs.difference_family(v,k,1,check=True)
Check a failing construction (trac ticket #17528):
sage: designs.difference_family(9,3) Traceback (most recent call last): ... NotImplementedError: No construction available for (9,3,1)-difference family
Todo
Implement recursive constructions from Buratti “Recursive for difference matrices and relative difference families” (1998) and Jungnickel “Composition theorems for difference families and regular planes” (1978)
-
sage.combinat.designs.difference_family.
group_law
(G)¶ Return a triple
(identity, operation, inverse)
that define the operations on the groupG
.EXAMPLES:
sage: from sage.combinat.designs.difference_family import group_law sage: group_law(Zmod(3)) (0, <built-in function add>, <built-in function neg>) sage: group_law(SymmetricGroup(5)) ((), <built-in function mul>, <built-in function inv>) sage: group_law(VectorSpace(QQ,3)) ((0, 0, 0), <built-in function add>, <built-in function neg>)
-
sage.combinat.designs.difference_family.
is_difference_family
(G, D, v=None, k=None, l=None, verbose=False)¶ Check wether
D
forms a difference family in the groupG
.INPUT:
G
- group of cardinalityv
D
- a set ofk
-subsets ofG
v
,k
andl
- optional parameters of the difference familyverbose
- whether to print additional information
See also
EXAMPLES:
sage: from sage.combinat.designs.difference_family import is_difference_family sage: G = Zmod(21) sage: D = [[0,1,4,14,16]] sage: is_difference_family(G, D, 21, 5) True sage: G = Zmod(41) sage: D = [[0,1,4,11,29],[0,2,8,17,21]] sage: is_difference_family(G, D, verbose=True) Too few: 5 is obtained 0 times in blocks [] 14 is obtained 0 times in blocks [] 27 is obtained 0 times in blocks [] 36 is obtained 0 times in blocks [] Too much: 4 is obtained 2 times in blocks [0, 1] 13 is obtained 2 times in blocks [0, 1] 28 is obtained 2 times in blocks [0, 1] 37 is obtained 2 times in blocks [0, 1] False sage: D = [[0,1,4,11,29],[0,2,8,17,22]] sage: is_difference_family(G, D) True sage: G = Zmod(61) sage: D = [[0,1,3,13,34],[0,4,9,23,45],[0,6,17,24,32]] sage: is_difference_family(G, D) True sage: G = AdditiveAbelianGroup([3]*4) sage: a,b,c,d = G.gens() sage: D = [[d, -a+d, -c+d, a-b-d, b+c+d], ....: [c, a+b-d, -b+c, a-b+d, a+b+c], ....: [-a-b+c+d, a-b-c-d, -a+c-d, b-c+d, a+b], ....: [-b-d, a+b+d, a-b+c-d, a-b+c, -b+c+d]] sage: is_difference_family(G, D) True
The following example has a third block with a non-trivial stabilizer:
sage: G = Zmod(15) sage: D = [[0,1,4],[0,2,9],[0,5,10]] sage: is_difference_family(G,D,verbose=True) It is a (15,3,1)-difference family True
The function also supports multiplicative groups (non necessarily Abelian):
sage: G = DihedralGroup(8) sage: x,y = G.gens() sage: i = G.one() sage: D1 = [[i,x,x^4], [i,x^2, y*x], [i,x^5,y], [i,x^6,y*x^2], [i,x^7,y*x^5]] sage: is_difference_family(G, D1, 16, 3, 2) True sage: from sage.combinat.designs.bibd import BIBD_from_difference_family sage: bibd = BIBD_from_difference_family(G,D1,lambd=2)
TESTS:
sage: K = GF(3^2,'z') sage: z = K.gen() sage: D = [[1,z+1,2]] sage: _ = is_difference_family(K, D, verbose=True) the number of differences (=6) must be a multiple of v-1=8 sage: _ False
-
sage.combinat.designs.difference_family.
one_cyclic_tiling
(A, n)¶ Given a subset
A
of the cyclic additive groupreturn another subset
so that
and
(i.e. any element of
is uniquely expressed as a sum
with
in
and
in
).
EXAMPLES:
sage: from sage.combinat.designs.difference_family import one_cyclic_tiling sage: tile = [0,2,4] sage: m = one_cyclic_tiling(tile,6); m [0, 3] sage: sorted((i+j)%6 for i in tile for j in m) [0, 1, 2, 3, 4, 5] sage: def print_tiling(tile, translat, n): ....: for x in translat: ....: print ''.join('X' if (i-x)%n in tile else '.' for i in range(n)) sage: tile = [0, 1, 2, 7] sage: m = one_cyclic_tiling(tile, 12) sage: print_tiling(tile, m, 12) XXX....X.... ....XXX....X ...X....XXX. sage: tile = [0, 1, 5] sage: m = one_cyclic_tiling(tile, 12) sage: print_tiling(tile, m, 12) XX...X...... ...XX...X... ......XX...X ..X......XX. sage: tile = [0, 2] sage: m = one_cyclic_tiling(tile, 8) sage: print_tiling(tile, m, 8) X.X..... ....X.X. .X.X.... .....X.X
ALGORITHM:
Uses dancing links
sage.combinat.dlx
-
sage.combinat.designs.difference_family.
one_radical_difference_family
(K, k)¶ Search for a radical difference family on
K
using dancing links algorithm.For the definition of radical difference family, see
radical_difference_family()
. Here, we consider only radical difference family with.
INPUT:
K
– a finite field of cardinality.
k
– a positive integer so thatdivides
.
OUTPUT:
Either a difference family or
None
if it does not exist.ALGORITHM:
The existence of a radical difference family is equivalent to a one dimensional tiling (or packing) problem in a cyclic group. This subsequent problem is solved by a call to the function
one_cyclic_tiling()
.Let
be the multiplicative group of the finite field
. A radical family has the form
, where
(for
odd) or
(for
even). Equivalently,
decomposes as:
We observe that
is a subgroup of the (cyclic) group
, that can thus be generated by some element
. Furthermore, we observe that
is always a union of cosets of
(which is twice larger than
).
where
Consequently,
is a radical difference family if and only if
is a partition of the cyclic group
.
EXAMPLES:
sage: from sage.combinat.designs.difference_family import ( ....: one_radical_difference_family, ....: is_difference_family) sage: one_radical_difference_family(GF(13),4) [[0, 1, 3, 9]]
The parameters that appear in [Bu95]:
sage: df = one_radical_difference_family(GF(449), 8); df [[0, 1, 18, 25, 176, 324, 359, 444], [0, 9, 88, 162, 222, 225, 237, 404], [0, 11, 140, 198, 275, 357, 394, 421], [0, 40, 102, 249, 271, 305, 388, 441], [0, 49, 80, 93, 161, 204, 327, 433], [0, 70, 99, 197, 230, 362, 403, 435], [0, 121, 141, 193, 293, 331, 335, 382], [0, 191, 285, 295, 321, 371, 390, 392]] sage: is_difference_family(GF(449), df, 449, 8, 1) True
-
sage.combinat.designs.difference_family.
radical_difference_family
(K, k, l=1, existence=False, check=True)¶ Return a
(v,k,l)
-radical difference family.Let fix an integer
and a prime power
. Let
be a field of cardinality
. A
-difference family is radical if its base blocks are either: a coset of the
-th root of unity for
odd or a coset of
-th root of unity and
if
is even (the number
is the number of blocks of that difference family).
The terminology comes from M. Buratti article [Bu95] but the first constructions go back to R. Wilson [Wi72].
INPUT:
K
- a finite fieldk
– positive integer, the size of the blocksl
– theparameter (default to
)
existence
– ifTrue
, then return eitherTrue
if Sage knows how to build such design,Unknown
if it does not andFalse
if it knows that the design does not exist.check
– boolean (default:True
). IfTrue
then the result of the computation is checked before being returned. This should not be needed but ensures that the output is correct.
EXAMPLES:
sage: from sage.combinat.designs.difference_family import radical_difference_family sage: radical_difference_family(GF(73),9) [[1, 2, 4, 8, 16, 32, 37, 55, 64]] sage: radical_difference_family(GF(281),5) [[1, 86, 90, 153, 232], [4, 50, 63, 79, 85], [5, 36, 149, 169, 203], [7, 40, 68, 219, 228], [9, 121, 212, 248, 253], [29, 81, 222, 246, 265], [31, 137, 167, 247, 261], [32, 70, 118, 119, 223], [39, 56, 66, 138, 263], [43, 45, 116, 141, 217], [98, 101, 109, 256, 279], [106, 124, 145, 201, 267], [111, 123, 155, 181, 273], [156, 209, 224, 264, 271]] sage: for k in range(5,10): ....: print "k = {}".format(k) ....: for q in range(k*(k-1)+1, 2000, k*(k-1)): ....: if is_prime_power(q): ....: K = GF(q,'a') ....: if radical_difference_family(K, k, existence=True): ....: print q, ....: _ = radical_difference_family(K,k) ....: print k = 5 41 61 81 241 281 401 421 601 641 661 701 761 821 881 1181 1201 1301 1321 1361 1381 1481 1601 1681 1801 1901 k = 6 181 211 241 631 691 1531 1831 1861 k = 7 337 421 463 883 1723 k = 8 449 1009 k = 9 73 1153 1873
-
sage.combinat.designs.difference_family.
radical_difference_set
(K, k, l=1, existence=False, check=True)¶ Return a difference set made of a cyclotomic coset in the finite field
K
and with paramtersk
andl
.Most of these difference sets appear in chapter VI.18.48 of the Handbook of combinatorial designs.
EXAMPLES:
sage: from sage.combinat.designs.difference_family import radical_difference_set sage: D = radical_difference_set(GF(7), 3, 1); D [[1, 2, 4]] sage: sorted(x-y for x in D[0] for y in D[0] if x != y) [1, 2, 3, 4, 5, 6] sage: D = radical_difference_set(GF(16,'a'), 6, 2) sage: sorted(x-y for x in D[0] for y in D[0] if x != y) [1, 1, a, a, a + 1, a + 1, a^2, a^2, ... a^3 + a^2 + a + 1, a^3 + a^2 + a + 1] sage: for k in range(2,50): ....: for l in reversed(divisors(k*(k-1))): ....: v = k*(k-1)//l + 1 ....: if is_prime_power(v) and radical_difference_set(GF(v,'a'),k,l,existence=True): ....: _ = radical_difference_set(GF(v,'a'),k,l) ....: print "{:3} {:3} {:3}".format(v,k,l) 3 2 1 4 3 2 7 3 1 5 4 3 7 4 2 13 4 1 11 5 2 7 6 5 11 6 3 16 6 2 8 7 6 9 8 7 19 9 4 37 9 2 73 9 1 11 10 9 19 10 5 23 11 5 13 12 11 23 12 6 27 13 6 27 14 7 16 15 14 31 15 7 ... 41 40 39 79 40 20 83 41 20 43 42 41 83 42 21 47 46 45 49 48 47 197 49 12
-
sage.combinat.designs.difference_family.
singer_difference_set
(q, d)¶ Return a difference set associated to the set of hyperplanes in a projective space of dimension
over
.
A Singer difference set has parameters:
The idea of the construction is as follows. One consider the finite field
as a vector space of dimension
over
. The set of
-lines in
is a projective plane and its set of hyperplanes form a balanced incomplete block design.
Now, considering a multiplicative generator
of
, we get a transitive action of a cyclic group on our projective plane from which it is possible to build a difference set.
The construction is given in details in [Stinson2004], section 3.3.
EXAMPLES:
sage: from sage.combinat.designs.difference_family import singer_difference_set, is_difference_family sage: G,D = singer_difference_set(3,2) sage: is_difference_family(G,D,verbose=True) It is a (13,4,1)-difference family True sage: G,D = singer_difference_set(4,2) sage: is_difference_family(G,D,verbose=True) It is a (21,5,1)-difference family True sage: G,D = singer_difference_set(3,3) sage: is_difference_family(G,D,verbose=True) It is a (40,13,4)-difference family True sage: G,D = singer_difference_set(9,3) sage: is_difference_family(G,D,verbose=True) It is a (820,91,10)-difference family True
-
sage.combinat.designs.difference_family.
twin_prime_powers_difference_set
(p, check=True)¶ Return a difference set on
.
The difference set is built from the following element of the cartesian product of finite fields
:
with any
with
and
squares
with
and
non-squares
For more information see Wikipedia article Difference_set.
INPUT:
check
– boolean (default:True
). IfTrue
then the result of the computation is checked before being returned. This should not be needed but ensures that the output is correct.
EXAMPLES:
sage: from sage.combinat.designs.difference_family import twin_prime_powers_difference_set sage: G,D = twin_prime_powers_difference_set(3) sage: G The cartesian product of (Finite Field of size 3, Finite Field of size 5) sage: D [[(1, 1), (1, 4), (2, 2), (2, 3), (0, 0), (1, 0), (2, 0)]]