小花梨的函数

题目： https://acm.ecnu.edu.cn/contest/173/statements/

题解：https://acm.ecnu.edu.cn/upload/mirror/cache/%E2%80%9C%E7%BE%8E%E7%99%BB%E6%9D%AF%E2%80%9D%E4%B8%8A%E6%B5%B7%E5%B8%82%E9%AB%98%E6%A0%A1%E5%A4%A7%E5%AD%A6%E7%94%9F%E7%A8%8B%E5%BA%8F%E8%AE%BE%E8%AE%A1%E9%82%80%E8%AF%B7%E8%B5%9B%EF%BC%88%E5%8D%8E%E4%B8%9C%E7%90%86%E5%B7%A5%E5%A4%A7%E5%AD%A6%EF%BC%89%E9%A2%98%E8%A7%A3.pdf

#include<bits/stdc++.h>
#define LL long long
#define ULL unsigned long long
#define rep(i,j,k) for(int i=j;i<=k;i++)
#define dep(i,j,k) for(int i=k;i>=j;i--)
#define INF 0x3f3f3f3f
#define mem(i,j) memset(i,j,sizeof(i))
#define make(i,j) make_pair(i,j)
#define pb push_back
using namespace std;
const int N = 2e7 + 10;
int f[N];
bool vis[N];
int pre[N];
int mu[N];
void init(int n) {
    int tot = 0;
    vis[1] = mu[1] = 1;
    rep(i, 2, n) {
        if(!vis[i]) pre[++tot] = i,mu[i] = -1;
        rep(j, 1, tot) {
            if(i * pre[j] > n) break;
            vis[i * pre[j]] = 1;
            if(i % pre[j]) mu[i * pre[j]] = -mu[i];
            if(i % pre[j] == 0) break;
        }
    }
    rep(i, 1, n) mu[i] += mu[i-1];
}
LL cal(int n,int m) {
    LL ans=0;
    for(LL l = 1, r; l <= n; l = r + 1) {
        r = min(n / (n / l),m / (m / l));
        ans += (mu[r] - mu[l - 1]) * (n / l) * (m / l);
    }
    return ans;
}
int main() {
    int n, m, a, p;
    scanf("%d %d %d %d",&n, &m, &a, &p);
    f[1] = f[2] = 1;
    int cnt;
    rep(i, 3, N-1) {
        f[i] =(f[i-1] + f[i-2]) % p;
        if(f[i] == 1 && f[i-1] == 1) {
            cnt = i - 2;
            break;
        }
    }
    if(n > m) swap(n, m);
    init(n);
    LL ans = 0; LL tmp = 1;
    rep(i, 1, n) {
        tmp = tmp * a % cnt;
        LL tot = cal(n / i, m / i) % p;///求gcd的 次数
        ans = (ans + tot * f[(tmp + cnt - 1) % cnt] % p) % p; ///求和；
    }
    printf("%lld
", ans);
    return 0;
}