Sumdiv
Time Limit: 1000MS | Memory Limit: 30000K | |
Total Submissions: 19576 | Accepted: 4930 |
Description
Consider two natural numbers A and B. Let S be the sum of all natural divisors of A^B. Determine S modulo 9901 (the rest of the division of S by 9901).
Input
The only line contains the two natural numbers A and B, (0 <= A,B <= 50000000)separated by blanks.
Output
The only line of the output will contain S modulo 9901.
Sample Input
2 3
Sample Output
15
Hint
2^3 = 8.
The natural divisors of 8 are: 1,2,4,8. Their sum is 15.
15 modulo 9901 is 15 (that should be output).
The natural divisors of 8 are: 1,2,4,8. Their sum is 15.
15 modulo 9901 is 15 (that should be output).
Source
首先因数和定理p^k=(1+p^1+p^2+...p^k);
之后构造矩阵求p^k%MOD
之前想过利用逆元来做,但是当(p-1)%MOD==0时不存在逆元,所以就没有继续想。但是逆元的思想比较直接简便。
#include<stdio.h>
#include<iostream>
#include<algorithm>
#include<string.h>
#include<string>
#include<math.h>
using namespace std;
#define ll long long
#define F(x,a,b) for(ll x=a;x<=b;x++)
#define me(x) memset(x,0,sizeof(x))
#define MOD 9901
#define _fast F(i,1,3)F(j,1,3)F(k,1,3) b[i][j]=(b[i][j]+a[i][k]*a[k][j])%MOD;
#define _reset F(i,1,3)F(j,1,3){a[i][j]=b[i][j]; b[i][j]=0;}
#define _orz F(i,1,3)F(j,1,3)F(k,1,3) b[i][j]=(b[i][j]+a[i][k]*c[k][j])%MOD;
ll p[100005],isp[100005],cnt,x,y;
ll b[10][10],a[10][10],c[10][10];
struct p{ll p,k;}rec[100000];
void _pr(ll x)
{
me(isp);me(p);
F(i,2,x)
{
if (!p[i]) isp[++p[0]]=i;
F(j,1,p[0])
{
if (isp[j]*i>x) break;
p[isp[j]*i]=1;
if (i%isp[j]==0) break;
}
}
}
void _s(ll x)
{
cnt=0;
me(rec);
int t=x;
F(i,1,p[0])
{
if (isp[i]*isp[i]>x) break;
if (t%isp[i]==0)
{
rec[++cnt].p=isp[i];
while (t%isp[i]==0)
{
t/=isp[i];
rec[cnt].k++;
}
}
}
if (t>1) {rec[++cnt].p=t;rec[cnt].k=1;}
}
void _f(ll kk)
{
if (kk==1) return;
if (kk&1){_f(kk-1); _orz _reset}
else {_f(kk/2); _fast _reset}
}
int main()
{
ll n,k;
_pr(100000);
cin>>n>>k;
if (k==0) {cout<<"1"<<endl;return 0;}
cnt=0;
_s(n);
ll ans=1; a[1][1]=1;a[1][3]=1;a[2][1]=1;
F(i,1,cnt)
{
me(a);me(b);me(c);
a[1][1]=1;a[1][3]=1;a[2][1]=1;
a[3][3]=rec[i].p;
F(i,1,3)F(j,1,3) c[i][j]=a[i][j];
_f(rec[i].k*k);
ans=ans*(a[1][1]*1+a[1][3]*rec[i].p)%MOD;
}
cout<<ans<<endl;
}