64bit-lcg-crack4.cpp

#include <cstdio>
#include <cstdint>
#include <cinttypes>
#include <cmath>
#include <algorithm>
#include <chrono>
#include <random>

typedef int64_t s64;
typedef uint64_t u64;
typedef int64_t s64;
typedef __int128 s128;

const u64 A = 0x5d588b656c078965ll;
const u64 B = 0x269ec3ll;
const s128 M = (s128)1 << 64;
const int BITS = 64;
const int P = 8;

class LCGOperator {
public:
    u64 a, b;
    LCGOperator(u64 aa, u64 bb) : a(aa), b(bb) {}

    LCGOperator o(LCGOperator y) {
        LCGOperator x = *this;
        LCGOperator res(x.a * y.a, x.a * y.b + x.b);
        return res;
    }

    LCGOperator pow(u64 n) {
        LCGOperator x = *this;
        if (n == 0) {
            return LCGOperator(1, 0);
        } else if (n % 2 == 1) {
            return x.o(x).pow(n / 2).o(x);
        } else {
            return x.o(x).pow(n / 2);
        }
    }

    u64 apply(u64 s) {
        return a * s + b;
    }
};

bool intersect(s128 a1, s128 b1, s128 a2, s128 b2) {
    if (b1 < a1 || b2 < a2) return false;
    if (a1 < a2) {
        return b1 < a2;
    } else {
        return a1 <= b2;
    }
}

void search(u64 s0, u64 sy, u64 sz) {
    printf("%lx %lx %lx\n", s0, sy, sz);
    LCGOperator op0(A, B);
    LCGOperator opy = op0.pow(74);
    LCGOperator opz = op0.pow(79);
    const u64 xstep = 0xc67e8b67;
    const u64 ystep = (xstep * opy.a) % M;
    const u64 zstep = (xstep * opz.a) % M;
    printf("%016lx %016lx\n", ystep, zstep);

    const s128 xstart = (s128)s0 << (BITS - P);
    const s128 xend = (s128)(s0 + 1) << (BITS - P);
    const s128 ystart = (s128)sy << (BITS - P);
    const s128 yend = (s128)(sy + 1) << (BITS - P);
    const s128 zstart = (s128)sz << (BITS - P);
    const s128 zend = (s128)(sz + 1) << (BITS - P);
    u64 found_count = 0;

    for (u64 i = 0; i < xstep; i ++) {
        s128 x = xstart + i;
        s128 y = opy.apply(x);
        s128 z = opz.apply(x);
        if (i%100000000==0) printf("i=%ld\n", i);

        while (x < xend) {
            while (ystart <= y && y < yend && zstart <= z && z < zend && x < xend) {
                if ((y % M) >> (BITS - P) == sy && (z % M) >> (BITS - P) == sz) {
                    //printf("found %016lx\n", (u64)x);
                    found_count ++;
                }
                x += xstep;
                y = (y + ystep) % M;
                z = (z + zstep) % M;
            }
            s128 ry = (yend - ystart) / 2 + 1;
            s128 rz = (zend - zstart) / 2 + 1;
            s128 y0 = (ystart + yend - 2 * y) / 2;
            s128 z0 = (zstart + zend - 2 * z) / 2;
            s128 ay = std::max((ystart - y) / ystep, (s128)0);
            s128 by = (yend - y) / ystep + ((yend - y) % ystep ? 1 : 0);
            s128 az = std::max((zstart - z) / zstep, (s128)0);
            s128 bz = (zend - z) / zstep + ((zend - z) % zstep ? 1 : 0);
            s128 s = std::max(ay, az);
            if (intersect(ay, by, az, bz)) {
                //printf("case 1 : s = %" PRId64 "\n", (s64)s);
            } else {
                s = (2 * M - (y0 + z0 + ry + rz)) / (ystep + zstep);
                //printf("case 2 : s = %" PRId64 "\n", (s64)s);
            }
            if (s <= 0) s = 1;
            //if (s < 1000) printf("s = %ld\n", (s64)s);
            x += xstep * s;
            y = (y + ystep * s) % M;
            z = (z + zstep * s) % M;
        }
    }
    printf("found_count=%ld\n", found_count);
}

int main() {
    std::mt19937 rnd;
    for (int i = 0; i < 1; i ++) {
        u64 seed = (u64)rnd() << 32 | rnd();
        printf("seed=%016lx\n", seed);
        LCGOperator op0(A, B);
        LCGOperator opy = op0.pow(74);
        LCGOperator opz = op0.pow(79);
        u64 s0 = seed >> (BITS - P);
        u64 sy = opy.apply(seed) >> (BITS - P);
        u64 sz = opz.apply(seed) >> (BITS - P);
        auto start = std::chrono::system_clock::now();
        search(s0, sy, sz);
        auto end = std::chrono::system_clock::now();
        printf("time=%ld\n", std::chrono::duration_cast<std::chrono::milliseconds>(end-start).count());
    }
    return 0;
}