Simplify •math.Comb implementation

This commit is contained in:
Marshall Lochbaum 2022-11-23 09:49:17 -05:00
parent 6eb21bd68c
commit 9a1268ea3a

View File

@ -360,8 +360,8 @@ AR_F_SCALAR("|", stile, pfmod(x.f, w.f))
AR_F_SCALAR("⋆⁼",log , log(x.f)/log(w.f)) AR_F_SCALAR("⋆⁼",log , log(x.f)/log(w.f))
#undef AR_F_SCALAR #undef AR_F_SCALAR
static f64 comb_nat(f64 k, f64 n, f64 j) { static f64 comb_nat(f64 k, f64 n) {
if (j < k) k = j; assert(k>=0 && n>=2*k);
if (k > 514) return INFINITY; if (k > 514) return INFINITY;
f64 p = 1; f64 p = 1;
for (usz i=0; i<(usz)k; i++) { for (usz i=0; i<(usz)k; i++) {
@ -375,15 +375,13 @@ static f64 comb(f64 k, f64 n) { // n choose k
bool jint = j == round(j); bool jint = j == round(j);
if (k == round(k)) { if (k == round(k)) {
if (jint) { if (jint) {
if (k<j) { f64 t=k; k=j; j=t; } // Now j<k
if (n >= 0) { if (n >= 0) {
if (!(k>=0 && j>=0)) return 0; // Negative phrasing to catch NaN return j<0? 0 : comb_nat(j, n);
return comb_nat(k, n, j);
} else { } else {
if (k<0) { if (k<0) return 0;
if (j<0) return 0; f64 l = -1-n; // l+k == -1-j
f64 t=k; k=j; j=t; // Swap so k is non-negative f64 r = comb_nat(k<l? k : l, -1-j);
}
f64 r = comb_nat(k, -1-j, -1-n);
return k<(1ull<<53) && ((i64)k&1)? -r : r; return k<(1ull<<53) && ((i64)k&1)? -r : r;
} }
} }