けんちょんの競プロ精進記録

競プロの精進記録や小ネタを書いていきます

New Year Contest 2015 L - 机のしみ

ガウス整数が使える問題と聞いて!!!

問題概要

二次元平面上に  N 個の格子点  (x_{i}, y_{i}) がある。これに最小個数の格子点を追加することで、格子点全体が「正方形グリッド」の模様を形作るようにしたい。追加すべき格子点の個数を求めよ。

制約

  •  1 \le N \le 10^{5}

ガウス整数

 x, y を整数として、 z = x + yi と表せる複素数ガウス整数と呼ぶ。ガウス整数の理論については以下の記事に書いた。

qiita.com

つまり、ガウス整数は

といった特徴を持っている。通常の整数と同様の理論が展開できるのだ。ただし、ユークリッドの互除法については注意が必要 (余りのノルムが商のノルムより小さくなるようにする必要がある)。ユークリッドの互除法における注意点は、maspy さんの記事にて。

maspypy.com

問題の解法

 z = x + yi とする。最終的に形作られる正方形グリッドの正方形の一辺を表すガウス整数を  d = a + bi とする。このとき、

  •  d は、 z_{2} - z_{1}, \dots, z_{N} - z_{1} の公約数

となっている。 d として最大公約数をとればよい。実際、 d が最大公約数でないとき、代わりに最大公約数  D をとれば、最終的な格子点数は  D/d のノルム倍になってしまう。

よって、ガウス整数に対するユークリッドの互除法を適用していけばよい。計算量は  O(N (\log x + \log y)) となる。

コード

#include <bits/stdc++.h>
using namespace std;
template<class T> inline bool chmax(T& a, T b) { if (a < b) { a = b; return 1; } return 0; }
template<class T> inline bool chmin(T& a, T b) { if (a > b) { a = b; return 1; } return 0; }

// Gauss Integer
struct gint {
    long long x, y;
    constexpr gint() : x(0), y(0) { }
    constexpr gint(long long x) : x(x), y(0) { }
    constexpr gint(long long x, long long y) : x(x), y(y) { }
    friend constexpr long long abs(const gint &r) noexcept {
        return r.x * r.x + r.y * r.y;
    }
    constexpr gint operator - () const noexcept {
        return gint(-x, -y);
    }
    constexpr gint operator + (const gint& r) const noexcept { return gint(*this) += r; }
    constexpr gint operator - (const gint& r) const noexcept { return gint(*this) -= r; }
    constexpr gint operator * (const gint& r) const noexcept { return gint(*this) *= r; }
    constexpr gint operator / (const gint& r) const noexcept { return gint(*this) /= r; }
    constexpr gint operator % (const gint& r) const noexcept { return gint(*this) %= r; }
    constexpr gint& operator += (const gint& r) noexcept {
        x += r.x, y += r.y;
        return *this;
    }
    constexpr gint& operator -= (const gint& r) noexcept {
        x -= r.x, y -= r.y;
        return *this;
    }
    constexpr gint& operator *= (const gint& r) noexcept {
        auto tx = x * r.x - y * r.y;
        auto ty = x * r.y + y * r.x;
        x = tx, y = ty;
        return *this;
    }
    constexpr gint& operator /= (const gint& r) noexcept {
        long long a = x, b = y, c = r.x, d = r.y;        
        assert(c != 0 || d != 0);
        long long tx = (a * c + b * d) / (c * c + d * d);
        long long ty = (b * c - a * d) / (c * c + d * d);
        gint res;
        long long dmin = -1;
        for (long long nx = tx - 1; nx <= tx + 1; ++nx) {
            for (long long ny = ty - 1; ny <= ty + 1; ++ny) {
                gint q(nx, ny);
                long long d = abs((*this) - q * r);
                if (dmin == -1 || dmin > d) {
                    dmin = d;
                    res = q;
                }
            }
        }
        return *this = res;
    }
    constexpr gint& operator %= (const gint& r) noexcept {
        gint q = (*this) / r;
        return (*this) -= q * r;
    }
    constexpr bool operator == (const gint& r) const noexcept {
        return x == r.x && y == r.y;
    }
    constexpr bool operator != (const gint& r) const noexcept {
        return x != r.x || y != r.y;
    }
    friend ostream& operator << (ostream &os, const gint& r) noexcept {
        if (r.x == 0 && r.y == 0) return os << "0";
        else if (r.x == 0) return os << r.y << "i";
        else os << r.x;
        if (r.y == 0) return os;
        else if (r.y > 0) return os << " + " << r.y << "i";
        else return os << " - " << (-r.y) << "i";
    }
    friend gint gcd(const gint&x, const gint&y) {
        if (y == 0) return x;
        else return gcd(y, x % y);
    }
};

const long long INF = 1LL<<60;
int main() {
    int N;
    cin >> N;
    vector<gint> A(N);
    for (int i = 0; i < N; ++i) cin >> A[i].x >> A[i].y;
    gint g = 0;
    for (int i = 1; i < N; ++i) g = gcd(g, A[i] - A[0]);

    long long xmin = INF, ymin = INF, xmax = -INF, ymax = -INF;
    for (int i = 0; i < N; ++i) {
        gint q = A[i] / g;
        chmin(xmin, q.x), chmin(ymin, q.y);
        chmax(xmax, q.x), chmax(ymax, q.y);
    }
    long long D = max(xmax - xmin, ymax - ymin) + 1;
    cout << D * D - N << endl;
}