IntroToCombinatorics
system:sage


<p>Trying to get oriented what is around in <strong>combinatorics</strong> ...</p>

{{{id=11|
#sage.combinat
///
}}}

<h2>Linear Diophantine equations</h2>

<p><strong>Problem: </strong>A postal clerk has only 14 and 21-cent stamps. What combinations can be used to make up 3.50 Dollars worth of stamps?</p>

{{{id=21|
#sage.combinat.integer_vector_weighted??
///
}}}

{{{id=23|
W = WeightedIntegerVectors(350,[14,21]); W
///
Integer vectors of 350 weighted by [14, 21]
}}}

{{{id=25|
W.cardinality()
///
9
}}}

{{{id=24|
W.list()
///
[[1, 16], [4, 14], [7, 12], [10, 10], [13, 8], [16, 6], [19, 4], [22, 2], [25, 0]]
}}}

<h2>Fibonacci numbers</h2>
<p>Let us experiment a little bit with the <strong>Fibonacci</strong> <strong>numbers</strong>. Recall that they are defined by the inital condition f_1=f_2=1 and the recursion f_{n+2} = f_{n+1} + f_n</p>

{{{id=26|
F = [fibonacci(n) for n in range(1,15)]; F
///
[1, 1, 2, 3, 5, 8, 13, 21, 34, 55, 89, 144, 233, 377]
}}}

{{{id=28|
[m for m in F if is_even(m)]
///
[2, 8, 34, 144]
}}}

{{{id=29|
sum(fibonacci(i) for i in range(1,5))
///
7
}}}

<p>Defining our first function:</p>

{{{id=30|
def sum_fib(n):
    return sum(fibonacci(i) for i in range(1,n+1))
///
}}}

{{{id=31|
[sum_fib(n) for n in range(1,13)]
///
[1, 2, 4, 7, 12, 20, 33, 54, 88, 143, 232, 376]
}}}

<p>Now we might conjecture that the sum of the first n Fibonacci numbers is equal to f_{n+2}-1. Let us test this conjecture!</p>

{{{id=32|
for n in range(1,100):
    assert( sum_fib(n) == fibonacci(n+2)-1)
print 'Done!'
///
Done!
}}}

<p>What would have happened if we had tried a wrong conjecture?</p>

{{{id=36|
for n in range(1,100):
    assert( sum_fib(n) == fibonacci(n+2))
print 'Done!'
///
Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
  File "_sage_input_115.py", line 10, in <module>
    exec compile(u'open("___code___.py","w").write("# -*- coding: utf-8 -*-\\n" + _support_.preparse_worksheet_cell(base64.b64decode("Zm9yIG4gaW4gcmFuZ2UoMSwxMDApOgogICAgYXNzZXJ0KCBzdW1fZmliKG4pID09IGZpYm9uYWNjaShuKzIpKQpwcmludCAnRG9uZSEn"),globals())+"\\n"); execfile(os.path.abspath("___code___.py"))' + '\n', '', 'single')
  File "", line 1, in <module>
    
  File "/private/var/folders/nR/nRzFEby2FV8wW6eFUOQcgU+++TI/-Tmp-/tmpRuobEL/___code___.py", line 4, in <module>
    assert( sum_fib(n) == fibonacci(n+_sage_const_2 ))
AssertionError
}}}

<h2>Partially ordered sets</h2>
<p>The set of all divisors of a number n form a <strong>poset</strong> by defining a&lt;=b if a divides b.</p>

{{{id=39|
n = 30
n.divisors()
///
[1, 2, 3, 5, 6, 10, 15, 30]
}}}

{{{id=40|
def order_relation(a,b):
    return b%a == 0
///
}}}

{{{id=60|
order_relation(15,30)
///
True
}}}

{{{id=61|
order_relation(6,15)
///
False
}}}

{{{id=41|
n = 30
P = Poset((n.divisors(), order_relation), cover_relations = False)
P.plot()
///
<html><font color='black'><img src='cell://sage0.png'></font></html>
}}}

{{{id=5|
B = Posets.BooleanLattice(3)
B.plot()
///
<html><font color='black'><img src='cell://sage0.png'></font></html>
}}}

{{{id=12|
n = 60
P = Poset((n.divisors(), order_relation), cover_relations = False)
P.plot()
///
<html><font color='black'><img src='cell://sage0.png'></font></html>
}}}

<p>We can ask about some properties of the poset:</p>

{{{id=13|
P.is_graded()
///
True
}}}

{{{id=18|
P.is_meet_semilattice()
///
True
}}}

{{{id=42|
P.is_join_semilattice()
///
True
}}}

{{{id=14|
P.minimal_elements()
///
[1]
}}}

{{{id=19|
P.maximal_elements()
///
[60]
}}}

{{{id=20|
P.mobius_function(P(1),P(30))
///
-1
}}}

{{{id=59|
P.mobius_function(P(1),P(4))
///
0
}}}

{{{id=43|
sum([P.mobius_function(P(1),v) for v in P])
///
0
}}}

<h2>Mersenne numbers</h2>
<p>Let us test some results on <strong>Mersenne numbers:</strong></p>

{{{id=2|
def mersenne(p):
    return 2^p-1
///
}}}

<p>Is p always prime when the Mersenne number is a prime?</p>

{{{id=46|
[is_prime(p) for p in range(1000) if is_prime(mersenne(p))]
///
[True, True, True, True, True, True, True, True, True, True, True, True, True, True]
}}}

<p>Is the converse also true?</p>

{{{id=48|
all( is_prime(mersenne(p)) for p in range(1000) if is_prime(p) )
///
False
}}}

<p>Finding the smallest counterexample:</p>

{{{id=50|
exists( (p for p in range(1000) if is_prime(p)), lambda p: not is_prime(mersenne(p)) )
///
(True, 11)
}}}

<h2>Game of poker</h2>
<p>&nbsp;We have 52=4x13 cards with suits Club, Spade, Heart and Diamond and values 2,3,...,Ace</p>

{{{id=64|
Suits = Set( ['Club', 'Spade', 'Heart', 'Diamond'] )
Values = Set( [2, 3, 4, 5, 6, 7, 8, 9, 10, 'Jack', 'Queen', 'King', 'Ace'] )
///
}}}

<p>The cards are the Cartesian product:</p>
<p>Cards = Suits x Values = {(s,v) | s in Suits, v in Values}</p>

{{{id=65|
Cards = CartesianProduct(Suits, Values)
///
}}}

{{{id=67|
Suits.cardinality()
///
4
}}}

{{{id=68|
Values.cardinality()
///
13
}}}

{{{id=69|
Cards.cardinality()
///
52
}}}

<p>We can deal a random element:</p>

{{{id=71|
Cards.random_element()
///
['Heart', 5]
}}}

<p>Now we need a brief technical digression:</p>

{{{id=73|
type(Cards.random_element())
///
<type 'list'>
}}}

<p>Python lists are not immutable, hence they cannot be used to put elements in a set. This can be fixed by transforming the lists into tuples:</p>

{{{id=74|
Cards = CartesianProduct(Suits, Values).map(tuple)
Cards.random_element()
///
('Diamond', 'King')
}}}

{{{id=76|
Set( [ Cards.random_element(), Cards.random_element() ] )
///
{('Diamond', 10), ('Spade', 7)}
}}}

<p>A poker hand consists of 5 cards. All cards are distinguishable and the order of the cards does not matter for a hand. A hand is hence a 5-subset&nbsp; of the set of all cards.</p>

{{{id=77|
Hands = Subsets(Cards,5)
Hands.random_element()
///
{('Spade', 'Queen'), ('Club', 6), ('Club', 2), ('Diamond', 6), ('Spade', 6)}
}}}

<p>Total number of hands:</p>

{{{id=80|
binomial(52,5)
///
2598960
}}}

{{{id=81|
Hands.cardinality()
///
2598960
}}}

<p>Let us now calculate the probability for a <strong>flush,&nbsp;</strong> that is, all cards in a given hand have the same suit:</p>

{{{id=83|
Colors = CartesianProduct(Suits, Subsets(Values,5))
Colors.cardinality()
///
5148
}}}

{{{id=84|
100.0 * Colors.cardinality() / Hands.cardinality()
///
0.198079231692677
}}}

<p>We can test our reasoning by doing a little simulation: we produce 10000 random hands and count how many of them have the same color.</p>

{{{id=87|
def is_same_suit(hand):
        return len(set(suit for (suit,value) in hand)) == 1
///
}}}

{{{id=88|
n = 10000
nsuits = 0
for i in range(n):
    hand = Hands.random_element()
    if is_same_suit(hand):
        nsuits += 1
print n, nsuits, 100.0 * nsuits / n
///
10000 19 0.190000000000000
}}}

<h2>Crystal bases!!</h2>

{{{id=52|
C = CrystalOfTableaux(['A',2], shape=[2,1])
C.cardinality()
///
8
}}}

<p>For the next command, you need to install dot2tex:</p>
<p>sage -f http://sage.math.washington.edu/home/nthiery/dot2tex-2.8.7-2.spkg</p>

{{{id=54|
view(C, pdflatex = True)
///
<html><img src="cell://sage0.png"></html>
}}}

{{{id=93|

///
}}}

{{{id=94|

///
}}}