More templated version

[Craig Black]

Here' a version that uses more templates.  It performs the same as the first
version I submitted.

import std.c.stdio;
import std.c.stdlib;
import std.c.time;

const int LONG_TIME=4000;

byte[] p;
byte[] t;
int q;

int main(char[][] args)
{
 if (args.length == 2) {
  sscanf(&args[1][0],"%d",&q);
 } else {
  printf("Usage: pi [precision]\n");
  exit(55);
 }

 if (q < 0)
 {
  printf("Precision was too low, running with precision of 0.\n");
  q = 0;
 }

 if (q > LONG_TIME)
 {
     printf("Be prepared to wait a while...\n");
 }

 // Compute one more digit than we display to compensate for rounding
 q++;

 p.length = q + 1;
 t.length = q + 1;

 // Compute pi
 time_t startime, endtime;
 std.c.time.time(&startime);
 arctan!(2);
 arctan!(3);
 fastmul!(4);
 std.c.time.time(&endtime);

 // Return to the number of digits we want to display
 q--;

 // Print pi
 printf("pi = %d.",cast(int)(p[0]));
 for (int i = 1; i <= q; i++)
  printf("%d",cast(int)(p[i]));
 printf("\n");
 printf("%ld seconds to compute pi with a precision of %d
digits.\n",endtime-startime,q);

 return 0;
}

/* Template that determines if a number is a power of
 * 2 at compile time
 */
template isPow2(uint n) { const bool isPow2 = (n % 2) != 0 ? false :
isPow2!(n / 2); }
template isPow2(uint n : 0) { const bool isPow2 = false; }
template isPow2(uint n : 1) { const bool isPow2 = false; }
template isPow2(uint n : 2) { const bool isPow2 = true; }

/* Use a template to get good compile time optimizations for
 * division and multiplication
 */
template fastdiv(int divisor)
{
 void fastdiv()
 {
  int r; // remainder


  for (int i = 0; i <= q; i++) {
   int b = t[i] + r * 10;
   int q = b / divisor;
   static if(isPow2!(divisor))
    r = b % divisor; // compiler can optimize modulus for powers of 2
   else
    r = b - q * divisor;
   t[i] = q;
  }
 }
}

/* Template to compute the arctangent.
 * The constant input parameter allows the compiler to optimize.
 */
template arctan(int s)
{
 void arctan()
 {
  int n;

  t[0] = 1;
  fastdiv!(s);
  add();
  n = 1;
  do {
   mul(n);
   fastdiv!(s*s);
   div(n += 2);
   if ((((n-1) / 2) % 2) == 0)
    add();
   else
    sub();
  } while (!tiszero());
 }
}

void add()
{
 for (int j = q; j >= 0; j--)
 {
  if (t[j] + p[j] > 9) {
   p[j] += t[j] - 10;
   p[j-1]++;
  } else
   p[j] += t[j];
 }
}

void sub()
{
 for (int j = q; j >= 0; j--)
  if (p[j] < t[j]) {
   p[j] -= t[j] - 10;
   p[j-1]--;
  } else
   p[j] -= t[j];
}


void mul(int multiplier)
{
 int c;

 for (int i = q; i >= 0; i--) {
  int b = t[i] * multiplier + c;
  c = b / 10;
  t[i] = b - c * 10;
 }
}

/* Fast multiplication using template */
template fastmul(int multiplier)
{
 void fastmul()
 {
  int c;

  for (int i = q; i >= 0; i--) {
   int b = p[i] * multiplier + c;
   c = b / 10;
   p[i] = b - c * 10;
  }
 }
}


void div(int divisor)
{
 int r; // remainder

        /* Optimization: It's faster to do floatint point multiplication
than
         * integer division. This is because there is no integer division
instruction.
         * Under the hood, integer division is essentially floating-point
division.
         */
        double idiv = 1.0 / divisor;

 for (int i = 0; i <= q; i++) {
                int b = t[i] + r * 10;
  int quotient = cast(int)(cast(double)b * idiv);
  r = b - quotient * divisor;
  t[i] = quotient;
 }
}

int tiszero()
{
 for (int k = 0; k <= q; k++)
  if (t[k] != 0)
   return false;
 return true;
}