Previous Up Next

Chapitre 9  Exercices de combinatoire

9.1  Fonction partage ou nombre de partitions de n∈ ℕ

9.1.1  L’énoncé

Définition
Un partition de n∈ ℕ est une suite de nombres entiers : λ1≥ λ2≥ ...λm>0 tels que : ∑j=1mλj =n.
On note p(n) la fonction partage de n : c’est le nombre de partitions distinctes de n∈ ℕ et on convient que p(0)=1.
Exemples
p(5)=7 car les 7 partitions distinctes de 5 sont :
1+1+1+1+1=5
2+1+1+1=5
2+2+1=5
3+1+1=5
3+2=5
4+1=5
5=5

  1. Chercher à la main : p(1),p(2),p(3),p(4),p(5),p(6),p(7)
  2. Écrire un programme qui renvoie les valeurs de p(0),p(1),p(2),p(3)..p(n)
  3. Montrer ce que Euler a remarqué à savoir que le développement en séries tronqué à l’ordre n de : ∏j=1n1/1−xj vaut ∑j=1np(j)xj.
  4. Écrire un programme qui renvoie les coefficients du développement en séries tronqué à l’ordre n de : ∏j=1n1/1−xj

9.1.2  La solution

  1. On trouve :
    p(0)=1, p(1)=1, p(2)=2, p(3)=3, p(4)=5, p(5)=7, p(6)=11, p(7)=15
  2. Soit A une matrice triangulaire inférieure tel que :
    si kj, A[j,k] represente le nombre de partitions de j tel que λ1=k≥ λ2≥ ...λm>0. Sur l’exemple
    p(5)=7 car les 7 partitions distinctes de 5 sont :
    1+1+1+1+1=5 donc A[5,1]=1
    2+1+1+1=5
    2+2+1=5 donc A[5,2]=2
    3+1+1=5
    3+2=5 donc A[5,3]=2
    4+1=5 donc A[5,4]=1
    5=5 donc A[5,5]=1
    On a alors p(j)=∑k=0jA[j,k].
    On remarque donc :
    A[0,0]=1 et A[j,0]=0 (car λm>0 sauf pour j=0 puisque p(0)=1)
    A[j,j]=1 puisque j=j
    A[j,1]=A[j−1,1]=1 puisque j=j−1+1=j−2+1+1=1+1...+1
    A[j,2]=A[j−2,1]+A[j−2,2] puisque j=2+1+1+...=2+2+....
    A[j,3]=A[j−3,1]+A[j−3,2]+A[j−3,3] puisque j=3+1+1+...1=3+2+....=3+3+....
    etc...
    A[j,j−1]=A[j−1,j−1]=1 puisque j=(j−1)+1
    ou encore :
    les coefficients d’indice 0 de toutes les lignes sont nuls et
    pour k=1..n le coefficient d’indice k de la n-ième ligne est obtenu en faisant la somme des coefficients d’indice 0..k de la nk-ième ligne c’est à dire :
    A[n,k]=
    k
    j=1
    A[nk,j]  et  
    A[n,n]=1   et  A[n,0]=0  si  n>0
    Voici l’illustration pour n=5 avec
    k=2, on a 5=0+1+4 :
    k=3, on a 8=0+1+3+4 :
    k=4, on a 9=0+1+3+3+2 :
    ....


































     
    1
    0
    0
    0
    0
    0
    0
    0
    0
    0
    0
     
    0
    1
    0
    0
    0
    0
    0
    0
    0
    0
    |
     
    0
    1
    1
    0
    0
    0
    0
    0
    0
    |
    0
     
    0
    1
    1
    1
    0
    0
    0
    0
    |
    00
     
    0
    1
    2
    1
    1
    0
    0
    |
    000
     
    0
    1
    2
    2
    1
    1
    |
    0000
     
    0
    1
    3
    3
    2
    |
    11000
     
    0
    1
    3
    4
    |
    321100
     
    0
    1
    4
    |
    5532110
    0
    1
    |
    47653211
    01589753211


































    On va faire un programme qui va calculer les coefficients de la matrice A.

    On tape :

    partitions(n):={
      local A,j,k,m,S;
      A:=idn(n+1);
      pour j de 2 jusque n faire
        pour k de 1 jusque j-1 faire
          S:=0;
          pour m de 1 jusque k faire
            S:=S+A[j-k,m];
          fpour;
          A[j,k]:=S;
        fpour;
      fpour;
      retourne A;
    }:;
    

    On tape : A:=partitions(10)
    On obtient :















    1000000000
    0100000000
    0110000000
    0111000000
    0121100000
    0122110000
    0133211000
    0134321100
    0145532110
    0147653211
    01589753211














    On tape : sum(A[k])$(k=0..10)
    On obtient : 1,1,2,3,5,7,11,15,22,30,42

    On peut faire un autre programme qui renverra la matrice B qui sera tel que : si kj, A[j,k] represente le nombre de partitions de j tel que k>=λ1≥ λ2≥ ...λm>0. Autrement dit les lignes de B sont les sommes partielles des lignes de A : B[j,k]=∑m=0kA[j,m] donc

    B[j,k]=B[jk,k]+B[j,k−1]

    On a donc p(j)=B[j,j].
    On tape :

    partition(n):={
      local B,j,k,m,S;
      B:=idn(n+1);
      pour k de 1 jusque n faire
        B[0,k]:=1;
      fpour;
      pour j de 1 jusque n faire
        pour k de 1 jusque j faire
          B[j,k]:=B[j-k,k]+B[j,k-1];
        fpour;
        pour k de j+1 jusque n faire
          B[j,k]:=B[j,j];
        fpour;
      fpour;
      retourne B;
    }
    :;
    

    On tape : B:=partition(10)
    On obtient :















    1111111111
    0111111111
    0122222222
    0123333333
    0134555555
    0135677777
    01479101111111111 
    014811131415151515 
    0151015182021222222 
    0151218232628293030 
    0161423303538404142














    On remarquera que l’on retrouve la matrice A dans les diagonales montantes de la matrice B.
    On tape : tran(B)[10]
    ou on tape : col(B,10)
    On obtient : [1,1,2,3,5,7,11,15,22,30,42] ou on tape : B[j,10]$(j=0..10)
    On obtient : 1,1,2,3,5,7,11,15,22,30,42 On tape : B:=partition(1000):;
    On obtient (Evaluation time: 972.05) enfin le résultat...
    On tape : B[200,1000]
    On obtient 3972999029388
    On tape : B[1000,1000]
    On obtient 24061467864032622473692149727991
    Remarque On peut aussi pour avoir une relation de récurrence considérer les les partitions P(n,p) de n en une somme d’exactement p éléments λ1≥ ...≥ λp>0. On dit alors :
    P(n,p) =(nombre de partage de n tels que λp=1)+(nombre de partage de n tels que λp>1).
    On a :
    (nombre de partage de n tels que λp=1)=(nombre de partage de n−1 en p−1 λj)=P(n−1,p−1)
    (nombre de partage de n tels que λp>1)=(nombre de partage de np en p λj)=P(np,p)

    P(n,p):=P(n−1,p−1)+P(np,p)

    avec P(n,0)=0, P(n,n)=1 et P(n,p)=0 si p>n
    On tape et on renvoie la somme des lignes de P :

    partage(n):={
      local P,j,k,m,S;
      P:=idn(n+1);
      pour k de 1 jusque n faire
        P[k,0]:=0;
      fpour;
      pour j de 1 jusque n faire
        pour k de 1 jusque j faire
          P[j,k]:=P[j-1,k-1]+P[j-k,k];
        fpour;
      fpour;
      S:=[];
      pour k de 1 jusque n faire
        S[k]:=sum(row(P,k));
      fpour;
      retourne S;
    }:;
    

    On tape : P:=partage(1000):;
    On obtient (Evaluation time: 455.81) enfin le résultat...mais c’est deux fois plus rapide qu’avec B:=partition(1000):;
    On tape : P[200]
    On obtient 3972999029388
    On tape : P[1000]
    On obtient 24061467864032622473692149727991
    relation entre P et A :
    On a :
    P[n,p]=P[n−1,p−1]+P[np,p]
    P[n−1,p−1]=P[n−2,p−2]+P[np,p−1]....
    donc

    P[n,p]= P[np,0]+
    p
    j=1
    P[np,j]

    avec P[n,0]=0, P[n,n]=1 et P[n,p]=0 si p>n
    Si on compare avec les relations obtenus pour A:

    A[n,k]=
    k
    j=1
    A[nk,j]  et  
    A[n,n]=1   et  A[n,0]=0  si  n>0

    On a les mêmes relations de récurrence, donc :

    P[n,p]=A[n,p]

    Représentation de Young
    On peut voir facilement cette égalité, grâce à la repésentation de Young qui repésente par exemple le partage :
    20=7+3+3+2+1+1+1+1+1
    par le tableau formé par 20 carrés disposés selon 9 lignes :
    7 carrés sur la 1-ière ligne,
    3 carrés sur la 2-ième ligne,
    3 carrés sur la 3-ième ligne,
    2 carrés sur la 4-ième ligne,
    1 carré sur de la 5-ième jusqu’à la 9-ième ligne :












    □ 
        
        
         
          
          
          
          
         











    Le nombre p de lignes correspond a un partage en une somme de p éléments et le plus grand élément m de ce partage est le nombre d’éléments de la première ligne.
    Si maintenant on échange les lignes et les colonnes (on prend le transposé de ce tableau) on obtient une transformation qui au partage :
    20=7+3+3+2+1+1+1+1+1 en une somme de 9 termes ayant comme plus grand élément 7 fait correspondre le partage :
    20=9+4+3+1+1+1+1 en une somme de 7 termes ayant comme plus grand élément 9 :









         
          
            
            
            
           








  3. Notations
    1 sera noté 11, 1+1 sera noté 21, 1+1+1 sera noté 31 etc...
    2 sera noté 22, 2+2 sera noté 42, 2+2+2 sera noté 62 etc...
    Le développement en série de ;
    1/1−x1 sera noté 1+x11+x21+...+xn1+...
    1/1−x2 sera noté 1+x22+x42+...+x2n2+...
    1/1−x3 sera noté 1+x33+x63+...+x3n3+...
    .....
    1/1−xp sera noté 1+xpp+x2pp+...+xpnp+...
    Lorsque l’on fait le produit de ces développements en série le terme en xj a pour coefficient p(j) car :
    xj=xj1=x(j−2)1x22=x(j−4)1x42 ...
    xj=x(j−2)j−2x21=x(j−2)j−2x22=x(j−1)j−1x11=xjj

    On tape : P:=truncate(product([(1-x^j)$(j=1..20)]),20)
    On obtient : -x^15-x^12+x^7+x^5-x^2-x+1
    On tape : A:=truncate(series(1/P,x=0,20),20)
    On obtient :
    627*x^20+490*x^19+385*x^18+297*x^17+231*x^16+176*x^15+ 135*x^14+101*x^13+77*x^12+56*x^11+42*x^10+30*x^9+22*x^8+ 15*x^7+11*x^6+7*x^5+5*x^4+3*x^3+2*x^2+x+1
    On tape :symb2poly(A)
    ou on tape : coeff(A,x)
    On obtient la liste p(20),p(19),..,p(1),p(0):
    [627,490,385,297,231,176,135,101,77,56,42,30,22,15,11, 7,5,3,2,1,1]
    On tape :lcoeff(A)
    On obtient p(20): 627

  4. On tape mais cela n’est pas très efficace lorsque n est grand :
    partitioner(n):={
    local A,P;
    P:=truncate(product([(1-x^j)$(j=1..n)]),n);
    A:=truncate(series(1/P,x=0,n),n);
    return revlist(symb2poly(A));
    }:;
    
    On tape : partitioner(20)
    On obtient la liste p(20),p(19),..,p(1),p(0) :
    [627,490,385,297,231,176,135,101,77,56,42,30,22,15,11, 7,5,3,2,1,1]

9.1.3  Une méthode plus rapide

On utilise le théorème du nombre pentagonal d’Euler. Ce théorème donne une relation de récurrence entre p(n) et p(j) pour j<n avec la convention que p(j)=0 lorsque j<0, p(0)=1. Cette relation est :

p(n)=p(n−1)+p(n−2)−p(n−5)−p(n−7)+p(n−12)+p(n−15)−p(n−22)−p(n−26)....
p(n)=
 
m>0
 (−1)m+1p(nm(3m−1)/2)+(−1)m+1p(nm(3m+1)/2)

Les nombres m(3m−1)/2 et m(3m+1)/2 sont les nombres pentagonaux généralisés. On tape :

partition_euler(n):={
local sg,s1,s2,m,k,j,p;
p:=[1];
pour j de 1 jusque n faire
m:=1;
sg:=1;
s1:=0;
k:=j-m*(3*m-1)/2;
tantque k>=0 faire 
s1:=s1+sg*p[k]
sg:=-sg;
m:=m+1;
k:=j-m*(3*m-1)/2;
ftantque;
m:=1;
sg:=1;
s2:=0;
k:=j-m*(3*m+1)/2;
tantque k>=0 faire 
s2:=s2+sg*p[k]
sg:=-sg;
m:=m+1;
k:=j-m*(3*m+1)/2;
ftantque;
p[j]:=s1+s2;
fpour;
retourne p;
}:;

On tape :
partition_euler(20)
On obtient :
[1,1,2,3,5,7,11,15,22,30,42,56,77,101,135,176,231,297,385,490,627]
On tape :
P:=partition_euler(1000)
On obtient la liste après 2.77s
On tape : P[200]
On obtient 3972999029388
On tape : P[1000]
On obtient 24061467864032622473692149727991

9.1.4  Estimation asymptotique

Ramanjan et Hardy ont trouvé une approximation de p(n) qui est :
p(n)≃ eπ2n/3/4n3 quand n tend vers +∞.
On tape :
P(n):=exp(pi*sqrt(2n/3))/(4n*sqrt(3))
floor(evalf(P(1000),31))
On obtient :
24401996316802476288263414942904 au lieu de
24061467864032622473692149727991

9.2  Un exercice de combinatoire et son programme

9.2.1  L’énoncé

Un jury est composé de p personnes que l’on choisit parmi n personnes : h hommes et f femmes (h+f=n).
Application numérique : p=10, n=17, h=9, f=8.

  1. Écrire un programme qui renvoie la liste contenant les différents jurys possibles.
  2. On impose comme contrainte que Monsieur Y et Madame X ne doivent pas se trouver ensemble.
    Écrire un programme qui renvoie la liste contenant tous les jurys respectant cette contrainte.
  3. On veut que le jury respecte la parité "homme-femme".
    Pour cela, si p=2*k ou p=2*k+1, on suppose fk et h=nfpk.
    Écrire un programme qui renvoie la liste contenant les jurys ayant pk hommes et k femmes.
  4. Avec l’hypothèse précédente, écrire un programme qui renvoie la liste contenant les jurys respectant la parité "homme-femme" et respectant la contrainte Monsieur Y et Madame X ne doivent pas se trouver ensemble.
  5. Rappels Avec Xcas, la commande comb(n,p) renvoie le nombre de combinaisons de p objets pris parmi n.
    Avec Xcas, la commande convert(k,base,2) renvoie la liste K de longueur l vérifiant k=∑j=0l−1K[j]*2j. revlist(K) renvoie alors l’écriture en base 2 de k.

9.2.2  Les programmes

  1. On peut avoir comb(n,p) jurys différents.
    On tape : comb(17,10)
    On obtient : 19448
    Supposons que n cases C numérotées de 0 à n−1 représentent les membres du jury.
    Dans ces cases on met 1 si la personne fait partie du jury et 0 sinon. Un jury est donc une liste de 0 et de 1 et cela fait penser à l’écriture en base 2 d’un entier.
    Avec Xcas, la commande convert(k,base,2) renvoie la liste K de longueur l vérifiant k=∑j=0l−1K[j]*2j On va donc considérer un jury comme un entier dont l’écriture en base 2 a au plus n chiffres et comporte exactement p 1. Pour faire le programme on va repésenter chaque jury par un nombre dont l’écriture en base 2 comporte au plus n chiffres dont p 1 et pas plus de np zéros.
    Le plus petit nombre correspondant à un jury est :
    m=sum(2^j,j,0,p-1)=2^p-1
    et le plus grand nombre correspondant à un jury est :
    M=sum(2^j,j,n-p,n-1)=2^n-2^(n-p)
    car les n cases sont numerotés de 0 à n−1.
    On peut initialiser la liste L soit avec L:=makelist(0,1,1); soit avec L:=makelist(0,1,comb(n,p)); puisque l’on connait sa longueur. On écrit le programme jurys(n,p) :
    jurys(n,p):={
    local j,k,L,m,M;
    //L:=makelist(0,1,comb(n,p));
    L:=makelist(0,1,1);
    k:=0;
    m:=2^p-1;
    M:=2^n-2^(n-p);
    for (j:=m;j<=M;j++){
    if (sum(convert(j,base,2))==p){
    L[k]=<j;
    k:=k+1;
    };
    }
    return L;
    }
    :;
    
    On tape : J:=jurys(17,10):;
    On obtient, au bout d’environ 11s, une liste J de taille 19448 et contenant des entiers qui lorsqu’on les décompose en base 2 donne la composition du jury.
    On tape : jurytire0:=J[rand(19448)]
    On obtient par exemple : 127255
    On tape : convert(127255,base,2)
    On obtient par exemple : [1,1,1,0,1,0,0,0,1,0,0,0,1,1,1,1,1]
    Donc, les personnes de numéros : 0,1,2,4,8,12,13,14,15,16 ont été tirées au sort pour faire partie du jury.
    On peut bien sûr taper directement :
    jurytire0:=convert(J[rand(19448)],base,2)
    On obtient par exemple : [1,0,0,0,1,1,1,0,1,1,1,0,0,1,0,1,1]
    Donc, les personnes de numéros : 0,4,5,6,8,9,10,13,15,16 ont été tirées au sort pour faire partie du jury.
  2. On peut supposer que la case 0 représente Monsieur Y et que la case n−1 représente Madame X. Donc le jury repésenté par le nombre j est valable si :
    irem(j,2)==0 ou si j<2^(n-1).
    On peut réecrire un programme jurysXY(n,p) qui utilise jurys ou un programme juryXY(n,p) qui n’utilise pas jurys ou utiliser la liste précédente J et rejeter les mauvais tirages.
  3. On peut supposer que les cases de 0 à h−1 représentent les hommes et que les cases de h à n−1 représentent les femmes. On écrit un nouveau programme jury55(h,f,p) qui utilise le programme jurys(n,p) écrit précédemment ou on se sert de la liste J trouvée précedemment.
  4. On peut supposer que la case 0 représente Monsieur Y et que la case n−1 représente Madame X. Donc le jury repésenté par le nombre j est valable si :
    irem(j,2)==0 ou si j<2^(n-1).
    On peut réecrire un programme jury55XY(h,f,p) qui utilise jury55, ou un programme jurys55XY(n,p) qui n’utilise pas jury55 mais jurys, ou utiliser la liste précédente J55 et rejeter les mauvais tirages, ou encore utiliser la liste J du début et rejeter les mauvais tirages.

9.3  Visualistion des combinaisons modulo 2

9.3.1  L’énoncé

On veut visualiser les points k,j (jk) tel que comb(k,j)=0 mod2

9.3.2  Le programme

On tape :

dessincomb(n):={
  local j,k,L;
  L:=NULL;
  pour j de 1 jusque n faire
    pour k de j jusque n faire
      si irem(comb(k,j),2) alors 
        L:=L,point(k,j,affichage=point_point);
      fsi;
    fpour;
  fpour;
  retourne L;
}:;

On tape :
dessincomb(255)
On obtient :

9.4  Un exercice

On tire 100 fois de suite au hasard de façon équiprobable un entier p parmi 1,2..100 et on note np le nombre de fois où il a été tiré.
Écrire un programme qui représente les points p,np On tape :

tirage():={
local p,np,j,L,P;
L:=0$(j=1..100)
pour j de 0 jusque 99 faire
p:=rand(100);
L[p]:=L[p]+1;
fpour;
P:=NULL;
pour p de 0 jusque 99 faire
np:=L[p];
si np!=0 alors P:=P,point(p+1+i*np); fsi; 
fpour;
print(max(L));
retourne P;
}:;

On tape tirage()
On obtient :

9.5  Valeur de e et le hasard

9.5.1  L’énoncé

On tire des nombres au hasard dans [0;1[ et on les ajoute. On s’arrête quand la somme obtenue est strictement plus grande que 1. On peut montrer que le nombres de tirages nécessaires est en moyenne égal à e. Écrire un programme vérifiant expérimentalement ce résultat.

9.5.2  La solution avec Xcas

Dans Xcas, on peut soit utiliser rand(0,1), soit utiliser la fonction f:=rand(0..1) pour obtenir un nombre au hasard entre 0 et 1.
On tape :rand(0.1)
On obtient : 0.482860322576
On tape : f:=rand(0..1)
puis : f() On obtient : 0.8627617261373
puis : f() On obtient : 0.3095522336662 etc....
On écrit le programme :

approxe(n):={
local j,S,k,N;
N:=0;
pour k de 1 jusque n faire
S:=0;
j:=0;
tantque S<=1 faire
S:=S+rand(0,1);
j:=j+1;
ftantque;
N:=N+j;
fpour;
retourne evalf(N/n);
}:;

On tape : approxe(100000)
On obtient : 2.71750000000
alors que evalf(e)=2.71828182846

9.6  Distance moyenne entre de 2 points

9.6.1  L’énoncé

On veut évaluer la distance moyenne entre 2 points choisis au hasard : dans un segment S de longueur 1, dans le carré unité C ou dans le cube unité K.

  1. Écrire un programme qui calcule cette distance moyenne lorsque l’on choisit au hasard n fois de suite 2 points dans S,dans C ou dans K, puis généraliser.
  2. Faire une représentation graphique permettant de visualiser la répartition des distances dans le cube unité K.

9.6.2  La solution avec Xcas

  1. distS(n):={
    local xa,xb,j,d,D;
    D:=0;
    pour j de 1 jusque n faire
    xa:=rand(0,1);
    xb:=rand(0,1);
    d:=abs(xa-xb);
    D:=D+d;
    fpour;
    retourne evalf(D/n);
    }:;
    
    distC(n):={
    local xa,ya,xb,yb,j,d,D;
    D:=0;
    pour j de 1 jusque n faire
    xa:=rand(0,1);ya:=rand(0,1);
    xb:=rand(0,1);yb:=rand(0,1);
    d:=distance([xa,ya],[xb,yb]);
    D:=D+d;
    fpour;
    retourne evalf(D/n);
    }:;
    
    distK(n):={
    local xa,ya,za,xb,yb,zb,j,d,D;
    D:=0;
    pour j de 1 jusque n faire
    xa:=rand(0,1);ya:=rand(0,1);za:=rand(0,1);
    xb:=rand(0,1);yb:=rand(0,1);zb:=rand(0,1);
    d:=distance([xa,ya,za],[xb,yb,zb]);
    D:=D+d;
    fpour;
    retourne evalf(D/n);
    }:;
    
    distp(p,n):={
    local a,b,j,d,D;
    D:=0;
    pour j de 1 jusque n faire
    a:=op(randMat(1,p,0..1));
    b:=op(randMat(1,p,0..1));
    d:=sqrt(sum((a-b)^2));
    D:=D+d;
    fpour;
    retourne evalf(D/n);
    }:;
    
    
    On tape : distS(100000)
    On obtient : 0.333424530881
    On tape : distC(100000)
    On obtient : 0.522140638365
    On tape : distK(100000)
    On obtient : 0.6627691590747
    On tape : distp(4,100000)
    On obtient : 0.777235961163
    On tape : distp(5,100000)
    On obtient : 0.877837288394
  2. graphdist(n):={
    local xa,ya,za,xb,yb,zb,j,d,L;
    L:=NULL;
    pour j de 1 jusque n faire
    xa:=rand(0,1);ya:=rand(0,1);za:=rand(0,1);
    xb:=rand(0,1);yb:=rand(0,1);zb:=rand(0,1);
    d:=distance([xa,ya,za],[xb,yb,zb]);
    L:=L,point(j+i*d);
    fpour;
    retourne L;
    }:;
    

On tape : graphdist(10000)
On obtient :


Previous Up Next