#include <assert.h>
#include <stdio.h>
#include <stdlib.h>
typedef long unsigned int uL;
typedef struct{uL q; uL r;} res;
res divq(res, uL);
res mulq(uL, uL);
res addq(res, res);
int min(int a, int b){return a<b?a:b;}
int max(int a, int b){return a<b?b:a;}
// llls long long left shift
res llls(res v, int n){return n<64?(res){(v.q<<n)|(n?(v.r>>64-n):0),v.r<<n}
     :(res){v.r<<(n-64), 0};}
int lleq(res a, res b){return a.q==b.q && a.r==b.r;}
res coer(uL a){return (res){0, a};}

int b; //calculated eventually ...
int const rs = 10; // number of bits in reciprocal approximation
int const lts = 9, ts = 1<<lts, bi = 8;
// reciprocal size and log table size.
short tab[ts];
// dv returns the quotient using the proposed scheme.
// the following variables would neither global nor mutable
// if clang had nested fuctions!!
uL delta;

uL rd(uL nnp, uL quo, int qb){
   if(!qb) return nnp+delta>=(1L<<(52+rs))?quo+1:quo;
   b = min(b, qb); // crude and safe estimate
   int splitpoint = 52+rs-b;
   int td = nnp>>splitpoint;
   uL nx = nnp - ((uL)td<<splitpoint); // split nnp.
   // below: if we were dividing nnp by nd would we be getting closer?
return rd((nx<<b)+td*delta, (quo<<b)+td, qb-b);}

uL dv(uL const num, uL const den){
  int big = num>=den; // not used soon; part of result
  b = bi;
  int N = tab[(den>>(64-(lts+1)))-ts];
  uL nn = N*(num>>(64-53)), nd = N*(den>>(64-53));
  if(0) printf("nd=%016lx\n", nd); assert(nd< 0x3fffffffffffffffL);
  delta = (1L<<(52+rs))-nd;
  return rd(nn, 0, 53-big);}

uL dvg(uL const num, uL const den){
  return divq(llls(coer(num>>11), 52+(num<den)), den>>11).q;}
void tsd(uL n, uL d){
  uL qg= dvg(n, d), qt = dv(n, d);
  assert(qt == qg);}

uL r52(){return (random()^((((uL)random())<<21) | (1L<<52)))<<11;}
int main (){
{int const mg = lts+rs; int j=ts; int y = 1<<mg, z=-1<<mg;
  while(j--) {int r = (1<<mg)/(ts+1+j);
    y = min(y, (ts+j)*r); z = max(z, (ts+1+j)*r);
    tab[j] = r; if(1) printf("%4d %4d, %5x\n", j, r, r*(ts+j));}
  printf("epsilon = %d = 0x%x; mx = %d = 0x%x\n", y, y, z, z);}
  if(1) {for(int i=0; i<ts; ++i) for(int j=0; j<ts; ++j) {
    uL n = (uL)i<<(64-(lts+1)), k=1L<<63,
       d = (uL)j<<(64-(lts+1)), m = (1L<<(64-(lts+1)))-(1<<11);
    tsd(k+n, k+d);   tsd(k+n+m, k+d);
    tsd(k+n, k+d+m); tsd(k+n+m, k+d+m);}
    printf("Good Patterns\n");}
{int j = 10000000; while(j--) tsd(r52(), r52());}
printf("Good Random\n"); return 0;}