/* alea2comb.c
 * combalea.c
 *
 * Recherche aléatoire d'une ou de plusieurs combinaisons de p éléments de
 * l'ensemble [1 n] = {1, 2, 3, 4, 5, ..., n}
 * 
 * Copyright (C) 2006 Jean-Paul Davalan <jpdvl@wanadoo.fr>
 *
 * This program is free software; you can redistribute it and/or
 * modify it under the terms of the GNU General Public License
 * as published by the Free Software Foundation; either version 2
 * of the License, or (at your option) any later version.
 *
 * This program is distributed in the hope that it will be useful,
 * but WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 * GNU General Public License for more details.
 *
 * You should have received a copy of the GNU General Public License
 * along with this program; if not, write to the Free Software
 * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301, USA.
 * 
 * compilation : 
 * -------------
 * fonction non récursive
 * gcc -O2 -o alea2comb alea2comb.c 
 * fonction recursive
 * gcc -DRECURSE -O2 -o alea2comb alea2comb.c
 * usage : alea2comb n p [k]
 *
 * quelques tests :
 * alea2comb 10 4 100|gawk '{for(i=1;i<=NF;i++)t[$i]++;};END{for (x in t) print x" "t[x];}'
 */
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <time.h>

int     alea(int n, int p);
void    combinaisonAleatoire(int n, int p);

int
main(int argc, char *argv[])
{
  int     n, p, k = 1, i, j, cnp, *L, *t;
  srandom(time(NULL));
  if (argc < 3) {
    printf("Combinaisons aléatoires\n"
       "Sous-ensembles au hasard de p éléments de [1 n]\n"
       "usage : %s n p [k]\n", argv[0]);
    exit(1);
  }
  n = atoi(argv[1]);		// lecture des paramètres
  p = atoi(argv[2]);
  if (argc > 3) {
    k = atoi(argv[3]);
  }
  // pas de combinaison
  if (n < 0 || p < 0 || p > n || k <= 0)
    return 0;
  // des combinaisons
  for (i = 0; i < k; i++) {
    combinaisonAleatoire(n, p);
    printf("\n");
  }
  return 0;
}

/*     Relation souvent utilisée pour construire le triangle de Pascal
 *     binom(n, p) = binom(n-1, p-1) + binom(n-1, p)
 *     On peut vérifier les deux égalités ci-dessous
 *     		binom(n-1, p-1)/binom(n, p) = p/n
 *     		binom(n-1, p)/binom(n, p) = (n-p)/n
 *
 *     Plusieurs exemples de constructions de combinaisons au hasard
 *     les probabilités sont proportionnelles aux Nbres de combinaisons
 *     et elles sont aussi proportionnelles à p et à n-p
 *     		de 3 éléments parmi 5
 *              1                    *                     *
 *              1  1                 .  1                  .  1
 *              1  2  1              .  *  .               .  .  2
 *              1  3  3  1           .  *  .  .            .  .  .  3
 *              1  4  6  4  1        .  .  4  .  .         .  .  .  *  .
 *              1  5 10 10  5 1      .  .  .  5  .  .      .  .  .  *  .  .
 *              
 *		de 4 éléments parmi 5
 *              *                    si on monte en diagonale, écrire le n° de la ligne
 *              .  1                 si on monte verticalement, ne rien écrire (ici une *)
 *     		.  .  2
 *     		.  .  *  .           1 2 4 5
 *     		.  .  .  4  .
 *     		.  .  .  .  5  .
 *     		
 *     algorithme pour construire une combinaison au hasard
 *     si p<0 ou n<p ou n<0 alors 
 *     		combinaison(n, p) ne retourne rien (PAS de COMBINAISON, binom(n,p)=0)
 *     sinon si p=0
 *     		combinaison(n, 0) = ENSEMBLE VIDE = ligne vide (UNE COMBIN. binom(n,p)=1)
 *     sinon
 *     		combinaison(n, p) = 
 *         		soit avec la probabilité p/n      =  combinaison(n-1, p-1) + " n"
 *         		soit avec la probabilité (n-p)/n  =  combinaison(n-1, p)
 *
 */

// retourne 1 avec la probabilité p/n 
// et retourne 0 le reste du temps (avec la proba contraire (n-p)/n )
int
alea(n, p)
{
  return ((double) n * random() / (RAND_MAX + 1.0) < (double) p) ? 1 : 0;
}

#ifdef RECURSE
void
combinaisonAleatoire(n, p)
{
  if (n < p || p < 0) {
    return;
  } else if (alea(n, p)) {
    combinaisonAleatoire(n - 1, p - 1);
    if (p > 1) {
      printf(" ");
    }
    printf("%d", n);
  } else {
    combinaisonAleatoire(n - 1, p);
  }
}
#else
void
combinaisonAleatoire(n, p)
{
  for (; n >= p && p >= 0; n--) {
    if (alea(n, p)) {
      printf("%d ", n);
      p--;
    }
  }
}
#endif