From 36fa19a38cf9a357f04d4ed76f25b1cbf44deedb Mon Sep 17 00:00:00 2001 From: Lucas Schwebler Date: Tue, 10 Sep 2024 21:50:42 +0200 Subject: new linear recurrence kthTerm code --- content/math/linearRecurence.cpp | 33 ----------------- content/math/linearRecurrence.cpp | 30 +++++++++++++++ content/math/linearRecurrenceOld.cpp | 33 +++++++++++++++++ content/math/math.tex | 5 ++- tcr.pdf | Bin 696418 -> 696332 bytes test/math/linearRecurence.cpp | 54 --------------------------- test/math/linearRecurenceOld.cpp | 54 +++++++++++++++++++++++++++ test/math/linearRecurrence.cpp | 57 +++++++++++++++++++++++++++++ test/math/linearRecurrenceSlowMul.cpp | 67 ++++++++++++++++++++++++++++++++++ 9 files changed, 244 insertions(+), 89 deletions(-) delete mode 100644 content/math/linearRecurence.cpp create mode 100644 content/math/linearRecurrence.cpp create mode 100644 content/math/linearRecurrenceOld.cpp delete mode 100644 test/math/linearRecurence.cpp create mode 100644 test/math/linearRecurenceOld.cpp create mode 100644 test/math/linearRecurrence.cpp create mode 100644 test/math/linearRecurrenceSlowMul.cpp diff --git a/content/math/linearRecurence.cpp b/content/math/linearRecurence.cpp deleted file mode 100644 index 2501e64..0000000 --- a/content/math/linearRecurence.cpp +++ /dev/null @@ -1,33 +0,0 @@ -constexpr ll mod = 1'000'000'007; -vector modMul(const vector& a, const vector& b, - const vector& c) { - ll n = sz(c); - vector res(n * 2 + 1); - for (int i = 0; i <= n; i++) { //a*b - for (int j = 0; j <= n; j++) { - res[i + j] += a[i] * b[j]; - res[i + j] %= mod; - }} - for (int i = 2 * n; i > n; i--) { //res%c - for (int j = 0; j < n; j++) { - res[i - 1 - j] += res[i] * c[j]; - res[i - 1 - j] %= mod; - }} - res.resize(n + 1); - return res; -} - -ll kthTerm(const vector& f, const vector& c, ll k) { - assert(sz(f) == sz(c)); - vector tmp(sz(c) + 1), a(sz(c) + 1); - tmp[0] = a[1] = 1; //tmp = (x^k) % c - - for (k++; k > 0; k /= 2) { - if (k & 1) tmp = modMul(tmp, a, c); - a = modMul(a, a, c); - } - - ll res = 0; - for (int i = 0; i < sz(c); i++) res += (tmp[i+1] * f[i]) % mod; - return res % mod; -} diff --git a/content/math/linearRecurrence.cpp b/content/math/linearRecurrence.cpp new file mode 100644 index 0000000..c15c25c --- /dev/null +++ b/content/math/linearRecurrence.cpp @@ -0,0 +1,30 @@ +// constexpr ll mod = 998244353; +// vector mul(const vector &a, const vector &b){ +// vector c(sz(a) + sz(b) - 1); +// for(int i = 0; i < sz(a); i++){ +// for(int j = 0; j < sz(b); j++){ +// c[i+j] += a[i]*b[j] % mod; +// } +// } +// for(ll &x : c) x %= mod; +// return c; +// } + +ll kthTerm(const vector& f, const vector& c, ll k){ + int n = sz(c); + vector q(n+1, 1); + for(int i = 1; i <= n; i++) q[i] = (mod-c[i-1])%mod; + vector p = mul(f, q); + p.resize(n); + p.push_back(0); + do{ + vector q2 = q; + for(int i = 1; i <= n; i += 2) q2[i] = (mod - q2[i]) % mod; + vector x = mul(p, q2), y = mul(q, q2); + for(int i = 0; i <= n; i++){ + p[i] = i == n ? 0 : x[2*i + (k&1)]; + q[i] = y[2*i]; + } + }while(k /= 2); + return p[0]; +} \ No newline at end of file diff --git a/content/math/linearRecurrenceOld.cpp b/content/math/linearRecurrenceOld.cpp new file mode 100644 index 0000000..2501e64 --- /dev/null +++ b/content/math/linearRecurrenceOld.cpp @@ -0,0 +1,33 @@ +constexpr ll mod = 1'000'000'007; +vector modMul(const vector& a, const vector& b, + const vector& c) { + ll n = sz(c); + vector res(n * 2 + 1); + for (int i = 0; i <= n; i++) { //a*b + for (int j = 0; j <= n; j++) { + res[i + j] += a[i] * b[j]; + res[i + j] %= mod; + }} + for (int i = 2 * n; i > n; i--) { //res%c + for (int j = 0; j < n; j++) { + res[i - 1 - j] += res[i] * c[j]; + res[i - 1 - j] %= mod; + }} + res.resize(n + 1); + return res; +} + +ll kthTerm(const vector& f, const vector& c, ll k) { + assert(sz(f) == sz(c)); + vector tmp(sz(c) + 1), a(sz(c) + 1); + tmp[0] = a[1] = 1; //tmp = (x^k) % c + + for (k++; k > 0; k /= 2) { + if (k & 1) tmp = modMul(tmp, a, c); + a = modMul(a, a, c); + } + + ll res = 0; + for (int i = 0; i < sz(c); i++) res += (tmp[i+1] * f[i]) % mod; + return res % mod; +} diff --git a/content/math/math.tex b/content/math/math.tex index dd88a5b..fb66110 100644 --- a/content/math/math.tex +++ b/content/math/math.tex @@ -136,9 +136,10 @@ sich alle Lösungen von $x^2-ny^2=c$ berechnen durch: Sei $f(n)=c_{0}f(n-1)+c_{1}f(n-2)+\dots + c_{n-1}f(0)$ eine lineare Rekurrenz. \begin{methods} - \method{kthTerm}{Berechnet $k$-ten Term einer Rekurrenz $n$-ter Ordnung}{\log(k)\cdot n^2} + \method{kthTerm}{Berechnet $k$-ten Term einer Rekurrenz $n$-ter Ordnung}{\log(k)\cdot \text{mul}(n)} \end{methods} - \sourcecode{math/linearRecurence.cpp} + Die Polynom-Multiplikation kann auch mit NTT gemacht werden! + \sourcecode{math/linearRecurrence.cpp} Alternativ kann der \mbox{$k$-te} Term in \runtime{n^3\log(k)} berechnet werden: $$\renewcommand\arraystretch{1.5} \setlength\arraycolsep{3pt} diff --git a/tcr.pdf b/tcr.pdf index 9ed6eae..1bf4c26 100644 Binary files a/tcr.pdf and b/tcr.pdf differ diff --git a/test/math/linearRecurence.cpp b/test/math/linearRecurence.cpp deleted file mode 100644 index a5290e5..0000000 --- a/test/math/linearRecurence.cpp +++ /dev/null @@ -1,54 +0,0 @@ -#include "../util.h" -#include - -struct RandomRecurence { - vector f, c, cache; - RandomRecurence(int n) : f(Random::integers(n, 0, mod)), c(Random::integers(n, 0, mod)), cache(f) {} - - ll operator()(ll k){ - while (sz(cache) <= k) { - ll cur = 0; - for (ll i = 0; i < sz(c); i++) { - cur += (c[i] * cache[sz(cache) - i - 1]) % mod; - } - cur %= mod; - cache.push_back(cur); - } - return cache[k]; - } -}; - -void stress_test() { - int queries = 0; - for (int i = 0; i < 10'000; i++) { - int n = Random::integer(1, 10); - RandomRecurence f(n); - for (int j = 0; j < 100; j++) { - ll k = Random::integer(0, 1000); - - ll got = kthTerm(f.f, f.c, k); - ll expected = f(k); - - if (got != expected) cerr << "got: " << got << ", expected: " << expected << FAIL; - queries++; - } - } - cerr << "tested random queries: " << queries << endl; -} - -constexpr int N = 1'000; -void performance_test() { - timer t; - RandomRecurence f(N); - t.start(); - hash_t hash = kthTerm(f.f, f.c, 1e18); - t.stop(); - if (t.time > 500) cerr << "too slow: " << t.time << FAIL; - cerr << "tested performance: " << t.time << "ms (hash: " << hash << ")" << endl; -} - -int main() { - stress_test(); - performance_test(); -} - diff --git a/test/math/linearRecurenceOld.cpp b/test/math/linearRecurenceOld.cpp new file mode 100644 index 0000000..dab2256 --- /dev/null +++ b/test/math/linearRecurenceOld.cpp @@ -0,0 +1,54 @@ +#include "../util.h" +#include + +struct RandomRecurence { + vector f, c, cache; + RandomRecurence(int n) : f(Random::integers(n, 0, mod)), c(Random::integers(n, 0, mod)), cache(f) {} + + ll operator()(ll k){ + while (sz(cache) <= k) { + ll cur = 0; + for (ll i = 0; i < sz(c); i++) { + cur += (c[i] * cache[sz(cache) - i - 1]) % mod; + } + cur %= mod; + cache.push_back(cur); + } + return cache[k]; + } +}; + +void stress_test() { + int queries = 0; + for (int i = 0; i < 10'000; i++) { + int n = Random::integer(1, 10); + RandomRecurence f(n); + for (int j = 0; j < 100; j++) { + ll k = Random::integer(0, 1000); + + ll got = kthTerm(f.f, f.c, k); + ll expected = f(k); + + if (got != expected) cerr << "got: " << got << ", expected: " << expected << FAIL; + queries++; + } + } + cerr << "tested random queries: " << queries << endl; +} + +constexpr int N = 1'000; +void performance_test() { + timer t; + RandomRecurence f(N); + t.start(); + hash_t hash = kthTerm(f.f, f.c, 1e18); + t.stop(); + if (t.time > 500) cerr << "too slow: " << t.time << FAIL; + cerr << "tested performance: " << t.time << "ms (hash: " << hash << ")" << endl; +} + +int main() { + stress_test(); + performance_test(); +} + diff --git a/test/math/linearRecurrence.cpp b/test/math/linearRecurrence.cpp new file mode 100644 index 0000000..50e98a0 --- /dev/null +++ b/test/math/linearRecurrence.cpp @@ -0,0 +1,57 @@ +#include "../util.h" +#include +#include +#include +#include + +struct RandomRecurence { + vector f, c, cache; + RandomRecurence(int n) : f(Random::integers(n, 0, mod)), c(Random::integers(n, 0, mod)), cache(f) {} + + ll operator()(ll k){ + while (sz(cache) <= k) { + ll cur = 0; + for (ll i = 0; i < sz(c); i++) { + cur += (c[i] * cache[sz(cache) - i - 1]) % mod; + } + cur %= mod; + cache.push_back(cur); + } + return cache[k]; + } +}; + +void stress_test() { + int queries = 0; + for (int i = 0; i < 1'000; i++) { + int n = Random::integer(1, 10); + RandomRecurence f(n); + for (int j = 0; j < 100; j++) { + ll k = Random::integer(0, 1000); + + ll got = kthTerm(f.f, f.c, k); + ll expected = f(k); + + if (got != expected) cerr << "got: " << got << ", expected: " << expected << FAIL; + queries++; + } + } + cerr << "tested random queries: " << queries << endl; +} + +constexpr int N = 100'000; +void performance_test() { + timer t; + RandomRecurence f(N); + t.start(); + hash_t hash = kthTerm(f.f, f.c, 1e18); + t.stop(); + if (t.time > 8000) cerr << "too slow: " << t.time << FAIL; + cerr << "tested performance: " << t.time << "ms (hash: " << hash << ")" << endl; +} + +int main() { + stress_test(); + performance_test(); +} + diff --git a/test/math/linearRecurrenceSlowMul.cpp b/test/math/linearRecurrenceSlowMul.cpp new file mode 100644 index 0000000..205e584 --- /dev/null +++ b/test/math/linearRecurrenceSlowMul.cpp @@ -0,0 +1,67 @@ +#include "../util.h" + +constexpr ll mod = 998244353; +vector mul(const vector &a, const vector &b){ + vector c(sz(a) + sz(b) - 1); + for(int i = 0; i < sz(a); i++){ + for(int j = 0; j < sz(b); j++){ + c[i+j] += a[i]*b[j] % mod; + } + } + for(ll &x : c) x %= mod; + return c; +} + +#include + +struct RandomRecurence { + vector f, c, cache; + RandomRecurence(int n) : f(Random::integers(n, 0, mod)), c(Random::integers(n, 0, mod)), cache(f) {} + + ll operator()(ll k){ + while (sz(cache) <= k) { + ll cur = 0; + for (ll i = 0; i < sz(c); i++) { + cur += (c[i] * cache[sz(cache) - i - 1]) % mod; + } + cur %= mod; + cache.push_back(cur); + } + return cache[k]; + } +}; + +void stress_test() { + int queries = 0; + for (int i = 0; i < 10'000; i++) { + int n = Random::integer(1, 10); + RandomRecurence f(n); + for (int j = 0; j < 100; j++) { + ll k = Random::integer(0, 1000); + + ll got = kthTerm(f.f, f.c, k); + ll expected = f(k); + + if (got != expected) cerr << "got: " << got << ", expected: " << expected << FAIL; + queries++; + } + } + cerr << "tested random queries: " << queries << endl; +} + +constexpr int N = 1'000; +void performance_test() { + timer t; + RandomRecurence f(N); + t.start(); + hash_t hash = kthTerm(f.f, f.c, 1e18); + t.stop(); + if (t.time > 500) cerr << "too slow: " << t.time << FAIL; + cerr << "tested performance: " << t.time << "ms (hash: " << hash << ")" << endl; +} + +int main() { + stress_test(); + performance_test(); +} + -- cgit v1.2.3