紫书 第7章 暴力求解法 习题7-7 Egyptian Fractions (HARD version) UVA - 12558 解题笔记

Egyptian Fractions (HARD version) - UVA 12558 - Virtual Judge

https://onlinejudge.org/external/125/12558.pdf

这道题是紫书中的例题里提到的埃及分数。在例题里出现过，那么很自然的就想到这道题求解的方法是必然是迭代加深搜素（IDA*）了。

但是实际做起来的时候，还是发现这道题有点恐怖。最后写出来代码没多少，但是需要考虑的edge case多的要命。导致我一开始想到的很多条件在调试的时候发现并不正确。（看来还是太菜了）

题外话：看到题目标题发现这道题竟然就是紫数作者出的，怪不得会拿来当例题：p

题目大意

给定一个分数 a/b （b%a != 0 且 b > a)，要求把这个分数展开成 a/b = 1/c1 + 1/c2 + 1/c3… 的形式。其中 c1 < c2 < c3 ….

要求展开的项数尽量小，项数相等的情况下，最小的项数要尽量小，如果有多解且最小项相等，那么次小的项数要尽量小。此外还会给出一些不能采用的分母，如果碰到就直接跳过。

题目分析

因为题目要求项数最小，用迭代加深搜索的话，我们就可以把深度当成项数。每次枚举可以拆分的下一个项。

比如说，假设我们现在的最大深度是maxd， 当前深度是d [0,maxd]，我们可以在d  = 0的时候把 a/b拆成 a/b = 1/c + an/bn，那么我们要做的就是枚举c的值，然后记录当前c，之后用回溯法递归an/bn (d=1)的值，直到d=maxd。

下一项a和b的计算（an和bn）

假设现在我们已经枚举了c， 我们就有 a/b =  1/c + an/bn，由此得到

an/bn = a/b - 1/c

通分后得到 an/bn = (a*c - b)/ b*c

因此：

an = ac - b

bn = b*c

坑：因为这里用了乘法，如果b和c的值够大用int的话bn可能会溢出！所以数据类型用long long会比较安全

Maxd的起始值

因为至少要有两项，我们可以把Maxd的起始值设置成2，然后从d = 1开始循环，我的代码里我用了maxd = 1，d = 0开始，其实也是一样的。

结束条件

当我们递归到d == maxd的时候就应该结束了。

如果这个时候我们能满足

1. 最后一项b%a == 0，也就是说所有项的分子都是1. 

2. 当前解能满足题目的条件：最小的项数要尽量小，如果有多解且最小项相等，那么次小的项数要尽量小

我们就把当前解当作最终答案。

注意 b%a只可能在最后一项发生。这是因为我们在做深度迭代，如果在之前的项发生了，那就应该在之前的深度就停止迭代了。

C的范围

这道题最难的就在这里。错误的范围会直接导致死循环或者超时。

C的下限

因为每一项都是由前边的a和b展开，所以我们可以得到a/b > 1/c， 由此可知 c > b/a, 所以我们可以写成 c = b/a + 1。

坑：这个条件其实是不充分的。因为当前项的分母必然要大于前一项的分母！所以当有前一项的时候，还得考虑前一项的值。

C的上限

这个条件很难被想出来！

因为之前说到的，假设现在我们的深度是d，最大深度是maxd，现在的展开式是 a/b = 1/c + an/bn，由此可知在此之后还会有maxd - d项。

我们又知道之后所有的项数都比 1/c小，也就是说， an/bn 必然要小于 (maxd - d)/c。

代入之前的算式，得到：

a/b < 1/c + (maxd - d)/c 

=> a * c < (1 + maxd - d) * b

所以当 a * c >= (1 + maxd - d) * b的时候，我们就要终止循环。

优化

另一个可以加的条件是如果现在已经有解，我们可以检查现在的c会不会比结果的最后一项大，如果比最后一项大的话，就算现在的c有解最终的解也不会被采用，所以就可以直接跳过了。

解答树的大小

限于自己水平不够，这道题解答树实在是难以估计。给定不同的数据可能会好解也可能会超时。 连题目中也自己表明了这个点，所以对于这道题，我没有浪费时间在计算解答树的大小上，而是直接上IDA* （毕竟这题已经很明显要用这个算法了）

AC 代码

#include <bits/stdc++.h>

using namespace std;

typedef long long ll;
int A, B, K;
ll cur[30], ans[30];
set<int> blocked;
bool found;

bool isSmaller(int d) {
    if (!ans[0]) return true;
    for (int i = d; i >= 0; i--) {
        if (ans[i] > cur[i]) return true;
        if (ans[i] < cur[i]) return false;
    }
    return false;
}

void bt(ll a, ll b, int d, int maxd) {
    ll c = b / a;  
    if (d == maxd) {
        if (b % a != 0 || blocked.count(c)) return;
        cur[maxd] = c;
        if (isSmaller(maxd)) memcpy(ans, cur, sizeof(cur));
        found = true;
        return;
    } 
    
    if (d > 0) c = max(cur[d - 1], c);
    for (;;) {
        c++;
        if (c * a >= b * (maxd - d + 1)) break;
        if (ans[maxd] && ans[maxd] <= c) break;
        if (blocked.count(c)) continue; 
        ll na = a*c - b, nb = b*c; cur[d] = c;
        bt(na, nb, d + 1, maxd);
    }
}

int main() {
    int t, kase = 1; cin >> t;
    while(t--) {
        cin >> A >> B >> K;
        blocked.clear();
        memset(ans, 0, sizeof(ans));
        for (int i = 0; i < K ; i++) {int tmp; cin >> tmp, blocked.insert(tmp);}
        for (int maxd = 1;; maxd++) {
            found = false;
            bt(A, B, 0, maxd);
            if (found) {
                printf("Case %d: %d/%d=", kase++, A, B);
                for (int i = 0; i <= maxd; i++) 
                    cout << "1/" << ans[i] << (i == maxd ? '\n' : '+');
                break;
            }
        }
    }
    return 0;
}