Mister雑記

競プロ記事がメインになりますが、内容は気分次第です。

AGC003-D 「Anticube」 (1100)

一応自力AC。一応。

atcoder.jp

概要

 N個の整数 \{s_i\}が与えられる。ここから以下の条件を満たすようにいくつかの数を選ぶ。

  • 任意の異なる2要素について、その積は立方数でない。

選べる個数の最大値を求めよ。

制約

  •  1 \leq N \leq 10^5
  •  1 \leq s_i \leq 10^{10}

考察

言い換えをすると、「掛け合わせると立方数になる整数間に辺が張られたグラフについて、最大独立集合のサイズを求めよ」となる。 しかし辺は密になる可能性があるのと、最大独立集合は計算量が指数オーダーなので普通にやると間に合わない。

グラフを書きながら色々考えるうちに、以下の事実に気づく。

  •  A B A' B A B'の積が立方数のとき、
    •  A' B'の積も立方数となる。
    • 素因数分解における指数を3で割った余りに変えたとき、 A A' B B'は一致する。

素因数分解における指数を3で割った余りに変える」というのは、例えば 2160 = 2^4 \cdot 3^3 \cdot 5^1 2^1 \cdot 3^0 \cdot 5^1 = 10に変換するということである。以降、この変換が施された後の値を「標準形」と呼ぶことにする。

また、ある標準形の値について、掛け合わせると立方数となる標準形の値を「補数」と呼ぶことにする。補数は素因数分解の指数の 1 2を交換することで得られる。

 \{s_i\}を全て標準形に変換したときに同じ値のものを1つの頂点とする。するとこの頂点から張られる辺は、補数からなる頂点のみとなる。補数からなる頂点も同様なので、2つのうちどちらか1つは選ぶことができる。

よって標準形が同じものを一まとまりにし、補数の集合と要素数が多い方を貪欲に選べばそれが最適解となる。

ただし標準形が 1の集合だけは例外で、ここからは1つだけ選ぶことができる。

高速化

しかし、これを愚直に実装すると間に合わない。というのも、標準形に変換するには素因数分解をする必要があり、素因数分解の計算量が O(\sqrt{s_i})であることから、全体の計算量は O(N \sqrt{s_i})になってしまうからだ。そこで何らかの高速化をする必要がある。

まず標準化だが、これは以下のような高速化をした。

  • 素数の立方であり、 10^{10}以下であるものを前計算で列挙する。
  • この全要素について割れるだけ割る。

立方数は 325個なので、これで間に合う。

あとは補数を高速に求める必要があるのだが、この方法が思いつかなかった。とりあえず 10^5以下の素数を前計算で列挙したものをダメ元で出してみたら、なんと通ってしまった。

実装例

以下のコードは2582msで動く(制約は5000ms)。想定解ではないけど、本番もこれくらいの実行時間の人がそこそこいたのでセーフってことで。

#include <iostream>
#include <map>
#include <vector>

using namespace std;
using ll = long long;

template <typename T, typename U>
T mypow(T b, U n) {
    if (n == 0) return 1;
    if (n == 1) return b;
    if (n % 2 == 0) {
        return mypow(b * b, n / 2);
    } else {
        return mypow(b, n - 1) * b;
    }
}

vector<ll> primes, cubes;

// 10^5以下の素数と、10^10以下の立方数を列挙する
void precalc() {
    bool isp[100010];
    fill(isp, isp + 100010, true);

    for (ll i = 2; i * i < 1e10; ++i) {
        if (!isp[i]) continue;
        primes.push_back(i);
        for (ll j = 2; i * j <= 100000; ++j) isp[i * j] = false;
    }

    for (ll p : primes) {
        if (p * p * p > 1e10) break;
        cubes.push_back(p * p * p);
    }
}

// 標準形に変換
ll delcube(ll n) {
    for (ll c : cubes) {
        while (n % c == 0) n /= c;
    }
    return n;
}

// 標準形の補数を求める
ll anticube(ll n) {
    ll ret = 1;

    for (ll p : primes) {
        if (p * p > n) break;
        if (n % p > 0) continue;

        int c = 0;
        while (n % p == 0) {
            n /= p;
            ++c;
        }
        ret *= mypow(p, 3 - c);
    }

    // 最後に残ったnは素数か1
    // 素数なら指数は1
    ret *= mypow(n, 2);
    return ret;
}

int main() {
    precalc();

    int N;
    cin >> N;

    // 標準形をキーに、その個数を保持する
    map<ll, int> cnt;
    for (int i = 0; i < N; ++i) {
        ll s;
        cin >> s;
        s = delcube(s);
        if (cnt.count(s) == 0) cnt[s] = 0;
        ++cnt[s];
    }

    int ans = 0;
    for (auto p : cnt) {
        if (p.second == 0) continue;
        if (p.first == 1) {  // 立方数の集合は例外
            ++ans;
            continue;
        }

        ll opp = anticube(p.first);

        // 補数集合と、大きい方を選ぶ
        if (cnt.count(opp)) {
            ans += max(cnt[p.first], cnt[opp]);
            cnt[p.first] = cnt[opp] = 0;
        } else {
            ans += cnt[p.first];
        }
    }

    cout << ans << endl;
    return 0;
}

感想

  • ガチで高速化したら1100点なんだろうなぁ。
    • 上の解法は900点くらいな気がする。
  • 考察パートだけで部分点がほしいね。
  • もっと考察が早くなりたい。
    • 累計4時間くらいの考察でようやく解法に至った。