Papan peringkat

 User            Language      Score
 =========================================
 Ell             C++11         293,619,555
 feersum         C++11         100,993,667
 Ell             C++11          78,824,732
 Geobits         Java           27,817,255
 Ell             Python         27,797,402
 Peter Taylor    Java                2,468
 <reference>     Julia                 530

Latar Belakang

Saat bekerja pada kisi-kisi integer 2-D, Anda terkadang ingin tahu apakah dua vektor (dengan komponen bilangan bulat) memiliki besaran yang sama. Tentu saja, dalam geometri Euclidean besarnya vektor (x,y)diberikan oleh

√(x² + y²)

Jadi implementasi yang naif mungkin menghitung nilai ini untuk kedua vektor dan membandingkan hasilnya. Tidak hanya itu menimbulkan perhitungan akar kuadrat yang tidak perlu, itu juga menyebabkan masalah dengan ketidakakuratan titik mengambang, yang mungkin menghasilkan positif palsu: vektor yang besarnya berbeda, tetapi di mana angka signifikan dalam representasi titik mengambang semuanya identik.

Untuk keperluan tantangan ini, kami mendefinisikan false positive sebagai sepasang pasangan koordinat (a,b)dan (c,d)yang:

• Besarnya kuadrat berbeda ketika diwakili sebagai integer 64-bit unsigned.
• Besarnya identik ketika direpresentasikan sebagai angka titik mengambang biner 64-bit dan dihitung melalui akar kuadrat 64-bit (sesuai IEEE 754 ).

Sebagai contoh, menggunakan representasi 16-bit (bukan 64), ¹ pasangan vektor terkecil yang menghasilkan false positive adalah

(25,20) and (32,0)

Besarnya kuadrat kuadrat adalah 1025dan 1024. Mengambil hasil akar kuadrat

32.01562118716424 and 32.0

Tetapi dalam 16-bit mengapung keduanya terpotong 32.0 .

Demikian juga, pasangan ² terkecil yang menghasilkan false positive untuk representasi 32-bit adalah

(1659,1220) and (1951,659)

_{¹ "Terkecil" yang diukur dengan besarnya floating point 16-bit.
² "Terkecil" yang diukur dengan besarnya floating point 32-bit.}

Terakhir, berikut adalah beberapa kasus 64-bit yang valid:

 (51594363,51594339) and (54792160,48184783)
 (54356775,54353746) and (54620742,54088476)
 (54197313,46971217) and (51758889,49645356)
 (67102042,  956863) and (67108864,       6) *

* Kasus terakhir adalah salah satu dari beberapa dengan besaran terkecil yang mungkin untuk false positive 64-bit.

Tantangan

Dalam kurang dari 10.000 byte kode, menggunakan satu utas, Anda akan menemukan banyak false positive untuk angka floating point 64-bit (biner) dalam rentang koordinat 0 ≤ y ≤ x(yaitu, hanya dalam oktan pertama bidang Euclidian) sedemikian rupa sehingga dalam 10 menitx² + y² ≤ 253 . Jika dua pengajuan mengikat untuk jumlah pasangan yang sama, pemecah dasi adalah waktu yang sebenarnya diambil untuk menemukan pasangan yang terakhir.

Program Anda tidak boleh menggunakan lebih dari 4 GB memori kapan saja (untuk alasan praktis).

Harus dimungkinkan untuk menjalankan program Anda dalam dua mode: satu yang menghasilkan setiap pasangan saat menemukannya, dan satu yang hanya menampilkan jumlah pasangan yang ditemukan di akhir. Yang pertama akan digunakan untuk memverifikasi validitas pasangan Anda (dengan melihat beberapa sampel output) dan yang kedua akan digunakan untuk menentukan waktu pengiriman Anda. Perhatikan bahwa pencetakan harus menjadi satu - satunya perbedaan. Secara khusus, program penghitungan tidak dapat meng-hardcode jumlah pasangan yang dapat ditemukannya. Itu masih harus melakukan loop yang sama yang akan digunakan untuk mencetak semua angka, dan hanya menghilangkan pencetakan itu sendiri!

Saya akan menguji semua pengiriman pada laptop Windows 8 saya, jadi silakan tanyakan di komentar jika Anda ingin menggunakan bahasa yang tidak terlalu umum.

Perhatikan bahwa pasangan tidak boleh boleh dihitung dua kali dengan beralih dari pasangan koordinat pertama dan kedua.

Perhatikan juga bahwa saya akan menjalankan proses Anda melalui pengontrol Ruby, yang akan mematikan proses Anda jika belum selesai setelah 10 menit. Pastikan untuk menampilkan jumlah pasangan yang ditemukan saat itu. Anda dapat melacak waktu sendiri dan mencetak hasilnya sesaat sebelum 10 menit berlalu, atau Anda hanya menampilkan jumlah pasangan yang ditemukan secara sporadis, dan saya akan mengambil angka terakhir sebagai skor Anda.

C ++, 275.000.000+

Kami akan merujuk pada pasangan yang besarnya dapat diwakili secara akurat, seperti (x, 0) , sebagai pasangan yang jujur dan untuk semua pasangan lainnya sebagai pasangan yang tidak jujur dari besarnya m , di mana m adalah besarnya pasangan yang dilaporkan salah. Program pertama pada posting sebelumnya menggunakan satu set pasangan yang terkait erat pasangan jujur ​​dan tidak jujur:
(x, 0) dan (x, 1) , masing-masing, untuk x yang cukup besar. Program kedua menggunakan pasangan pasangan yang tidak jujur ​​yang sama tetapi memperluas pasangan pasangan yang jujur ​​dengan mencari semua pasangan jujur ​​yang besarnya tidak terpisahkan. Program tidak berhenti dalam sepuluh menit, tetapi menemukan sebagian besar hasilnya sangat awal, yang berarti bahwa sebagian besar runtime menjadi sia-sia. Alih-alih terus mencari pasangan jujur ​​yang semakin jarang, program ini menggunakan waktu luang untuk melakukan hal logis berikutnya: memperluas perangkat tidak jujur pasangan .

Dari posting sebelumnya kita tahu bahwa untuk semua bilangan bulat besar cukup r , sqrt (r ² + 1) = r , di mana sqrt adalah floating-point fungsi akar kuadrat. Rencana kami dari serangan adalah untuk menemukan pasangan P = (x, y) sehingga x ² + y ² = r ² + 1 untuk beberapa besar cukup bilangan bulat r . Itu cukup sederhana untuk dilakukan, tetapi secara naif mencari pasangan seperti itu terlalu lambat untuk menarik. Kami ingin menemukan pasangan ini secara massal, seperti yang kami lakukan untuk pasangan jujur ​​di program sebelumnya.

Biarkan { v , w } menjadi pasangan vektor ortonormal. Untuk semua skalar nyata r , || r v + w || ² = r ² + 1 . Dalam ℝ ² , ini adalah akibat langsung dari teorema Pythagoras:

Kami sedang mencari vektor v dan w sedemikian rupa sehingga ada bilangan r yang x dan y juga bilangan bulat. Sebagai catatan tambahan, perhatikan bahwa himpunan pasangan tidak jujur ​​yang kami gunakan dalam dua program sebelumnya hanyalah kasus khusus ini, di mana { v , w } adalah basis standar dari ℝ ² ; kali ini kami ingin menemukan solusi yang lebih umum. Di sinilah kembar tiga Pythagoras (bilangan bulat bilangan bulat (a, b, c) memuaskan a ² + b ² = c ², yang kami gunakan dalam program sebelumnya) membuat comeback mereka.

Biarkan (a, b, c) menjadi triplet Pythagoras. Vektor v = (b / c, a / c) dan w = (-a / c, b / c) (dan juga
w = (a / c, -b / c) ) adalah ortonormal, karena mudah untuk memverifikasi . Ternyata, untuk sembarang pilihan triplet Pythagoras, ada bilangan r sehingga x dan y adalah bilangan bulat. Untuk membuktikan ini, dan untuk secara efektif menemukan r dan P , kita memerlukan sedikit teori / kelompok; Saya akan menyimpan rinciannya. Either way, misalkan kita memiliki r integral , x dan y . Kami masih kekurangan beberapa hal: kami perlu rcukup besar dan kami ingin metode cepat untuk mendapatkan lebih banyak pasangan serupa dari yang satu ini. Untungnya, ada cara sederhana untuk mencapai ini.

Perhatikan bahwa proyeksi P ke v adalah r v , maka r = P · v = (x, y) · (b / c, a / c) = xb / c + ya / c , semua ini untuk mengatakan bahwa xb + ya = rc . Akibatnya, untuk semua bilangan bulat n , (x + bn) ² + (y + an) ² = (x ² + y ² ) + 2 (xb + ya) n + (a ² + b ² ) n ² = ( r ² + 1) + 2 (rc) n + (c ² ) n ² = (r + cn) ² + 1. Dengan kata lain, besarnya kuadrat pasangan bentuk
(x + bn, y + an) adalah (r + cn) ² + 1 , yang merupakan jenis pasangan yang kami cari! Untuk n yang cukup besar , ini adalah pasangan magnitudo r + cn yang tidak jujur .

Selalu menyenangkan untuk melihat contoh konkret. Jika kita mengambil triplet Pythagoras (3, 4, 5) , maka pada r = 2 kita memiliki P = (1, 2) (Anda dapat memeriksa bahwa (1, 2) · (4/5, 3/5) = 2 dan, jelas, 1 ² + 2 ² = 2 ² + 1. ) Menambahkan 5 ke r dan (4, 3) ke P membawa kita ke r '= 2 + 5 = 7 dan P' = (1 + 4, 2 + 3) = (5, 5) . Lihatlah, 5 ² + 5 ² = 7 ² + 1. Koordinat berikutnya adalahr '' = 12 dan P '' = (9, 8) , dan lagi, 9 ² + 8 ² = 12 ² +1 , dan seterusnya, dan seterusnya ...

Setelah r cukup besar, kita mulai mendapatkan pasangan yang tidak jujur ​​dengan kenaikan 5 . Itu kira-kira 27.797.402 / 5 pasangan tidak jujur.

Jadi sekarang kita memiliki banyak pasangan tidak jujur ​​yang besarnya tidak terpisahkan. Kita dapat dengan mudah memasangkan mereka dengan pasangan jujur ​​dari program pertama untuk membentuk false-positive, dan dengan hati-hati kita juga dapat menggunakan pasangan jujur ​​dari program kedua. Inilah yang pada dasarnya dilakukan oleh program ini. Seperti program sebelumnya, ia juga menemukan sebagian besar hasilnya sangat awal --- ia mendapatkan 200.000.000 positif palsu dalam beberapa detik --- dan kemudian melambat jauh.

Kompilasi dengan g++ flspos.cpp -oflspos -std=c++11 -msse2 -mfpmath=sse -O3. Untuk memverifikasi hasil, tambahkan -DVERIFY(ini akan lebih lambat.)

Jalankan dengan flspos. Argumen baris perintah apa pun untuk mode verbose.

#include <cstdio>
#define _USE_MATH_DEFINES
#undef __STRICT_ANSI__
#include <cmath>
#include <cfloat>
#include <vector>
#include <iterator>
#include <algorithm>

using namespace std;

/* Make sure we actually work with 64-bit precision */
#if defined(VERIFY) && FLT_EVAL_METHOD != 0 && FLT_EVAL_METHOD != 1
#   error "invalid FLT_EVAL_METHOD (did you forget `-msse2 -mfpmath=sse'?)"
#endif

template <typename T> struct widen;
template <> struct widen<int> { typedef long long type; };

template <typename T>
inline typename widen<T>::type mul(T x, T y) {
    return typename widen<T>::type(x) * typename widen<T>::type(y);
}
template <typename T> inline T div_ceil(T a, T b) { return (a + b - 1) / b; }
template <typename T> inline typename widen<T>::type sq(T x) { return mul(x, x); }
template <typename T>
T gcd(T a, T b) { while (b) { T t = a; a = b; b = t % b; } return a; }
template <typename T>
inline typename widen<T>::type lcm(T a, T b) { return mul(a, b) / gcd(a, b); }
template <typename T>
T div_mod_n(T a, T b, T n) {
    if (b == 0) return a == 0 ? 0 : -1;
    const T n_over_b = n / b, n_mod_b = n % b;
    for (T m = 0; m < n; m += n_over_b + 1) {
        if (a % b == 0) return m + a / b;
        a -= b - n_mod_b;
        if (a < 0) a += n;
    }
    return -1;
}

template <typename T> struct pythagorean_triplet { T a, b, c; };
template <typename T>
struct pythagorean_triplet_generator {
    typedef pythagorean_triplet<T> result_type;
private:
    typedef typename widen<T>::type WT;
    result_type p_triplet;
    WT p_c2b2;
public:
    pythagorean_triplet_generator(const result_type& triplet = {3, 4, 5}) :
        p_triplet(triplet), p_c2b2(sq(triplet.c) - sq(triplet.b))
    {}
    const result_type& operator*() const { return p_triplet; }
    const result_type* operator->() const { return &p_triplet; }
    pythagorean_triplet_generator& operator++() {
        do {
            if (++p_triplet.b == p_triplet.c) {
                ++p_triplet.c;
                p_triplet.b = ceil(p_triplet.c * M_SQRT1_2);
                p_c2b2 = sq(p_triplet.c) - sq(p_triplet.b);
            } else
                p_c2b2 -= 2 * p_triplet.b - 1;
            p_triplet.a = sqrt(p_c2b2);
        } while (sq(p_triplet.a) != p_c2b2 || gcd(p_triplet.b, p_triplet.a) != 1);
        return *this;
    }
    result_type operator()() { result_type t = **this; ++*this; return t; }
};

int main(int argc, const char* argv[]) {
    const bool verbose = argc > 1;

    const int min = 1 << 26;
    const int max = sqrt(1ll << 53);

    const size_t small_triplet_count = 1000;
    vector<pythagorean_triplet<int>> small_triplets;
    small_triplets.reserve(small_triplet_count);
    generate_n(
        back_inserter(small_triplets),
        small_triplet_count,
        pythagorean_triplet_generator<int>()
    );

    int found = 0;
    auto add = [&] (int x1, int y1, int x2, int y2) {
#ifdef VERIFY
        auto n1 = sq(x1) + sq(y1), n2 = sq(x2) + sq(y2);
        if (x1 < y1 || x2 < y2 || x1 > max || x2 > max ||
            n1 == n2 || sqrt(n1) != sqrt(n2)
        ) {
            fprintf(stderr, "Wrong false-positive: (%d, %d) (%d, %d)\n",
                    x1, y1, x2, y2);
            return;
        }
#endif
        if (verbose) printf("(%d, %d) (%d, %d)\n", x1, y1, x2, y2);
        ++found;
    };

    int output_counter = 0;
    for (int x = min; x <= max; ++x) add(x, 0,    x, 1);
    for (pythagorean_triplet_generator<int> i; i->c <= max; ++i) {
        const auto& t1 = *i;

        for (int n = div_ceil(min, t1.c); n <= max / t1.c; ++n)
            add(n * t1.b, n * t1.a,    n * t1.c, 1);

        auto find_false_positives = [&] (int r, int x, int y) {
            {
                int n = div_ceil(min - r, t1.c);
                int min_r = r + n * t1.c;
                int max_n = n + (max - min_r) / t1.c;
                for (; n <= max_n; ++n)
                    add(r + n * t1.c, 0,    x + n * t1.b, y + n * t1.a);
            }
            for (const auto t2 : small_triplets) {
                int m = div_mod_n((t2.c - r % t2.c) % t2.c, t1.c % t2.c, t2.c);
                if (m < 0) continue;
                int sr = r + m * t1.c;
                int c = lcm(t1.c, t2.c);
                int min_n = div_ceil(min - sr, c);
                int min_r = sr + min_n * c;
                if (min_r > max) continue;
                int x1 = x + m * t1.b, y1 = y + m * t1.a;
                int x2 = t2.b * (sr / t2.c), y2 = t2.a * (sr / t2.c);
                int a1 = t1.a * (c / t1.c), b1 = t1.b * (c / t1.c);
                int a2 = t2.a * (c / t2.c), b2 = t2.b * (c / t2.c);
                int max_n = min_n + (max - min_r) / c;
                int max_r = sr + max_n * c;
                for (int n = min_n; n <= max_n; ++n) {
                    add(
                        x2 + n * b2, y2 + n * a2,
                        x1 + n * b1, y1 + n * a1
                    );
                }
            }
        };
        {
            int m = div_mod_n((t1.a - t1.c % t1.a) % t1.a, t1.b % t1.a, t1.a);
            find_false_positives(
                /* r = */ (mul(m, t1.c) + t1.b) / t1.a,
                /* x = */ (mul(m, t1.b) + t1.c) / t1.a,
                /* y = */ m
            );
        } {
            int m = div_mod_n((t1.b - t1.c % t1.b) % t1.b, t1.a, t1.b);
            find_false_positives(
                /* r = */ (mul(m, t1.c) + t1.a) / t1.b,
                /* x = */ m,
                /* y = */ (mul(m, t1.a) + t1.c) / t1.b
            );
        }

        if (output_counter++ % 50 == 0)
            printf("%d\n", found), fflush(stdout);
    }
    printf("%d\n", found);
}

Python, 27.797.402

Hanya untuk mengatur bar sedikit lebih tinggi ...

from sys import argv
verbose = len(argv) > 1
found = 0
for x in xrange(67108864, 94906266):
    found += 1
    if verbose:
        print "(%d, 0) (%d, 1)" % (x, x)
print found

Sangat mudah untuk memverifikasi bahwa untuk semua 67.108.864 <= x <= 94.906.265 = lantai (sqrt (2 ⁵³ )) pasangan (x, 0) dan (x, 1) adalah positif palsu.

Mengapa ini berhasil : 67.108.864 = 2 ²⁶ . Karenanya, semua angka x dalam rentang di atas adalah dari bentuk 2 ²⁶ + x ' untuk beberapa 0 <= x' <2 ²⁶ . Untuk semua e positif , (x + e) ² = x ² + 2xe + e ² = x ² + 2 ²⁷ e + 2x'e + e ² . Jika kita ingin memiliki
(x + e) ² = x ² + 1 kita membutuhkan setidaknya 2 ²⁷ e <= 1 , yaitu, e <= 2 ^-27 Namun, karena mantissa angka floating point presisi ganda adalah lebar 52-bit, e terkecil sehingga x + e> x . Karena mode pembulatan IEEE-754 default adalah bulat-ke-terdekat; ties-to-even, itu akan selalu membulatkan ke x dan tidak pernah ke x + 2 ^-26 (di mana tie-break benar-benar hanya relevan untuk x = 67.108.864is e = 2 ^{26 - 52} = 2 ^-26 . Dengan kata lain, angka keterwakilan terkecil yang lebih besar dari x adalah x + 2 ^-26 sedangkan hasil sqrt (x ² + 1) paling banyak adalah x + 2 ^-27, jika sama sekali. Angka yang lebih besar akan membulatkan ke x terlepas).

C ++, 75.000.000+

Ingat bahwa 3 ² + 4 ² = 5 ² . Ini artinya bahwa titik (4, 3) terletak pada lingkaran jari-jari 5 yang berpusat di sekitar titik asal. Sebenarnya, untuk semua bilangan bulat n , (4n, 3n) terletak pada lingkaran jari-jari 5n . Untuk n yang cukup besar (yaitu, sehingga 5n> = 2 ²⁶ ), kita sudah tahu false-positive untuk semua poin pada lingkaran ini: (5n, 1) . Bagus! Itu yang lain 27.797.402 / 5 pasangan positif-palsu gratis di sana! Tetapi mengapa berhenti di sini? (3, 4, 5) bukan satu-satunya triplet tersebut.

Program ini terlihat untuk semua kembar tiga bilangan bulat positif (a, b, c) sehingga suatu ² + b ² = c ² , dan jumlah false-positif dalam mode ini. Ia mencapai 70.000.000 positif palsu cukup cepat tetapi kemudian melambat ketika jumlahnya bertambah.

Kompilasi dengan g++ flspos.cpp -oflspos -std=c++11 -msse2 -mfpmath=sse -O3. Untuk memverifikasi hasil, tambahkan -DVERIFY(ini akan lebih lambat.)

Jalankan dengan flspos. Argumen baris perintah apa pun untuk mode verbose.

#include <cstdio>
#include <cmath>
#include <cfloat>

using namespace std;

/* Make sure we actually work with 64-bit precision */
#if defined(VERIFY) && FLT_EVAL_METHOD != 0 && FLT_EVAL_METHOD != 1
#   error "invalid FLT_EVAL_METHOD (did you forget `-msse2 -mfpmath=sse'?)"
#endif

template <typename T> inline long long sqr(T x) { return 1ll * x * x; }
template <typename T>
T gcd(T a, T b) { while (b) { T t = a; a = b; b = t % b; } return a; }

int main(int argc, const char* argv[]) {
    const bool verbose = argc > 1;

    const int min = 1 << 26;
    const int max = sqrt(1ll << 53);

    int found = 0;
    auto add = [=, &found] (int x1, int y1, int x2, int y2) {
#ifdef VERIFY
        auto n1 = sqr(x1) + sqr(y1), n2 = sqr(x2) + sqr(y2);
        if (n1 == n2 || sqrt(n1) != sqrt(n2)) {
            fprintf(stderr, "Wrong false-positive: (%d, %d) (%d, %d)\n",
                    x1, y1, x2, y2);
            return;
        }
#endif
        if (verbose) printf("(%d, %d) (%d, %d)\n", x1, x2, y1, y2);
        ++found;
    };

    for (int x = min; x <= max; ++x) add(x, 0,    x, 1);

    for (int a = 1; a < max; ++a) {
        auto a2b2 = sqr(a) + 1;
        for (int b = 1; b <= a; a2b2 += 2 * b + 1, ++b) {
            int c = sqrt(a2b2);
            if (a2b2 == sqr(c) && gcd(a, b) == 1) {
                int max_c = max / c;
                for (int n = (min + c - 1) / c; n <= max_c; ++n)
                    add(n * a, n * b,    n * c, 1);
            }
        }

        if (a % 512 == 0) printf("%d\n", found), fflush(stdout);
    }

    printf("%d\n", found);
}

C ++ 11 - 100.993.667

EDIT: Program baru.

Yang lama menggunakan terlalu banyak memori. Ini membagi dua penggunaan memori dengan menggunakan array vektor raksasa, bukan tabel hash. Itu juga menghilangkan cruft utas acak.

   /* by feersum  2014/9
   http://codegolf.stackexchange.com/questions/37627/false-positives-on-an-integer-lattice */
#include <iostream>
#include <cmath>
#include <cstdlib>
#include <cstring>
#include <functional>
#include <vector>
using namespace std;
#define ul unsigned long long

#define K const



#define INS(A)   { bool already = false; \
    for(auto e = res[A.p[0][0]].end(), it = res[A.p[0][0]].begin(); it != e; ++it) \
        if(A.p[0][1] == it->y1 && A.p[1][0] == it->x2 && A.p[1][1] == it->y2) { \
            already = true; \
            break; } \
    if(!already) res[A.p[0][0]].push_back( {A.p[0][1], A.p[1][0], A.p[1][1]} ), ++n; }

#define XMAXMIN (1<<26)

struct ints3 {
    int y1, x2, y2;
};


struct pparm {
    int a,b,c,d;
    int E[4];
    pparm(int A,int B,int C, int D):
        E{B*B+D*D,A*B+C*D,A*A+C*C+2*(B+D),A+C}
    {
        a=A;b=B;c=C;d=D;
    }

};

struct ans {
    int p[2][2];

};
ostream&operator<<(ostream&o, ans&a)
{
    o<<'('<<a.p[0][0]<<','<<a.p[0][1]<<"),("<<a.p[1][0]<<','<<a.p[1][1]<<')'<<endl;
    return o;
}



vector<ints3> res[XMAXMIN];

bool print;
int n;

void gen(K pparm&p1, K pparm&p2)
{
#ifdef DBUG
    for(int i=0;i<2;i++){
    K pparm&p=i?p2:p1;
    cout<<' '<<p.a<<' '<<p.b<<' '<<p.c<<' '<<p.d<<' ';}
    cout<<endl;
#endif

    for(ul x = 0; ; ++x) {
        ans a;
        ul s[2];
        for(int i = 0; i < 2; i++) {
            K pparm &p = i?p2:p1;
            int *pt = a.p[i];
            pt[0] = p.b+x*(p.a+x);
            pt[1] = p.d+x*(p.c+x);
            s[i] = (ul)pt[0]*pt[0] + (ul)pt[1]*pt[1];
        }
        if(*s >> 53)
            break;
if(s[1] - s[0] != 1)
exit(4);

        if(sqrt(s[0]) == sqrt(s[1])) {
             for(int i = 0; i < 2; i++)
                if(a.p[i][0] > a.p[i][1])
                    swap(a.p[i][0], a.p[i][1]);
            if(a.p[0][0] > a.p[0][1])
                for(int i = 0; i < 2; i++)
                    swap(a.p[0][i], a.p[1][i]);
            INS(a)
        }
    }
}



int main(int ac, char**av)
{
    for(int i = 1; i < ac; i++) {
        print |= !strcmp(av[1], "-P");
    }


    #define JMAX 43000000
    for(ul j = 0; j < JMAX; j++) {
        pparm p1(-~j,j,~-j,0),p2(j,1,j,j);
        gen(p1,p2);
        if(!print && !(j%1024))
#ifdef DBUG
            cout<<j<<' ',
#endif
            cout<<n<<endl;

    }
    if(print) 
        for(vector<ints3>& v: res)
            for(ints3& i: v)
                printf("(%d,%d),(%d,%d)\n", &v - res, i.y1, i.x2, i.y2);

    return 0;
}

Jalankan dengan a -P argumen untuk mencetak poin, bukan nomornya.

Bagi saya dibutuhkan kurang dari 2 menit dalam mode penghitungan dan sekitar 5 menit dengan pencetakan diarahkan ke file (~ 4 GB), sehingga tidak menjadi I / O terbatas.

Program asli saya rapi, tetapi saya membuang sebagian besar karena hanya bisa menghasilkan urutan 10 ^ 5 hasil. Apa yang dilakukannya adalah mencari parameterisasi bentuk (x ^ 2 + Ax + B, x ^ 2 + Cx + D), (x ^ 2 + ax + b, x ^ 2 + cx + d) sedemikian rupa sehingga untuk setiap x, (x ^ 2 + Kapak + B) ^ 2 + (x ^ 2 + Cx + D) ^ 2 = (x ^ 2 + kapak + b) ^ 2 + (x ^ 2 + cx + d) ^ 2 + 1. Ketika ia menemukan seperangkat parameter {a, b, c, d, A, B, C, D} itu dilanjutkan untuk memeriksa semua nilai x di bawah maksimum. Sambil melihat output debug saya dari program ini, saya melihat parameterisasi tertentu dari parameterisasi yang memungkinkan saya menghasilkan banyak angka dengan mudah. Saya memilih untuk tidak mencetak nomor Ell karena saya punya banyak nomor. Mudah-mudahan sekarang seseorang tidak akan mencetak kedua set angka kami dan mengklaim sebagai pemenang :)

 /* by feersum  2014/9
   http://codegolf.stackexchange.com/questions/37627/false-positives-on-an-integer-lattice */
    #include <iostream>
    #include <cmath>
    #include <cstdlib>
    #include <cstring>
    #include <functional>
    #include <unordered_set>
    #include <thread>
using namespace std;
#define ul unsigned long long

#define h(S) unordered_##S##set
#define P 2977953206964783763LL
#define K const

#define EQ(T, F)bool operator==(K T&o)K{return!memcmp(F,o.F,sizeof(F));}

struct pparm {
    int a,b,c,d;
    int E[4];
    pparm(int A,int B,int C, int D):
        E{B*B+D*D,A*B+C*D,A*A+C*C+2*(B+D),A+C}
    {
        a=A;b=B;c=C;d=D;
    }
    EQ(pparm,E)
};

struct ans {
    int p[2][2];
    EQ(ans,p)
};
ostream&operator<<(ostream&o, ans&a)
{
    o<<'('<<a.p[0][0]<<','<<a.p[0][1]<<"),("<<a.p[1][0]<<','<<a.p[1][1]<<')'<<endl;
    return o;
}

#define HASH(N,T,F) \
struct N { \
    size_t operator() (K T&p) K { \
        size_t h = 0; \
        for(int i = 4; i--; ) \
            h=h*P+((int*)p.F)[i]; \
        return h; \
    }};

#define INS(r, a) { \
    bool new1 = r.insert(a).second; \
    n += new1; \
    if(print && new1) \
        cout<<a; }

HASH(HA,ans,p)

bool print;
int n;

void gen(h()<ans,HA>&r, K pparm&p1, K pparm&p2)
{
#ifdef DBUG
    for(int i=0;i<2;i++){
    K pparm&p=i?p2:p1;
    cout<<' '<<p.a<<' '<<p.b<<' '<<p.c<<' '<<p.d<<' ';}
    cout<<endl;
#endif

    for(ul x = 0; ; ++x) {
        ans a;
        ul s[2];
        for(int i = 0; i < 2; i++) {
            K pparm &p = i?p2:p1;
            int *pt = a.p[i];
            pt[0] = p.b+x*(p.a+x);
            pt[1] = p.d+x*(p.c+x);
            s[i] = (ul)pt[0]*pt[0] + (ul)pt[1]*pt[1];
        }
        if(*s >> 53)
            break;
if(s[1] - s[0] != 1)
exit(4);

        if(sqrt(s[0]) == sqrt(s[1])) {
             for(int i = 0; i < 2; i++)
                if(a.p[i][0] > a.p[i][1])
                    swap(a.p[i][0], a.p[i][1]);
            INS(r,a)
        }
    }
    //if(!print) cout<<n<<endl;
}

void endit()
{
    this_thread::sleep_for(chrono::seconds(599));
    exit(0);
}

int main(int ac, char**av)
{
    bool kill = false;
    for(int i = 1; i < ac; i++) {
        print |= ac>1 && !stricmp(av[1], "-P");
        kill |= !stricmp(av[i], "-K");
    }

    thread KILLER;
    if(kill)
        KILLER = thread(endit);

    h()<ans, HA> res;
    res.reserve(1<<27);

    #define JMAX 43000000
    for(ul j = 0; j < JMAX; j++) {
        pparm p1(-~j,j,~-j,0),p2(j,1,j,j);
        gen(res,p1,p2);
        if(!print && !(j%1024))
#ifdef DBUG
            cout<<j<<' ',
#endif
            cout<<n<<endl;

    }
    exit(0);
}

Pemindaian lingkaran Java, Bresenham-esque

Secara heuristik saya berharap untuk mendapatkan lebih banyak tabrakan dengan mulai dari ujung annulus yang lebih luas. Saya berharap mendapatkan beberapa peningkatan dengan melakukan satu pemindaian untuk setiap tabrakan, mencatat nilai yang surplusada di antaranya 0dan r2max - r2inklusif, tetapi dalam pengujian saya yang terbukti lebih lambat daripada versi ini. Demikian pula upaya untuk menggunakan int[]buffer tunggal daripada banyak pembuatan array dan daftar dua elemen. Optimalisasi kinerja memang binatang yang aneh.

Jalankan dengan argumen baris perintah untuk output dari pasangan, dan tanpa untuk penghitungan sederhana.

import java.util.*;

public class CodeGolf37627 {
    public static void main(String[] args) {
        final int M = 144;
        boolean[] possible = new boolean[M];
        for (int i = 0; i <= M/2; i++) {
            for (int j = 0; j <= M/2; j++) {
                possible[(i*i+j*j)%M] = true;
            }
        }

        long count = 0;
        double sqrt = 0;
        long r2max = 0;
        List<int[]> previousPoints = null;
        for (long r2 = 1L << 53; ; r2--) {
            if (!possible[(int)(r2 % M)]) continue;

            double r = Math.sqrt(r2);
            if (r != sqrt) {
                sqrt = r;
                r2max = r2;
                previousPoints = null;
            }
            else {
                if (previousPoints == null) previousPoints = findLatticePointsBresenham(r2max, (int)r);

                if (previousPoints.size() == 0) {
                    r2max = r2;
                    previousPoints = null;
                }
                else {
                    List<int[]> points = findLatticePointsBresenham(r2, (int)r);
                    for (int[] p1 : points) {
                        for (int[] p2 : previousPoints) {
                            if (args.length > 0) System.out.format("(%d, %d) (%d, %d)\n", p1[0], p1[1], p2[0], p2[1]);
                            count++;
                        }
                    }
                    previousPoints.addAll(points);
                    System.out.println(count);
                }
            }
        }
    }

    // Surprisingly, this seems to be faster than doing one scan for all two or three r2s.
    private static List<int[]> findLatticePointsBresenham(long r2, long r) {
        List<int[]> rv = new ArrayList<int[]>();
        // Require 0 = y = x
        long x = r, y = 0, surplus = r2 - r * r;
        while (y <= x) {
            if (surplus == 0) rv.add(new int[]{(int)x, (int)y});

            // Invariant: surplus = r2 - x*x - y*y >= 0
            y++;
            surplus -= 2*y - 1;
            if (surplus < 0) {
                x--;
                surplus += 2*x + 1;
            }
        }

        return rv;
    }
}

Jawa - 27.817.255

Sebagian besar sama dengan apa yang ditunjukkan Ell , dan sisanya didasarkan pada (j,0) (k,l). Untuk masing-masing j, saya berjalan mundur beberapa kotak dan memeriksa apakah sisanya memberikan hasil positif palsu. Ini pada dasarnya mengambil seluruh waktu dengan hanya 25rb (sekitar 0,1%) lebih dari hanya (j,0) (j,1), tetapi keuntungan adalah keuntungan.

Ini akan selesai dalam waktu kurang dari sepuluh menit pada mesin saya, tetapi saya tidak tahu apa yang Anda miliki. Karena alasan, jika tidak selesai sebelum waktu habis, skornya akan jauh lebih buruk. Dalam hal ini, Anda dapat mengubah pembagi pada baris 8 sehingga selesai tepat waktu (ini hanya menentukan seberapa jauh berjalannya untuk masing-masing j). Untuk beberapa pembagi yang berbeda, nilainya adalah:

11    27817255 (best on OPs machine)
10    27818200
8     27820719
7     27822419 (best on my machine)

Untuk mengaktifkan output untuk setiap pertandingan (dan, Tuhan, lambat jika Anda melakukannya), cukup batalkan komentar pada baris 10 dan 19.

public class FalsePositive {
    public static void main(String[] args){
        long j = 67108864;
        long start = System.currentTimeMillis();
        long matches=0;
        while(j < 94906265 && System.currentTimeMillis()-start < 599900){
            long jSq = j*j;
            long limit = (long)Math.sqrt(j)/11; // <- tweak to fit inside 10 minutes for best results
            matches++; // count an automatic one for (j,0)(j,1)
            //System.out.println("("+j+",0) ("+j+",1)");        
            for(int i=1;i<limit;i++){
                long k = j-i;
                long kSq = k*k;
                long l = (long)Math.sqrt(jSq-kSq);
                long lSq = l*l;
                if(kSq+lSq != jSq){
                    if(Math.sqrt(kSq+lSq)==Math.sqrt(jSq)){
                        matches++;
                        //System.out.println("("+j+",0) ("+k+","+l+")");        
                    }
                }
            }
            j++;
        }
        System.out.println("\n"+matches+" Total matches, got to j="+j);
    }
}

Untuk referensi, 20 output pertama yang diberikan (untuk pembagi = 7, tidak termasuk (j,0)(j,1)jenis) adalah:

(67110083,0) (67109538,270462)
(67110675,0) (67109990,303218)
(67111251,0) (67110710,269470)
(67111569,0) (67110668,347756)
(67112019,0) (67111274,316222)
(67112787,0) (67111762,370918)
(67115571,0) (67115518,84346)
(67117699,0) (67117698,11586)
(67117971,0) (67117958,41774)
(67120545,0) (67120040,260368)
(67121043,0) (67120118,352382)
(67122345,0) (67122320,57932)
(67122449,0) (67122444,25908)
(67122633,0) (67122328,202348)
(67122729,0) (67121972,318784)
(67122849,0) (67122568,194224)
(67124195,0) (67123818,224970)
(67125201,0) (67125172,62396)
(67125705,0) (67124632,379540)
(67126195,0) (67125882,204990)

Julia, 530 positif salah

Berikut ini adalah pencarian brute force yang sangat naif, yang dapat Anda lihat sebagai implementasi referensi.

num = 0
for i = 60000000:-1:0
    for j = i:-1:ifloor(0.99*i)
        s = i*i + j*j
        for x = ifloor(sqrt(s/2)):ifloor(sqrt(s))
            min_y = ifloor(sqrt(s - x*x))
            max_y = min_y+1
            for y = min_y:max_y
                r = x*x + y*y
                if r != s && sqrt(r) == sqrt(s)
                    num += 1
                    if num % 10 == 0
                        println("Found $num pairs")
                    end
                    #@printf("(i,j) = (%d,%d); (x,y) = (%d,%d); s = %d, r = %d\n", i,j,x,y,s,r)
                end
            end
        end
    end
end

Anda dapat mencetak pasangan (dan besarnya kuadrat persisnya) dengan menghapus komentar pada @printfgaris.

Pada dasarnya ini memulai pencarian pada x = y = 6e7pasangan koordinat pertama dan memindai sekitar 1% dari jalan ke sumbu x sebelum mengurangi x. Kemudian untuk setiap pasangan koordinat seperti itu, ia memeriksa seluruh busur dengan besaran yang sama (membulatkan ke atas dan ke bawah) untuk tabrakan.

Kode mengasumsikan bahwa itu dijalankan pada sistem 64-bit, sehingga tipe integer dan floating-point default adalah 64-bit (jika tidak, Anda dapat membuatnya dengan int64()dan float64()konstruktor).

Itu menghasilkan sedikit hasil 530.