태그 보관물: arithmetic

arithmetic

정수 격자에 대한 오탐 C++11

리더 보드

 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

배경

정수 좌표의 2 차원 그리드에서 작업 할 때 정수 성분이있는 두 벡터의 크기가 같은지 여부를 알고 싶은 경우가 있습니다. 물론, 유클리드 기하학에서 벡터의 크기가 (x,y)주어진다

√(x² + y²)

따라서 순진한 구현은 두 벡터 모두에 대해이 값을 계산하고 결과를 비교할 수 있습니다. 이로 인해 불필요한 제곱근 계산이 발생할뿐만 아니라 부동 소수점 부정확성에 문제가 발생하여 오 탐지가 발생할 수 있습니다.

이 과제의 목적을 위해, 우리는 정의 위양성 좌표 쌍의 한 쌍으로 (a,b)하고 (c,d)있는을 :

  • 64 비트 부호없는 정수로 표시 될 때 제곱 크기가 다릅니다.
  • 크기는 64 비트 이진 부동 소수점 숫자로 표시되고 64 비트 제곱근을 통해 계산 될 때 동일합니다 ( IEEE 754에 따름 ).

예를 들어, 64 비트 대신 16 비트 표현 을 사용하면 오 탐지를 생성하는 가장 작은 1 쌍의 벡터가

(25,20) and (32,0)

그들의 제곱 제곱 크기는 10251024입니다. 제곱근 수확량을 복용

32.01562118716424 and 32.0

그러나 16 비트 플로트에서는 둘 다 잘립니다 32.0.

마찬가지로 32 비트 표현에 대해 오 탐지가 발생하는 가장 작은 2 쌍은 다음과 같습니다.

(1659,1220) and (1951,659)


1 16 비트 부동 소수점 크기로 측정 한 “가장 작은”. 32 비트 부동 소수점 크기로 측정 된 “가장 작은”
2 개 .

마지막으로 유효한 64 비트 사례는 다음과 같습니다.

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

* 마지막 경우는 64 비트 오 탐지에 대해 가능한 가장 작은 크기를 가진 몇 가지 중 하나입니다.

도전

단일 스레드를 사용하는 10,000 바이트 미만의 코드에서는 좌표 범위 (즉, Euclidian 평면의 첫 번째 영역 내에서만 )에있는 64 비트 (이진) 부동 소수점 숫자 에 대해 많은 오 탐지를 찾을 수 있습니다. 0 ≤ y ≤ x그런 그 안에 십분 . 두 개의 제출이 동일한 수의 쌍으로 묶인 경우, 동점 차단기는 해당 쌍의 마지막을 찾는 데 걸리는 실제 시간입니다.x² + y² ≤ 253

프로그램은 실제적인 이유로 언제든지 4GB 이상의 메모리를 사용해서는 안됩니다.

프로그램은 두 가지 모드로 실행할 수 있어야합니다. 하나는 찾은대로 모든 쌍을 출력하고 다른 하나는 끝에 발견 된 쌍의 수만 출력합니다. 첫 번째는 출력의 일부 샘플을보고 쌍의 유효성을 확인하는 데 사용되고 두 번째는 실제로 제출 타이밍을 지정하는 데 사용됩니다. 인쇄는 유일한 차이점 이어야합니다 . 특히, 카운팅 프로그램은 찾을 있는 쌍 수를 하드 코딩하지 않을 수 있습니다. 모든 숫자를 인쇄하는 데 사용되는 것과 동일한 루프를 수행해야하며 인쇄 자체 만 생략해야합니다!

Windows 8 랩톱에서 모든 제출물을 테스트 할 것이므로 자주 사용하지 않는 언어를 사용하려면 의견을 보내주십시오.

첫 번째와 두 번째 좌표 쌍을 전환 할 때 쌍을 두 번 계산 해서는 안됩니다 .

또한 루비 컨트롤러를 통해 프로세스를 실행합니다. 10 분 후에 프로세스가 완료되지 않으면 프로세스가 종료됩니다. 그때 찾은 쌍 수를 출력하십시오. 시간을 직접 추적하고 10 분이 경과하기 직전에 결과를 인쇄하거나 산발적으로 발견 된 쌍 수를 출력하면 점수와 같은 마지막 숫자를 사용합니다.



답변

C ++, 275,000,000+

(x, 0) 과 같은 크기 를 정직한 쌍으로 , 다른 모든 쌍을 부정직 한 크기의 m 으로 표현할 수있는 쌍을 참조합니다 . 여기서 m 은 쌍의 잘못보고 된 크기입니다. 이전 게시물 의 첫 번째 프로그램은 충분히 큰 x에 대해
(x, 0)(x, 1) 과 밀접하게 관련된 정직하고 부정직 한 쌍을 사용했습니다.. 두 번째 프로그램은 동일한 부정직 한 쌍을 사용했지만 모든 정직한 크기의 정직한 쌍을 찾아서 정직한 쌍을 확장했습니다. 이 프로그램은 10 분 안에 종료되지 않지만 대부분의 결과는 매우 초기에 발견되므로 대부분의 런타임이 낭비됩니다. 자주없는 정직한 쌍을 계속 찾는 대신이 프로그램은 여가 시간을 사용하여 다음 논리적 인 일을 수행합니다 . 부정직 한 쌍 집합 확장 .

이전 게시물에서 우리는 모든 큰 정수 r , sqrt (r 2 + 1) = r 에서 sqrt 는 부동 소수점 제곱근 함수 라는 것을 알고 있습니다. 공격 우리 계획을 찾을 것이다 쌍 P = (x, y)로 되도록 X 2 + y 2 = r에 2 + 1 일부 충분히 큰 정수 대한 R . 그것은 충분히 간단하지만, 순전히 그러한 쌍을 찾는 것은 너무 느리기 때문에 흥미 롭습니다. 이전 프로그램에서 정직한 쌍을 위해했던 것처럼이 쌍을 대량으로 찾고 싶습니다.

{ v , w }를 정규 직교 벡터 쌍으로 하자 . 모든 실제 스칼라 r , || r v + w || 2 = r 2 + 1 입니다. 에서는 2 , 이것은 피타고라스의 정리의 직접적인 결과이다 :

우리는 벡터를 찾고 V 거기에 존재하도록 정수 r에 있는 X를 하고 y는 또한 정수입니다. 참고로, 앞의 두 프로그램에서 사용한 부정직 한 쌍의 집합은이 경우의 특별한 경우였습니다. 여기서 { v , w }2 의 표준 기준입니다 . 이번에는보다 일반적인 해결책을 찾고자합니다. 피타고라스의 삼중 항 (a, b, c)이 a 2 + b 2 = c 2를 만족시키는 곳입니다 .(이전 프로그램에서 사용한)는 다시 돌아옵니다.

(a, b, c) 를 피타고라스 삼중 항으로 하자 . 벡터 v = (b / c, a / c)w = (-a / c, b / c) (및 또한
w = (a / c, -b / c) )는 검증하기 쉬운 직교 정규화입니다. . 그것이 나오는 것에 따라, 피타고라스의 삼중의 선택의 여지를 들어, 정수가 존재 r은 그런 Xy는 정수입니다. 이를 증명하고 효과적으로 rP를 찾으려면 약간의 숫자 / 그룹 이론이 필요합니다. 세부 사항을 아끼지 않겠습니다. 어느 쪽이든 정수 r , xy 가 있다고 가정하십시오 . 우리는 여전히 몇 가지 부족합니다 .r이 필요합니다.충분히 커야하고 우리는 이것으로부터 더 많은 유사한 쌍을 도출하는 빠른 방법을 원합니다. 다행히도이를 달성하는 간단한 방법이 있습니다.

참고 투영 것을 P, V는R의 V , 따라서 R = P · V = (X, Y) · (B / C, A / C) = XB / C + 나중에 / C , 모든 말 것을 XB의 + ya = rc . 결과적으로 모든 정수 n 에 대해 (x + bn) 2 + (y + an) 2 = (x 2 + y 2 ) + 2 (xb + ya) n + (a 2 + b 2 ) n 2 = ( r 2 + 1) + 2 (rc) n + (c 2 ) n 2 = (r + cn) 2 + 1. 다시 말해,
(x + bn, y + an) 형식의 쌍의 제곱 크기 는 (r + cn) 2 + 1입니다 . 이것은 우리가 찾고있는 쌍의 종류입니다! 충분히 큰 n에 대해 , 이들은 부정직 한 크기의 r + cn 쌍입니다 .

구체적인 예를 보는 것이 항상 좋습니다. 우리는 피타고라스의 삼중 걸릴 경우 (3, 4, 5) , 그 다음에 R = 2 우리가 P는 = (1, 2) (당신이 확인할 수 (1, 2) · (4/5, 3/5) = 2 그리고, 명확하게, 1 2 + 2 2 = 2 2 + 1 ). 추가 5R(4,3)P는 로 데려 간다 R ‘= 2 + 5 = 7P’= (+ 4 1 2 + 3) = (5, 5) 입니다. 보라, 5 2 + 5 2 = 7 2 + 1. 다음 좌표는 r ”= 12P ”= (9, 8) 이며, 9 2 + 8 2 = 12 2 + 1 등입니다.

일단 r은 충분히 큰이며, 우리의 크기 증가와 부정직 한 쌍을 받기 시작 5 . 약 27,797,402 / 5 부정직 한 쌍입니다.

이제 우리는 정수 크기의 부정직 한 쌍을 많이 가지고 있습니다. 우리는 그것들을 첫 프로그램의 정직한 쌍과 쉽게 결합하여 오탐 (false-positive)을 형성 할 수 있으며,주의를 기울이면 두 번째 프로그램의 정직한 쌍을 사용할 수도 있습니다. 이것이 기본적으로이 프로그램이하는 일입니다. 이전 프로그램과 마찬가지로이 프로그램은 결과의 대부분을 매우 일찍 발견합니다. 몇 초 내에 200,000,000 개의 오 탐지에 도달 한 후 상당히 느려집니다.

로 컴파일하십시오 g++ flspos.cpp -oflspos -std=c++11 -msse2 -mfpmath=sse -O3. 결과를 확인하려면 추가하십시오 -DVERIFY(특히 느려질 수 있습니다).

로 실행하십시오 flspos. 상세 모드에 대한 모든 명령 줄 인수

#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);
}

답변

파이썬, 27,797,402

바를 조금 더 높이려면 …

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

모든 67,108,864 <= x <= 94,906,265 = floor (sqrt (2 53 ))에 대해 쌍 (x, 0)(x, 1) 이 오 탐지인지 쉽게 확인할 수 있습니다.

작동 이유 : 67,108,864 = 2 26 . 따라서, 상기 범위의 모든 숫자 x0 <= x ‘<2 26에 대해 2 26 + x’ 형식 입니다. 모든 양의 e 에 대해 (x + e) 2 = x 2 + 2xe + e 2 = x 2 + 2 27 e + 2x’e + e 2 입니다. 우리가 (x + e) 2 = x 2 + 1 을 원한다면 적어도 2 27 e <= 1 , 즉 e <= 2 -27이 필요합니다

그러나 배정 밀도 부동 소수점 수의 가수가 52 비트 폭이므로, 작은 전자 되도록 X + E> X는 이고 E = 2 26-52 = 2 -26 . 즉, x 보다 큰 표현 가능한 가장 작은 숫자 는 x + 2-26 이고 sqrt (x 2 + 1) 의 결과는 최대 x + 2 -27 입니다. 기본 IEEE-754 반올림 모드는 반올림입니다. 동점-짝수-항상 x로 반올림하고 x + 2 -26으로 반올림 하지 않습니다 (동점 구분은 실제로 x = 67,108,864 에만 해당됩니다)라면 더 큰 숫자는 x에 관계없이 반올림합니다 ).


C ++, 75,000,000+

그 리콜 3 2 + 4 2 = 5 2 . 이것이 의미하는 것은 점 (4, 3) 이 원점을 중심으로하는 반경 5 의 원에 있다는 것입니다. 실제로, 모든 정수 n에 대해 , (4n, 3n) 은 반지름 5n의 원에 있습니다. 충분히 큰 n (즉, 5n> = 2 26 )의 경우이 원의 모든 점에 대해 (5n, 1) 의 위양성을 이미 알고 있습니다. 큰! 그것은 또 다른 27,797,402 / 5 무료 오 탐지 쌍입니다! 근데 왜 여기서 멈춰? (3, 4, 5) 만이 그러한 삼중 항이 아닙니다.

모든 양의 정수 세 쌍둥이를위한이 프로그램은 외모 (A, B, C) 등이 2 + B 2 = C 2 이런 식으로, 그리고 카운트 가양. 그것은에 도달 70,000,000 꽤 빨리 가양하지만 숫자가 증가 상당히으로 속도가 느려집니다.

로 컴파일하십시오 g++ flspos.cpp -oflspos -std=c++11 -msse2 -mfpmath=sse -O3. 결과를 확인하려면 추가하십시오 -DVERIFY(특히 느려질 수 있습니다).

로 실행하십시오 flspos. 상세 모드에 대한 모든 명령 줄 인수

#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

편집 : 새로운 프로그램.

오래된 메모리는 너무 많은 메모리를 사용했습니다. 이것은 해시 테이블 대신 거대한 벡터 배열을 사용하여 메모리 사용량을 절반으로 줄입니다. 또한 임의의 스레드 부스러기를 제거합니다.

   /* 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;
}

-P인수 대신 숫자로 포인트를 인쇄 하십시오 .

필자의 경우 카운팅 모드에서 2 분, 파일로 인쇄하는 데 약 5 분 (~ 4GB)이 걸리므로 I / O 제한이되지 않았습니다.

내 원래 프로그램은 깔끔했지만 10 ^ 5 정도의 결과 만 생성 할 수 있기 때문에 대부분을 삭제했습니다. 그것이 한 것은 (x ^ 2 + Ax + B, x ^ 2 + Cx + D), (x ^ 2 + ax + b, x ^ 2 + cx + d) 형식의 매개 변수를 찾는 것입니다. x, (x ^ 2 + Ax + B) ^ 2 + (x ^ 2 + Cx + D) ^ 2 = (x ^ 2 + ax + b) ^ 2 + (x ^ 2 + cx + d) ^ 2 + 1. 그러한 파라미터 세트 {a, b, c, d, A, B, C, D}를 찾았을 때, 최대 값 이하의 모든 x 값을 확인하는 과정을 진행했다. 이 프로그램의 디버그 출력을 보면서 많은 숫자를 쉽게 생성 할 수있는 매개 변수화의 매개 변수화에 대한 특정 매개 변수화를 발견했습니다. 나는 내 자신의 많은 것을 가지고 있기 때문에 엘의 숫자를 인쇄하지 않기로 선택했습니다. 바라건대 이제 누군가 우리의 숫자 세트를 모두 인쇄하지 않고 승자라고 주장합니다. 🙂

 /* 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);
}

답변

자바, Bresenham-esque 서클 스캔

경험적으로 나는 고리의 더 넓은 끝에서 시작하여 더 많은 충돌을 기대합니다. 나는 각 충돌에 대해 하나의 스캔을 수행하여 surplus사이 0r2max - r2포함 된 값을 기록 하지만이 버전보다 느린 테스트에서 약간의 개선을 기대했습니다 . 마찬가지로 int[]두 요소로 구성된 배열과 목록을 많이 만들지 않고 단일 버퍼 를 사용하려고 합니다. 성능 최적화는 실제로 이상한 짐승입니다.

쌍의 출력을 위해 명령 행 인수로 실행하고 간단한 계수를 사용하지 마십시오.

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;
    }
}

답변

자바-27,817,255

이것들의 대부분은 Ell이 보여주는 것과 동일 하고 나머지는에 기초합니다 (j,0) (k,l). 각각에 대해 j, 나는 사각형을 뒤로 걸어 나머지가 거짓 긍정인지 확인합니다. 이것은 기본적으로 전체에 25k (약 0.1 %)의 이득만으로 전체 시간을 차지 (j,0) (j,1)하지만 이득은 이득입니다.

내 컴퓨터에서 10 분 안에 완료되지만, 당신이 무엇을하는지 모르겠습니다. 이 경우 다음과 같은 이유 때문에 하지 않는 시간이 다 떨어지기 전에 완료, 그것은 크게 나쁜 점수를해야합니다. 이 경우 8 번 줄에서 제수를 조정하여 시간이 지날 수 있습니다 (이것은 단순히 각각에 대해 얼마나 멀리 뒤로 이동하는지 결정합니다 j). 다양한 제수의 점수는 다음과 같습니다.

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

각 경기에 대한 출력을 켜려면 (그리고, 신이 있다면 느립니다) 10 행과 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);
    }
}

참고로 처음 20 개의 출력 (제수 = 7의 경우 (j,0)(j,1)유형 제외 )은 다음과 같습니다.

(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)

답변

줄리아, 530 개의 오탐

다음은 참조 구현으로 볼 수있는 매우 순진한 무차별 검색입니다.

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

@printf줄 을 주석 해제하여 쌍 (및 정확한 제곱 크기)을 인쇄 할 수 있습니다 .

기본적으로 이것은 x = y = 6e7첫 번째 좌표 쌍 에 대한 검색을 시작하고 x를 감소시키기 전에 x 축으로가는 길의 약 1 %를 스캔합니다. 그런 다음 이러한 각 좌표 쌍에 대해 충돌에 대해 동일한 크기 (반올림 및 내림)의 전체 호를 확인합니다.

이 코드는 64 비트 시스템에서 실행되는 것으로 가정하므로 기본 정수 및 부동 소수점 유형은 64 비트 유형입니다 (그렇지 않은 경우 int64()float64()생성자를 사용하여 작성할 수 있음 ).

그 결과 빈약 한 530 결과가 나옵니다.


답변