#include <stdio.h>
#include <stdint.h>
#include <stdlib.h>

/* http://marc-b-reynolds.github.io/math/2017/09/18/ModInverse.html */
static uint32_t modinv(uint32_t a)
{
    uint32_t x;
    x = a;      //  3 bits
    x *= 2 - a*x;   //  6
    x *= 2 - a*x;   // 12
    x *= 2 - a*x;   // 24
    x *= 2 - a*x;   // 48 / retaining bottom 32
    return x;
}

/* http://number-none.com/blow/blog/programming/2016/07/08/fabian-on-lcg-fast-forward.html */
uint32_t skip(uint32_t x, uint32_t a, uint32_t c, unsigned n)
{
    while (n) {
            if (n & 1)
                x = a * x + c;
        c *= a + 1;
        a *= a;
        n >>= 1;
    }
    return x;
}

int main(int argc, char **argv)
{
    srandom(time(NULL));
    unsigned i, s = 1, a = random() | 1UL, c =  random() | 1UL;

    /* move forward one step at a time */
    for (i = 0; i < 100; i++)
        s = s * a + c;
    printf("%08x\n", s);
    /* skip forwaward */
    printf("%08x\n", skip(1, a, c, 100));
    /* you could loop through the full LCG cycle */
    printf("%08x\n", skip(s, a, c, -100));
    /* or inverse of LCG is still LCG */
    a = modinv(a);
    c = -c * a;
    /* skip backward */
    printf("%08x\n", skip(s, a, c, 100));
    /* move backward one step at a time */
    for (i = 0; i < 100; i++)
        s = s * a + c;
    printf("%08x\n", s);
}