Skip to content

bpo-37295: Use constant-time comb() for larger n depending on k #30305

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
serhiy-storchaka merged 7 commits into from
Jan 9, 2022
Merged

bpo-37295: Use constant-time comb() for larger n depending on k #30305

Changes from 1 commit
Commits
Show all changes
7 commits
Select commit Hold shift + click to select a range
f9804ac
bpo-37295: Use constant-time comb() for larger n depending on k
serhiy-storchaka Dec 30, 2021
0029e7e
Optimize also perm().
serhiy-storchaka Dec 30, 2021
4cbeee6
Silence complier warnings.
serhiy-storchaka Dec 31, 2021
f9b0f28
Merge branch 'main' into fast_comb
serhiy-storchaka Dec 31, 2021
d57c36b
Optimize limit tables and tiny refactoring.
serhiy-storchaka Jan 1, 2022
0d1ed7c
Refactor tests.
serhiy-storchaka Jan 1, 2022
78efc24
Improve comments.
serhiy-storchaka Jan 1, 2022
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
Prev Previous commit
Improve comments.
  • Loading branch information
@serhiy-storchaka
serhiy-storchaka committed Jan 1, 2022
commit 78efc2457d5fc89d451e0df3e94de2f0053ea57e
23 changes: 17 additions & 6 deletions Modules/mathmodule.c
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -3358,6 +3358,9 @@ perm_comb_small(unsigned long long n, unsigned long long k, int iscomb)
/* For small enough n and k the result fits in the 64-bit range and can
* be calculated without allocating intermediate PyLong objects. */
if (iscomb) {
/* Maps k to the maximal n so that 2*k-1 <= n <= 127 and C(n, k)
* fits into a uint64_t. Exclude k = 1, because the second fast
* path is faster for this case.*/
static const unsigned char fast_comb_limits1[] = {
0, 0, 127, 127, 127, 127, 127, 127, // 0-7
127, 127, 127, 127, 127, 127, 127, 127, // 8-15
Expand All @@ -3367,8 +3370,7 @@ perm_comb_small(unsigned long long n, unsigned long long k, int iscomb)
};
if (k < Py_ARRAY_LENGTH(fast_comb_limits1) && n <= fast_comb_limits1[k]) {
/*
For 0 <= k <= n <= fast_comb_limits1[k], comb(n, k) always fits
into a uint64_t. We compute it as
comb(n, k) fits into a uint64_t. We compute it as

comb_odd_part << shift

Expand All @@ -3386,13 +3388,15 @@ perm_comb_small(unsigned long long n, unsigned long long k, int iscomb)
return PyLong_FromUnsignedLongLong(comb_odd_part << shift);
}

/* Only contains items larger than in fast_comb_limits1. */
/* long long is at least 64 bit. */
/* Maps k to the maximal n so that 2*k-1 <= n <= 127 and C(n, k)*k
* fits into a long long (which is at least 64 bit). Only contains
* items larger than in fast_comb_limits1. */
static const unsigned long long fast_comb_limits2[] = {
Comment on lines +3391 to 3394
Copy link
Copy Markdown
Member

@arhadthedev arhadthedev Jan 1, 2022

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

which is at least 64 bit

C99 provides types like uint_least64_t that express such comments explicitly. So I think it will be better to #include <stdint.h> and use it:

Suggested change
/* Maps k to the maximal n so that 2*k-1 <= n <= 127 and C(n, k)*k
* fits into a long long (which is at least 64 bit). Only contains
* items larger than in fast_comb_limits1. */
static const unsigned long long fast_comb_limits2[] = {
/* Maps k to the maximal n so that 2*k-1 <= n <= 127 and C(n, k)*k
* fits into a long long. Only contains
* items larger than in fast_comb_limits1. */
static const uint_least64_t fast_comb_limits2[] = {

Copy link
Copy Markdown
Member Author

@serhiy-storchaka serhiy-storchaka Jan 1, 2022

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

It should contain at least LLONG_MAX. There is no guarantee that uint_least64_t is at least so large as positive long long.

0, ULLONG_MAX, 4294967296ULL, 3329022, 102570, 13467, 3612, 1449, // 0-7
746, 453, 308, 227, 178, 147, // 8-13
};
if (k < Py_ARRAY_LENGTH(fast_comb_limits2) && n <= fast_comb_limits2[k]) {
/* C(n, k) = C(n, k-1) * (n-k+1) / k */
unsigned long long result = n;
for (unsigned long long i = 1; i < k;) {
result *= --n;
Expand All @@ -3402,20 +3406,24 @@ perm_comb_small(unsigned long long n, unsigned long long k, int iscomb)
}
}
else {
/* Maps k to the maximal n so that k <= n and P(n, k)
* fits into a long long (which is at least 64 bit). */
static const unsigned long long fast_perm_limits[] = {
0, ULLONG_MAX, 4294967296ULL, 2642246, 65537, 7133, 1627, 568, // 0-7
259, 142, 88, 61, 45, 36, 30, 26, // 8-15
24, 22, 21, 20, 20, // 16-20
};
if (k < Py_ARRAY_LENGTH(fast_perm_limits) && n <= fast_perm_limits[k]) {
if (n <= 127) {
/* P(n, k) fits into a uint64_t. */
uint64_t perm_odd_part = reduced_factorial_odd_part[n]
* inverted_factorial_odd_part[n - k];
int shift = factorial_trailing_zeros[n]
- factorial_trailing_zeros[n - k];
return PyLong_FromUnsignedLongLong(perm_odd_part << shift);
}

/* P(n, k) = P(n, k-1) * (n-k+1) */
unsigned long long result = n;
for (unsigned long long i = 1; i < k;) {
result *= --n;
Expand All @@ -3425,8 +3433,11 @@ perm_comb_small(unsigned long long n, unsigned long long k, int iscomb)
}
}

/* For larger n use recursive formula. */
/* C(n, k) = C(n, j) * C(n-j, k-j) // C(k, j) */
/* For larger n use recursive formulas:
*
* P(n, k) = P(n, j) * P(n-j, k-j)
* C(n, k) = C(n, j) * C(n-j, k-j) // C(k, j)
*/
unsigned long long j = k / 2;
PyObject *a, *b;
a = perm_comb_small(n, j, iscomb);
Expand Down