ProjectEuler习题

<翻译+题解。>
既然要做英文题，那么就好好做，认真写结题报告。

Problem 401

题意：求$\sum_{i=1}^{n}g(i)$其中$g(n) = \sum_{d|n}d^2$。
$n \leq 10^{15}$，答案对$10^9$取模。

---

题解：由反演的一些经验，我们知道可以枚举因数$d$.
答案变成了：

$$\sum_{d = 1}^{n}\lfloor \frac {n}{d}\rfloor *d^2$$
注意到$\lfloor \frac {n}{d}\rfloor$的取值是根号级别的。
分块暴力，统计平方和。
但是怎么算平方和呢？
有公式：$\sum_{i=1}^{n}i^2 = \frac {(n+1)*n*(2n+1)}{6}$
看似问题是解决了，但是我们注意到6是没有逆元的，我们需要知道一个式子：

$$a/b\ mod\ m \equiv a \ mod \ (bm) / b$$
那么我们只需要对 $2*10^9$取模，然后求出 $3$的逆元即可。
最后的答案记得除以2，并且要开unsigned long long.

#include <bits/stdc++.h>
#define Rep(i,n) for(int i = 1;i <= n;++ i)
#define rep(i,a,b) for(int i = a;i <= b;++ i)
#define Dwn(i,n) for(int i = n;i ;i --)
#define dwn(i,a,b) for(int i = a;i >= b;-- i)
using namespace std;
const int Mod = 2e9;
typedef unsigned long long ull;
ull rev = 666666667;
ull Get(ull n)
{
    n %= Mod;
    return (n + 1) * n % Mod * (2ull * n + 1) % Mod * rev % Mod;
}
ull calc(long long i,long long j,long long d)
{
    return d % Mod * ((Get(j) - Get(i - 1) + Mod) % Mod) % Mod;
}
int main()
{
    ull x = 1e15;
    ull ans = 0;
    ull i = 1,n = x;
    for(ull j;i <= n;)
    {
        j = n / (n / i);
        ans = (ans + calc(i,j,n / i)) % Mod;
        i = j + 1;
    }
    printf("%llu\n",ans / 2);
    return 0;
}

知识点：$n/d$的取值，逆元的转化，[1,n]平方和公式。

Problem 464

先拆前缀和，如果不满足的话，要么N大要么P大。
把式子化出来的话：

1. $$99N(a,b) > 100P(a,b)$$ $$99N(b) - 100P(b) > 99N(a-1)-100P(a-1)$$
2. $$99P(a,b) > 100N(a,b)$$ $$99P(b) - 100N(b) > 99P(a-1)-100N(a-1)$$

然后这样我们用个set查一下就没了。
……set不滋磁这个操作……我写了个TreapQAQ

#include <bits/stdc++.h>
#define Rep(i,n) for(int i = 1;i <= n;++ i)
#define v edge[i].to
#define RepG(i,x) for(int i = head[x];~ i; i = edge[i].next)
#define fl edge[i].f
#define vfl edge[i ^ 1].f

#define u t[x]
#define o t[y]
#define lc ch[0]
#define rc ch[1]
#define tc ch[ty]
#define vc ch[!ty]
using namespace std;
const int N = 20000005;
typedef long long LL;
class HH
{
    public:
    int rt,tot;
    struct Node{int ch[2],sz,fix,cnt;LL val;void Set(LL _){val = _;fix = rand() + 1 + rand() * 4;cnt = sz = 1;ch[0] = ch[1] = 0;}}t[N];
    void Upd(int x){u.sz = t[u.lc].sz + t[u.rc].sz + u.cnt;}
    void Rot(int &x,bool ty){int y = u.tc;u.tc = o.vc,o.vc = x,Upd(x),Upd(x = y);}
    void ins(int &x,LL _)
    {
        if(!x){x = ++ tot,u.Set(_);return;}
        if(u.val == _){u.cnt ++;u.sz ++;return;}
        bool ty = u.val < _;
        ins(u.tc,_);
        (u.fix > t[u.tc].fix) ? Upd(x) : Rot(x,ty);
    }
    int Qry(int x,LL _)
    {
        if(!x)return 0;
        if(u.val < _)return t[u.lc].sz + u.cnt + Qry(u.rc,_);
        else return Qry(u.lc,_);
    }
}a,b;
int n,m,p[N],mu[N],cnt,sn[N],sp[N];
bool isp[N];

int main()
{
    n = 20000000;
    LL ans = 0;
    mu[1] = 1;
    for(int i = 1;i <= n;++ i)ans += (n - i + 1);
    for(int i = 2;i <= n;++ i)
    {
        if(!isp[i])mu[i] = -1,p[++ cnt] = i;
        for(int j = 1;j <= cnt && i * p[j] <= n;++ j)
        {
            isp[i * p[j]] = 1;
            if(i % p[j] == 0){mu[i * p[j]] = 0;break;}
            else mu[i * p[j]] = -mu[i];
        }
    }
    for(int i = 1;i <= n;++ i)sn[i] = sn[i - 1] + (mu[i] == -1),sp[i] = sp[i - 1] + (mu[i] == 1);
    for(int i = 1;i <= n;++ i)
    {
        a.ins(a.rt,99ll * sn[i - 1] - 100ll * sp[i - 1]);
        ans -= a.Qry(a.rt,99ll * sn[i] - 100ll * sp[i]); 

        b.ins(b.rt,99ll * sp[i - 1] - 100ll * sn[i - 1]);
        ans -= b.Qry(b.rt,99ll * sp[i] - 100ll * sn[i]);
    }
    printf("%lld\n",ans);
    return 0;
}

知识点：平衡树，前缀和拆分，二维数点。

Problem 565

C++11大法吼！
首先我们注意到，$\omega(p^k) = 1+p+p^2+p^3...$
2017是个质数。
我们先考虑筛出根号$n$这个范围内的所有素数。
之后找所有$p$满足$p = k*2017-1$
这个咋找呢？
我们可以用小于$\sqrt{N/2017}$的数字去筛这些数字对吧。
瓶颈在筛这个素数上。
如果P为合数，那么有$d*P/d = 2016(mod\ 2017),d<\sqrt n$
那么有$k*2017- 1 \mod d = 0$。
$k=2017^{-1}mod d$
那么我们就可以筛$k$，看这个东西是否满足。
$k = td + (2017^{-1} \mod d)$
枚举$t$,用$d$去筛，这一步复杂度是$loglogn$。
然后我们就要考虑重合的部分，但是其他地方其实都显然了……
因为要去吃饭，所以到时候再补。

#include <bits/stdc++.h>
#define Rep(i,n) for(int i = 1;i <= n;++ i)
#define v edge[i].to
#define RepG(i,x) for(int i = head[x];~ i; i = edge[i].next)
#define fl edge[i].f
#define vfl edge[i ^ 1].f

using namespace std;
typedef long long LL;
const int N = 1000000;
const int NR = N + 5;
LL n,arr[(LL)(1e11+1) / 2017],m = 2017;
vector<int>pr;
vector<pair<LL,LL> >re;
bool isp[NR * 50];

void init(int ZZ)
{
    for(int i = 2;i * i <= ZZ;++ i)
        if(!isp[i])
            for(int j = i * i;j <= ZZ;j += i)isp[j] = 1;
    for(int i = 2;i <= ZZ;++ i)if(!isp[i])pr.push_back(i);
}

int Pw(int a,int b)
{
    int cur = 1;
    int p = b - 2;
    while(p)
    {
        if(p & 1)cur = 1ll * cur * a % b;
        p >>= 1;
        a = 1ll * a * a % b;
    }
    return cur;
}
int main()
{
    n = 1e11;
    LL L = n / m;
    LL SQ = (int)sqrt(n);
    init(SQ);
    for(LL i = m - 1;i <= n;i += m)arr[i / m] = i;
    int LEN = (LL)(1e11+1) / 2017;
    for(int i : pr)
        if(i != m)
        {
            int rev = Pw(m,i);
            for(int j = ((arr[rev - 1] == i) ? rev + i : rev) - 1;j < LEN;j += i)
                arr[j] = 0;
        }
    for(int z = 0;z < LEN;++ z)
    {
        LL i = arr[z];
        if(i)
        {
            if(i <= SQ)re.push_back(make_pair(i,i*i));
            else re.push_back(make_pair(i,0ll));
        }
    }
    for(LL i : pr)
    {
        for(LL j = 1 + i + i * i,a = i * i;a <= n;a *= i,j += a)
        {
            if(j % m == 0)
            {
                if(i * a <= n)re.push_back(make_pair(a,i * a));
                else re.push_back(make_pair(a,0ll));
            }
        }
    }
    sort(re.begin(),re.end());
    LL ans = 0;
    for(int i = 0;i < re.size();++ i)
    {
        LL a = re[i].first,b = re[i].second;
        ans += a * (n / a) * (n / a + 1) / 2;
        if(b)ans -= b * (n / b) * (n / b + 1) / 2;
    }
    for (int i = 0; i < re.size(); i ++)
        for (int j = i + 1; j < re.size(); j ++) 
        {
            LL temp = re[i].first * re[j].first;
            if (temp > n)
            {
                if (j == i + 1)i = re.size();
                break;
            }
            for (LL k = temp; k <= n; k += temp)
                ans -= k;
        }
    cout << ans << endl;
    return 0;
}