【SDOI2017】数字表格 [莫比乌斯反演]

  1 #include<iostream>  
  2 #include<string>  
  3 #include<algorithm>  
  4 #include<cstdio>  
  5 #include<cstring>  
  6 #include<cstdlib>  
  7 #include<cmath>
  8 using namespace std; 
  9 typedef long long s64;
 10 
 11 const int ONE = 1e6+2;
 12 const int MOD = 1e9+7;
 13 const int PHI = 1e9+6;
 14 
 15 int T;
 16 int n,m;
 17 int prime[ONE],p_num,miu[ONE];
 18 int Niyu[ONE];
 19 int F[ONE];
 20 bool isp[ONE];
 21 int Ans;
 22 
 23 namespace input
 24 {
 25         const int BufferSize = 1 << 16 | 1;
 26         
 27         char buffer[BufferSize];
 28         char *head = buffer + BufferSize;
 29         const char *tail = head;
 30         
 31         inline char nextChar()
 32         {
 33             if (head == tail)
 34             {
 35                 fread(buffer, 1, BufferSize, stdin);
 36                 head = buffer;
 37             }
 38             return *head++;
 39         }
 40         
 41         inline int get()
 42         {
 43             static char c;
 44             while ((c = nextChar()) < ‘0‘ || c > ‘9‘);
 45         
 46             int res = c - ‘0‘;
 47             while ((c = nextChar()) >= ‘0‘ && c <= ‘9‘)
 48                 res = res * 10 + c - ‘0‘;
 49             return res;
 50         }
 51 }
 52 using input::get;
 53 
 54 inline int Quickpow(int a,int b)
 55 {
 56         int res = 1;
 57         while(b)
 58         {
 59             if(b&1) res = (s64)res*a%MOD;
 60             a = (s64)a*a%MOD;
 61             b>>=1;
 62         }
 63         return res;
 64 }
 65 
 66 inline void Deal_first(int MaxN)
 67 {
 68         F[0]=0;    F[1]=1;
 69         int val=1;
 70         for(int i=2; i<=MaxN; i++)
 71         {
 72             F[i] = F[i-1]+F[i-2];
 73             if(F[i] >= MOD) F[i] -= MOD;
 74             val=(s64)val*F[i]%MOD;
 75         }
 76         Niyu[MaxN] = Quickpow(val, MOD-2);
 77         for(int i=MaxN-1;i>=1;i--) Niyu[i] = (s64)Niyu[i+1] * F[i+1] % MOD;
 78         Niyu[0] = Niyu[1]; F[0]=1;
 79         for(int i=1; i<=MaxN; i++) F[i] = (s64)F[i]*F[i-1] % MOD;
 80         
 81         miu[1] = 1;
 82         for(int i=2; i<=MaxN; i++)
 83         {
 84             if(!isp[i])
 85                 prime[++p_num] = i, miu[i] = -1;
 86             for(int j=1; j<=p_num, i*prime[j]<=MaxN; j++)
 87             {
 88                 isp[i * prime[j]] = 1;
 89                 if(i % prime[j] == 0)
 90                 {
 91                     miu[i * prime[j]] = 0;
 92                     break;
 93                 }
 94                 miu[i * prime[j]] = -miu[i];
 95             }
 96             miu[i] += miu[i-1];
 97         }
 98 }
 99 
100 inline int f(int n,int m)
101 {
102         if(n > m) swap(n,m);
103         int Ans = 0;
104         for(int i=1,j=0; i<=n; i=j+1)
105         {
106             int x=n/i, y=m/i;
107             j = min(n/x, m/y);
108             Ans = ((s64)Ans + (s64)x * y%PHI * ((s64)miu[j] - miu[i-1] + PHI) )%PHI;
109         }
110         return Ans;
111 }
112 
113 inline void Solve()
114 {
115         n=get();    m=get();
116         if(n > m) swap(n,m);
117         Ans = 1;
118         for(int i=1,j=0; i<=n; i=j+1)
119         {
120             int x=n/i, y=m/i;
121             j = min(n/x, m/y);
122             Ans = (s64)Ans * Quickpow( (s64)F[j] * Niyu[i-1] % MOD , f(x,y) ) % MOD;
123         }
124         printf("%d\n",Ans);
125 }
126 
127 int main()
128 {
129         Deal_first(ONE-1);
130         T = get();
131         while(T--)
132             Solve();
133         return 0;
134 }