An incidence structure is specified by a list of points, blocks, and an incidence matrix ([1], [2]).
Classes:
IncidenceStructure
This software is released under the terms of the GNU General Public License, version 2 or above (your choice). For details on licensing, see the accompanying documentation.
REFERENCES:
This is a significantly modified form of part of the module block_design.py (version 0.6) written by Peter Dobcsanyi peter@designtheory.org.
Copyright 2007-2008 by David Joyner wdjoyner@gmail.com, Peter Dobcsanyi peter@designtheory.org.
Bases: object
This the base class for block designs.
Returns the subgroup of the automorphism group of the incidence graph which respects the P B partition. This is (isomorphic to) the automorphism group of the block design, although the degrees differ.
EXAMPLES:
sage: from sage.combinat.designs.block_design import BlockDesign
sage: BD = BlockDesign(7,[[0,1,2],[0,3,4],[0,5,6],[1,3,5],[1,4,6],[2,3,6],[2,4,5]])
sage: G = BD.automorphism_group(); G
Permutation Group with generators [(4,5)(6,7), (4,6)(5,7), (2,3)(6,7), (2,4)(3,5), (1,2)(5,6)]
sage: BD = BlockDesign(4,[[0],[0,1],[1,2],[3,3]],test=False)
sage: G = BD.automorphism_group(); G
Permutation Group with generators [()]
sage: BD = BlockDesign(7,[[0,1,2],[0,3,4],[0,5,6],[1,3,5],[1,4,6],[2,3,6],[2,4,5]])
sage: G = BD.automorphism_group(); G
Permutation Group with generators [(4,5)(6,7), (4,6)(5,7), (2,3)(6,7), (2,4)(3,5), (1,2)(5,6)]
This is not a wrapper for GAP Design’s IsBlockDesign. The GAP Design function IsBlockDesign http://www.gap-system.org/Manuals/pkg/design/htm/CHAP004.htmSSEC001.1 apparently simply checks the record structure and no mathematical properties. Instead, the function below checks some necessary (but not sufficient) “easy” identities arising from the identity.
INPUT:
WARNING: This is very fast but can return false positives.
EXAMPLES:
sage: from sage.combinat.designs.block_design import BlockDesign
sage: BD = BlockDesign(7,[[0,1,2],[0,3,4],[0,5,6],[1,3,5],[1,4,6],[2,3,6],[2,4,5]])
sage: BD.parameters()
(2, 7, 3, 1)
sage: BD.block_design_checker(2, 7, 3, 1)
True
sage: BD.block_design_checker(2, 7, 3, 1,"binary")
True
sage: BD.block_design_checker(2, 7, 3, 1,"connected")
True
sage: BD.block_design_checker(2, 7, 3, 1,"simple")
True
Return a list of block’s sizes.
EXAMPLES:
sage: from sage.combinat.designs.block_design import BlockDesign
sage: BD = BlockDesign(7,[[0,1,2],[0,3,4],[0,5,6],[1,3,5],[1,4,6],[2,3,6],[2,4,5]])
sage: BD.block_sizes()
[3, 3, 3, 3, 3, 3, 3]
Return the list of blocks.
EXAMPLES:
sage: from sage.combinat.designs.block_design import BlockDesign
sage: BD = BlockDesign(7,[[0,1,2],[0,3,4],[0,5,6],[1,3,5],[1,4,6],[2,3,6],[2,4,5]])
sage: BD.blocks()
[[0, 1, 2], [0, 3, 4], [0, 5, 6], [1, 3, 5], [1, 4, 6], [2, 3, 6], [2, 4, 5]]
Wraps GAP Design’s DualBlockDesign (see [1]). The dual of a block design may not be a block design.
Also can be called with dual_design.
REQUIRES: method=”gap” option requires GAP’s Design package. method=None option does not require GAP’s Design.
EXAMPLES:
sage: from sage.combinat.designs.block_design import BlockDesign
sage: D = BlockDesign(4, [[0,2],[1,2,3],[2,3]], test=False)
sage: D
Incidence structure with 4 points and 3 blocks
sage: D.dual_design()
Incidence structure with 3 points and 4 blocks
sage: print D.dual_design(method="gap") # optional - gap_design
IncidenceStructure<points=[0, 1, 2], blocks=[[0], [0, 1, 2], [1], [1, 2]]>
sage: BD = IncidenceStructure(range(7),[[0,1,2],[0,3,4],[0,5,6],[1,3,5],[1,4,6],[2,3,6],[2,4,5]], name="FanoPlane")
sage: BD
Incidence structure with 7 points and 7 blocks
sage: print BD.dual_design(method="gap") # optional - gap_design
IncidenceStructure<points=[0, 1, 2, 3, 4, 5, 6], blocks=[[0, 1, 2], [0, 3, 4], [0, 5, 6], [1, 3, 5], [1, 4, 6], [2, 3, 6], [2, 4, 5]]>
sage: BD.dual_incidence_structure()
Incidence structure with 7 points and 7 blocks
REFERENCE:
Wraps GAP Design’s DualBlockDesign (see [1]). The dual of a block design may not be a block design.
Also can be called with dual_design.
REQUIRES: method=”gap” option requires GAP’s Design package. method=None option does not require GAP’s Design.
EXAMPLES:
sage: from sage.combinat.designs.block_design import BlockDesign
sage: D = BlockDesign(4, [[0,2],[1,2,3],[2,3]], test=False)
sage: D
Incidence structure with 4 points and 3 blocks
sage: D.dual_design()
Incidence structure with 3 points and 4 blocks
sage: print D.dual_design(method="gap") # optional - gap_design
IncidenceStructure<points=[0, 1, 2], blocks=[[0], [0, 1, 2], [1], [1, 2]]>
sage: BD = IncidenceStructure(range(7),[[0,1,2],[0,3,4],[0,5,6],[1,3,5],[1,4,6],[2,3,6],[2,4,5]], name="FanoPlane")
sage: BD
Incidence structure with 7 points and 7 blocks
sage: print BD.dual_design(method="gap") # optional - gap_design
IncidenceStructure<points=[0, 1, 2, 3, 4, 5, 6], blocks=[[0, 1, 2], [0, 3, 4], [0, 5, 6], [1, 3, 5], [1, 4, 6], [2, 3, 6], [2, 4, 5]]>
sage: BD.dual_incidence_structure()
Incidence structure with 7 points and 7 blocks
REFERENCE:
Returns the incidence graph of the design, where the incidence matrix of the design is the adjacency matrix of the graph.
EXAMPLE:
sage: BD = IncidenceStructure(range(7),[[0,1,2],[0,3,4],[0,5,6],[1,3,5],[1,4,6],[2,3,6],[2,4,5]])
sage: BD.incidence_graph()
Bipartite graph on 14 vertices
sage: A = BD.incidence_matrix()
sage: Graph(block_matrix([A*0,A,A.transpose(),A*0])) == BD.incidence_graph()
True
REFERENCE:
Return the incidence matrix A of the design. A is a (v x b) matrix defined by: A[i,j] = 1 if i is in block B_j 0 otherwise
EXAMPLES:
sage: BD = IncidenceStructure(range(7),[[0,1,2],[0,3,4],[0,5,6],[1,3,5],[1,4,6],[2,3,6],[2,4,5]])
sage: BD.block_sizes()
[3, 3, 3, 3, 3, 3, 3]
sage: BD.incidence_matrix()
[1 1 1 0 0 0 0]
[1 0 0 1 1 0 0]
[1 0 0 0 0 1 1]
[0 1 0 1 0 1 0]
[0 1 0 0 1 0 1]
[0 0 1 1 0 0 1]
[0 0 1 0 1 1 0]
Returns a pair True, pars if the incidence structure is a t-design, for some t, where pars is the list of parameters [t, v, k, lmbda]. The largest possible t is returned, provided t=10.
EXAMPLES:
sage: BD = IncidenceStructure(range(7),[[0,1,2],[0,3,4],[0,5,6],[1,3,5],[1,4,6],[2,3,6],[2,4,5]])
sage: BD.is_block_design()
(True, [2, 7, 3, 1])
sage: BD.block_design_checker(2, 7, 3, 1)
True
sage: BD = WittDesign(9) # requires optional gap package 'design'
sage: BD.is_block_design() # requires optional gap package 'design'
(True, [2, 9, 3, 1])
sage: BD = WittDesign(12) # requires optional gap package 'design'
sage: BD.is_block_design() # requires optional gap package 'design'
(True, [5, 12, 6, 1])
sage: BD = AffineGeometryDesign(3, 1, GF(2))
sage: BD.is_block_design()
(True, [2, 8, 2, 2])
Returns (t,v,k,lambda). Does not check if the input is a block design. Uses t=2 by default.
EXAMPLES:
sage: from sage.combinat.designs.block_design import BlockDesign
sage: BD = BlockDesign(7,[[0,1,2],[0,3,4],[0,5,6],[1,3,5],[1,4,6],[2,3,6],[2,4,5]], name="FanoPlane")
sage: BD.parameters()
(2, 7, 3, 1)
sage: BD.parameters(t=3)
(3, 7, 3, 0)
Returns the list of points.
EXAMPLES:
sage: from sage.combinat.designs.block_design import BlockDesign
sage: BD = BlockDesign(7,[[0,1,2],[0,3,4],[0,5,6],[1,3,5],[1,4,6],[2,3,6],[2,4,5]])
sage: BD.points()
[0, 1, 2, 3, 4, 5, 6]
Literally pushes this block design over to GAP and returns the points of that. Other than debugging, usefulness is unclear.
REQUIRES: GAP’s Design package.
EXAMPLES:
sage: from sage.combinat.designs.block_design import BlockDesign
sage: BD = BlockDesign(7,[[0,1,2],[0,3,4],[0,5,6],[1,3,5],[1,4,6],[2,3,6],[2,4,5]])
sage: BD.points_from_gap() # requires optional gap package 'design'
[1, 2, 3, 4, 5, 6, 7]
M must be a (0,1)-matrix. Creates a set of “points” from the rows and a set of “blocks” from the columns.
EXAMPLES:
sage: from sage.combinat.designs.block_design import BlockDesign
sage: BD1 = BlockDesign(7,[[0,1,2],[0,3,4],[0,5,6],[1,3,5],[1,4,6],[2,3,6],[2,4,5]])
sage: M = BD1.incidence_matrix()
sage: BD2 = IncidenceStructureFromMatrix(M)
sage: BD1 == BD2
True
L is a list of n-vectors or lists all of length n with a common parent. This returns the vector whose i-th coordinate is the product of the i-th coordinates of the vectors.
EXAMPLES:
sage: from sage.combinat.designs.incidence_structures import coordinatewise_product
sage: L = [[1,2,3],[-1,-1,-1],[5,7,11]]
sage: coordinatewise_product(L)
[-5, -14, -33]