From ddbb892748230a601d2365168be9e7138867f214 Mon Sep 17 00:00:00 2001 From: filifa Date: Tue, 16 Sep 2025 22:23:39 -0400 Subject: [PATCH] massively improve performance with dynamic programming --- cmd/stirling.go | 9 ++--- internal/lib/stirling.go | 72 ++++++++++++++++++++++------------------ 2 files changed, 45 insertions(+), 36 deletions(-) diff --git a/cmd/stirling.go b/cmd/stirling.go index 91dfec5..27eccb8 100644 --- a/cmd/stirling.go +++ b/cmd/stirling.go @@ -19,6 +19,7 @@ package cmd import ( "fmt" "math/big" + "strconv" "github.com/spf13/cobra" "scm.dairydemon.net/filifa/mathtools/internal/lib" @@ -31,13 +32,13 @@ var stirlingBottom string var stirlingUnsigned bool func stirling(cmd *cobra.Command, args []string) { - n, ok := new(big.Int).SetString(stirlingTop, 10) - if !ok { + n, err := strconv.Atoi(stirlingTop) + if err != nil { cobra.CheckErr("invalid input " + stirlingTop) } - k, ok := new(big.Int).SetString(stirlingBottom, 10) - if !ok { + k, err := strconv.Atoi(stirlingBottom) + if err != nil { cobra.CheckErr("invalid input " + stirlingBottom) } diff --git a/internal/lib/stirling.go b/internal/lib/stirling.go index b3f99ff..55ced30 100644 --- a/internal/lib/stirling.go +++ b/internal/lib/stirling.go @@ -20,49 +20,57 @@ import ( "math/big" ) -// TODO: implement both of these with dynamic programming - -func Stirling1(n, k *big.Int) *big.Int { - if n.Cmp(big.NewInt(0)) == 0 && k.Cmp(big.NewInt(0)) == 0 { - return big.NewInt(1) +// given a slice where vals[i] = Stirling1(i+k-1, k-1) for a given k, nextStirling1 updates the slice so vals[i] = Stirling1(i+k, k) using the property that Stirling1(n, k) = -(n-1)*Stirling1(n-1, k) + Stirling1(n-1, k-1) +func nextStirling1(k int, vals []*big.Int) { + for i := 1; i < len(vals); i++ { + n := int64(k + i - 1) + v := big.NewInt(-n) + v.Mul(v, vals[i-1]) + vals[i].Add(vals[i], v) } +} - if n.Cmp(big.NewInt(0)) == 0 || k.Cmp(big.NewInt(0)) == 0 { +func Stirling1(n, k int) *big.Int { + if k > n { return big.NewInt(0) } - newN := new(big.Int).Set(n) - newN.Sub(newN, big.NewInt(1)) + vals := make([]*big.Int, n-k+1) + for i := range vals { + vals[i] = big.NewInt(0) + } + vals[0] = big.NewInt(1) - result := Stirling1(newN, k) - result.Mul(result, newN) - result.Neg(result) - - newK := new(big.Int).Set(k) - newK.Sub(newK, big.NewInt(1)) - - result.Add(result, Stirling1(newN, newK)) - return result -} - -func Stirling2(n, k *big.Int) *big.Int { - if n.Cmp(k) == 0 { - return big.NewInt(1) + for i := 1; i <= k; i++ { + nextStirling1(i, vals) } - if n.Cmp(big.NewInt(0)) == 0 || k.Cmp(big.NewInt(0)) == 0 { + return vals[n-k] +} + +// given a slice where vals[i] = Stirling2(i+k-1, k-1) for a given k, nextStirling2 updates the slice so vals[i] = Stirling2(i+k, k) using the property that Stirling2(n, k) = k*Stirling2(n-1, k) + Stirling2(n-1, k-1) +func nextStirling2(k int64, vals []*big.Int) { + for i := 1; i < len(vals); i++ { + v := big.NewInt(k) + v.Mul(v, vals[i-1]) + vals[i].Add(vals[i], v) + } +} + +func Stirling2(n, k int) *big.Int { + if k > n { return big.NewInt(0) } - newN := new(big.Int).Set(n) - newN.Sub(newN, big.NewInt(1)) + vals := make([]*big.Int, n-k+1) + for i := range vals { + vals[i] = big.NewInt(0) + } + vals[0] = big.NewInt(1) - result := Stirling2(newN, k) - result.Mul(result, k) + for i := 1; i <= k; i++ { + nextStirling2(int64(i), vals) + } - newK := new(big.Int).Set(k) - newK.Sub(newK, big.NewInt(1)) - - result.Add(result, Stirling2(newN, newK)) - return result + return vals[n-k] }