← Back to branch summary

~ubuntu-branches/ubuntu/quantal/genometools/quantal-backports

« back to all changes in this revision

Viewing changes to src/core/combinatorics.h

  • Committer: Package Import Robot
  • Author(s): Sascha Steinbiss
  • Date: 2012-07-09 14:10:23 UTC
  • Revision ID: package-import@ubuntu.com-20120709141023-juuu4spm6chqsf9o
Tags: upstream-1.4.1
Import upstream version 1.4.1

Show diffs side-by-side

added added

removed removed

Lines of Context:
src/core/combinatorics.h
 
1
/*
 
2
  Copyright (c) 2007 Thomas Jahns <Thomas.Jahns@gmx.net>
 
3
 
 
4
  Permission to use, copy, modify, and distribute this software for any
 
5
  purpose with or without fee is hereby granted, provided that the above
 
6
  copyright notice and this permission notice appear in all copies.
 
7
 
 
8
  THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL WARRANTIES
 
9
  WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED WARRANTIES OF
 
10
  MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE AUTHOR BE LIABLE FOR
 
11
  ANY SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL DAMAGES OR ANY DAMAGES
 
12
  WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR PROFITS, WHETHER IN AN
 
13
  ACTION OF CONTRACT, NEGLIGENCE OR OTHER TORTIOUS ACTION, ARISING OUT OF
 
14
  OR IN CONNECTION WITH THE USE OR PERFORMANCE OF THIS SOFTWARE.
 
15
*/
 
16
 
 
17
#ifndef COMBINATORICS_H
 
18
#define COMBINATORICS_H
 
19
 
 
20
#include "core/assert_api.h"
 
21
#include "core/minmax.h"
 
22
 
 
23
/**
 
24
 * \file combinatorics.h
 
25
 * \brief Simple routines for distribution computations relevant to
 
26
 * combinatorial problems.
 
27
 * \author Thomas Jahns <Thomas.Jahns@gmx.net>
 
28
 */
 
29
 
 
30
/**
 
31
 * Computes \f$n! = 1 \cdot 2 \cdot \dots \cdot (n - 1) \cdot n\f$
 
32
 * @param n number of which to compute factorial
 
33
 * @return \f$n!\f$
 
34
 */
 
35
static inline unsigned long
 
36
factorial(int n)
 
37
{
 
38
  unsigned long k = 1;
 
39
  while (n > 1)
 
40
    k *= n--;
 
41
  return k;
 
42
}
 
43
 
 
44
/**
 
45
 * \brief Compute binomial coefficient \f[{n\choose k} = \frac{n!}{k!\cdot (n
 
46
 * - k)!}\f]
 
47
 * @param n
 
48
 * @param k
 
49
 * @return \f$n\choose{k}\f$
 
50
 */
 
51
static inline unsigned long
 
52
binomialCoeff(unsigned long n, unsigned long k)
 
53
{
 
54
  unsigned long accum;
 
55
  gt_assert(k <= n);
 
56
  if (k == 0 || k == n)
 
57
    return 1;
 
58
  else if (k < n/2)
 
59
    return binomialCoeff(n, n - k);
 
60
  else
 
61
  {
 
62
    unsigned long i = k;
 
63
    accum = ++i;
 
64
    while (i < n)
 
65
      accum *= ++i;
 
66
    return accum /= factorial(n - k);
 
67
  }
 
68
}
 
69
 
 
70
/**
 
71
 * \brief Compute multinomial coefficient
 
72
 * \f[{n\choose k_1, k_2,\dots,k_m} = \frac{n!}{k_1!\cdot
 
73
 * k_2!\cdot\dots\cdot k_m!}\f] where \f$m\f$ equals \link numBins\endlink
 
74
 * @param n
 
75
 * @param numBins
 
76
 * @param binSizes points to array containing \link numBins\endlink values which
 
77
 * represent the \f$k_i\f$
 
78
 * @return \f$n\choose{k_1, k_2,\dots,k_m}\f$
 
79
 */
 
80
static inline unsigned long
 
81
multinomialCoeff(unsigned n, size_t numBins, const unsigned binSizes[])
 
82
{
 
83
  unsigned long accum = 1, nfac;
 
84
  size_t i, maxBin = 0, maxBinSize = 0;
 
85
#ifndef NDEBUG
 
86
  unsigned long binSum = 0;
 
87
#endif
 
88
  gt_assert(n > 0 && numBins > 0 && binSizes);
 
89
  for (i = 0; i < numBins; ++i)
 
90
  {
 
91
#ifndef NDEBUG
 
92
    binSum += binSizes[i];
 
93
#endif
 
94
    if (binSizes[i] > maxBinSize)
 
95
    {
 
96
      maxBinSize = binSizes[i];
 
97
      maxBin = i;
 
98
    }
 
99
  }
 
100
  gt_assert(binSum <= n);
 
101
  for (nfac = maxBinSize + 1; nfac <= n; ++nfac)
 
102
    accum *= nfac;
 
103
  for (i = 0; i < numBins; ++i)
 
104
    if (i != maxBin)
 
105
      accum /= factorial(binSizes[i]);
 
106
  return accum;
 
107
}
 
108
 
 
109
static inline unsigned long long
 
110
iPow(unsigned long long x, unsigned i)
 
111
{
 
112
   unsigned long long result = 1;
 
113
   while (i)
 
114
   {
 
115
     if (i & 1)
 
116
       result *= x;
 
117
     x *= x;
 
118
     i >>= 1;
 
119
   }
 
120
   return result;
 
121
}
 
122
 
 
123
#endif