use strict;
use warnings;
use integer;

# PRNG parameters
my $m  = 0x343fd;
my $b  = 0x269ec3;
my $m1 = 0x39b33155; # multiplicative inverse of $m mod 2**31

# build the lookup tables
my (@diff_tbl, @prod_tbl);
for my $x (0 .. 0x2bc) {  # 0x2bc is a magic number... see below
  my $p = ($x * $m1) & 0x7fffffff;
  my $e = ($p >> 22) & 0x1ff;
  for ($e, ($e-1)&0x1ff, ($e-2)&0x1ff) {
    push @{$diff_tbl[$_]}, $x;
    push @{$prod_tbl[$_]}, $p;
  }
}

# "Unknown" values: internal states of the PRNG
my $X0 = int(rand 2**31);
my $X1 = ($X0 * $m + $b) & 0x7fffffff;
my $X2 = ($X1 * $m + $b) & 0x7fffffff;
my $X3 = ($X2 * $m + $b) & 0x7fffffff;

# "Known" values: a few bytes of output from the PRNG
my $x0 = $X0 & 0x7f800000;
my $x1 = $X1 & 0x7f800000;
my $x2 = $X2 & 0x7f800000;
my $x3 = $X3 & 0x7f800000;

# find the lowest nonnegative $d such that
#     $x0 == ((($x1 + $d - $b) * $m1) & 0x7f800000)
my $prev = (($x1 - $b) * $m1) & 0x7fffffff;
my $ent = (($x0 - $prev) >> 22) & 0x1ff;
my ($d, $p);
for (0 .. $#{$prod_tbl[$ent]}) {
  $p = $prev + $prod_tbl[$ent][$_];
  if (($p & 0x7f800000) == $x0) { $d = $diff_tbl[$ent][$_]; last }
}
die "can't happen\n" unless defined $d;
$p &= 0x7fffff;
my $q = ($x1 + $d) * $m + $b;

# modest (about 2**15 steps) brute-force search for the right $d
while ($d < 0x800000) {
  # invariant: $p == ((($x1 + $d - $b) * $m1) & 0x7fffff)
  # invariant: $q == ($x1 + $d) * $m + $b

  # is this a solution?
  if (($q & 0x7f800000) == $x2) {
    my $r = $q * $m + $b;
    if (($r & 0x7f800000) == $x3) {
      printf "guessed x0 = %08x\n", (($x1 + $d - $b) * $m1) & 0x7fffffff;
    }
  }

  # find the next admissible value of $d
  # there is much magic here... see below for a hint
  if    ($p <  0x67ceeb) { $d += 0x0c1; $p += 0x183115; $q += 0x27641bd }
  elsif ($p >= 0x6a1b54) { $d += 0x1fc; $p -= 0x6a1b54; $q += 0x67aea0c }
  else                   { $d += 0x2bd; $p -= 0x51ea3f; $q += 0x8f12bc9 }
}

printf "actual  X0 = %08x\n", $X0;

# calculate the magic numbers
my $low = 0;
my $high = 0x800000;
my $x = 0;
while ($low < $high) {
  $x++;
  my $y = ($x * $m1) & 0x7fffffff;
  my $yhi = ($y >> 23) & 0xff;
  my $yneg = (-$y) & 0x7fffff;
  if ($yhi == 0 && $yneg > $low) {
    printf "if y <  0x%x then x += 0x%x, y += 0x%x\n", $yneg, $x, $y;
    $low = $yneg;
  }
  elsif ($yhi == 0xff && $yneg < $high) {
    printf "if y >= 0x%x then x += 0x%x, y -= 0x%x\n", $yneg, $x, $yneg;
    $high = $yneg;
  }
}