30 分题解

30 分的话，可以当做一个大模拟+基础计算几何的题目来做，需要选手具备初中平面几何相关的知识点。

我们对平面上的点设置一个基础的类，并且重载基础的一些四则运算，以及如何计算点的距离这些的，面向对象的基础写法会方便我们后续的书写。

struct pt
{
    double x, y;
    pt(double xx = 0, double yy = 0):x(xx), y(yy) {}
    // standard
    pt operator+(const pt& o) const { return pt(x + o.x, y + o.y); }
    pt operator-(const pt& o) const { return pt(x - o.x, y - o.y); }
    pt operator*(double d) const { return pt(x * d, y * d); }
    pt operator/(double d) const { return pt(x / d, y / d); }
    pt operator-() const { return pt() - (*this); }
    double dis2(const pt& o) const { return (x - o.x) * (x - o.x) + (y - o.y) * (y - o.y); }
    void operator=(const pt& o) { x = o.x, y = o.y; }
}

操作 1：平移

这个最简单，直接调用上面的加法运算就行。

操作 2：旋转

首先计算当前点与轴点距离的向量，将向量逆时针旋转 $\theta$ 之后，加上轴点的坐标即为新的位置。

$(x,y)$ 绕原点逆时针旋转 $\theta$ 角之后答案是 $(x\cos\theta - y\sin \theta ,x\sin \theta +y\cos \theta)$。对应操作也就有了。

    // oper 2 : rotate
    void rotate(const pt& pivot, double theta)
    {
        pt vec = (*this) - pivot;
        double _cos = cos(theta), _sin = sin(theta);
        vec = pt(_cos * vec.x - _sin * vec.y, _sin * vec.x + _cos * vec.y);
        (*this) = (vec + pivot);
    }

操作 3：缩放

这个比 2 简单，计算距离向量之后，直接乘以缩放比例，再加上轴点就得到了。

    // oper 3 : scale
    void scale(const pt& pivot, double lambda)
    {
        pt vec = (*this) - pivot;
        (*this) = pivot + (vec * lambda);
    }

操作 5：投影

显然，操作 4 和 5 显然是先计算 5 更合适。因为计算出 5 的点之后，操作 4 的结果就是原来的点以操作 5 的结果为轴点做轴对称就行了。

此处给出一个不依赖 $\tan\theta$ 存在的解法（也就是说平行于 $y$ 轴的直线也能做）。

我们把这个直线改为斜截式 $ax+by=c$，将原本的 $(x_1,y_1)$ 向 $ax+by=c$ 做投影，则伸过去的直线方程即为 $b(x-x_1)=a(y-y_1)$，两条直线解方程不难得到结果是 $(\dfrac{b^2x_1-aby_1+ac}{a^2+b^2},\dfrac{-abx_1+a^2y_1+bc}{a^2+b^2})$。

然后将原本直线带入进来，$a=-\tan\theta,b=1,c=y_0$，进行运算就行了。

操作 4：轴对称

设原本点是 $(x_1,y_1)$，操作 5 得到的结果是 $(x_2,y_2)$，则操作 4 的结果就是 $(2x_2-x_1,2y_2-y_1)$。

操作 6 和 7

这个不赘述了，有上面的基础操作就可以做。

然后我们在区间操作的时候使用 switch-case 判断使用哪种操作，不同模块封装一下即可。

时间复杂度为 $O(nq)$，空间复杂度 $O(n)$。由于官方数据的区间操作纯随机，所以每次生成区间长度的期望为 $n/3$，这导致 30 分做法在 5 秒左右就能跑完，可以直接获得 100 分。在复刻数据当中，我们加入了加强数据，重新卡回了 30 分。

其中 FastIO 可以做到等价于关闭同步流的效果，将输入输出分离之后方便对拍。

#include <stdio.h>
#include <math.h>
struct fastIO
{
    static const int BUFF_SZ = 1 << 18;
    char inbuf[BUFF_SZ], outbuf[BUFF_SZ];
    fastIO()
    {
        setvbuf(stdin, inbuf, _IOFBF, BUFF_SZ);
        setvbuf(stdout, outbuf, _IOFBF, BUFF_SZ);
    }
} IO;
const int N = 500010;
int n, q, op, l, r;
double a, b, theta, lambda, _y0;
struct pt
{
    double x, y;
    pt(double xx = 0, double yy = 0):x(xx), y(yy) {}
    // standard
    pt operator+(const pt& o) const { return pt(x + o.x, y + o.y); }
    pt operator-(const pt& o) const { return pt(x - o.x, y - o.y); }
    pt operator*(double d) const { return pt(x * d, y * d); }
    pt operator/(double d) const { return pt(x / d, y / d); }
    pt operator-() const { return pt() - (*this); }
    double dis2(const pt& o) const { return (x - o.x) * (x - o.x) + (y - o.y) * (y - o.y); }
    void operator=(const pt& o) { x = o.x, y = o.y; }

    // oper 2 : rotate
    void rotate(const pt& pivot, double theta)
    {
        pt vec = (*this) - pivot;
        double _cos = cos(theta), _sin = sin(theta);
        vec = pt(_cos * vec.x - _sin * vec.y, _sin * vec.x + _cos * vec.y);
        (*this) = (vec + pivot);
    }
    // oper 3 : scale
    void scale(const pt& pivot, double lambda)
    {
        pt vec = (*this) - pivot;
        (*this) = pivot + (vec * lambda);
    }
    // oper 5 : project
    // y = (tan\theta)x + y0 -> ax + by = c
    // a = -tan\theta b = 1 c = y0
    void project(double theta, double y0)
    {
        double a = -tan(theta), b = 1, c = y0;
        (*this) = pt(b * b * x - a * b * y + a * c, -a * b * x + a * a * y + b * c) / (a * a + b * b);
    }
    // oper 4 : mirror
    void mirror(double theta, double y0)
    {
        pt c = (*this);
        c.project(theta, y0);
        (*this) = (c * 2.0) - (*this);
    }
} p[N], tmp;
void op1()
{
    scanf("%d%d%lf%lf", &l, &r, &a, &b), tmp = pt(a, b);
    for (int i = l; i <= r; ++i) p[i] = p[i] + tmp;
}
void op2()
{
    scanf("%d%d%lf%lf%lf", &l, &r, &a, &b, &theta), tmp = pt(a, b);
    for (int i = l; i <= r; ++i) p[i].rotate(tmp, theta);
}
void op3()
{
    scanf("%d%d%lf%lf%lf", &l, &r, &a, &b, &lambda), tmp = pt(a, b);
    for (int i = l; i <= r; ++i) p[i].scale(tmp, lambda);
}
void op4()
{
    scanf("%d%d%lf%lf", &l, &r, &theta, &_y0);
    for (int i = l; i <= r; ++i) p[i].mirror(theta, _y0);
}
void op5()
{
    scanf("%d%d%lf%lf", &l, &r, &theta, &_y0);
    for (int i = l; i <= r; ++i) p[i].project(theta, _y0);
}
void op6()
{
    scanf("%d%d", &l, &r), tmp = pt();
    for (int i = l; i <= r; ++i) tmp = tmp + p[i];
    tmp = tmp / (r - l + 1);
    printf("%lf %lf\n", tmp.x, tmp.y);
}
void op7()
{
    double res = 0;
    scanf("%d%d%lf%lf", &l, &r, &a, &b), tmp = pt(a, b);
    for (int i = l; i <= r; ++i) res += p[i].dis2(tmp);
    printf("%lf\n", res);
}
int main()
{
    scanf("%d%d", &n, &q);
    for (int i = 1; i <= n; ++i) scanf("%lf%lf", &p[i].x, &p[i].y);
    while (q--)
    {
        scanf("%d", &op);
        switch (op)
        {
        case 1: op1(); break;
        case 2: op2(); break;
        case 3: op3(); break;
        case 4: op4(); break;
        case 5: op5(); break;
        case 6: op6(); break;
        case 7: op7(); break;
        default: break;
        }
    }
}

100 分题解

显然，我们需要一个类似 $O(q\log n)$ 的做法才能完成，不难想到可以利用线段树的方式维护。且在这种区间操作较为复杂的情况下，我们可以考虑用矩阵乘法的方式来进行维护。

不考虑维护平方时的操作

由于所有操作均为旋转、平移、投影操作组合而成，我们可以在线段树上维护 $[x,y,1]^T$ 这样的向量。以下具体操作我们不再重复讲解，只给出对应的操作矩阵，表示如何从 $[x_1,y_1,1]^T$ 变换到 $[x_2,y_2,1]^T$。

操作 1：

$$\left[\begin{matrix} x_2\\y_2\\1\end{matrix}\right]=\left[\begin{matrix} 1&0&a\\0&1&b\\0&0&1\end{matrix}\right]\left[\begin{matrix} x_1\\y_1\\1\end{matrix}\right]\\ $$

操作 2：

$$\left[\begin{matrix} x_2\\y_2\\1\end{matrix}\right]=\left[\begin{matrix} 1 & 0 &a\\0&1&b\\0&0&1\end{matrix}\right]\left[\begin{matrix} \cos\theta & -\sin \theta &0\\\sin \theta&\cos\theta&0\\0&0&1\end{matrix}\right]\left[\begin{matrix} 1 & 0 &-a\\0&1&-b\\0&0&1\end{matrix}\right]\left[\begin{matrix} x_1\\y_1\\1\end{matrix}\right]\\=\left[\begin{matrix} \cos\theta & -\sin \theta & a(1-\cos\theta)+b\sin \theta\\\sin \theta&\cos\theta&b(1-\cos\theta)-a\sin\theta \\0&0&1\end{matrix}\right]\left[\begin{matrix} x_1\\y_1\\1\end{matrix}\right] \\ $$

操作 3：

$$ \left[\begin{matrix} x_2\\y_2\\1\end{matrix}\right]=\left[\begin{matrix} 1 & 0 &a\\0&1&b\\0&0&1\end{matrix}\right]\left[\begin{matrix} \lambda & 0 &0\\0&\lambda&0\\0&0&1\end{matrix}\right]\left[\begin{matrix} 1 & 0 &-a\\0&1&-b\\0&0&1\end{matrix}\right]\left[\begin{matrix} x_1\\y_1\\1\end{matrix}\right]\\=\left[\begin{matrix} \lambda & 0 & (1-\lambda)a\\0&\lambda&(1-\lambda)b\\0&0&1\end{matrix}\right]\left[\begin{matrix} x_1\\y_1\\1\end{matrix}\right] \\ $$

操作 5：

$$\left[\begin{matrix} x_2\\y_2\\1\end{matrix}\right]=\left[\begin{matrix} \dfrac{b^2}{a^2+b^2} & -\dfrac{ab}{a^2+b^2} & \dfrac{ac}{a^2+b^2}\\-\dfrac{ab}{a^2+b^2}&\dfrac{a^2}{a^2+b^2}&\dfrac{bc}{a^2+b^2}\\0&0&1\end{matrix}\right] \left[\begin{matrix} x_1\\y_1\\1\end{matrix}\right]\\ $$

其中 $a=-\tan\theta,b=1,c=y_0$。

操作 4：

$$ \left[\begin{matrix} x_2\\y_2\\1\end{matrix}\right] =\left[\begin{matrix} \dfrac{2b^2}{a^2+b^2}-1 & -\dfrac{2ab}{a^2+b^2} & \dfrac{2ac}{a^2+b^2}\\-\dfrac{2ab}{a^2+b^2}&\dfrac{2a^2}{a^2+b^2}-1&\dfrac{2bc}{a^2+b^2}\\0&0&1\end{matrix}\right]\left[\begin{matrix} x_1\\y_1\\1\end{matrix}\right] \\ $$

其中 $a=-\tan\theta,b=1,c=y_0$。

至此，我们可以维护区间的向量和以及矩阵懒标记，懒标记下放采用矩阵和向量的乘法操作，区间和用向量加和维护即可。

考虑维护平方和的操作

我们会发现操作 7 要维护的距离平方这个东西其实也很难维护。这个时候我们可以推算一下 $\sum x^2$ 和 $\sum y^2$ 需要怎么通过 $\sum x$ 和 $\sum y$ 维护，不难发现我们还需要维护 $\sum xy$ 一项才能操作（也就是 $xy$ 乘积的区间和）。

上述 5 种操作的变换矩阵都满足 $\left[\begin{matrix} F_1 & F_2 & F_3\\F_4&F_5&F_6\\0&0&1\end{matrix}\right]$ 的形式，因此我们只需要做一次推算即可扩展到所有的运算。

当 $x'=F_1x+F_2y+F_3,~y'=F_4x+F_5y+F_6$ 时，我们有：

• $x'^2=(F_1x+F_2y+F_3)^2=F_1^2x^2+F_2^2y^2+2F_1F_2xy+2F_1F_3x+2F_2F_3y+F_3^2$
• $y'^2=(F_4x+F_5y+F_6)^2=F_4^2x^2+F_5^2y^2+2F_4F_5xy+2F_4F_6x+2F_5F_6y+F_6^2$
• $x'y'=F_1F_4x^2+F_2F_5y^2+(F_1F_5+F_2F_4)xy+(F_1F_6+F_3F_4)x+(F_2F_6+F_3F_5)y+F_3F_6$

所以我们最终要维护的向量以及变换矩阵为：

$$\left[\begin{matrix} x'\\y'\\x'^2\\y'^2\\x'y'\\1\end{matrix}\right]=\left[\begin{matrix} F_1 & F_2&0&0&0 & F_3\\F_4&F_5&0&0&0&F_6\\2F_1F_3&2F_2F_3&F_1^2&F_2^2&2F_1F_2&F_3^2\\2F_4F_6&2F_5F_6&F_4^2&F_5^2&2F_4F_5&F_6^2\\F_1F_6+F_3F_4&F_2F_6+F_3F_5&F_1F_4&F_2F_5&F_1F_5+F_2F_4&F_3F_6\\0&0&0&0&0&1\end{matrix}\right] \left[\begin{matrix} x\\y\\x^2\\y^2\\xy\\1\end{matrix}\right]\\ $$

对于前面 5 种操作，把他们对应的 F1 到 F6 传入，然后转成这个通用的 $6\times 6$ 的变换矩阵就行了。

对于操作 6，向量中自带 $\sum x$ 和 $\sum y$，维护出来然后除以区间长度 $(r-l+1)$ 就行了。

对于操作 7，$\sum(x-a)^2+(y-b)^2=\sum x^2+\sum y^2-2a\sum x-2b\sum y+(r-l+1)(a^2+b^2)$ 因此在得到区间向量和之后，带入对应公式计算就行了。

综上，时间复杂度为 $O(k^3q\log n)$，空间复杂度 $O(k^2n)$，其中 $k=6$，代表矩阵的存储以及矩阵乘法的运算。使用动态开点的线段树建树方式可以大幅降低空间常数，顺带优化时间常数，期望得分 100。

#include <stdio.h>
#include <string.h>
#include <math.h>
struct fastIO
{
    static const int BUFF_SZ = 1 << 18;
    char inbuf[BUFF_SZ], outbuf[BUFF_SZ];
    fastIO()
    {
        setvbuf(stdin, inbuf, _IOFBF, BUFF_SZ);
        setvbuf(stdout, outbuf, _IOFBF, BUFF_SZ);
    }
} IO;
const double eps = 3e-4;
struct matrix
{
    double a[6][6];
    matrix(bool identity = 1)
    {
        memset(a, 0, sizeof(a));
        if (identity)
            for (int i = 0; i < 6; ++i)
                a[i][i] = 1.0;
    }
    matrix operator*(const matrix& o) const
    {
        matrix ret = matrix(0);
        for (int i = 0; i < 6; ++i)
            for (int j = 0; j < 6; ++j)
                for (int k = 0; k < 6; ++k)
                    ret.a[i][j] += a[i][k] * o.a[k][j];
        return ret;
    }  
};
struct _vector
{
    double a[6]; // [\sum x, \sum y, \sum x^2, \sum y^2, \sum xy, 1]^T
    _vector() { memset(a, 0, sizeof(a)); }
    _vector(double x, double y)
    {
        a[0] = x, a[1] = y, a[2] = x * x, a[3] = y * y, a[4] = x * y, a[5] = 1.0;
    }
    _vector operator+(const _vector& o) const
    {
        _vector ret;
        for (int i = 0; i < 6; ++i)
            ret.a[i] = a[i] + o.a[i];
        return ret;
    }
    friend _vector operator* (const matrix& l, const _vector& r)
    {
        _vector ret;
        for (int i = 0; i < 6; ++i)
            for (int j = 0; j < 6; ++j)
                ret.a[i] += l.a[i][j] * r.a[j];
        return ret;
    }
};
matrix op;
void set_op(double F1, double F2, double F3, double F4, double F5, double F6)
{
    op.a[0][0] = F1, op.a[0][1] = F2, op.a[0][5] = F3;
    op.a[1][0] = F4, op.a[1][1] = F5, op.a[1][5] = F6;
    op.a[2][0] = 2 * F1 * F3, op.a[2][1] = 2 * F2 * F3, op.a[2][2] = F1 * F1;
    op.a[2][3] = F2 * F2, op.a[2][4] = 2 * F1 * F2, op.a[2][5] = F3 * F3;
    op.a[3][0] = 2 * F4 * F6, op.a[3][1] = 2 * F5 * F6, op.a[3][2] = F4 * F4;
    op.a[3][3] = F5 * F5, op.a[3][4] = 2 * F4 * F5, op.a[3][5] = F6 * F6;
    op.a[4][0] = F1 * F6 + F3 * F4, op.a[4][1] = F2 * F6 + F3 * F5, op.a[4][2] = F1 * F4;
    op.a[4][3] = F2 * F5, op.a[4][4] = F1 * F5 + F2 * F4, op.a[4][5] = F3 * F6;
}

const int N = 500010;
int n, q;
int opnum, l, r;

double a, b, theta, lambda, yzero;


double x[N], y[N];

struct node
{
    int l, r;
    int lc, rc;
    _vector sum;
    bool tag_mul;
    matrix tag_mat;
    void mul(const matrix& m)
    {
        tag_mul = 1;
        tag_mat = m * tag_mat, sum = m * sum;
    }
} tr[N << 1];
int cnt, rt;
void pushup(int u) { tr[u].sum = tr[tr[u].lc].sum + tr[tr[u].rc].sum; }
void pushdown(int u)
{
    if (tr[u].tag_mul)
    {
        tr[tr[u].lc].mul(tr[u].tag_mat), tr[tr[u].rc].mul(tr[u].tag_mat);
        tr[u].tag_mul = 0, tr[u].tag_mat = matrix();
    }
}
int _build(int l, int r)
{
    int u = ++cnt;
    tr[u].l = l, tr[u].r = r;
    if (l == r)
        return tr[u].sum = _vector(x[l], y[l]), u;
    int m = (l + r) >> 1;
    tr[u].lc = _build(l, m), tr[u].rc = _build(m + 1, r);
    pushup(u);
    return u;
}
void build() { rt = _build(1, n); }
void _modify(int u, int l, int r, const matrix& m)
{
    if (l > tr[u].r || r < tr[u].l)
        return;
    if (l <= tr[u].l && tr[u].r <= r)
    {
        tr[u].mul(m);
        return;
    }
    pushdown(u);
    _modify(tr[u].lc, l, r, m), _modify(tr[u].rc, l, r, m);
    pushup(u);
}
void modify(int l, int r, const matrix& m) { _modify(rt, l, r, m); }
_vector _query(int u, int l, int r)
{
    if (l > tr[u].r || r < tr[u].l)
        return _vector();
    if (l <= tr[u].l && tr[u].r <= r)
        return tr[u].sum;
    pushdown(u);
    return _query(tr[u].lc, l, r) + _query(tr[u].rc, l, r);
}
_vector query(int l, int r) { return _query(rt, l, r); }
void set_op1(double a, double b)
{
    set_op(1, 0, a, 0, 1, b);
}

void set_op2(double a, double b, double theta)
{
    set_op(cos(theta), -sin(theta), a * (1 - cos(theta)) + b * sin(theta), sin(theta), cos(theta), b * (1 - cos(theta)) - a * sin(theta));
}

void set_op3(double a, double b, double lambda)
{
    set_op(lambda, 0, (1 - lambda) * a, 0, lambda, (1 - lambda) * b);
}
// y = \tan\theta x + y0 -> ax + by = c
void set_op4(double a, double b, double c)
{
    double deno = a * a + b * b;
    set_op((2 * b * b) / deno - 1, (-2 * a * b) / deno, (2 * a * c) / deno, (-2 * a * b) / deno, (2 * a * a) / deno - 1, (2 * b * c) / deno);
}
void set_op5(double a, double b, double c)
{
    double deno = a * a + b * b;
    set_op((b * b) / deno, (-a * b) / deno, (a * c) / deno, (-a * b) / deno, (a * a) / deno, (b * c) / deno);
}
int main()
{
    scanf("%d%d", &n, &q);
    for (int i = 1; i <= n; ++i)
        scanf("%lf%lf", &x[i], &y[i]);
    build();
    while (q--)
    {
        scanf("%d%d%d", &opnum, &l, &r);
        switch (opnum)
        {
        case 1:
            scanf("%lf%lf", &a, &b);
            set_op1(a, b);
            break;
        case 2:
            scanf("%lf%lf%lf", &a, &b, &theta);
            set_op2(a, b, theta);
            break;
        case 3:
            scanf("%lf%lf%lf", &a, &b, &lambda);
            set_op3(a, b, lambda);
            break;
        case 4:
            scanf("%lf%lf", &theta, &yzero);
            set_op4(-tan(theta), 1, yzero);
            break;
        case 5:
            scanf("%lf%lf", &theta, &yzero);
            set_op5(-tan(theta), 1, yzero);
            break;
        case 7:
            scanf("%lf%lf", &a, &b);
            break;
        default:
            break;
        }
        if (opnum <= 5)
            modify(l, r, op);
        else if (opnum == 6)
        {
            _vector tmp = query(l, r);
            printf("%lf %lf\n", (tmp.a[0] + eps) / (r - l + 1), (tmp.a[1]  + eps) / (r - l + 1));
        }
        else
        {
            _vector tmp = query(l, r);
            double x = tmp.a[2] - 2 * tmp.a[0] * a + (r - l + 1) * a * a;
            double y = tmp.a[3] - 2 * tmp.a[1] * b + (r - l + 1) * b * b;
            printf("%lf\n", x + y + eps);
            
        }
    }
}