/* It turned up up just at the right moment: 
http://www.cs.mtu.edu/~shene/COURSES/cs201/NOTES/chap04/gcd.html 
 
PROGRAM  GreatestCommonDivisor
   IMPLICIT  NONE
 
   INTEGER   :: a, b, c
 
   WRITE(*,*) 'Two positive integers please --> '
   READ(*,*)  a, b
   IF (a < b) THEN       ! since a >= b must be true, they
      c = a              ! are swapped if a < b
      a = b
      b = c
   END IF
 
   DO                    ! now we have a <= b
      c = MOD(a, b)      !    compute c, the reminder
      IF (c == 0) EXIT   !    if c is zero, we are done.  GCD = b
      a = b              !    otherwise, b becomes a
      b = c              !    and c becomes b
   END DO                !    go back
 
   WRITE(*,*) 'The GCD is ', b
 
END PROGRAM  GreatestCommonDivisor
 */
int GreatestCommonDivisor(int a, int b)
{
    int c;
 
    if (a < 0) a = -a;
    if (b < 0) b = -b;
 
    if (a < b)
    {
        c = a; a = b; b = c;
    }
    if (b == 0)
       return a;
    for (;;)
   {
        c = a % b;
        if (c == 0)
            break;
        a = b;
        b = c;
    }
    return b;
}
 
char *double2fraction(char* buf, int bufsize, double x, int maxdigits)
{
    double intpart;
    double f;
    int negative = 0;
    int i, gcd, numerator, divisor;
    char temp[64];
    int templen;
 
    if (!buf || bufsize < 2)
        return 0;
    else if (x == 0.0)
    {
        buf[0] = '0'; buf[1] = '\0';
    }
 
    if (x < 0.0)
    {
        x = -x; negative = 1;
    }
 
    if (maxdigits <= 0)
        maxdigits = 5;
    else if (maxdigits > 9)
        maxdigits = 9;
 
    divisor = 1;
    for (i = 0; i < maxdigits; ++i)
        divisor *= 10;
    f = modf(x,&intpart);
 
    numerator = (int)(f*divisor+0.5);
 
    if (numerator != 0)
    {
        while ((gcd = GreatestCommonDivisor(divisor,numerator)) > 1)
        {
            divisor   /= gcd;
            numerator /= gcd;
        }
        sprintf(temp,"%s%.0f/%d",(negative?"-":""),intpart*divisor+numerator,divisor);
    }
    else
        sprintf(temp,"%.0f",intpart);
    templen = strlen(temp);
    if (templen+1 < bufsize)
        strcpy(buf,temp);
    else
    {
        buf[0] = '*'; buf[1] = '\0';
    }
    return buf;
}

int greatestCommonDivisor(int a, int b)
{
    int t;
 
    if (a < 0) a = -a;
    if (b < 0) b = -b;
    if (a < b)
    {
        t = a; a = b; b = t;
    }
   while (b)
   {
        t = b;
        b = a % b;
        a = t;
    }
    return a;
}
 
char *double2fraction(char* buf, int bufsize, double x, int maxdigits)
{
    double intpart;
    double f;
    int negative = 0;
    int i, gcd, numerator, divisor;
    char temp[64];
    int templen;
 
    if (!buf || bufsize < 2)
        return 0;
 
    if (x == 0.0)
   {
        buf[0] = '0'; buf[1] = '\0';
        return buf;
    }
 
    if (x < 0.0)
    {
        x = -x; negative = 1;
    }
 
    if (maxdigits <= 0)
        maxdigits = 5;
    else if (maxdigits > 9)
        maxdigits = 9;
 
    divisor = 1;
    for (i = 0; i < maxdigits; ++i)
        divisor *= 10;
 
    f = modf(x,&intpart);
    numerator = (int)(f*divisor+0.5);
 
    if (numerator != 0)
    {
        while ((gcd = greatestCommonDivisor(divisor,numerator)) > 1)
        {
            divisor   /= gcd;
            numerator /= gcd;
        }
        sprintf(temp,"%s%.0f/%d",(negative?"-":""),
        intpart*divisor+numerator,divisor);
    }
    else
        sprintf(temp,"%.0f",intpart);
 
    templen = strlen(temp);
    if (templen+1 < bufsize)
        strcpy(buf,temp);
    else
    {
        buf[0] = '*'; buf[1] = '\0';
    }
    return buf;
}