#include <iostream>
#include <cstdio>
#include <cstring>
#include <cmath>
#include <algorithm>
#include <vector>
#include <stdlib.h>

using namespace std;
int n;
const int N = 1010;
const double eps = 1e-7;

struct node {
    int x, y ,r;
    bool operator < (const node &a) const {
        if(x - r != a.x - a.r) return x - r < a.x - a.r;
        return r < a.r;
    }
}cc[N];

inline int cmp(double x, double y) {
    if(fabs(x - y) < eps) return 0;
    return x < y ? -1 : 1;
}

using pdd = pair<double, double>;
pdd q[N];
vector<pdd>segx;

inline double f(double _x) {

    int cnt = 0;
    for(int i = 0; i < n; ++i)
    {
        double xx = fabs(_x - cc[i].x),
               r = cc[i].r;

        if(cmp(xx, r) < 0) {
            double yy = sqrtl(r*r - xx*xx);
            q[cnt++] = {cc[i].y - yy, cc[i].y + yy};
        }
    }
    if(!cnt) return 0;

    sort(q, q+cnt);
    double res = 0;
    double _y1 = q[0].first,
           _y2 = q[0].second;

    for(int i = 1; i < cnt; ++i)
    {
        if(_y2 >= q[i].first)
            _y2 = max(_y2, q[i].second);
        else
            res += _y2 - _y1, 
            _y1 = q[i].first,
            _y2 = q[i].second;
    }
    return res + _y2 - _y1;

}

inline double simpson(double l, double r) {
    double mid = l / 2.0 + r / 2.0;                 // 积分部分
    return (r - l) * (f(l) + f(r) + 4 * f(mid)) / 6;
}

double recursive_calc(double l, double r, double s) {
    double mid = (r + l) / 2.00;
    double _L = simpson(l, mid), 
           _R = simpson(mid, r);
    if(!cmp(_L + _R, s))
        return _L + _R;
    return recursive_calc(l, mid, _L) + 
           recursive_calc(mid, r, _R);              // 递归调用部分
}

const int delta = 100;

void merge() {
    sort(segx.begin(), segx.end());
    vector<pdd> tmp;

    double st = -1e20,
           ed = -1e20;
    for(auto p : segx)
    {
        if(cmp(p.first, ed) <= 0) ed = max(ed, p.second);
        else {
            if(cmp(st, -1e20)) tmp.push_back({st, ed});
            st = p.first, ed = p.second;
        }
    }
    tmp.push_back({st, ed});
    segx = tmp;
}


int main() {
    scanf("%d", &n);

    for(int i = 0; i < n; ++i)
        scanf("%d %d %d", &cc[i].x, &cc[i].y, &cc[i].r),
        segx.push_back({cc[i].x - cc[i].r, cc[i].x + cc[i].r});

    merge();
    double res = 0;
    for(auto p: segx)
        res += recursive_calc(p.first, p.second,
                              simpson(p.first, p.second));

    printf("%.3lf", res);
    return 0;
}

// 好的数据结构，不应该写着写着就忘记最初要干什么
#include <iostream>
#include <cstring>
#include <cstdio>
#include <algorithm>
#include <cmath>

using namespace std;
using ll = long long;

const int N = 5E5+10, inf = 1e9;
int n, m;

// 常数有点大
struct splay {      // 需要父节点得信息可以被两个孩子算出来
    int son[2];     // 左右儿
    int val;        // 当前节点的值
    int size;       // 节点大小
    int fa;         // 父节点

    int rev;        // 整个区间是否翻转
    int same;       // 整个区间是否要变成相同的数
    // 两个的操作顺序对其它值得影响要考虑好，不然会很惨
    // 此处记录得是操作之后的值


    // 类似线段树的操作
    int f_sum;      // 最大字段和

    // 用来处理跨越区间的情况
    int pre_sum,    // 最大前缀和
        suf_sum;    // 最大后缀和
    int sum;        // 区间和

    void init(int _fa, int _val) {
        son[0] = son[1] = 0;
        size = 1;   fa = _fa;
        val = _val;
        rev = same = 0;

        sum = f_sum = val;
        pre_sum = suf_sum = max(0, val);
    }
}tr[N];

int root, nodes[N], tt;
int w[N];

// 每个点维护的是一个序列
// 为了找到第 k 个数，还需要 size

// 首先需要把序列变成 splay
// 比较平衡的做法是通过递归建成二叉树

// 内存回收记得做一做
// 把点放到垃圾桶里，用的时候直接从垃圾桶要。

/**
 *      套路：每次要操作区间的时候就把区间的前一个数和
 *            后一个数记录下来，然后把前一个数转到根，
 *            后一个数转到右儿子，那么待操作区间就在右儿子的左儿子上。
 *            o
 *             \
 *              o
 *             /
 *           <>
 *         <>ST<>
 *        <><><><>
**/

void push_up(int x) {
    auto &u = tr[x];
    auto &l = tr[u.son[0]],
         &r = tr[u.son[1]];
    u.size = l.size + r.size + 1;       // 包含自己
    u.sum = u.val + l.sum + r.sum;      // 所有值等于左 + 右 + 自己

    // 按照原先构想更新标记
    u.pre_sum = max(l.pre_sum, l.sum + u.val + r.pre_sum);
    u.suf_sum = max(r.suf_sum, r.sum + u.val + l.suf_sum);
    u.f_sum = max({l.f_sum, r.f_sum, 
                   l.suf_sum + u.val + r.pre_sum});
}   


inline void push_down(int x) {
    auto &u = tr[x];
    auto &l = tr[u.son[0]],
         &r = tr[u.son[1]];
    if(u.same) {
        u.same = u.rev = 0;
        if(u.son[0]) l.same = 1, l.val = u.val, 
                     l.sum = l.val * l.size;
        if(u.son[1]) r.same = 1, r.val = u.val, 
                     r.sum = r.val * r.size;

        if(u.val > 0)
        {
            if(u.son[0]) l.f_sum = l.pre_sum = l.suf_sum = l.sum;
            if(u.son[1]) r.f_sum = r.pre_sum = r.suf_sum = r.sum;
        }
        else {
            if(u.son[0]) l.f_sum = l.val, l.pre_sum = l.suf_sum = 0;
            if(u.son[1]) r.f_sum = r.val, r.pre_sum = r.suf_sum = 0;

        }
    }
    else if(u.rev) {
        u.rev = 0;
        l.rev ^= 1, r.rev ^= 1;
        swap(l.pre_sum, l.suf_sum);
        swap(r.pre_sum, r.suf_sum);
        swap(l.son[1], l.son[0]);
        swap(r.son[1], r.son[0]);
    }
}

inline void rotate(int x) {        // 模板
    int y = tr[x].fa;
    int z = tr[y].fa;
    int k = tr[y].son[1] == x;

    tr[z].son[tr[z].son[1] == y] = x, 
    tr[x].fa = z;

    tr[y].son[k] = tr[x].son[k ^ 1], 
    tr[tr[x].son[k ^ 1]].fa =y;

    tr[x].son[k ^ 1] = y;
    tr[y].fa = x;
    push_up(y), push_up(x);
}


void splay(int x, int k) {          // 模板
    while(tr[x].fa ^ k) {
        int y = tr[x].fa,
            z = tr[y].fa;
        if(z ^ k)
            if((tr[y].son[1] == x) ^ (tr[z].son[1] == y) )
                rotate(x);
            else rotate(y);
        rotate(x);
    }
    if(!k) root = x;
}

int get_k(int k) {                  // 模板
    int ru = root;
    while(ru) {
        push_down(ru);
        if(tr[tr[ru].son[0]].size >= k) ru = tr[ru].son[0];
        else if(tr[tr[ru].son[0]].size + 1 == k) return ru;
        else k -= tr[tr[ru].son[0]].size + 1,
             ru = tr[ru].son[1];
    }
}


int build (int l, int r, int fa)
{
    int mid = l + r >> 1;
    int ru = nodes[tt--];       // 从垃圾回收站里拿一个节点出来
    tr[ru].init(fa, w[mid]);    // 初始化节点

    if(l < mid) tr[ru].son[0] = build(l, mid - 1, ru);
    if(mid < r) tr[ru].son[1] = build(mid + 1, r, ru);
    push_up(ru);
    return ru;

}

void meme_dfs(int u) {          // 垃圾回收，从原有的树上面薅下来
    if(tr[u].son[0]) meme_dfs(tr[u].son[0]);
    if(tr[u].son[1]) meme_dfs(tr[u].son[1]);
    nodes[++tt] = u;
}


int main() {
    for(int i = 1; i < N; i++)
        nodes[++tt] = i;
    scanf("%d %d", &n, &m);

    tr[0].f_sum = -inf, w[0] = w[n + 1] = -inf;         // 初值还是要的
    for(int i = 1; i <= n; i++) scanf("%d", w + i);
    root = build(0, n + 1, 0);                          // 然后建树

    char op[32];
    while(m--)
    {
        scanf("%s", op);
        if(!strcmp(op, "INSERT")) {
            int posi, tot;
            scanf("%d %d", &posi, &tot);
            for(int i = 0; i < tot; ++i)
                scanf("%d", w + i);
            int l = get_k(posi + 1),        
                r = get_k(posi + 2);        // 因为是插入到后面
            splay(l, 0), splay(r, l);       // 转成上面讲的方式



            int u = build(0, tot - 1, r);

            tr[r].son[0] = u;
            push_up(r), push_up(l);
        }
        else if(!strcmp(op, "DELETE")) {
            int posi, tot;
            scanf("%d %d", &posi, &tot);

            int l = get_k(posi),
                r = get_k(posi + tot + 1);

            splay(l, 0), splay(r, l);

            meme_dfs(tr[r].son[0]);         // 转成上面讲的方式之后，直接回收孩子
            tr[r].son[0] = 0;               // 简单置零就好
            push_up(r), push_up(l);         // 然后将数据上传

        }
        else if(!strcmp(op, "MAKE-SAME")) {
            int posi, tot, c;
            scanf("%d %d %d", &posi, &tot, &c);

            int l = get_k(posi),
                r = get_k(posi + tot + 1);
            splay(l, 0), splay(r, l);

            auto &son = tr[tr[r].son[0]];                   // 左子树

            son.same = 1, son.val = c, son.sum = c * son.size;
            if(c > 0) son.f_sum = son.pre_sum = son.suf_sum = son.sum;
            else son.f_sum = c, son.pre_sum = son.suf_sum = 0;

            push_up(r), push_up(l);
        }
        else if(!strcmp(op, "REVERSE")) {
            int posi, tot;
            scanf("%d %d %d", &posi, &tot);

            int l = get_k(posi),
                r = get_k(posi + tot + 1);
            splay(l, 0), splay(r, l);
            auto &son = tr[tr[r].son[0]];

            son.rev ^= 1;
            swap(son.pre_sum, son.suf_sum);
            swap(son.son[0], son.son[1]);
            push_up(r), push_up(l);
        }
        else if(!strcmp(op, "GET-SUM")) {
            int posi, tot;
            scanf("%d %d %d", &posi, &tot);

            int l = get_k(posi),
                r = get_k(posi + tot + 1);
            splay(l, 0), splay(r, l);
            printf("%d\n", tr[tr[r].son[0]].sum);

        }
        else printf("%d\n",tr[root].f_sum);
    }   
    return 0l;
}

#include <iostream>
#include <cstring>
#include <cstdio>
#include <cmath>
#include <algorithm>

using namespace std;
const int N = 1E5 + 10;

int a[N];
int F, n;
/* 表述得不是太好 */
// 应该说只围成一块，围起来的地要大于等于 F 块
// 那么我们计算上就可以使用双指针来做完这件事情
double pre[N];

bool check(double x)        // 判断大多数方案是否能被满足
{
    // 如何判断测算的平均值是可以的？首先拿所有值减去 x 然后算前缀和

    for(int i = 1; i <= n; ++i)
        pre[i] = pre[i - 1] + a[i] - x;
    double minx = 2e9+10;
    for(int i = 0, j = F; j <=n; j++, i++)
    {
        // 如何保证大于 F 块？
        // 很简单，只要统计顺过去的最小值就可以保证
        minx = min(minx, pre[i]);
        if(pre[j] >= minx) return 1;    // 找到
    }
    return 0;               // 找不到
}

int main() {
    // N 块地划分为至少 F 块地时，每块地包含的牛的平均值最大值 
    scanf("%d %d", &n, &F);

    double l = 0, r = 0;
    for(int i = 1; i <= n; ++i)
        scanf("%d", a + i), r = max(r, (double)a[i]);

    while(fabs(l - r) > 1e-6)
    {
        double mid = (l + r) / 2.0;
        if(check(mid)) l = mid;     // 仍有进步空间
        else r = mid;               // “进步过头”
    }
    printf("%lld", (long long)(r * 1000));
    return 0;
}

还注意到有一种实现二分的方式。利用了分得越小，剩得越多得原理。

bool check(double x)        // 判断大多数方案是否能被满足
{
    for(int i = 1; i <= n; ++i) pre[i] = pre[i - 1] + a[i] - x; 
    // 如何判断测算的平均值是可以的？首先拿所有值减去 x 然后算前缀和
    double minx = 2e9, maxn = -2e9;
    for(int i = F; i <= n; ++i)
    {
        minx = min(minx, (double)pre[i - F]);
        maxn = max(maxn, (double)pre[i] - minx);
    }
    if(maxn >= 0) return 1;
    return 0;               // 找不到
}

也完全可以不用二分去做。$O(n)$ 就办完了。

#include <iostream>
#include <cstring>
#include <cstdio>
#include <cmath>
#include <algorithm>

using namespace std;
const int N = 1E5 + 10;

int a[N];
int F, n;
/* 表述得不是太好 */
// 应该说只围成一块，围起来的地要大于等于 F 块
// 那么我们计算上就可以使用双指针来做完这件事情
double pre[N];
using ld = long double;
using ll = long long;

inline ld k(int l, int r) {return (ld)(pre[r] - pre[l]) / (r - l);}

int q[N << 1];

int main() {
    scanf("%d %d", &n, &F);
    for(int i = 1; i <= n; ++i)
        scanf("%d", a + i),
        pre[i] = pre[i-1] + a[i];

    /** 什么是平均数？
     * (pre[i] - pre[i - F]) / F 就相当于是平均数，而且不难发现，分得越小对平均数越好。
     * 那么，假定就一直在 F 取到。
     * 如果是为了找最大值而特意去算这个式子，也就需要维护一个上凸包？
     * 那么出列的条件依次是 队尾比上一个更劣，队头比下一个更劣。
     * 可以公式化它们的关系。
    **/

    int head = 1, tail = 0;

    ld res = 0.0L;
    for(int i = F; i <= n; ++i)
    {
        int delta = i - F;
        while (head < tail and k(q[tail], delta) <= k(q[tail - 1], delta))
            tail--;
        q[++tail] = delta;

        while(head < tail and k(q[head], i) <= k(q[head + 1], i))
            head++;
        res = max(res, k(q[head], i));
    }
    printf("%lld", (ll)(res * 1000));
    return 0;
}

时间上：

• $O(n\log n)$：133ms。
• $O(n)$：64ms。

对初学的人来说，理解还是要花很大的功夫……

#include <iostream>
#include <cstring>
#include <cstdio>
#include <algorithm>
#include <string>

/****************************************************************************************************
 *   后缀数组 用倍增 nlogn 排完所有的后缀
 *   需要：
 *      排名为 i 的后缀的在原字符串中的次序                             @sa[i]
 *      第 i 个后缀的排名                                               @rk[i]
 *      排名为 i 和 i - 1 的最长公共前缀的值                            @height[i]
 *      第 i 个后缀和排名在第 i 个后缀前面一位后缀的最长公共前缀长度    @h[i] = height[rk[i]]
 *
 *      对于字符串 aanklllv 而言
 *      后缀有    v                     次序是      8 
 *                lv                                7
 *                llv                               6
 *                lllv                              5
 *                klllv                             4
 *                nklllv                            3
 *                anklllv                           2
 *                aanklllv                          1                       
 ****************************************************************************************************/

const int N = 1E6+10;
char str[N];
int sa[N], x[N], y[N], cnt[N];
int rk[N], height[N];
int n, m;

/* 映射关系有点复杂 */
void get_sa() {

    /* 第一次基数排序 */
    for(int i = 1; i <= n; ++i) cnt[x[i] = str[i]] ++;                  // 统计各字符的个数 
                                                                        // 并计算各个位置第一关键字离散化后是什么值
    for(int i = 2; i <= m; ++i) cnt[i] += cnt[i - 1];                   // 计算**字典序**下的，小于等于当前关键字的数量
    for(int i = n; i; --i) sa[cnt[x[i]] --] = i;                        // 从后往前确定排序后排名在原字符串中的次序
                                                                        // 并将关键字数量 --

    // 做完基数排序后，视字符串后缀为有序后的比较。

    for(int k = 1; k <= n; k <<= 1)                                     // 倍增
    {
        int num = 0;

        for(int i = n - k + 1; i <= n; ++i) y[++num] = i;               // 第二关键字排名最小的部分（因为没有）   
        for(int i = 1 ; i <= n; ++i)
            if(sa[i] > k)                                               // 某个后缀的第二关键字是另一个后缀的第一关键字
                y[++num] = sa[i] - k;                                   // 要求某个后缀的位置要 > k
        /* 计算第二关键字并按第二关键字排序 */


        /* 第二次基数排序 */
        for(int i = 1; i <= m; ++i) cnt[i] = 0;                     
        for(int i = 1; i <= n; ++i) cnt[x[i]] ++;                       
        for(int i = 2; i <= m; ++i) cnt[i] += cnt[i-1];
        for(int i = n; i ; i--) sa[cnt[x[y[i]]]--] = y[i], y[i] = 0;    // 第二关键字来排，排完以后 y 没有用了

        std::swap(x, y);                                                // 把 x 的信息存到 y 里
        x[sa[1]] = 1, num = 1l;

        for(int i = 2; i <= n; ++i)                                     // 离散化前 2k 个字符
            x[sa[i]] = (y[sa[i - 1] + k] == y[sa[i] + k] &&
                        y[sa[i]] == y[sa[i - 1]]) ?
                        num : ++num;

        if(num == n) break;                                             // 完全区别开所有后缀
        m = num;                                                        // 不行就付给下一轮
    }
}

void get_height() {                                                     // 最长公共前缀
    for(int i = 1; i <= n; i ++) rk[sa[i]] = i;
    for(int i = 1, k = 0; i <= n; ++i) {
        if(rk[i] == 1) continue;
        if(k) k--;
        int j = sa[rk[i] - 1];
        while(i + k <= n and j + k <= n and str[i + k] == str[j + k]) k++;
        height[rk[i]] = k;
    }
}

int main (){

    scanf("%s", str + 1);
    n = strlen(str + 1);
    m = 122;

    get_sa();
    get_height();
    for(int i = 1l; i <= n; ++i) printf("%d ", sa[i]);
    puts("");
    for(int i = 1l; i <= n; ++i) printf("%d ", height[i]);
    return 0;
} 

#include <iostream>
#include <cstdio>
#include <cstring>
#include <cmath>
#include <algorithm>
using namespace std;
using db = long double;
const db pi = (db)acos(-1);
const int M = 3e5 + 10;
char sa[M], sb[M];
class cp
{
public:
    db x, y;
    cp(db xx = 0, db yy = 0) : x(xx), y(yy) {}
    cp operator+(const cp &a) { return cp(a.x + x, y + a.y); }
    cp operator-(const cp &a) { return cp(x - a.x, y - a.y); }
    cp operator*(const cp &a) { return cp(x * a.x - y * a.y, x * a.y + y * a.x); }
} a[M];
int lena, lenb, tot, bits;
int rev[M];
int res[M];
void fft(cp a[], int len, int sig = 1)
{
    for (int i = 0; i < len; ++i)
        if (i < rev[i])
            swap(a[i], a[rev[i]]);
    for (int mid = 1; mid < len; mid <<= 1)
    {
        cp w(cos(pi / mid), sig * sin(pi / mid));
        for (int i = 0; i < len; i += mid << 1)
        {
            cp w0(1, 0);
            for (int j = 0; j < mid; j++, w0 = w0 * w)
            {
                auto x = a[i + j], y = w0 * a[i + j + mid];
                a[i + j] = x + y, a[i + j + mid] = x - y;
            }
        }
    }
}
int main()
{
    scanf("%s%s", sa, sb);
    lena = strlen(sa), lenb = strlen(sb);
    int La = lena - 1, Lb = lenb - 1;
    for (int i = 0; i < lena; ++i)
        a[i].x = sa[La - i] - 48;
    for (int i = 0; i < lenb; ++i)
        a[i].y = sb[Lb - i] - 48;
    while ((1 << bits) < lena + lenb)
        bits++;
    tot = 1 << bits;
    for (int i = 0; i < tot; ++i)
        rev[i] = (rev[i >> 1] >> 1) | ((i & 1) << (bits - 1));
    fft(a, tot);
    for (int i = 0; i < tot; ++i)
        a[i] = a[i] * a[i];
    fft(a, tot, -1);
    int t = 0, k = 0;
    for (int i = 0; i < tot || t; ++i)
    {
        a[i].y /= 2.0 * tot; // 这里最需要去掉前导零。
        t += a[i].y + 0.5;
        res[k++] = t % 10;
        t /= 10;
    }
    while (k - 1 > 0 and !res[k - 1])
        k--;
    for (int i = k - 1; ~i; --i)
        printf("%d", res[i]);
    return 0;
}

上高中以后一直好奇最速下降曲线是怎么算出来的。到了大学怎么还能留着问题呢？今天就用漏洞百出的证法把它扬了。

前提

假定无摩擦，起点为 $(0,0)$，终点是 $(x_0, y_0)$，整套系统只有重力场和最速降线，并建立 $y$ 轴竖直向下，$x$ 轴水平向右的左手坐标系。

分析过程

一

对质点而言，由经典力学动能定理有，$\dfrac{1}{2}mv^2=mgy$。消去质量，并变形得到 $v^2=2gy$。

同时又由于质点在最速降线上运动，此时也有 $\dfrac{\sqrt{(\mathrm{d}x)^2 + (\mathrm{d}y)^2}}{\mathrm{d} t} = v $。进一步变形得到 $\dfrac{ \sqrt{1 + (y’)^2}\mathrm{d}x}{\mathrm{d} t} = v $。

那么对于时间的微分而言，有 $\mathrm{d} t=\sqrt{\dfrac{ 1 + (y’)^2}{2gy}}\mathrm{d}x$。

积分回来就有

$$t=\int_0^{x_0}\sqrt{\dfrac{ 1 + (y’)^2}{2gy}}\mathrm{d}x$$

问题便是选取合适的 $y$ 使得 $t$ 在积分之后取到最小值。

不难发现此时的函数涉及到了三个变量：自变量 $x$，因变量 $y$ 以及因变量的导数 $y’$。

除此之外还能注意到 $y(0) = 0, y(x_0) = y_0$。然后就似乎没有额外的限制条件了。

二

我们从物理常识可以知道待求的曲线一定存在、存在极值、连续、可导且二阶可导。

而为了研究如何找到这样的函数，我们可以很自然地类比微分寻找整个函数的极值点一样 引入函数的变动量（或者称之为函数的变分） $\delta(x)$，在两点间 所有可能 的函数为 $\bar{y}$，最优解为 $y$。

那么，$\bar{y}=y+\delta$。从而有 $\bar{y}’=y’+\delta’$。相应地，要求 $\delta(0)=\delta(x_0) = 0$ 以固定起始点和终点。

而为了使用到极限和微分求极值的理论，我们还需要对能式子中任意性最大的 $\delta(x)$ 做变换 $\delta(x) = \xi\eta (x)$。其中 $\xi$ 是常数。

具体而言，我们考虑三个变量所构成的多元函数 $\bar{F}(x,y+\xi\eta ,y’ + \xi\eta’ )=\bar{F}(x,\bar{y},\bar{y}’)$ 以及形成最优解函数的多元函数 $F(x,y,y’)$，其中最优解函数存在于 $\lim\limits_{\xi \to 0}\bar{F}(x,\bar{y},\bar{y}’)$ 所确定的函数族之中。

同时，$\bar{F}(x,\bar{y},\bar{y}’)$ 也是一个 关于 $\xi$ 的函数。

$$I(\xi) = t = \int_0^{x_0}\bar{F}(x,\bar{y},\bar{y}’)\mathrm d x$$

最优解可通过 $\xi$ 逐渐趋近以至于成为零得到，且这个过程隐含一个条件：当 $\xi$ 趋近于零时，$I(\xi)$ 整体对 $\xi$ 的微分一定是零。这样才能保证 $I(\xi)$ 取到极值。

那么整个过程可以被表示为
$$\lim_{\xi \to 0}\dfrac{\mathrm d I(\xi)}{\mathrm d \xi} = \lim_{\xi \to 0}\Big(\dfrac{\mathrm d }{\mathrm d \xi}\int_0^{x_0}\bar{F}(x,\bar{y},\bar{y}’)\mathrm d x\Big) = \lim_{\xi \to 0}\int_0^{x_0}\dfrac{\mathrm d \bar{F}}{\mathrm d \xi}\mathrm d x = 0$$

根据函数性质，以及链式法则，可以得到

$$\dfrac{\mathrm d \bar{F}}{\mathrm d \xi} = \dfrac{\partial \bar{F}}{\partial x}\cdot\dfrac{\partial x}{\partial \xi} + \dfrac{\partial \bar{F}}{\partial \bar{y}}\cdot\dfrac{\partial \bar{y}}{\partial \xi} + \dfrac{\partial \bar{F}}{\partial \bar{y}’}\cdot\dfrac{\partial \bar{y}’}{\partial \xi} = 0 + \eta\dfrac{\partial \bar{F}}{\partial \bar{y}} + \eta’\dfrac{\partial \bar{F}}{\partial \bar{y}’}$$

进一步有

$$\begin{aligned}\lim_{\xi \to 0}\dfrac{\mathrm d I(\xi)}{\mathrm d \xi} &= \lim_{\xi \to 0}\int_0^{x_0} \eta\dfrac{\partial \bar{F}}{\partial \bar{y}} + \eta’\dfrac{\partial \bar{F}}{\partial \bar{y}’}\mathrm d x \\
&= \int_0^{x_0} \lim_{\xi \to 0}\Big(\eta\dfrac{\partial \bar{F}}{\partial (y + \xi\eta)} + \eta’\dfrac{\partial \bar{F}}{\partial (y’ + \xi\eta’)}\Big)\mathrm d x \\ \\
&= \int_0^{x_0}\Big( \eta\dfrac{\partial F}{\partial y} + \eta’\dfrac{\partial F}{\partial y’}\Big)\mathrm d x = \int_0^{x_0}\eta\dfrac{\partial F}{\partial y} \mathrm d x + \left. \eta\dfrac{\partial F}{\partial y’}\right|_0^{x_0} - \int_0^{x_0}\eta \mathrm d \dfrac{\partial F}{\partial y’} \\
&=\int_0^{x_0}\eta\Big(\dfrac{\partial F}{\partial y} - \dfrac{\mathrm d}{\mathrm d x }\dfrac{\partial F}{\partial y’}\Big) \mathrm d x + \color{red}{0} = 0\end{aligned}$$

标红是因为代入上下限始终有 $\eta(0) = \eta(x_0) = 0$。

然后我们考虑一个问题：当满足 $\eta(x_1) = \eta(x_2) = 0$ 时，积分 $\displaystyle\int_{x_1}^{x_2}\eta(x)f(x)\mathrm d x=0$ 会导出什么结论。

由于 $\eta(x)$ 随意，不妨构造为 $\eta(x) = -f(x)(x - x_1)(x - x_2)$ 以满足题目条件。此时可以有 $\displaystyle \int_{x_1}^{x_2} -[f(x)]^2(x - x_1)(x - x_2)\mathrm d x$。

因为 $[f(x)]^2 \geq 0;\forall x \in [x_1,x_2],-(x-x_1)(x- x_2) \geq 0$，如果存在值使得 $[f(x)]^2$ 大于零，则会直接导致算出的积分值大于零，进一步与已知的函数值产生矛盾，所以 唯一的可能 只有 $[f(x)]^2 = f(x) = 0$。

那么回归到刚刚的讨论，此时就可以有 $\dfrac{\partial F}{\partial y} - \dfrac{\mathrm d}{\mathrm d x }\dfrac{\partial F}{\partial y’} = 0$。

三

再来看一看 最优解所构成的多元函数 对于 $x$ 的全微分。

$$\dfrac{\mathrm d F}{\mathrm d x} = \dfrac{\partial F}{\partial x} + \dfrac{\partial F}{\partial y}\dfrac{\partial y}{\partial x} + \dfrac{\partial F}{\partial y’}\dfrac{\partial y’}{\partial x} = \dfrac{\partial F}{\partial x} + \dfrac{\partial F}{\partial y}y’ + \dfrac{\partial F}{\partial y’}y’'$$

注意到 $\dfrac{\mathrm d}{\mathrm d x}\Big( y’\dfrac{\partial F}{\partial y’}\Big) = \color{red}{y’'\dfrac{\partial F}{\partial y’}} + y’\dfrac{\mathrm d}{\mathrm d x}\dfrac{\partial F}{\partial y’}$，那么上式减左式得到

$$\dfrac{\mathrm d}{\mathrm d x}\Big(F - y’\dfrac{\partial F}{\partial y’}\Big) = \dfrac{\partial F}{\partial x} + y’\Big( \dfrac{\partial F}{\partial y} - \dfrac{\mathrm d}{\mathrm d x}\dfrac{\partial F}{\partial y’}\Big)$$

根据第二节的结论，此时有 $\dfrac{\mathrm d}{\mathrm d x}\Big(F - y’\dfrac{\partial F}{\partial y’}\Big) = \dfrac{\partial F}{\partial x}$。

当整个函数没有出现 $x$ 的时候，$\dfrac{\mathrm d}{\mathrm d x}\Big(F - y’\dfrac{\partial F}{\partial y’}\Big) = 0$。这个时候就有 $F - y’\dfrac{\partial F}{\partial y’} = C$。

四

将上一节算得的结论运用到 $t=\displaystyle\int_0^{x_0}\sqrt{\dfrac{ 1 + (y’)^2}{2gy}}\mathrm{d}x$ 中的 $\sqrt{\dfrac{ 1 + (y’)^2}{2gy}}$ 就有

$$
\sqrt{\dfrac{ 1 + (y’)^2}{2gy}} - y’\cdot\dfrac{2y’}{\sqrt{ 2gy}}\cdot\dfrac{1}{2\sqrt{ 1 + (y’)^2}}=\dfrac{1+(y’)^2 - (y’)^2}{\sqrt{ 2gy( 1 + (y’)^2)}} = C
$$

于是我们得到微分方程 $y( 1 + (y’)^2) = \dfrac{C}{2g}$。将 $\dfrac{C}{2g}$ 设为 $k$，移项开方整理得到

$$\dfrac{\mathrm d y}{\sqrt{\dfrac{k-y}{y}}} = \mathrm d x$$

做变量代换 $y = \dfrac{k}{2}(1-\cos\theta)$ 得到

$$\dfrac{k}{2}\sqrt{\dfrac{1-\cos\theta}{1+\cos\theta}}\sin\theta\mathrm d \theta = \dfrac{k}{2}(1-\cos\theta)\mathrm d \theta = \mathrm d x$$

进而有 $x = \dfrac{k}{2}(\theta-\sin\theta) + C$。

因此 $\begin{cases}x = \dfrac{k}{2}(\theta-\sin\theta) + C \\ y = \dfrac{k}{2}(1-\cos\theta)\end{cases}$，其中更具体的表达式需要由 $(0,0),(x_0,y_0)$ 确定。

对于 $(0,0)$，此时有 $\begin{cases}x = \dfrac{k}{2}(\theta-\sin\theta) \\ y = \dfrac{k}{2}(1-\cos\theta)\end{cases}$。而对于另一个坐标而言，既然已经找到了表达式，就不必多言了。

悬链线

在悬挂问题中，悬挂点是固定的两个点，绳子可看成经过这两个点且总长为 $L$ 的曲线。

根据能量最低原理，绳子对应曲线的形状将会使得绳子的质心最低，而曲线可用一个函数表示。

那么问题就变成求解函数的表示。

我们以地面水平方向为 $x$ 轴，两悬挂点水平方向中垂线为 $y$ 轴，则两悬挂坐标分别为 $-\dfrac{a}{2}$ 和 $\dfrac{a}{2}$ ，绳子对应的曲线函数为 $f(x)$，绳长 $L$，质量为 $m$，则其质心的 $y$ 坐标值 $y_\text{center} $ 不难根据质心方程列出¹

$$
y_\text{center} = \dfrac{\displaystyle\int_{0}^{L}f(x)\cdot\dfrac{m}{L} \mathrm{d}{s}}{m} = \dfrac{1}{L}\int_{-\frac{a}{2}}^{\frac{a}{2}} f(x)\sqrt{1+[f’(x)]^2} \mathrm{d}{x}
$$

沿用先前算得的结论，即可得到以下的微分方程。

$$f(x)\sqrt{1+[f’(x)]^2} - f(x)\dfrac{[f’(x)]^2}{\sqrt{1+[f’(x)]^2}}=\dfrac{f(x)}{\sqrt{1+[f’(x)]^2}}=C$$

进一步有

$$\dfrac{\mathrm d f(x)}{\mathrm d x}=\sqrt{\dfrac{[f(x)]^2 - C}{C}}$$

从而

$$\dfrac{\mathrm d f(x)}{\sqrt{[f(x)]^2 - C}}=\dfrac{\mathrm d x}{\sqrt{C}}$$

后面就不多解释了。

参考资料

1. 泛函（最简泛函）导数解法
2. 分析力学—最速降线与变分法

前言

这个月没什么好分享的了，干脆就把自己通过计算几何的计算过程来学 markdown 和 $\LaTeX$ 时写的文章从洛谷 第二次 转投这里。

第一次是在能让人绝望的 CSDN 里投的，现在想想完全不如在 AcWing 分享好。所以已经删掉了。虽然这里也可能还是有借鉴不带出处的情况出现，但反正无所谓，源和思路都开了，现在就任其恣意 抄袭 吧。

然而在此之前，我还有一件事要说——

什么杂七杂八的网站我就懒得管了，但在第一次公开发布文章的洛谷看见了，我还是非常不想吐槽洛谷的管理员把我自己修了两次的文章批掉，另写了一篇题解，却完全 不引用 的做法。即使不是大学生，也多少要有点著作权意识吧。哪天被查抄的是你以后自己“写好”的论文，我看你以后怎么办？

对了，也正是 因为不好维权，所以抄袭成风。这就进一步成为我们的 AI 做不过 ChatGPT、我们的开源生态总是失败的又一个因素。

哦还有，这场 Codeforces 之所以 unrate 也是因为出题人“借鉴了”一道题目的思路而被外国人扒了出来。

这里附上自己所看见的 “影子”。

还有管理员的文章。我别的不要，其他人也不需要主动对线浪费自己的时间和精力，如果你在某种情况下又再次看到了这里记得来道歉:) 看不到的话以后总有人会拿这篇文章说事的，即使 我强调过不要主动对线。

具体使用到的细节转换了讲法。但是作为偶然或者必然间看到的材料，我认为，完全有必要指出出处。因为这是 别人先想到的思路而不是你先想到的。

你只是通过各种各样的途径看到然后拿来主义罢了。单凭没有说明清楚就写了文章的这点，我只能表示怀疑以及非常的遗憾。希望下次自己自私一点，不要再写这种无利可图的文章了。

下面是我第一次投稿的时间。

活着很累的，能写一篇文章说说内心的感受或者思考问题的心路历程都是一种 奢侈的浪费。

临近期末外加决定是否继续学业的年纪，头痛到不管哪一门语言的代码都不想碰。不过还是要感谢这个学期，让我学会了盲打键盘并见到了一堆极其难绷的事情。反正，剩下的内容就是正文了。总之我还要留一个:)。

题面翻译

大意

以逆时针顺序给出笛卡尔平面上 $n (4 \leqslant n \leqslant 5000)$ 个点 $A_1,A_2,\ldots,A_n$ 的坐标，从而确定一凸多边形。该凸多边形的性质有如下保证：
1. 所有点的坐标均为整数且 $-10^{9} \leqslant x,y \leqslant 10^{9}$ 。
2. 对于该凸多边形的任意内角 $\alpha$ ，均有 $90^{\circ} \leqslant \alpha < 180^{\circ}$ 。

现从凸多边形的任意两边上分别取出两点，将所有由两点所确定的、长度不超过 $1$ 的弦所经过的凸多边形内部区域标红。求标红区域的面积。

输入格式

第一行输入 $n$ 。

之后 $n$ 行每行给出 $1$ 个点的坐标。

输出格式

只输出标红区域的面积，保留 $11$ 位小数。误差不超过 $10^{-4}$ 即视为正确。

没错，我就是译者。

正文

首先要说明的是这道题，代码实现并不是太难，难的是推式子还有调试。我会试着把过程讲解得尽可能细致。当然如有问题还请不吝指出，我会尽力和大伙一起解决的。

大致所需知识

• 曲线上点坐标的参数方程表示。
• 函数图像变换。
• 三角函数积化和差与和差化积公式。
• 定积分。
• 二分求根。

程序构建脉络

题目的描述已经很清楚了，对于这类考数学推演运算而不是数据结构或者纯算法的题目，我个人一开始的想法就是研究样例是怎样的。

先优先考虑凸多边形各边长度 均大于 $1$ 的情况。注意到弦长为 $1$ 的弦只会出现在凸多边形的两邻边之间。这启发我们：算出每一个角然后用角度计算相应面积即可。

我们先来看看单个符合要求 $(90^{\circ}\leqslant \alpha <180^{\circ})$ 的角 $\alpha$。

按照题目限制下的逻辑，可以想象一个长度为 $1$ 的刚性杆从边 $\text{BD}$ 顺畅地滑到边 $\text{BE}$。整个过程中，杆所经过的区域即为我们的目标。

显然由滑动的过程以及圆自身的性质可知，曲线绝不会是圆弧。那么此刻，问题转化为求出杆运动所生成包络线的函数表达式及其积分。

由于我们只在逻辑上知道这条线性质如何如何，而不知道这条线到底可以由哪一个初等的函数表出，为此，我们引入参数方程来代替相对更直观的函数表达式。

按上图，包络线明显与 $\alpha$ 有关。又注意到，在杆的运动过程中， $\theta$ 也是一直在变化的。那么我们理应使用 $\alpha$ 和 $\theta$ 作为参数方程的参数。

统一记包络线上的点为 $G$，其坐标均为 $(x,y)$ 的抽象形式。则点 $G$ 可由两参数 $\theta$ 与 $\alpha$ 表示为 $\begin{matrix}x(\theta,\alpha)\ y(\theta,\alpha)\end{matrix}$。同时，直线 $l_{HI}$ 的方程可经以下推导得出。

在 $\triangle \text{HIB}$ 中，由正弦定理得到 $\dfrac{x_\text{I}}{\sin(\pi-\theta-\alpha)}=\dfrac{\text{HI}}{\sin\alpha}$，那么，进一步就可以有 $\dfrac{x_\text{I}}{\sin(\theta+\alpha)}=\dfrac{1}{\sin\alpha}$。

从而，我们可以得到以下结果。

$$ \begin{aligned}
&\Rightarrow x_\text{I}=\dfrac{\sin(\theta+\alpha)}{\sin\alpha},y_\text{J}=\dfrac{\sin(\theta+\alpha)}{ \sin\alpha}\tan\theta
\\& \therefore \ l_\text{HI}:\dfrac{\sin\alpha \cdot x}{\sin(\theta+\alpha)}+\dfrac{\sin\alpha \cdot y}{\tan\theta\sin(\theta+\alpha)}=1
\end{aligned}$$

我们还可从包络线切线的性质出发，得到 $l_\text{HI}$ 的斜率。注意，按照图上的标法，斜率需要加负号。这里再提供另一种不需要管斜率，算起来相对更简单的推导过程。

此处不以题解中显示的锐角 $\theta$ 作为参数，而是以其补角作为 以下的参数 $\theta$。

参数方程部分：
$$\begin{aligned}&\dfrac{x_\text{I}}{\sin(\theta-\alpha)}=\dfrac{1}{\sin\alpha},y_\text{J}=-\tan\theta \cdot x_\text{I}\\
\\ \Rightarrow & x_\text{I}=\dfrac{\sin(\theta-\alpha)}{\sin\alpha}
, y_\text{J}=-\dfrac{\sin(\theta-\alpha)}{\sin\alpha}\tan\theta
\\
\\l:&\dfrac{x(\theta,\alpha)}{x_\text{I}}+\dfrac{y(\theta,\alpha)}{y_\text{J}}=1,\ \dfrac{\mathrm{d}y}{\mathrm{d}x}=\dfrac{y_\theta}{x_\theta}=\tan\theta
\\ \Rightarrow & -\tan\theta\cdot x(\theta,\alpha)+y(\theta,\alpha)=-\dfrac{\sin(\theta-\alpha)}{\sin\alpha}\tan\theta
,\ y_\theta=\tan\theta \cdot x_\theta
\\
\therefore &-\dfrac{x(\theta,\alpha)}{\cos^2\theta}-x_\theta\tan\theta+y_\theta=\csc\alpha(-\cos(\theta-\alpha)\tan\theta-\sin(\theta-\alpha)\sec^2\theta)\\
\\\Rightarrow &-\dfrac{x(\theta,\alpha)}{\cos^2\theta}=
-\dfrac{\sin(\theta-\alpha)+\sin\theta\cos\theta\cos(\theta-\alpha)}{\sin\alpha\cos^2\theta}
\\\Rightarrow &x(\theta,\alpha)=\dfrac{\sin(\theta-\alpha)+\sin\theta\cos\theta\cos(\theta-\alpha)}{\sin\alpha},y(\theta,\alpha)=\dfrac{\sin^2\theta\cos(\theta-\alpha)}{\sin\alpha}
\end{aligned}
$$
积分部分：
$$\begin{aligned}
x_\theta&=\csc\alpha[\cos(\theta-\alpha)+\cos2\theta\cos(\theta-\alpha)-\dfrac{1}{2}\sin2\theta\sin(\theta-\alpha)]
\\ \therefore
\text{S}&=\int_\alpha^\pi y(\theta,\alpha)\cdot x_\theta\mathrm{d}\theta
\\
\Rightarrow \text{S}&=\csc^2\alpha\int_\alpha^\pi\sin^2\theta\cos(\theta-\alpha)[(1+\cos2\theta)\cos(\theta-\alpha)-\dfrac{1}{2}\sin2\theta\sin(\theta-\alpha)]\mathrm{d}\theta
\\&=
\dfrac{\csc^2\alpha}{4}\int_\alpha^\pi 8\cos^2\theta\sin^2\theta\cos^2(\theta-\alpha)-\sin^2\theta\sin2\theta\sin(2\theta-2\alpha)\mathrm{d}\theta
\\&=
\dfrac{\csc^2\alpha}{4}\int_\alpha^\pi 2\sin^22\theta\cos^2(\theta-\alpha)-\sin^2\theta\sin2\theta\sin(2\theta-2\alpha)\mathrm{d}\theta
\\&=
\dfrac{\csc^2\alpha}{4}\int_\alpha^\pi \sin^22\theta+\sin^22\theta\cos(2\theta-2\alpha)-\sin^2\theta\sin2\theta\sin(2\theta-2\alpha)\mathrm{d}\theta
\\&=
\dfrac{\csc^2\alpha}{8}\int_\alpha^\pi 1-\cos4\theta+
(1-\cos4\theta)\cos(2\theta-2\alpha)-(1-\cos2\theta)\sin2\theta\sin(2\theta-2\alpha)\mathrm{d}\theta
\\&=
\dfrac{\csc^2\alpha}{16}\int_\alpha^\pi 2-2\cos4\theta+
2\cos(2\theta-2\alpha)-2\cos4\theta\cos(2\theta-2\alpha)
-(1-\cos2\theta)(\cos(-2\alpha)-\cos(4\theta-2\alpha))
\mathrm{d}\theta
\\&=
\dfrac{\csc^2\alpha}{16}\int_\alpha^\pi 2-2\cos4\theta+
2\cos(2\theta-2\alpha)-\cos(6\theta-2\alpha)-\cos (2\alpha+2\theta)
-\cos2\alpha+\cos(4\theta-2\alpha)
+\cos2\theta\cos2\alpha-\cos2\theta\cos(4\theta-2\alpha)
\mathrm{d}\theta
\\&=
\dfrac{\csc^2\alpha}{32}\int_\alpha^\pi 4-4\cos4\theta+
4\cos(2\theta-2\alpha)-2\cos(6\theta-2\alpha)-2\cos (2\alpha+2\theta)
-2\cos2\alpha+2\cos(4\theta-2\alpha)
+2\cos2\alpha\cos2\theta-\cos(6\theta-2\alpha)-\cos(2\theta-2\alpha)
\mathrm{d}\theta
\\&=
\dfrac{\csc^2\alpha}{32}\int_\alpha^\pi 4-4\cos4\theta+
4\cos(2\theta-2\alpha)-3\cos(6\theta-2\alpha)-\cos (2\alpha+2\theta)
-2\cos2\alpha+2\cos(4\theta-2\alpha)
\mathrm{d}\theta
\\&=
\dfrac{\csc^2\alpha}{64}\int_\alpha^\pi 8-8\cos4\theta+
8\cos(2\theta-2\alpha)-6\cos(6\theta-2\alpha)-2\cos (2\alpha+2\theta)
-4\cos2\alpha+4\cos(4\theta-2\alpha)
\mathrm{d}\theta
\\
\Rightarrow \text{S}&=\left. \dfrac{\csc^2\alpha}{64}(8\theta-2\sin4\theta+4\sin(2\theta-2\alpha)-\sin(6\theta-2\alpha)-\sin(2\theta+2\alpha)
-4\cos2\alpha\cdot\theta+\sin(4\theta-2\alpha))\right|_{\alpha}^{\pi}
\end{aligned} $$

注记：
1. 对于不包含交点的面积，需要再加上 $\dfrac{1}{2}\sin\alpha\cos\alpha$。
2. 对于含有交点的部分，只需换 $\alpha$ 为二分套二分求得的 $\theta_i$ 即可。

$$-\tan\theta=\dfrac{\mathrm{d}y}{\mathrm{d}x}=\dfrac{\mathrm{d}y(\theta,\alpha)}{\mathrm{d}x(\theta,\alpha)}$$

还观察到当两向量给定时，其夹角 $\alpha$ 也便随之确定。所以此时的斜率可认为只与 $\theta$ 有关。那么：

$$-\tan\theta=\dfrac{y_\theta}{x_\theta}\Rightarrow y_\theta=-\tan\theta \cdot x_\theta$$

其中 $y_\theta,x_\theta$ 分别表示 $y(\theta,\alpha),x(\theta,\alpha)$ 对 $\theta$ 求偏导。联立直线方程与斜率导出式有如下方程组。

$$
\left\{\begin{matrix}
\dfrac{\sin\alpha \cdot x(\theta,\alpha)}{\sin(\theta+\alpha)}+\dfrac{\sin\alpha \cdot y(\theta,\alpha)}{\tan\theta\sin(\theta+\alpha)}=1
\\ \\ y_\theta=-\tan\theta \cdot x_\theta
\end{matrix}\right.
$$

对第一式化简整理，等式两边同时对 $\theta$ 求偏导，与第二式的条件一起，即可求出参数方程。具体步骤如下。
$$
\begin{aligned}
&\tan\theta \cdot x(\theta,\alpha)+y(\theta,\alpha)=\dfrac{\tan\theta\sin(\theta+\alpha)}{\sin\alpha}
\\&\Rightarrow \dfrac{x(\theta,\alpha)}{\cos^2\theta}+\tan\theta \cdot x_\theta+y_\theta=\frac{\tan\theta\cos(\theta+\alpha)}{\sin\alpha}+\dfrac{\sin(\theta+\alpha)}{\sin\alpha\cos^2\theta}
\\&\Rightarrow x(\theta,\alpha)+\overbrace{\tan\theta \cdot x_\theta+y_\theta}^0=\dfrac{\sin(\theta+\alpha)+\sin\theta \cos\theta \cos(\theta+\alpha)}{\sin\alpha}
\\&\Rightarrow y(\theta,\alpha)=\dfrac{\tan\theta\sin(\theta+\alpha)}{\sin\alpha}-\dfrac{\tan\theta\sin(\theta+\alpha)+\sin^2\theta\cos(\theta+\alpha)}{\sin\alpha}
\\
\\
& \therefore
\left\{\begin{matrix}
\begin{aligned}
&x(\theta,\alpha)=\dfrac{\sin2\theta\cdot \cos(\theta+\alpha)+2\sin(\theta+\alpha)}{2\sin\alpha}
\\&y(\theta,\alpha)= -\dfrac{\sin^2\theta\cdot \cos(\theta+\alpha)}{\sin\alpha}
\end{aligned}
\end{matrix}\right.
\end{aligned}
$$
然后重点来了，我们要求的是面积，而得到的只是参数方程。接下来怎么办？其实也还好，我们带进去硬解就完了。要求面积，就必须先知道参数的范围。我们可以由图上的关系，轻松得出 $\theta \in [0,\pi-\alpha]$。既然上下限已知，那我们可以动手了。对于参数方程，求面积一般有以下形式。

$$
\text{S}=\int y(t_1,t_2,\ldots,t_n)\mathrm{d}[x(t_1,t_2,\ldots,t_n)]
$$
其中 $t_1,t_2,\ldots,t_n$ 代表参数。当只有一个参数变化，其它参数保持不变时，即可使用链式法则积分之。对本道题而言，有下面的形式。
$$
\text{S}=\int y(\theta,\alpha)\mathrm{d}[x(\theta,\alpha)]=\int y(\theta,\alpha)\frac{\partial x(\theta,\alpha)}{\partial\theta}\mathrm{d}\theta
$$
对于合乎几何逻辑的上下限，图上已经很明确了，因此可以列得下式。
$$
\text{S}=\int^{0}_{\pi-\alpha} y(\theta,\alpha)\frac{\partial x(\theta,\alpha)}{\partial\theta}\mathrm{d}\theta
$$
上式才能保证从左往右积。抽时间，或者直接找 WolframAlpha 计算整理。具体说，先输入 $x$ 的参数方程求解其关于 $\theta$ 的偏微分。得到答案后，再计算其与 $y$ 参数方程的乘积，最后再算二者乘积的不定积分即可。但总之，经由暴力计算，可以有下式。
$$
\begin{aligned}
\text{S}&=\left. \dfrac{1}{64\sin^2\alpha}
[
(8-4\cos2\alpha)\theta-2\sin4\theta+4\sin(2\theta+2\alpha)
-\sin(6\theta+2\alpha)
+\sin(2\alpha-2\theta)
+\sin(4\theta+2\alpha)]
\right|_{0}^{\pi-\alpha}
\end{aligned}
$$
同时，作为第二次提交失败而得到的反馈，这里我将给出我个人的证明。
$$
\begin{aligned}
\text{S}&=\int_{\pi-\alpha}^{0}y(\theta,\alpha)\mathrm{d}[x(\theta,\alpha)]
\\ &=\int_{0}^{\pi-\alpha}\dfrac{\sin^2\theta\cdot \cos(\theta+\alpha)}{\sin\alpha}\mathrm{d}[\dfrac{\sin2\theta\cdot \cos(\theta+\alpha)+2\sin(\theta+\alpha)}{2\sin\alpha}]
\\ &=\dfrac{1}{2\sin^2\alpha}\int_{0}^{\pi-\alpha}
\sin^2\theta\cdot \cos[\theta+\alpha](2\cos2\theta\cdot\cos(\theta+\alpha)-\sin2\theta\cdot\sin(\theta+\alpha)+2\cos(\theta+\alpha))\mathrm{d\theta}
\\ &=\dfrac{1}{2\sin^2\alpha}
\int_{0}^{\pi-\alpha}
\sin^2\theta [2\cos^2(\theta+\alpha)(1+\cos2\theta)
-\sin2\theta\sin(\theta+\alpha)\cos(\theta+\alpha)]
\mathrm{d\theta}
\\ &=\dfrac{1}{4\sin^2\alpha}
\int_{0}^{\pi-\alpha}
(1-\cos2\theta) [(1+\cos(2\theta+2\alpha))(1+\cos2\theta)
-\sin2\theta\sin(\theta+\alpha)\cos(\theta+\alpha)]
\mathrm{d\theta}
\\ &=\dfrac{1}{4\sin^2\alpha}
\int_{0}^{\pi-\alpha}
(1+\cos(2\theta+2\alpha))(1-\cos^22\theta)
-\sin2\theta\sin(\theta+\alpha)\cos(\theta+\alpha)
+\sin2\theta\cos2\theta\sin(\theta+\alpha)\cos(\theta+\alpha)
\mathrm{d\theta}
\\ &=\dfrac{1}{16\sin^2\alpha}
\int_{0}^{\pi-\alpha}
2(1+\cos(2\theta+2\alpha))(1-\cos4\theta)
-2\sin2\theta\sin(2\theta+2\alpha)
+\sin4\theta\sin(2\theta+2\alpha)
\mathrm{d\theta}
\\ &=\dfrac{1}{16\sin^2\alpha}
\int_{0}^{\pi-\alpha}
2-2\cos4\theta+2\cos(2\theta+2\alpha)-2\cos(2\theta+2\alpha)\cos4\theta
-\cos2\alpha+\cos(4\theta+2\alpha)
+\dfrac{1}{2}\cos(2\alpha-2\theta)-\dfrac{1}{2}\cos(2\alpha+6\theta)
\mathrm{d\theta}
\\ &=\dfrac{1}{32\sin^2\alpha}
\int_{0}^{\pi-\alpha}
4-4\cos4\theta+4\cos(2\theta+2\alpha)
-2\cos(6\theta+2\alpha)-2\cos(2\alpha-2\theta)
-2\cos2\alpha+2\cos(4\theta+2\alpha)
+\cos(2\alpha-2\theta)-\cos(2\alpha+6\theta)
\mathrm{d\theta}
\\ &=\dfrac{1}{32\sin^2\alpha}
\int_{0}^{\pi-\alpha}
4-4\cos4\theta+4\cos(2\theta+2\alpha)
-3\cos(6\theta+2\alpha)-\cos(2\alpha-2\theta)
-2\cos2\alpha+2\cos(4\theta+2\alpha)
\mathrm{d\theta}
\\ &=\left. \dfrac{1}{32\sin^2\alpha}
[
4\theta-\sin4\theta+2\sin(2\theta+2\alpha)-\dfrac{1}{2}\sin(6\theta+2\alpha)+\dfrac{1}{2}\sin(2\alpha-2\theta)
-2\cos2\alpha\cdot\theta+\dfrac{1}{2}\sin(4\theta+2\alpha)]
\right|_
{0}^{\pi-\alpha}
\\ &=\left. \dfrac{1}{64\sin^2\alpha}
[
(8-4\cos2\alpha)\theta-2\sin4\theta+4\sin(2\theta+2\alpha)
-\sin(6\theta+2\alpha)
+\sin(2\alpha-2\theta)
+\sin(4\theta+2\alpha)]
\right|_{0}^{\pi-\alpha}
\end{aligned}
$$

到这一步，兴许有人会问：“为什么不直接写结果呢，直接算到最后不行吗？”对此，有两点原因。一是后面也会用到这个式子，直接写成代入式会很方便。二是计算几何中喜闻乐见的精度问题，直接换成结果有时是不会 AC 的。

总之，它代表的是如图所示的面积。而为了计算的正确性，我们需要减掉直线下方那块多余的三角形面积。

由 $\sin(\pi-\alpha)=\sin\alpha,\cos(\pi-\alpha)=-\cos\alpha$ 恒等式得到三角形的面积计算式 $\text{S}_\triangle =-\dfrac{1}{2}\sin\alpha\cos\alpha$。所以，我们减去的时候就相当于加上 $\dfrac{1}{2}\sin\alpha\cos\alpha$了。

到这一步，按照上述逻辑编出来的正确代码已经可以通过样例一了。但第二个问题随之而来：为什么第二个样例跑出来的结果和给的答案对不上？

那是因为我们还没有处理重叠的部分。容易验证，当给出所有点的坐标分量都是长整数时，重叠面积只会在某边边长为 $1$ 或 $\sqrt{2}$ 时出现。这里由以下三张图辅助感性理解。

边长为 $1$ 时。

边长为 $\sqrt{2}$ 时。

边长为 $\color{red}{2}$ 时。而往后由于长度的增加，两条线一定不会相交。

而对上述结论做较为理性严谨的定性分析，可以按以下方式来看。

毕竟这道题是通过杆滑动的过程受到启发而建构起来的。我们又知道杆最终是会停在 $x$ 轴上的。加上杆的长度是 $1$，如果边的长度不是 $2$ 又严格小于 $2$，最终两根杆停下来，一定是会交叠在一起的。但从现实来看，显然不可能有一个叠着一个的情况，运动途中两杆必然发生了碰撞。那么对应回本道题，当边长严格小于 $2$ 时，便一定有相交的点。

回到相交的讨论，这里假设绿线上的坐标为 $\big(x(\theta_1,\alpha_1),y(\theta_1,\alpha_1)\big)$，变换前橙线上的坐标为 $\big(x(\theta_2,\alpha_2),y(\theta_2,\alpha_2)\big)$。第一个问题是：变换后橙线的参数方程怎么求，求完了以后怎么求出参数？

根据高中所学知识，橙线上点的坐标理应为 $\big(\text{Len}-x(\theta_2,\alpha_2),y(\theta_2,\alpha_2)\big)$。其中 $\text{Len}$ 代表凸多边形上两邻近顶点的边长。对于这一点的证明，具体来说有如下操作。

右侧的确定角度的直线可以理解为左侧的直线先平移，再关于 $y$ 轴做对称变换得到。那么对应地，曲线也是执行了上述的操作（也即平移到 $x=x_0$ 的位置后再关于 $x=x_0$ 对称）。只不过，这样容易忽略 $x$ 轴上平移的“左加右减”含义是在不加负号的前提下成立的事实，从而容易出错。

以线性代数的角度理解，上述的操作等价于对二维向量 $\boldsymbol{a}=(x,y)^T$ 施行变换 $f(\boldsymbol{a})=\boldsymbol{A}(\boldsymbol{a}+\boldsymbol{b})$，其中 $\boldsymbol{A}=\begin{pmatrix}
-1 &0 \\
0 &1
\end{pmatrix},\boldsymbol{b}=(-\operatorname{Len},0)^T$。所以有结论$\big(\text{Len}-x(\theta_2,\alpha_2),y(\theta_2,\alpha_2)\big)$。

为了求出交点对应的两个参数 $\theta_1,\theta_2$，我们发扬人类智慧，根据并利用两条曲线的单调性，从两条曲线纵坐标值的初始上界 $\sin[\min(\alpha_1,\alpha_2)]$ 与下界 $0$ 开始，以上下界所确定的中心线与曲线的交点为基准，逐步二分到与交点坐标的误差不超过 $10^{-13}$ 的位置；同时在二分缩小纵坐标误差的过程中，每次对属于 $[0,\pi-\alpha]$ 区间内的 $\theta_u$ 二分，以求得当前纵坐标为 $y_u$ 时两条曲线的对应参数。
所以，当发现边长不够 $2$ 及以上时，可以用前文得到的积分式，把积分上界换成两条曲线经由二分求出交点所对应的参数 $\theta_{u,1},\theta_{u,2}$，再减掉即可。
$$
\begin{aligned}
\text{S}&=\left. \dfrac{1}{64\sin^2\alpha}
[
(8-4\cos2\alpha)\theta-2\sin4\theta+4\sin(2\theta+2\alpha)
-\sin(6\theta+2\alpha)
+\sin(2\alpha-2\theta)
+\sin(4\theta+2\alpha)]
\right|_{0}^{\theta_{u,i}}
\end{aligned}
$$
这个时候你可能会问：“如果图形是下面这样的，我们求交点然后减掉面积，这样不会多减吗？”

对于这个担忧，由于曲线上点的坐标是由 $\theta,\alpha$ 所确定的，而对于同一个角 $\alpha_k$，若在上一次计算参数时是通过形如 $\text{Len}-x(\theta_k,\alpha_k)$ 的方式得到的话，在下一次计算时就必须得通过 $x(\theta_k,\alpha_k)$ 的方式得到，才保证不会重复。而程序按循环执行下去，只要测试数据没有锅、取模写得对，这是一定不会重复计算的。

到这里，我们仍需要由第一个样例的暗示注意到，某些矩形的宽为 $1$。这使得整个图形的面积都是标红的。好在题目的条件 $90^{\circ} \leq \alpha<180^{\circ}$ 限制了更为复杂的情况出现，所以我们只需为顶点数为 $4$ 的图形做一次特判就行了。

最后，还有两个得加以关注的点是，由于 $-10^{9} \leq x,y\leq10^{9}$, 当两个点的坐标分量 $x,y$ 都差得非常大时，用 long double 算向量模长，精度是不够的。得用 long long 才能跑过最毒瘤的那几个大样例。而这个，是我前前后后花了一个月时间才弄清楚的。

第二个是夹角的计算。我们反应最快的公式是 $\alpha=\arccos\dfrac{\boldsymbol{a}\cdot \boldsymbol{b}}{|\boldsymbol{a}|\cdot|\boldsymbol{b}|}$，但若是追求答案的精度，这也是不够用的。为了尽量降低浮点误差的累积，以 C++ 编写时，推荐使用 atan2(y,x) 函数。相应地，这会使得叉乘的应用成为必然，因为确定正切不只是需要余弦值，还需要正弦值。

代码实现

全都被 H4ck 掉了。想跑对自己想办法更正吧。

C++ 代码。

#include <iostream>
#include <cstdio>
#include <cstring>
#include <algorithm>
#include <cmath>
using namespace std;
using ll = long long;
using db = long double;
const int N = 5001;
const db eps = 1e-6, pi = acos((db)-1); 

class vec{
public:
    ll x, y;
    vec(ll xx = 0, ll yy = 0) : x(xx), y(yy) {}
    vec(const vec &a) : x(a.x), y(a.y) {}
    vec operator-(const vec &a) const { return vec(x - a.x, y - a.y); }
    db len() { return (db)sqrt(x * x + y * y); }
    vec rev() { return vec(-x, -y); }
} dot[N], edge[N];

inline db angle(vec a, vec b) {return acos(a.x * b.x / a.len() / b.len() 
                                         + a.y * b.y / a.len() / b.len());} 
                                      // atan2(a.x * b.y - a.y * b.x, a.x * b.x + a.y * b.y)

inline db interg(db theta, db alpha){ // 上下界代入式
    return ((8.0 - 4.0 * cos(2.0 * alpha)) * theta - 2.0 * sin(4 * theta) + 4.0 * sin(2.0 * (theta + alpha)) 
          - sin(2.0 * (8.0 * theta + alpha)) + sin(2.0 * (alpha - theta)) + sin(2.0 * (2.0 * theta + 3.0 * alpha)))
          / sin(alpha) / sin(alpha) / 64.0;
}

inline db overlap(db mida, db midb, db ang1, db ang2){return interg(midb, ang1) + interg(mida, ang2) 
                                                           - interg(0, ang1) - interg(0, ang2);}

void calc(db &h, db &l, db ang, db Y){ // 求两条曲线交点的参数theta
    while (fabs(h - l) > eps){
        db mid = h / 2.0 + l / 2.0;
        db yy = -pow(sin(mid), 2) / sin(ang) * cos(mid + ang);
        if (yy < Y) h = mid;
        if (yy > Y) l = mid;
    }
    return;
}
db getarea(int _n, int n, int nt){
    auto now = edge[n], lst = edge[_n], nxt = edge[nt];
    db res = 0.0, len = now.len(),
       alpha1 = angle(now, lst.rev()), alpha2 = angle(nxt, now.rev());
    if (len < 2.0){
        db la = 0.0, lb = 0.0,
           ra = pi - alpha1, rb = pi - alpha2,
           cy = sin(min(alpha1, alpha2)), fy = 0, mida, midb;

        while (fabs(cy - fy) > eps){   // 二分套二分
            db midy = cy / 2.0 + fy / 2.0;
            la = 0.0, lb = 0.0, ra = pi - alpha1, rb = pi - alpha2;
            calc(ra, la, alpha1, midy), calc(rb, lb, alpha2, midy);

            mida = la / 2.0 + ra / 2.0, midb = lb / 2.0 + rb / 2.0;
            db x1 = (sin(2.0 * mida) * cos(alpha1 + mida) + 2 * sin(mida + alpha1)) / 2.0 / sin(alpha1),
               x2 = len - (sin(2.0 * midb) * cos(alpha2 + midb) + 2 * sin(midb + alpha2)) / 2.0 / sin(alpha2);

            if (x1 < x2) cy = midy;
            if (x1 > x2) fy = midy;
        }

        res -= overlap(mida, midb, alpha1, alpha2);
    }
    res += interg(pi - alpha1, alpha1) - interg(0, alpha1) + 0.25 * sin(alpha1) * cos(alpha1); 
                                     // 推荐积分求解，不推荐写成 3*(pi-alpha)/16 + 1/16/tan(alpha) + (pi-alpha)/16/tan(alpha)/tan(alpha)
    return res;
}
int main(){
    int n;
    scanf("%d", &n);
    for (int i = 0; i < n; ++i)
        scanf("%ld %ld", &dot[i].x, &dot[i].y);
    for (int i = 0; i < n; ++i)
        edge[i] = vec(dot[(i + 1) % n] - dot[i]);
    if (n == 4){                    // 特判4个点时的情况
        db len1 = edge[0].len(), len2 = edge[1].len(),
           len3 = edge[2].len(), len4 = edge[3].len();

        if (len1 == len3 and len2 == len4 and
            (len1 == 1 or len2 == 1) and
            fabs(angle(edge[0], edge[1]) - pi / 2.0) < eps)
            printf("%.11lf", len1 * len2), exit(0);
    }
    db res = 0.0;
    for (int i = 0; i < n; ++i)
        res += getarea((i - 1 + n) % n, i, (i + 1) % n);
    printf("%.11lf", res);
    return 0;
}

Python 3 的代码。

import math
pi = 3.14159265358979323846264338327950288419716939937510
eps, sq2 = 1e-6, pow(2, 0.5)
x, y, n = [], [], 0


def fabs(xx):
    return math.fabs(xx)


def Sin(xx):
    return math.sin(xx)


def Cos(xx):
    return math.cos(xx)


def Acos(xx):
    return math.acos(xx)


def get_len(ux, uy):
    return pow((pow(ux, 2) + pow(uy, 2)), 0.5)


def angel_(ux, uy, vx, vy, u, v):
    return Acos((ux * vx + uy * vy) / u / v)


def small(_a: float, _b: float) -> bool:
    return fabs(_a - _b) < eps


def big(_a: float, _b: float) -> bool:
    return fabs(_a - _b) > eps


def inter_grate(theta, alpha):
    ans = ((8 - 4 * Cos(2 * alpha)) * theta - 2 * Sin(4 * theta) + 4 * Sin(2 * theta + 2 * alpha)
           - Sin(6 * theta + 2 * alpha) + Sin(2 * alpha - 2 * theta)
           + Sin(4 * theta + 2 * alpha)) / (64 * pow(Sin(alpha), 2))
    return ans


def calc(la, ra, mid_y, alpha_1):
    while big(ra, la):
        mid_a = ra / 2.0 + la / 2.0
        yy = - pow(Sin(mid_a), 2) / Sin(alpha_1) * Cos(alpha_1 + mid_a)
        if yy > mid_y:
            la = mid_a
        if yy < mid_y:
            ra = mid_a
    return la, ra


def formulas_x(la, ra, alpha_1):
    return (2.0 * Sin(la / 2.0 + ra / 2.0 + alpha_1) + Sin(la + ra)
            * Cos(la / 2.0 + ra / 2.0 + alpha_1)) / (2 * Sin(alpha_1))


def binary_find(la, lb, ra, rb, cy, fy, alpha_1, alpha_2, ab):
    while big(cy, fy):
        mid_y = cy / 2.0 + fy / 2.0
        la, lb, ra, rb = 0.0, 0.0, pi - alpha_1, pi - alpha_2

        la, ra = calc(la, ra, mid_y, alpha_1)
        lb, rb = calc(lb, rb, mid_y, alpha_2)

        x1, x2 = formulas_x(la, ra, alpha_1), ab - formulas_x(lb, rb, alpha_2)

        if x1 < x2:
            cy = mid_y
        if x1 > x2:
            fy = mid_y
    return la, lb, ra, rb, cy, fy


def get_area(_i, ni, i_, i_2):
    ans = 0.00
    ux, uy, vx, vy = x[_i] - x[ni], y[_i] - y[ni], x[i_] - x[ni], y[i_] - y[ni]

    ab, ad = get_len(ux, vy), get_len(vx, uy)
    alpha_1 = angel_(ux, vy, vx, uy, ab, ad)

    if small(ab, sq2) or small(ab, 1.00):
        wx, wy = x[i_2] - x[i_], y[i_2] - y[i_]

        bc = get_len(wx, wy)
        alpha_2 = angel_(wx, wy, vx, vy, ab, bc)

        la, lb, ra, rb, cy, fy = 0.0, 0.0, pi - alpha_1, pi - alpha_2, Sin(min(alpha_1, alpha_2)), 0.0000
        la, lb, ra, rb, cy, fy = binary_find(la, lb, ra, rb, cy, fy, alpha_1, alpha_2, ab)

        ans -= inter_grate(la / 2.0 + ra / 2.0, alpha_1) - inter_grate(0, alpha_1)
        ans -= inter_grate(lb / 2.0 + rb / 2.0, alpha_2) - inter_grate(0, alpha_2)

    ans += inter_grate(pi - alpha_1, alpha_1) - inter_grate(0, alpha_1) + 0.5 * Sin(alpha_1) * Cos(alpha_1)
    return ans


if __name__ == "__main__":
    n = eval(input())
    for i in range(1, n + 1):
        a, b = input().split()
        a, b = int(eval(a)), int(eval(b))
        x.append(a), y.append(b)
    if n == 4:
        Ax, Ay, Bx, By, Cx, Cy, Dx, Dy = x[0], y[0], x[1], y[1], x[2], y[2], x[3], y[3]
        ABx, ABy, ADx, ADy = Bx - Ax, By - Ay, Dx - Ax, Dy - Ay
        CBx, CBy, CDx, CDy = Bx - Cx, By - Cy, Dx - Cx, Dy - Cy

        LenAB, LenAD, LenCB, LenCD = get_len(ABx, ABy), get_len(ADx, ADy), get_len(CBx, CBy), get_len(CDx, CDy)

        a, b = angel_(ABx, ABy, ADx, ADy, LenAB, LenAD), angel_(ABx, ABy, CBx, CBy, LenAB, LenCB)
        c, d = angel_(CDx, CDy, CBx, CBy, LenCD, LenCB), angel_(CDx, CDy, ADx, ADy, LenCD, LenAD)

        if small(a, pi / 2.0) and small(b, pi / 2.0) and small(c, pi / 2.0) and small(d, pi / 2.0) and \
         ((small(LenAB, 1.00) and small(LenAB, LenCD)) or (small(LenCB, 1.00) and small(LenCB, LenAD))):
            print('%.11Lf' % (LenAB * LenCB)), exit(0)

    res = 0.0000
    for i in range(1, n + 1):
        res += get_area((i - 1 + n) % n, i % n, (i + 1) % n, (i + 2) % n)

    print('%.11Lf' % res)

问题拓展

如果完全抛开凸多边形的鬼畜计算和直角钝角的性质不谈，只考虑单个锐角不取到 $0^\circ,180^\circ$ 以外的情况，又会是个什么情况？

首先可以想象出下面的样子。

可观察到杆 $AC$ 位于 $x$ 轴的端点一直从点 $C$ 运动到点 $G$ 后收回到 $F$。而另一个端点则是先上升到点 $D$ 后一直沿着夹角 $\alpha$ 的一边往 $C$ 走。

同时也不难想象围成图形的面积由两部分构成，即三角形面积加上包络线的积分值。也即 $\dfrac{1}{2\tan\alpha}+\displaystyle \int y\mathrm{d}x$。

而包络线的计算方式可以继续沿用之前的思想。

假如让 $\angle \text{DEC}$ 作为 $\theta$，那么包络线上的点 $(x,y)$ 必定经过直线 $l:y=-\tan\theta(x-x_0)$，且存在关系

$$\begin{aligned}&\dfrac{x_0}{\sin(\pi-\alpha-\theta)} = \dfrac{1}{\sin\alpha}\Rightarrow x_0 = \dfrac{\sin(\alpha+\theta)}{\sin\alpha} \\
&\dfrac{\mathrm d y}{\mathrm d x} = \dfrac{y_\theta}{x_\theta} = -\tan\theta\Rightarrow y_\theta = -\tan\theta x_\theta\end{aligned}$$

变形整理得到
$$\begin{cases}x = \dfrac{2\sin(\alpha + \theta) + \sin2\theta\cos(\alpha + \theta)}{2\sin\alpha} \\ y = -\dfrac{\sin^2\theta\cos(\alpha + \theta)}{2\sin\alpha} \end{cases}$$

结果是完全一致的。我们再来考虑积分的上下限。

由图片不难看到，最开始的下限是 $\dfrac{\pi}{2}$，上限是 $\dfrac{\pi}{2}-\alpha$。在经过上文的积分后即可，无需减掉三角形面积，保留即可得到正确答案。

所以，对于单个锐角而言，其覆盖的面积为

$$\begin{aligned}\text{S}(\alpha) &= \dfrac{1}{2\tan\alpha} + \Big(\left. \dfrac{1}{64\sin^2\alpha}
[
(8-4\cos2\alpha)\theta-2\sin4\theta+4\sin(2\theta+2\alpha)
-\sin(6\theta+2\alpha)
+\sin(2\alpha-2\theta)
+\sin(4\theta+2\alpha)]
\right|_{\frac{\pi}{2} - \alpha}^{\frac{\pi}{2}}\Big).
\end{aligned}$$

参考网页与在线工具

1. https://bilibili.com/video/av370957522

2. https://www.bilibili.com/video/BV13s4y1X7yb

3. https://codeforces.com/blog/entry/106675

4. https://www.wolframalpha.com/

5. https://www.latexlive.com/

6. https://www.geogebra.org/calculator

最后，欢迎纠错。

#include <iostream>
#include <cstdio>
#include <cstring>
#include <algorithm>
#include <cmath>

using namespace std;
const int N = 5E4 + 10,
          M = 2E5 + 10,
          S = 1E6 + 10;
int n, m, len = 1,
    a[N], res[M];

int cnt[S];                         // 用于计数用

inline int get(int i) {return i / len;}
inline void add(int x, int &ans) {if(++cnt[x] == 1) ans ++;}
inline void sub(int x, int &ans) {if(!--cnt[x])     ans --;}

struct node {int id, l , r;}q[M];


int main() {

    scanf("%d", &n);
    for(int _ = 1; _ <= n; ++_)
        scanf("%d", a + _);
    scanf("%d", &m);

    len = max(1, (int)sqrtf((double)n * n / m));
    int l, r;
    for(int i = 0; i < m; ++i)
    {
        scanf("%d %d", &l, &r);
        q[i] = {i, l, r};
    }

    sort(q, q + m, [&](node &a, node &b) {
        int al = get(a.l), bl = get(b.l);
        if(al ^ bl) return al < bl;
        return a.r > b.r;
    });                             // 横跨整段的优先

    for(int k = 0, i = 1, j = 0, ans = 0; k < m; ++k)
    {
        int id = q[k].id, l = q[k].l, r = q[k].r;

        while(j < r) add(a[++j], ans);
        while(j > r) sub(a[j--], ans);
        while(i < l) sub(a[i++], ans);
        while(i > l) add(a[--i], ans);

        res[id] = ans;
    }
    for(int _ = 0; _ < m; ++_)
        printf("%d\n", res[_]);
    return 0;
}

好久没能写算法题了。

今天下午五点，老师布置给我的这个实验就结束了。

事实上，这个汇编实验很有意思，需要 IDA 和 gdb 配合着来做才会做得比较顺畅。单纯 gdb 看汇编代码实在有点辛苦了。但是今天看见隔壁班学委完完全全拿 gdb 肝过去了，不妨在此 %%% 一下他 。

里头有 $7$ 个 phase（一个隐藏，其余的正常显示），而学校专门为班里的同学生成了不同的 $32$ 位和 $64$ 位程序用于读汇编代码。下面我只将我做了 $64$ 位的部分讲解清楚。

打开 IDA 导入 $64$ 位的 bomb，然后忽略其它部分，按 F5 主动查看文件的反汇编伪代码（注意，有时看了不如直接看汇编代码，但作为经验先看能简化），然后左键双击——

phase_1：字符串比较 难度：简单

只要点进去后这题就结束了。

所以答案是 I was trying to give Tina Fey more material.。做完后回退到 main 函数中。

对于 gdb，在 bomb 所在目录下打开终端，使用下列指令查看程序的汇编代码。

gdb -q bomb # -q 表示忽略载入信息。
disassem phase_1 

然后需要按照 gdb 给出的注释来检查内存情况。这样子也能得到答案。

后面找数组、节点、二叉树节点等内容也可以按照上面的方式做。检查的参数视情况而定，感兴趣的读者还请自己参阅 gdb 检查内存内容的一节。

phase_2：等比数列 难度：还好

在汇编的世界中，指针或变量极其灵活。

这题直接看汇编代码会有障碍，但是 IDA 就没问题，基本也是一眼搞完。

首先看到 v5 = __readfsqword(0x28u); 和最后的 return v5 - __readfsqword(0x28u); 这个只是反栈溢出的函数，静态分析可以暂时不用管它。因为程序最重要的逻辑都与 v3 数组有关。

可以很明显地看到程序调用了一个『读六个数』的函数。其中传入的参数是 v3 数组和传入的字符串指针 a1。跟进审阅可以得到

在计算机中，一个 int 占 $4$ 个字节。所以上面的偏移不过是需要依次输够六个数到 v3 数组之中（或者更多，但前提是不影响下一题的输入），不然就爆炸罢了。

然后 phase_2 中要求 v3[0] = 1，不然爆炸拆弹失败。其次可以看见 v1 取了 v3 的地址之后就进了一个 do{...}while(v1!=&v4); 循环。

将指针的表示转换一下，就有要求 v1[1] == 2 * v1[0]。也就是说，后一个数总是前一个数的 $2$ 倍。

而边界则是通过 phase_2 中的变量声明来看的。

int  v3[5]; // rsp+0h = rsp
char v4; // rsp+14h = rsp + 20

此时我们要关心的是注释而非变量本身。可以看到从栈指针所指的 v3 开始，到向 上 偏移 $20$ 个字节后指到 v4 的地址，这段区域实际上就是读六个数之后所存放的位置。

到此，不难看出 while(v1!=&v4); 的意义以及最终的答案 1 2 4 8 16 32。

phase_3：神秘的分支选择 难度：中等

可以明确，程序要求我们输入两个数，第一个数要小于等于 $7$。然后发现往下 IDA 就解析不出来了。汇编代码的意义是跳到寄存器所指的地址，然而问题是我们不知道 rax 里头的值。

通常来讲一般是开了优化来优化 switch 才会出现这种状况的。恢复跳转表可以参考 这篇 文章

所以这个时候最好在 gdb 里动态调试了。首先使用到断点，在确保 bomb.c 和 bomb 文件位于同一目录后，终端里进入 gdb 后敲入 l 然后连续回车到对应的位置。

进而可以确定断点的位置是在第 $89$ 行。于是可以在终端中打入

b 89
r  # 键入前两题的正确答案。 
disassem phase_3  # 等执行到断点的位置反汇编 phase_3。
# 然后根据反汇编情况多次打入 si 单次执行汇编语句。

当然也可以一步到位，直接 b phase_3。

具体而言

最后调试得到

所以一个有效答案是 1 -582。

phase_4：斐波那契前缀和 难度：还好

这个时候 IDA 就简单了。

不难将 func4 转写为

def func4(a1, a2):
    if not a1:
        return 0
    elif a1 == 1:
        return a2
    return func(a1-1,a2)+func(a1-2,a2)+a2;

那么

$$\begin{aligned}
\operatorname{func}(6,a_2)&=\operatorname{func}(5,a_2)+\operatorname{func}(4,a_2)+a_2
\\&=2\operatorname{func}(4,a_2)+\operatorname{func}(3,a_2)+2a_2
\\&=3\operatorname{func}(3,a_2)+2\operatorname{func}(2,a_2)+4a_2
\\&=5\operatorname{func}(2,a_2)+3\operatorname{func}(1,a_2)+7a_2
\\&=8\operatorname{func}(1,a_2)+5\operatorname{func}(0,a_2)+12a_2
\\&=20a_2\end{aligned}
$$
所以只需要输一个解：60 3 即可。

回过头来看，拿斐波那契数列的前缀和数列来看就有

$\begin{aligned}&1,1,2,3,5,8,\ldots\\&1,2,4,7,12,20\ldots\end{aligned}$

phase_5： 模拟 难度：还好

老样子，继续用 IDA 静态分析。

跟进分析 array_0 数组，有

通过 array_0 数组可以 pwn 掉密码。从 IDA 中提取数据的过程是：选中数据块 -> 选择菜单栏上的『编辑』-> 『导出数据』->如图所示的选项。然后复制到代码编辑器里就好。

// windows 下跑的程序
#include<cstring>
#include<cstdio>
#include<iostream>

using namespace std;

int array_0[16] ={ 10, 2, 14, 7, 8, 12, 15, 11, 0, 4, 1, 13, 3, 9, 6, 5 };

int main()
{
    int v2 = 0, i = 0, 
        v1 = 0, v3 = 0;
    for(;i < 16; ++i)
    {
        v2 = v3 = 0;
        v1 = i;
        do
        {
            ++ v3;
            v1 = array_0[v1];
            v2 += v1;
        } while (v1 != 15);
        printf("I = %d\t v2 = %d\t v3 = %d\n",i, v2, v3);
    }
    system("pause");
    return 0;
}

答案就是下面显示的 5 115。

phase_6： 链表数据段排序 难度：稍难

一开始我看了头都有点昏。反汇编代码里充满了指针。但是认真思考，又没那么难。

跟进分析 node1 得到暂且不成看的内容。

在 IDA 中为每个数据选中起始位置后，使用两次 ‘D’ 键使单字节的 db 合并为四字节 dd，从而有规整的形式。

然后只需要读出数据位数据的大小关系就可以了。我得到的关系是 $[2] > [6] > [1] > [4] > [5] > [3]$，再拿 $7$ 去减有答案 5 1 6 3 2 4。

secret_phase： 二叉树 + 递归 难度：一般难

从逻辑来说，应该是用 gdb 读源代码时会在最后面看到有关 secret_phase 的暗示。然而我自己使用了 IDA 的缘故，很早的时候就在 IDA 的函数栏中看到 secret_phase 了。跟进分析，有

再转到 fun7 和 n1 中查看情况。

然后就不用多说了。现在第二个问题是，如何找到调用 secret_phase 函数的函数？

很简单，右击函数名，选中『跳转到交叉引用』然后点击『ok』。

然后就发现只在 phase_defused 中调用了它。要求是调用时 num_input_strings == 6 而且先输两个随便的数之后再输入 DrEvil。

但是如果只是猜在解第六题的环节这么做，事实上经过验证是调用不到的。而自己最近也是因为事情太多没能深挖什么时候 num_input_strings == 6。但是据经验来看，应该是前面多打几轮 DrEvil 就好了。而这样，整个实验就搞完了。

最终结果：

前言

三个流程，让班上三十几个学生集体 “坐牢” 一个半小时。

学校开设的这门课，给我的感觉是：

• 浪费青春；
• 很有用，但很没用，下次不要再来了；
• 混沌。实验步骤等有用信息还得自己思考或筛选得到，理解成本很高。

大体流程

1. 用示波器抓取 STC 单片计算机板发送的 UART 信号以分析该信号的 波特率；
2. 根据识别与读取到的结果更改 Linux 程序，只编译执行一个 main 程序后（操作上面完全没有写），在终端中读取串口收到的数据；
3. 中止程序，继续在终端中输入形如 curl "ip.ip.ip.ip:port-number/checkBaud?id=iiiiiiiiiiii&v=112233445566778899AABB" 的指令提交读取到的数据至测评服务器，完成实验并得到指导老师的确认验收。其中：
   • ip.ip.ip.ip:port-number 表示服务器的 IP 地址与接收的端口号（出于保护需要，此处不提供实际的形式）；
   • iiiiiiiiiiii 是需要换成的学号；
   • 112233445566778899AABB 需要换成读取到的序列号。

准备工具

两人小组需要取 两块不同的 STC 单片计算机板、两台笔记本电脑【一台最好直接装有 Linux 系统（我自己是直接把电脑配成 win10 + kali 双系统的。当时实验被虚拟机坑到吐血，灵机一动，想到直接用自己的 kali 直连单片计算机板试一试，然后实验就好了）或者安装 VMWare 的 ubuntu 虚拟机，按照操作文档来 除了恶心也没别的不好。另一台笔记本只需在 windows 系统中运行即可】、示波器、数据线与探针。

关于如何在一台电脑上装配双系统，大致流程是先下载 kali 的镜像文件然后用 rufus 配自己的一个 U 盘成启动盘，在 win10 文件资源管理器中的左侧栏中 右击 『此电脑』，选中『管理』后再选中『磁盘管理』，为 kali 分配足够大小的磁盘空间。

然后重启进入 BIOS 系统勾销安全启动项后，于再次重启的过程中插入启动盘，进入引导安装界面配置 kali 系统。

这里再提供一个避坑方案：选择语言的时候不要提前选中文，优先选英文。

一是方便后面配置回中文，二是用久了中不中文其实无所谓，用英文还能锻炼英文能力。

上述过程会遇到很多技术细节，在此不做过多阐述，还请自己查阅相关资料。充分准备后再谨慎地安装。

具体步骤

烧录与测量

首先是波特率。所谓波特率，就是 $1\ \operatorname{s}$ 内传输的 bits 总数。更具体地说，是在配置好示波器属性（这里我个人隐约记得是触发方式选用下降沿触发，耦合设为直流）测量工具连接好单片计算机板后，通过测量单个脉冲 1 的时间间隔后拿 $1\ \text{s}$ 去除即可。

也即下图的两个纵向光标 a,b 的位置与时间间隔 $560 \ \mu\operatorname{s}$ 所示。此时计算得到 $\dfrac{1}{560}\times10^6\approx 1786 \approx 1800$。所以之后的程序里需要改参数为 1800（这个是照着提供的表改的，后面会更具体地说）。

但在此之前，单片计算机板上啥也没有，于是就需要烧录。学校提供的烧录软件是 stc -isp-v6.88F.exe，配合驱动程序 ch341ser.exe 以及烧录文件 dut7. hex 和单片计算机板来使用。

在一台笔记本的 windows 系统下下载好三者，先后安装好 ch341ser.exe 与 stc -isp-v6.88F.exe 后，通过数据线连接单片计算机板并运行 stc -isp-v6.88F.exe。当 正常连接 时，此程序 会显示 USB-SERAL CH340 和对应串口号，否则会是其它项，遇到这种情况最好检查一下连接情况。

然后就是使用程序中的『打开程序文件』来导入待烧录的数据。

总之，安好并导入完 hex 文件后，界面会是下面这样的（这里偷懒直接用学习通上文档的示意图）。

接着点『下载/编程』就可以烧录了。中途会中断是因为我们需要自己按单片计算机上的复位按钮才能擦除原有数据。具体位置见下。

待显示出 “操作成功 ！” 时，这块板子就烧录好了。调节示波器后，图像就很好地显示出来了。

编程读出数据

学校为被迫做这个实验的学生提供了三个只能在 Linux 中编译的 c 程序，分别是 com(1).c、com.h、main (1).c。特别想吐槽命名，事实上不改命名，终端会报错。对于我，我在 kali 中把那两个带有 (1) 的文件名都去掉了 (1)。所以从此至下文起，一直用不带 (1) 的名称来分别指代三者。

其中在烧录与测量环节中得到的波特率需要更改到其中的 com.c。具体如下。

到三个文件所在的目录中打开终端，敲指令进去编辑。

gedit com.c

然后弹出来改掉要改的内容就是了。

int openSerial(char *cSerialName)
{
    int iFd;

    struct termios opt;  // 省略原有注释...

    iFd = open(cSerialName, O_RDWR | O_NOCTTY |O_NONBLOCK); //...
    //... //...
    if(iFd < 0) {
        perror(cSerialName);
        return -1;
    }

    tcgetattr(iFd, &opt); //... ...  


    // 对于这块板子，只需于下两行代码中修改 B4800 为 B1800 即可；
    // 算出其它的结果是根据测量结果的不同而不同，需要对照文档所提供的数字来写。

    cfsetispeed(&opt, B4800); //...
    cfsetospeed(&opt, B4800); //...
    /*省略后文*/

最后就只需要将单片计算机连接到 kali，编译运行 main.c 就能得到串口正在输出的数据序列了。相应在终端上的指令是

gcc main.c -o main
./main

注意，取的是稳定以后的序列（总共 $13$ 个十六进制数）。当时由于做出后太过激动，没能截取程序运行状况。但是可以得到一个答案：a0 1d d0 33 1a 5a 7b 07 69 ab 93 aa 55。队友由于本身的板子有误，也因此和我提交了一样的答案。可能会扣点分吧，就让他们扣吧，咱也不稀罕这点分，这点分代表不了什么。

提交答案

答案序列是连起来写的，中间没有空格。具体形式如下。

输出 OK 拿去给老师验收就做完这个实验了。如果是 ERROR 或者别的什么，（悲）继续做实验罢。

这次参考资料硬要说有的话，也许就是那份极其垃圾的说明文档吧。