🌑

小羊儿的心情天空

Camp Day3 Div2

Jan 22, 2019 由 小羊

Problem B - 集合

一个计算几何题目,如果两人线段与圆相交,那么作切线然后再绕圆弧走;如果两人线段与圆相离,那么两人汇合点与圆心连线是汇合点与两人之间的角平分线,二分即可。(目前不是AC代码)

#include <bits/stdc++.h>
using namespace std;
typedef long long lld;
const double eps = 1e-8;
const double PI = acos(-1.0);
inline int sgn(double t) { return fabs(t) < eps ? 0 : t > 0 ? 1 : -1; }

struct point
{
    double x, y;
    point operator+(const point &a) const { return point(x+a.x, y+a.y); }
    point operator-(const point &a) const { return point(x-a.x,y-a.y); }
    double operator^(const point &a) const { return x*a.y-y*a.x; }
    double operator*(const point &a) const { return x*a.x+y*a.y; }
    point operator/(double a) const { return point(x/a, y/a); }
    point operator*(double a) const { return point(x*a, y*a); }
    double length() const { return sqrt(x*x+y*y); }
    point(double _x, double _y) noexcept { x = _x; y = _y; }
    double rad(const point &a, const point &b) { point q(a-*this), r=(b-*this); return acos(q*r/q.length()/r.length()); }
    point regulize(double len) { return *this/(length()/len); }
    double arg() const { return atan2(y, x); }
};

struct line
{
    double A, B, C;
    const point &s, &e;
    line(const point &a, const point &b) noexcept : s(a), e(b) { A = B = C = 0; }
    void load() { A = e.y - s.y; B = s.x - e.x; C = s.y*(e.x-s.x) - (e.y-s.y)*s.x; }
    double dis(const point &a) const { return fabs(A*a.x + B*a.y + C) / sqrt(A*A + B*B); }
    double disP2seg(point p) { if (sgn((p-s)*(e-s))<0 || sgn((p-e)*(s-e))<0) return min((p-s).length(),(p-e).length()); else return dis(p); }
};

struct circle : point
{
    double r;
    circle() noexcept : point(0,0) { r = 0; }
    bool seg(line v) { double dst = v.disP2seg(*this); return sgn(dst-r)<0; }
};

point A1(0,0), A2(0,0);
line px(A1, A2);
circle C;

double ans_seg()
{
    // line corrupts circle
    double bigRad = C.rad(A1, A2);
    double d1 = (A1-C).length();
    double d2 = (A2-C).length();
    double sRad1 = acos(C.r / d1);
    double sRad2 = acos(C.r / d2);
    double dd1 = sqrt(d1*d1-C.r*C.r);
    double dd2 = sqrt(d2*d2-C.r*C.r);
    double ans = dd1 + dd2 + C.r*(bigRad-sRad1-sRad2);
    return ans;
}

double f(double alpha)
{
    point A3 = (A1-C)*alpha + (A2-C)*(1-alpha) + C;
    return A3.rad(A1, C) - A3.rad(A2, C);
}

double ans_bs()
{
    double l = 0.0, r = 1.0, m;
    double fl = f(l), fr = f(r), fm;

    while (r - l > eps)
    {
        m = (r+l)/2;
        fm = f(m);
        if (sgn(fl) * sgn(fm) > 0)
            l = m, fl = fm;
        else
            r = m, fr = fm;
    }

    // so now m is the correct alpha percent.
    point A3 = ((A1-C)*m + (A2-C)*(1-m)).regulize(C.r) + C;
    return (A1-A3).length() + (A2-A3).length();
}

void solve()
{
    scanf("%lf %lf %lf %lf", &A1.x, &A1.y, &A2.x, &A2.y);
    scanf("%lf %lf %lf", &C.x, &C.y, &C.r); px.load();
    printf("%.3lf\n", C.seg(px) ? ans_seg() : ans_bs());
}

int main()
{
    int T; scanf("%d", &T);
    while (T--) solve();
    return 0;
}

Problem F - 小清新数论

对于 1n1071 \le n \le 10^7 (for div2) 或 1n10111 \le n \le 10^{11}

ans=i=1nj=1nμ(gcd(i,j))ans = \sum_{i=1}^n\sum_{j=1}^n\mu(\gcd(i,j))

那么对于div2来说,可以暴力枚举 gcd(i,j)\gcd(i,j) 的值,也可以进一步继续推算,即

ans=d=1nμ(d)×i=1nj=1n[gcd(i,j)=d]=d=1nμ(d)×k=1n/dμ(k)ndk2=t=1nnt2dk=tμ(d)μ(k)\begin{aligned} ans & = \sum_{d=1}^n\mu(d)\times\sum_{i=1}^n\sum_{j=1}^n[\gcd(i,j)=d] \\ & = \sum_{d=1}^n\mu(d)\times\sum_{k=1}^{\lfloor n/d\rfloor}\mu(k)\left\lfloor\frac{n}{dk}\right\rfloor^2 \\ & = \sum_{t=1}^n\left\lfloor\frac{n}{t}\right\rfloor^2\sum_{dk=t}\mu(d)\mu(k) \end{aligned}

那么显然后面那个是 μμ\mu * \mu,也是个积性函数,可以用线性筛来完成。对于div1来说,这时候需要用到 O(n23)O(n^\frac23) 的杜教筛来完成工作。

f(n)=dk=nμ(d)μ(k)f(n)=\sum_{dk=n}\mu(d)\mu(k),现在欲求 F(n)=i=1nf(i)F(n)=\sum_{i=1}^nf(i)。考虑 M(n)=i=1nμ(i)M(n)=\sum_{i=1}^n\mu(i),以及 μμ1=μe=μ\mu * \mu * 1 = \mu * e = \mu,那么有

F(n)=M(n)d=2nFndF(n)=M(n)-\sum_{d=2}^nF\left\lfloor\frac{n}{d}\right\rfloor

由于杜教筛主要计算的是 F(n/i)F(n/i)M(n/i)M(n/i) 的值,并且存在记忆化,那么两个函数的计算实际上仍然是原来的复杂度,只要线性筛预处理前 O(n23)O(n^\frac23) 的前缀和,后面这个函数计算复杂度依然是 O(n23)O(n^\frac23)

然后在实际写代码的时候遇到了一个坑点:(n/k)*(n/k)神tm会爆long long我没注意到。。

#include <bits/stdc++.h>
typedef long long lld;
using namespace std;
const size_t MAXN = 1e7+5;
const int mod = 998244353;
char isnp[MAXN], miu[MAXN];
int pcn, pri[MAXN], mm[MAXN];
int _M1[MAXN], _F1[MAXN];
map<lld,lld> _M2, _F2;

void init_prime()
{
    miu[1] = mm[1] = _M1[1] = _F1[1] = 1;
    for (int i = 2; i < MAXN; i++)
    {
        if (!isnp[i])
        {
            miu[i] = -1;
            mm[i] = -2;
            pri[pcn++] = i;
        }

        for (int j = 0; j < pcn; j++)
        {
            int k = i * pri[j];
            if (k >= MAXN) break;
            isnp[k] = 1;
            if (i % pri[j] == 0)
            {
                miu[k] = 0;
                mm[k] = ((i/pri[j])%pri[j]==0) ? 0 : mm[i/pri[j]];
                break;
            }
            else
            {
                miu[k] = -miu[i];
                mm[k] = -2*mm[i];
            }
        }

        _M1[i] = (_M1[i-1] + miu[i]) % mod;
        _F1[i] = (_F1[i-1] + mm[i]) % mod;
    }
}

lld M(lld n)
{
    if (n < MAXN) return _M1[n];
    else if (_M2.count(n)) return _M2[n];

    lld ans = 1;
    for (lld k = 2, fk; k <= n; k = fk+1)
    {
        fk = n / (n/k);
        ans -= (fk-k+1)*M(n/k)%mod;
    }

    ans = (ans % mod + mod) % mod;
    return _M2[n] = ans;
}

lld F(lld n)
{
    if (n < MAXN) return _F1[n];
    else if (_F2.count(n)) return _F2[n];

    lld ans = M(n);
    for (lld k = 2, fk; k <= n; k = fk+1)
    {
        fk = n / (n/k);
        ans -= (fk-k+1)*F(n/k)%mod;
    }

    ans = (ans % mod + mod) % mod;
    return _F2[n] = ans;
}

int main()
{
    init_prime();
    lld n; scanf("%lld", &n);
    lld ans = 0;

    for (lld k = 1, fk; k <= n; k = fk+1)
    {
        fk = n / (n/k);
        ans += (((n/k)%mod)*((n/k)%mod))%mod * (F(fk)-F(k-1)) % mod;
    }

    ans = (ans % mod + mod) % mod;
    printf("%lld\n", ans);
    return 0;
}

Problem G - 排列

简单题,过。