Check-in [a6b1fd4985]
Not logged in
Overview
SHA1 Hash:a6b1fd4985ed49811e32f53c5e9adb82fbf2e507
Date: 2011-09-28 15:31:48
User: kinaba
Comment:518
Timelines: family | ancestors | descendants | both | trunk
Downloads: Tarball | ZIP archive
Other Links: files | file ages | manifest
Tags And Properties
Changes

Added SRM/518-U/1C.cpp version [1b6b6f1bf178804b]

> 1 #include <iostream> > 2 #include <sstream> > 3 #include <iomanip> > 4 #include <vector> > 5 #include <string> > 6 #include <map> > 7 #include <set> > 8 #include <algorithm> > 9 #include <numeric> > 10 #include <iterator> > 11 #include <functional> > 12 #include <complex> > 13 #include <queue> > 14 #include <valarray> > 15 #include <stack> > 16 #include <cmath> > 17 #include <cassert> > 18 #include <cstring> > 19 using namespace std; > 20 typedef long long LL; > 21 typedef complex<double> CMP; > 22 > 23 static const int MODVAL = 1000000007; > 24 struct mint > 25 { > 26 int val; > 27 mint():val(0){} > 28 mint(int x):val(x%MODVAL) {} > 29 mint(size_t x):val(x%MODVAL) {} > 30 mint(LL x):val(x%MODVAL) {} > 31 }; > 32 mint& operator+=(mint& x, mint y) { return x = x.val+y.val; } > 33 mint& operator-=(mint& x, mint y) { return x = x.val-y.val+MODVAL; } > 34 mint& operator*=(mint& x, mint y) { return x = LL(x.val)*y.val; } > 35 mint POW(mint x, LL e) { mint v=1; for(;e;x*=x,e>>=1) if(e&1) v*=x; return v; } > 36 mint& operator/=(mint& x, mint y) { return x *= POW(y, MODVAL-2); } > 37 mint operator+(mint x, mint y) { return x+=y; } > 38 mint operator-(mint x, mint y) { return x-=y; } > 39 mint operator*(mint x, mint y) { return x*=y; } > 40 mint operator/(mint x, mint y) { return x/=y; } > 41 mint inv2 = mint(1) / 2; > 42 > 43 class Nim { public: > 44 int count(int K, int L) > 45 { > 46 valarray<mint> v(65536); > 47 v[slice(2,L-1,1)] = 1; > 48 for(unsigned p=2; p<=L; ++p) > 49 if( v[p].val ) > 50 for(unsigned q=p*p; q<=L; q+=p) > 51 v[q] = 0; > 52 > 53 pre(v, 0, 65536); > 54 valarray<mint> a(1, 65536); > 55 for(int i=1; i<=K; i<<=1, v*=v) > 56 if( K & i ) > 57 a *= v; > 58 post(a, 0, 65536); > 59 return a[0].val; > 60 } > 61 > 62 void pre(valarray<mint>& v, int s, int e) > 63 { > 64 if( s+1 == e ) > 65 return; > 66 int m = (s+e)/2; > 67 pre(v, s, m); > 68 pre(v, m, e); > 69 for(int i=s,j=m; i<m; ++i,++j) { > 70 mint vi=v[i], vj=v[j]; > 71 v[i] = vi - vj; > 72 v[j] = vi + vj; > 73 } > 74 } > 75 > 76 void post(valarray<mint>& v, int s, int e) > 77 { > 78 if( s+1 == e ) > 79 return; > 80 int m = (s+e)/2; > 81 for(int i=s,j=m; i<m; ++i,++j) { > 82 mint dif=v[i], sum=v[j]; > 83 v[i] = (sum + dif) * inv2; > 84 v[j] = (sum - dif) * inv2; > 85 } > 86 post(v, s, m); > 87 post(v, m, e); > 88 } > 89 }; > 90 > 91 // BEGIN CUT HERE > 92 #include <ctime> > 93 double start_time; string timer() > 94 { ostringstream os; os << " (" << int((clock()-start_time)/CLOCKS_PER_SEC*1000) > 95 template<typename T> ostream& operator<<(ostream& os, const vector<T>& v) > 96 { os << "{ "; > 97 for(typename vector<T>::const_iterator it=v.begin(); it!=v.end(); ++it) > 98 os << '\"' << *it << '\"' << (it+1==v.end() ? "" : ", "); os << " }"; return > 99 void verify_case(const int& Expected, const int& Received) { > 100 bool ok = (Expected == Received); > 101 if(ok) cerr << "PASSED" << timer() << endl; else { cerr << "FAILED" << timer() > 102 cerr << "\to: \"" << Expected << '\"' << endl << "\tx: \"" << Received << '\"' > 103 #define CASE(N) {cerr << "Test Case #" << N << "..." << flush; start_time=clock( > 104 #define END verify_case(_, Nim().count(K, L));} > 105 int main(){ > 106 > 107 CASE(0) > 108 int K = 3; > 109 int L = 7; > 110 int _ = 6; > 111 END > 112 CASE(1) > 113 int K = 4; > 114 int L = 13; > 115 int _ = 120; > 116 END > 117 CASE(2) > 118 int K = 10; > 119 int L = 100; > 120 int _ = 294844622; > 121 END > 122 CASE(3) > 123 int K = 123456789; > 124 int L = 12345; > 125 int _ = 235511047; > 126 END > 127 CASE(4) > 128 int K = 1000000000; > 129 int L = 50000; > 130 int _ = 428193537; > 131 END > 132 /* > 133 CASE(5) > 134 int K = ; > 135 int L = ; > 136 int _ = ; > 137 END > 138 */ > 139 } > 140 // END CUT HERE

Added lib/doc/fft.txt version [893844b7e5ed0726]

> 1 > 2 http://apps.topcoder.com/wiki/display/tc/SRM+518 > 3 > 4 > 5 以下のような関数を考える。 > 6 > 7 int[] xor_count(int[] a, int[] b, int N) > 8 { > 9 int[] c = {0} * N; > 10 for(int i=0; i<N; ++i) > 11 for(int k=0; k<N; ++k) > 12 c[i ^ k] += a[i] * b[k]; > 13 return c; > 14 } > 15 > 16 高速化。 > 17 > 18 int[] xor_count_ultra(int[] a, int[] b, int N) > 19 { > 20 a = xor_pre(a); > 21 b = xor_pre(b); > 22 c = a * b; // element-wise product > 23 return xor_post(c); > 24 } > 25 > 26 このようになるマジカルな変換 xor_pre, xor_post が存在する。 > 27 > 28 N=1 の場合、pre : {x} --> {x} > 29 N=2 の場合、pre : {x,y} --> {x-y, x+y} > 30 xor_count( {x,y}, {u,v} ) = {xu+yv, yu+xv} > 31 と {x-y, x+y} * {u-v, u+v} = {xu+yv-xu-xv, xu+yv+yu+xv} > 32 から確認できる。 > 33 一般にNが2のベキの場合、xorの作用するビットのパターンが再帰的なので > 34 pre( [X,Y] ) = [pre(X)-pre(Y), pre(X)+pre(Y)] > 35 と再帰的にすれば成り立つ > 36 > 37 > 38 > 39 > 40 > 41 > 42 xorでなくて足し算の場合 > 43 > 44 int[] plus_count(int[] a, int[] b, int N) > 45 { > 46 int[] c = {0} * 2N; > 47 for(int i=0; i<N; ++i) > 48 for(int k=0; k<N; ++k) > 49 c[i + k] += a[i] * b[k]; > 50 return c; > 51 } > 52 > 53 も同様に。 > 54 > 55 int[] plus_count_ultra(int[] a, int[] b, int N) > 56 { > 57 a = plus_pre(a); > 58 b = plus_pre(b); > 59 c = a * b; // element-wise product > 60 return plus_post(c); > 61 } > 62 > 63 とできるのが FFT