// gcd128.c
// read and caclulate the GCD of two 128-bit unsigned integers
//      (<=38 decimal digits)
//Uses a C structure of an array of 4 ints. (assumes int is 32 bits long)
//

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

typedef struct {
    unsigned int i[4];
    } int128;

void ODD128 ( int128 *n);
int GC2int128 ( int128 *n, int128 *p);
void GCD2ODD128( int128 *n, int128 *p);
int128 * ptrmax128 ( int128 *x, int128 *y );
int128 * ptrmin128 ( int128 *x, int128 *y );
void atoi128( int128 *n, char *s);
void i128toa (char *s, int128 *n);

int main(int argc, char *args[])
{
    int128 *P, *Q;
    int  R;
    char Pstr[40], Qstr[40];

    P = (int128 * ) malloc( sizeof( int128 ) );
    Q = (int128 * ) malloc( sizeof( int128 ) );

    printf ("Enter two large integers (up to 38 digits in decimal\n");

    scanf("%s %s", Pstr, Qstr);
    atoi128( P, Pstr);
    atoi128( Q, Qstr);
    i128toa(Pstr,P);
    i128toa(Qstr,Q);
    printf(" The GCD of \n  %s  and\n  %s \n is calculated as follows:\n", Pstr, Qstr);

    R = GC2int128(P, Q);
    printf(" The greatest power of 2 dividing both numbers is 2^%d\n", R);

    ODD128(P);
    ODD128(Q);
    i128toa(Pstr, P);
    i128toa(Qstr, Q);
    printf(" The remaining odd part of P is %s\n", Pstr);
    printf(" The remaining odd part of Q is %s\n", Qstr);

    GCD2ODD128(P,Q);
    i128toa( Pstr, P);
    printf(" The GCD of these odd parts is %s \n", Pstr);
    printf(" The GCD of the original pair is 2^%d * %s \n", R, Pstr);

    return 0;
}
// GC2int128: repeatedly divides two int128's by 2 if BOTH are even
//  until at least one is odd.

int GC2int128 ( int128 *n, int128 *p)
{   int e=0;
    while ( !(n->i[0]&1) &&  !(p->i[0]&1)  ) {
        //n->i[0]>>=1;
        //p->i[0]>>=1;
        //e++;
        _asm{
            push esi
            push edi
            mov esi, n
            mov edi, p
            again:
                test [esi], 1
                jnz done
                test [edi], 1
                jnz done
                shr DWORD PTR 12[esi], 1
                ror DWORD PTR 8[esi], 1
                ror DWORD PTR 4[esi], 1
                ror DWORD PTR 0[esi], 1

                shr DWORD PTR 12[edi], 1
                ror DWORD PTR 8[edi], 1
                ror DWORD PTR 4[edi], 1
                ror DWORD PTR 0[edi], 1

                inc e
                jmp again

            done:
                pop edi
                pop esi

        }

    }
    return (e);
}
// ODD128: repeatedly divides an int128 by 2 until it is odd.
void ODD128 ( int128 *n)
{
    while( !(n->i[0]&1) ) n->i[0]>>=1;
    return;
}
//GCD2DOO128: repeatedly subtracts the smaller odd int128 from the larger one
//  and calls ODD128 to make them both odd. Stops when they are equal.
void GCD2ODD128( int128 *n, int128 *p)
{
    int128 *lg, *sm;
    while (n->i[0] != p->i[0])
    {
        lg =  ptrmax128(n,p);
        sm = ptrmin128(n,p);
         lg->i[0]  -= sm->i[0];
         p = lg;
         n = sm;
         ODD128 (p);
        }
    return;
}
//ptrmax128: given two pointers to int128's, returns the pointer to the larger int128.
int128 *ptrmax128(int128 *x, int128 *y )
{
    return (x->i[0]) >=( y->i[0]) ? x: y;
}

//ptrmin128: given two pointers to int128's, returns the pointer to the smaller int128.
int128 *ptrmin128 (int128 *x, int128 *y )
{
    return  (x->i[0] )<= ( y->i[0] )? x: y;
}
//atoi128: converts a numeric string to int128.
void atoi128( int128 *n, char *s)
{
    n->i[0] = atoi(s);
    n->i[1] = 0;
    n->i[2] = 0;
    n->i[3] = 0;
    return;
}
// itoa128: converts a 128-bit integer to a string
void i128toa (char *s, int128 * n)
{
    sprintf(s,"%u", n->i[0]);
    return;
}