// CDC 6600 and IBM Stretch multiply. <- multiply.html
#include <stdio.h>
#include <assert.h>
typedef long int L;
typedef unsigned int ui;
typedef struct{L hi; L lo;} rp;
#define C const
char p = 1;
long unsigned int dump;
void discard(ui x){assert(x<22); dump = (dump>>3) + ((L)x<<(64-5));}
rp cook(L x, L y, L z){return (rp){x&y|x&z|y&z, x^y^z};}
typedef unsigned long int ul;
typedef struct{ul q; ul r;} res;
res divq(res, ul);
res mulq(ul, ul);
res addq(res, res);
void prr(char * m, res x){
   printf("[%016lx %016lx] = %s\n", x.q, x.r, m);}
res llsrs(res v, int n){return n<64?(res){(L)v.q>>n, (n?(v.q<<(64-n)):0)|(v.r>>n)}
          :(res){(L)v.q>>63, (L)v.q>>(n-64)};}
int eqr(res a, res b){return a.q==b.q && a.r==b.r;}
res mul (L C a, L b){L C th = a+(a<<1); ui carry = 0; dump = 0;
  if(p) {res C pr = mulq(a, b);
    printf("%016lx*%016lx=%016lx%016lx\n", a, b, pr.q, pr.r);}
  L sel(ui n){assert(n<8); n += carry; carry = (n>3);
    return carry?~(n>6?(n==7?a:0):(n==5?th:(a<<(n==4?2:1)))):
                 (n>2?th:!n?0:(n==1?a:a<<1));}
  {rp sum = {0, 0}; int k; for(k=0; k<48; k+=3) {
   {L C z = sel(b&7); sum = cook((sum.hi<<1)|carry, sum.lo, z);}
   b>>=3;
   {discard(((7&sum.hi)<<1)+(7&sum.lo));
     sum = (rp){sum.hi>>3, sum.lo>>3};}}
   if(carry) sum = cook(sum.hi<<1, sum.lo, a);
   return addq(
         llsrs((res){(sum.hi<<1)+sum.lo, 0}, 16),
         llsrs((res){0, dump}, 14));}}

#include <stdlib.h>
void tst(L a, L b){res hp = mul(a, b); res qp = mulq(a, b);
  if(!eqr(qp, hp)) {printf("Bad!\n");
    prr("hp", hp); prr("qp", qp); exit(0);}}
L rL(){return ((L)random()<<17)^random();}
int main(){tst(0, 4);
  {int j=64; tst(0, 7237578L); while(j--) tst (1, j<<3);}
  {tst(0, 234587237578L); tst(1, 4); tst(0, 63);
    tst(1L<<47, 1L<<2); tst(1L<<4, 7L<<44); tst(1L<<47, 7L<<44);}
  {int j = 100000000; while(j--){
     p = !(((1<<18)-1)&random()); tst(rL(), rL());}}
  {int j = 48; while(j--) {int k = 48; while(k--) tst(1L<<j, 1L<<k);}}
printf("Good\n");
return 0;}