/* ANSI-C(?)によるπ計算   VERSION 2.0 */
/* Copyright(c) Daisuke Takahashi 1991 */

/* 任意の瞬間にキー入力によって計算を中断したいのですが */
/* High-Cで signal.h が使えないので、一定時間毎にセーブさせています */
/* その代わり、DOS.H未使用の為、殆どのＣでコンパイルできます(^_^) */

#include <stdio.h>
#include <stdlib.h>
#include <ctype.h>
#include <time.h>
#include <string.h>

/* ----------- definitions ------------ */
#define TEN8   100000000L  /* １億 10^8 */
#define MAX    4190000000L /* max of USL - TEN8 */
#define ON     (1)
#define OFF    (0)
#define SET    (2)
#define LAP    (1)
#define DONE   (0)
#define DATA   (1)
#define PLUS   (1)
#define MINUS  (-1)
#define NUL    (\0)
#define USL    unsigned long     /* typedef にすべきでしょうが(^_^;) */
#define RATE   (int)( 1000 * (long)ex.head / (long)ex.max ) 

/* ---------------- prototypes ------------------ */
  void set_ex(int,char**);              /* 引数→変数 */
  void get_memory(void);
  void back_memory(void);               /* ワーク返還 */
  void formal(void);                    /* 計算公式   */
  void arctan(int,int,int);             /* arctan設定 */
  long bench_mark(long);                /* 速度測定   */
  void message(int,int,int);            /* 計算量表示 */
  void c_std(int,USL,USL,long*,long*);  /* 41以下     */
  void c_low(int,USL,USL,long*,long*);  /* 6.4万以下  */
  void c_mid(int,USL,USL,long*,long*);  /* 200万以下  */
  void c_high(int,USL,USL,long*,long*); /* 256万以下  */
  void c_huge(int,USL,USL,long*,long*); /* 1600万以下 */
  void regular(int);                    /* 正規化     */
  void load(char *);                    /* 経過ロード */
  void save(void);                      /* 経過セーブ */
  void show_ans(void);                  /* 結果の表示 */
  void time_set(void);                  /* 開始時刻記録 */
  void time_end(void);                  /* 計算時間更新 */
  void usage(int);                      /* 使い方 */

/* ------------- structs --------- */
struct {
  /* ファイルセーブする変数 */
  int  ver;   /* Version  */
  long sec;   /* 計算時間 */
  long cent;  /* 1/100秒  */
  long keta;  /* 計算桁数 */
  int  max;   /* 配列数   */
  int  type;  /* 計算公式 */
  int  head;  /* 動的先頭 */
  int  tail;  /* 動的最終 */
  USL  save;  /* 退避flag */
  USL  con;   /*          */
  USL  var;   /*          */
  /* セーブ必要なし */
  long *arc;      /* 計算作業領域 */
  long *ans;      /* 計算結果領域 */
  int  load;      /* file読込済   */
  time_t  s_sec;  /* 計算開始時刻 */
  clock_t s_cent; /* 〃1/100 秒   */
} ex;

int main( int argc, char **argv ) {
  set_ex( argc, argv );  /* 引き数から変数設定 */

  formal();      /* 計算 */

  show_ans();    /* 結果と時間を表示 */
  back_memory(); /* ワーク返還 */
  return 0;      /* 正常終了   */
}

void set_ex( int argc, char **argv ) {
  char *buf;

  if (argc == 1)                       /* オプション無しなら */
    usage( !DATA );                  /* 普通説明表示       */
  if (strchr( argv[1], '.' ) != NULL)  /* ファイル名であれば */
    load( argv[1] );                 /* 途中結果を読み込む */
  else if ( !(isdigit( *argv[1] ) ) )  /* 無効オプションなら */
    usage( DATA );                   /* データ説明         */
  else {            /* 桁数指定オプション  */
    ex.ver  = 0;  /* On memory Version   */
    ex.sec  = 0;  /* 計算時間(秒単位)    */
    ex.cent = 0;  /* 計算時間(0.01秒単位 */
    ex.keta = atol( argv[1] );         /* 計算桁を代入 */
    ex.max  = (int)(ex.keta / 8 + 2);  /* 配列数 */
    ex.save = 0;  /* セーブモードでない */
    get_memory(); /* 必要メモリを確保   */
  }
  if (argc == 3) {
    buf = argv[2];
    while ( *buf ) {
      if ( isdigit(*buf) )                    /* 数字の場合は */
        ex.save = strtol( buf, &buf, 10 );  /* 分単位のセーブ時間 */
      else {
        if (strchr( "gmrsGMRS", *buf ) != NULL && ex.load == OFF)
          ex.type = *buf;  /* ロード後の公式変更は不可 */
        buf++;               /* 計算公式を代入 */
      }
    }
  }
}

void get_memory( void ) {      /* ワークを確保 */
  ex.arc = (long *)calloc( (size_t)ex.max, (size_t)sizeof(long) ); 
  ex.ans = (long *)calloc( (size_t)ex.max, (size_t)sizeof(long) );
  if (ex.ans == NULL) {
    fputs( "メモリーが足らんぜよ(;_;) \n", stderr );
    exit( 1 );
  }
}

void back_memory( void ) {   /* ワークを返還 */
  free( ex.arc );
  free( ex.ans );
}

void formal( void ) {    /* 計算公式の値をセット */
  switch ( tolower( ex.type ) ) {

    case 'g': /* ガウスの式 */
      /* 48*arctan(1/18) + 32*arctan(1/57) - 20*arctan(1/239) */
      arctan( PLUS, 48, 18 );
      arctan( PLUS, 32, 57 );
      arctan( MINUS, 20, 239 );
      break;
    case 'r': /* ラザフォードの式 */
      /* 16*arctan(1/5) - 4*arctan(1/70) + 4*arctan(1/99) */
      arctan( PLUS, 16, 5 );
      arctan( MINUS, 4, 70 );
      arctan( PLUS, 4, 99 );
      break;
    case 's': /* シュテルマーの式 */
      /* 24*arctan(1/8) + 8*arctan(1/57) + 4*arctan(1/239) */
      arctan( PLUS, 24, 8 );
      arctan( PLUS, 8, 57 );
      arctan( PLUS, 4, 239 );
      break;
    default:  /* デフォルトは最速マチン式 */
    case 'm': /* マチンの式 */
      /* 16*arctan(1/5) - 4*arctan(1/239) */
      arctan( PLUS, 16, 5 );
      arctan( MINUS, 4, 239 );
      break;
  }
  if (ex.load == ON) 
    fputs( "Not calculated !!!\n", stderr );
}

void arctan( int sign, int mag, int base ) {
  /* arctan(1/base) = 1/base - 1/(3*base^3) + 1/(5*base^5) - */
  time_t cur_time, pre_time;
  long min_one;     /* １分相当の var */

  message( SET, sign, SET );  /* 3番目はダミー：計算式の符号 */
  if (ex.load == OFF ) {    /* 新規計算開始 */
    *ex.arc = (long)base * (long)mag;  /* 初期値設定 */
    ex.con  = (long)base * (long)base;
    ex.head = 0;
    ex.var  = 1;
    if ( base == 5 || base == 8 )  /* 10進数で有限小数なら */
      ex.tail = 1;                 /* 最終位置は動的に変化 */
    else                 /* 10進数で循環小数なら最終位置は物理的最終 */
      ex.tail = ex.max;
  } else {                  /* 計算再開 */
    if (ex.con != (long)base * (long)base) {
      message( DONE, mag, base );
      return;  /* この項は計算済 */
    } else 
      ex.load = OFF;   /* この項から計算再開：ロードフラグクリア */
  }
  message( LAP, mag, base );

  time_set();   /* ここから計算時間測定開始 */
  time( &pre_time );  /* 表示間隔初期値 */
  min_one = bench_mark( 0 );  /* 時間測定初期化 */

  while ( (ex.head < ex.max) || *(ex.arc + ex.max) ) {
    USL    check;       /* どのcalc_*_? を使うかの判定 */

    if ( ex.var > 64000L && (TEN8 % ex.var) < (MAX / ex.var) ) { 
      check = 60000L;     /* enough under 64000L */
      /* ans_r * T8_v_r < 41億 なら calc_low() で計算 */
    } else 
      check = ex.var;
    if      (check < 42)
      c_std( sign, ex.con, ex.var, ex.arc, ex.ans ); 
    else if (check < 64000L) 
      c_low( sign, ex.con, ex.var, ex.arc, ex.ans );
    else if (check < 2030000L)
      c_mid( sign, ex.con, ex.var, ex.arc, ex.ans );
    else if (check < 2560000L)
      c_high( sign, ex.con, ex.var, ex.arc, ex.ans );
    else /*  check < 1600万 */
      c_huge( sign, ex.con, ex.var, ex.arc, ex.ans );

    sign *= MINUS; /* 符号逆転 */
      if ( ex.tail < ex.max && *(ex.arc + ex.tail) )
        ex.tail++;
      if ( !*(ex.arc + ex.head) && ex.head < ex.max ) 
        ex.head++;  /* 有効数字の先頭に移動 */
      ex.var += 2;  /* 次の展開項を計算     */

    if ( ex.var % 80 == 1 ) { /* ex.var 80  毎に */
      regular( OFF );       /* ans[]を正規表現化（高速モード） */

      if (min_one < 80)    /* １分相当値未計算なら */
        min_one = bench_mark( min_one );  /* 計算する */

      if ( (min_one >= 80) && ex.var % min_one == 1 ) {  /* １分経過なら */
        time( &cur_time );  /* 現在時刻の取得 */

        if (ex.save) {   /* セーブ指定有りで時間なら */
          if ( (long)difftime( cur_time, ex.s_sec ) > 60 * ex.save )
            save();      /* セーブするだよん♪(^0^) */
        }
        if ( (long)difftime( cur_time, pre_time ) <= 40 ) {
          min_one *= 2;   /* 表示間隔が短くなったら倍にする */
        }
        pre_time = cur_time;
        message( LAP, mag, base ); /* 表示時 Ctrl-C有効 */
      }
    }
  }
  regular( ON );
  time_end();    /* 次のarctan()の為 */
  message( DONE, mag, base );
  if (ex.save)  /* 一つのarctan(1/n)終了時も */
    save();     /* saveする */
}

long bench_mark( long min_one ) {   /* 1 分相当の ex.var 差を返す */
  time_t cur_time;
  long   sec;
  static USL st_var;

  if (min_one == 0) { /* 初回計測？   */
    st_var = ex.var;  /* 計測開始時刻 */
    return 1;         /* 初回時は計測 */
  }
  if (min_one >= 2)   /* まだまだ時間がある */
    return min_one - 1;

  time( &cur_time );
  sec = (long)difftime( cur_time, ex.s_sec );
  if (sec >= 60)
    return  ex.var - st_var; /* １分相当値：多桁のロスを防ぐ為80加算 */
  else if (sec >= 30)       /* 80加算の理由は上と同じ */
    return (ex.var - st_var) *60 /sec /80 *80;  /* 1分近似値 */
  else                      /* 30秒以上余裕がある場合は */
    return (60 - sec) / (sec + 1);  /* 残り回数の推算 */
}

void message( int mode, int mag, int base ) {
  static int sign;

  if (mode == SET)
    sign = mag;  /* 計算式全体のフラグ */

  if (mode == LAP)
    fprintf( stderr, "%darctan(1/%d) %2d.%d%%\r",
                      sign * mag, base, RATE/10, RATE%10 );
  if (mode == DONE)
    fprintf( stderr, "%darctan(1/%d) completed !!!\n", sign * mag, base );
}

void c_std( int sign, USL con, USL var, long *arc, long *ans ) {
  /* var must be under 42 */
  int i, j;
  USL co, va, buf, arc_r, ans_r;
  USL con_q = TEN8 / con;  /* arctan(1/5)以外 */
  USL con_r = TEN8 % con;  /* arctan(1/5)以外 */

  if (sign == PLUS) {
    if (con <= 40) {   /* arctan(1/5)のみ */
      for (co = con, va = var, arc_r = 0, ans_r = 0,
           i = ex.head, j = ex.tail; i <= j; i++) {
        buf      = *(arc+i) + arc_r * TEN8;
        *(arc+i) = buf / co;
        arc_r    = buf % co;
        buf       = *(arc+i) + ans_r * TEN8;
        *(ans+i) += buf / va;
        ans_r     = buf % va;
      }
      for (j = ex.max; i <= j; i++) {
        buf       = ans_r * TEN8;
        *(ans+i) += buf / va;
        ans_r     = buf % va;
      }
    } else {         /* arctan(1/5)以外 */
      for (co = con, va = var, arc_r = 0, ans_r = 0,
           i = ex.head, j = ex.tail; i <= j; i++) {
        buf      = *(arc+i) + arc_r * con_r;
        *(arc+i) = buf / co + arc_r * con_q;
        arc_r    = buf % co;
        buf       = *(arc+i) + ans_r * TEN8;
        *(ans+i) += buf / va;
        ans_r     = buf % va;
      }
      for (j = ex.max; i <= j; i++) {
        buf       = ans_r * TEN8;
        *(ans+i) += buf / va;
        ans_r     = buf % va;
      }
    }
  } else {    /* sign == MINUS */
    if (con <= 40) {  /* arctan(1/5)のみ */
      for (co = con, va = var, arc_r = 0, ans_r = 0,
           i = ex.head, j = ex.tail; i <= j; i++) {
        buf      = *(arc+i) + arc_r * TEN8;
        *(arc+i) = buf / co;
        arc_r    = buf % co;
        buf       = *(arc+i) + ans_r * TEN8;
        *(ans+i) -= buf / va;
        ans_r     = buf % va;
      }
      for (j = ex.max; i <= j; i++) {
        buf       = ans_r * TEN8;
        *(ans+i) -= buf / va;
        ans_r     = buf % va;
      }
    } else {         /* arctan(1/5)以外 */
      for (co = con, va = var, arc_r = 0, ans_r = 0,
           i = ex.head, j = ex.tail; i <= j; i++) {
        buf      = *(arc+i) + arc_r * con_r;
        *(arc+i) = buf / co + arc_r * con_q;
        arc_r    = buf % co;
        buf       = *(arc+i) + ans_r * TEN8;
        *(ans+i) -= buf / va;
        ans_r     = buf % va;
      }
      for (j = ex.max; i <= j; i++) {
        buf       = ans_r * TEN8;
        *(ans+i) -= buf / va;
        ans_r     = buf % va;
      }
    }
  }
}

void c_low( int sign, USL con, USL var, long *arc, long *ans ) {
  /* var must be under 64_000 */
  int i, j;
  USL co, va, buf, arc_r, ans_r;
  USL con_q = TEN8 / con;  /* arctan(1/5)以外 */
  USL con_r = TEN8 % con;  /* arctan(1/5)以外 */
  USL var_q = TEN8 / var;
  USL var_r = TEN8 % var;

  if (sign == PLUS) {
    if (con <= 40) {  /* arctan(1/5)のみ */
      for ( co = con, va = var, arc_r = 0, ans_r = 0,
            i = ex.head, j = ex.tail; i <= j; i++) {
        buf      = *(arc+i) + arc_r * TEN8;
        *(arc+i) = buf / co;
        arc_r    = buf % co;
        buf       = *(arc+i) + ans_r * var_r;
        *(ans+i) += buf / va + ans_r * var_q;
        ans_r     = buf % va;
      }
      for (j = ex.max; i <= j; i++ ) {
        buf       = ans_r * var_r;
        *(ans+i) += buf / va + ans_r * var_q;
        ans_r     = buf % va;
      }
    } else {         /* arctan(1/5)以外 */
      for (co = con, va = var, arc_r = 0, ans_r = 0,
           i = ex.head, j = ex.tail; i <= j; i++) {
        buf      = *(arc+i) + arc_r * con_r;
        *(arc+i) = buf / co + arc_r * con_q;
        arc_r    = buf % co;
        buf       = *(arc+i) + ans_r * var_r;
        *(ans+i) += buf / va + ans_r * var_q;
        ans_r     = buf % va;
      }
      for (j = ex.max; i <= j; i++ ) {
        buf       = ans_r * var_r;
        *(ans+i) += buf / va + ans_r * var_q;
        ans_r     = buf % va;
      }
    }
  } else {    /* sign == MINUS */
    if (con <= 40) {  /* arctan(1/5)のみ */
      for (co = con, va = var, arc_r = 0, ans_r = 0,
           i = ex.head, j = ex.tail; i <= j; i++) {
        buf      = *(arc+i) + arc_r * TEN8;
        *(arc+i) = buf / co;
        arc_r    = buf % co;
        buf       = *(arc+i) + ans_r * var_r;
        *(ans+i) -= buf / va + ans_r * var_q;
        ans_r     = buf % va;
      }
      for (j = ex.max; i <= j; i++ ) {
        buf       = ans_r * var_r;
        *(ans+i) -= buf / va + ans_r * var_q;
        ans_r     = buf % va;
      }
    } else {         /* arctan(1/5)以外 */
      for (co = con, va = var, arc_r = 0, ans_r = 0,
           i = ex.head, j = ex.tail; i <= j; i++) {
        buf      = *(arc+i) + arc_r * con_r;
        *(arc+i) = buf / co + arc_r * con_q;
        arc_r    = buf % co;
        buf       = *(arc+i) + ans_r * var_r;
        *(ans+i) -= buf / va + ans_r * var_q;
        ans_r     = buf % va;
      }
      for (j = ex.max; i <= j; i++ ) {
        buf       = ans_r * var_r;
        *(ans+i) -= buf / va + ans_r * var_q;
        ans_r     = buf % va;
      }
    }
  }
}

void c_mid( int sign, USL con, USL var, long *arc, long *ans ) {
  /* var must be under 2_030_000 */
  int i, j;
  USL co, va, buf, buf2, quot;
  USL arc_r, ans_r;
  USL con_q  = TEN8 / con;  /* arctan(1/5)以外 */
  USL con_r  = TEN8 % con;  /* arctan(1/5)以外 */
  USL var_q  = TEN8 / var;
  USL var_r1 = (TEN8 % var) / 1024;
  USL var_r2 = (TEN8 % var) % 1024;

  if (sign == PLUS) {
    if (con <= 40) {  /* arctan(1/5)のみ */
      for (co = con, va = var, arc_r = 0, ans_r = 0,
           i = ex.head, j = ex.tail; i <= j; i++) {
        buf      = *(arc+i) + arc_r * TEN8;
        arc_r    = buf % co;
        *(arc+i) = buf / co;
        buf    = ans_r * var_r1;
        quot   = ( (buf / va) << 10 ) + ans_r * var_q;
        buf2   = ( (buf % va) << 10 ) + ans_r * var_r2 + *(arc+i);
        *(ans+i) += buf2 / va + quot;
        ans_r     = buf2 % va;
      }
      for (j = ex.max; i <= j; i++) {
        buf    = ans_r * var_r1;
        quot   = ( (buf / va) << 10 ) + ans_r * var_q;
        buf2   = ( (buf % va) << 10 ) + ans_r * var_r2;
        *(ans+i) += buf2 / va + quot;
        ans_r     = buf2 % va;
      }
    } else {         /* arctan(1/5)以外 */
      for (co = con, va = var, arc_r = 0, ans_r = 0,
           i = ex.head, j = ex.tail; i <= j; i++) {
        buf      = *(arc+i) + arc_r * con_r;
        *(arc+i) = buf / co + arc_r * con_q;
        arc_r    = buf % co;
        buf    = ans_r * var_r1;
        quot   = ( (buf / va) << 10 ) + ans_r * var_q;
        buf2   = ( (buf % va) << 10 ) + ans_r * var_r2 + *(arc+i);
        *(ans+i) += buf2 / va + quot;
        ans_r     = buf2 % va;
      }
      for (j = ex.max; i <= j; i++) {
        buf    = ans_r * var_r1;
        quot   = ( (buf / va) << 10 ) + ans_r * var_q;
        buf2   = ( (buf % va) << 10 ) + ans_r * var_r2;
        *(ans+i) += buf2 / va + quot;
        ans_r     = buf2 % va;
      }
    }
  } else {    /* sign == MINUS */
    if (con <= 40) {  /* arctan(1/5)のみ */
      for (co = con, va = var, arc_r = 0, ans_r = 0,
           i = ex.head, j = ex.tail; i <= j; i++) {
        buf      = *(arc+i) + arc_r * TEN8;
        arc_r    = buf % co;
        *(arc+i) = buf / co;
        buf    = ans_r * var_r1;
        quot   = ( (buf / va) << 10 ) + ans_r * var_q;
        buf2   = ( (buf % va) << 10 ) + ans_r * var_r2 + *(arc+i);
        *(ans+i) -= buf2 / va + quot;
        ans_r     = buf2 % va;
      }
      for (j = ex.max; i <= j; i++) {
        buf    = ans_r * var_r1;
        quot   = ( (buf / va) << 10 ) + ans_r * var_q;
        buf2   = ( (buf % va) << 10 ) + ans_r * var_r2;
        *(ans+i) -= buf2 / va + quot;
        ans_r     = buf2 % va;
      }
    } else {         /* arctan(1/5)以外 */
      for (co = con, va = var, arc_r = 0, ans_r = 0,
           i = ex.head, j = ex.tail; i <= j; i++) {
        buf      = *(arc+i) + arc_r * con_r;
        *(arc+i) = buf / co + arc_r * con_q;
        arc_r    = buf % co;
        buf    = ans_r * var_r1;
        quot   = ( (buf / va) << 10 ) + ans_r * var_q;
        buf2   = ( (buf % va) << 10 ) + ans_r * var_r2 + *(arc+i);
        *(ans+i) -= buf2 / va + quot;
        ans_r     = buf2 % va;
      }
      for (j = ex.max; i <= j; i++) {
        buf    = ans_r * var_r1;
        quot   = ( (buf / va) << 10 ) + ans_r * var_q;
        buf2   = ( (buf % va) << 10 ) + ans_r * var_r2;
        *(ans+i) -= buf2 / va + quot;
        ans_r     = buf2 % va;
      }
    }
  }
}

void c_high( int sign, USL con, USL var, long *arc, long *ans ) {
  /* var must be under 2_560_000 */
  int i, j;
  USL co, va, buf, buf2, quot;
  USL arc_r, ans_r;
  USL con_q  = TEN8 / con;  /* arctan(1/5)以外 */
  USL con_r  = TEN8 % con;  /* arctan(1/5)以外 */
  USL var_q  = TEN8 / var;
  USL var_r1 = (TEN8 % var) / 1600;
  USL var_r2 = (TEN8 % var) % 1600;

  if (sign == PLUS) {
    if (con <= 40) {
      for (co = con, va = var, arc_r = 0, ans_r = 0,
           i = ex.head, j = ex.tail; i <= j; i++) {
        buf      = *(arc+i) + arc_r * TEN8;
        *(arc+i) = buf / co;
        arc_r    = buf % co;
        buf    = ans_r * var_r1;
        quot   = ans_r * var_q;
        quot  += 1600 * (buf / va);
        buf2   = 1600 * (buf % va) + ans_r * var_r2 + *(arc+i);
        *(ans+i) += buf2 / va + quot;
        ans_r     = buf2 % va;
      }
      for (j = ex.max; i <= j; i++) {
        buf    = ans_r * var_r1;
        quot   = ans_r * var_q;
        quot  += 1600 * (buf / va);
        buf2   = 1600 * (buf % va) + ans_r * var_r2;
        *(ans+i) += buf2 / va + quot;
        ans_r     = buf2 % va;
      }
    } else {       
      for (co = con, va = var, arc_r = 0, ans_r = 0,
           i = ex.head, j = ex.tail; i <= j; i++) {
        buf      = *(arc+i) + arc_r * con_r;
        *(arc+i) = buf / co + arc_r * con_q;
        arc_r    = buf % co;
        buf    = ans_r * var_r1;
        quot   = ans_r * var_q;
        quot  += 1600*(buf / va);
        buf2   = 1600*(buf % va) + ans_r * var_r2 + *(arc+i);
        *(ans+i) += buf2 / va + quot;
        ans_r     = buf2 % va;
      }
      for (j = ex.max; i <= j; i++) {
        buf    = ans_r * var_r1;
        quot   = ans_r * var_q;
        quot  += 1600 * (buf / va);
        buf2   = 1600 * (buf % va) + ans_r * var_r2;
        *(ans+i) += buf2 / va + quot;
        ans_r     = buf2 % va;
      }
    }
  } else {
    if (con <= 40) {
      for (co = con, va = var, arc_r = 0, ans_r = 0,
           i = ex.head, j = ex.tail; i <= j; i++) {
        buf      = *(arc+i) + arc_r * TEN8;
        *(arc+i) = buf / co;
        arc_r    = buf % co;
        buf    = ans_r * var_r1;
        quot   = ans_r * var_q;
        quot  += 1600 * (buf / va);
        buf2   = 1600 * (buf % va) + ans_r * var_r2 + *(arc+i);
        *(ans+i) -= buf2 / va + quot;
        ans_r     = buf2 % va;
      }
      for (j = ex.max; i <= j; i++) {
        buf    = ans_r * var_r1;
        quot   = ans_r * var_q;
        quot  += 1600 * (buf / va);
        buf2   = 1600 * (buf % va) + ans_r * var_r2;
        *(ans+i) -= buf2 / va + quot;
        ans_r     = buf2 % va;
      }
    } else {
      for (co = con, va = var, arc_r = 0, ans_r = 0,
           i = ex.head, j = ex.tail; i <= j; i++) {
        buf      = *(arc+i) + arc_r * con_r;
        *(arc+i) = buf / co + arc_r * con_q;
        arc_r    = buf % co;
        buf    = ans_r * var_r1;
        quot   = ans_r * var_q;
        quot  += 1600 * (buf / va);
        buf2   = 1600 * (buf % va) + ans_r * var_r2 + *(arc+ i);
        *(ans+i) -= buf2 / va + quot;
        ans_r     = buf2 % va;
      }
      for (j = ex.max; i <= j; i++) {
        buf    = ans_r * var_r1;
        quot   = ans_r * var_q;
        quot  += 1600 * (buf / va);
        buf2   = 1600 * (buf % va) + ans_r * var_r2;
        *(ans+i) -= buf2 / va + quot;
        ans_r     = buf2 % va;
      }
    }
  }
}

void c_huge( int sign, USL con, USL var, long *arc, long *ans ) {
  /* var must be under 15_625_000 (=250^3) */
  int i, j;
  USL co, va, buf, buf2, quot;
  USL arc_r, ans_r;
  USL con_q = TEN8 / con;  /* arctan(1/5)以外 */
  USL con_r = TEN8 % con;  /* arctan(1/5)以外 */
  USL var_q = TEN8 / var;
  USL var_r1 = (TEN8 % var) / 62500L;
  USL var_r2 = ( (TEN8 % var) % 62500L ) / 250;
  USL var_r3 = (TEN8 % var) % 250;

  if (sign == PLUS) {
    if (con <= 40) {     /* arctan(1/5)のみ */
      for (co = con, va = var, arc_r = 0, ans_r = 0,
           i = ex.head, j = ex.tail; i <= j; i++) {
        buf      = *(arc+i) + arc_r * TEN8;
        *(arc+i) = buf / co;
        arc_r    = buf % co;
        buf   = ans_r * var_r1;
        quot  = ans_r * var_q;
        quot += 62500L * (buf / va);
        buf2  = 250 * (buf % va) + ans_r * var_r2;
        quot += 250 * (buf2 / va);
        buf   = 250 * (buf2 % va) + ans_r * var_r3 + *(arc+i);
        *(ans+i) += buf / va + quot;
        ans_r     = buf % va;
      }
      for (j = ex.max; i <= j; i++) {
        buf   = ans_r * var_r1;
        quot  = ans_r * var_q;
        quot += 62500L * (buf / va);
        buf2  = 250 * (buf % va) + ans_r * var_r2;
        quot += 250 * (buf2 / va);
        buf   = 250 * (buf2 % va) + ans_r * var_r3;
        *(ans+i) += buf / va + quot;
        ans_r     = buf % va;
      }
    } else {             /* arctan(1/5)以外 */
      for (co = con, va = var, arc_r = 0, ans_r = 0,
           i = ex.head, j = ex.tail; i <= j; i++) {
        buf      = *(arc+i) + arc_r * con_r;
        *(arc+i) = buf / co + arc_r * con_q;
        arc_r    = buf % co;
        buf   = ans_r * var_r1;
        quot  = ans_r * var_q;
        quot += 62500L * (buf / va);
        buf2  = 250 * (buf % va) + ans_r * var_r2;
        quot += 250 * (buf2 / va);
        buf   = 250 * (buf2 % va) + ans_r * var_r3 + *(arc+i);
        *(ans+i) += buf / va + quot;
        ans_r     = buf % va;
      }
      for (j = ex.max; i <= j; i++) {
        buf   = ans_r * var_r1;
        quot  = ans_r * var_q;
        quot += 62500L * (buf / va);
        buf2  = 250 * (buf % va) + ans_r * var_r2;
        quot += 250 * (buf2 / va);
        buf   = 250 * (buf2 % va) + ans_r * var_r3;
        *(ans+i) += buf / va + quot;
        ans_r     = buf % va;
      }
    }
  } else {    /* sign == MINUS */
    if (con <= 40) {     /* arctan(1/5)のみ */
      for (co = con, va = var, arc_r = 0, ans_r = 0,
           i = ex.head, j = ex.tail; i <= j; i++) {
        buf      = *(arc+i) + arc_r * TEN8;
        *(arc+i) = buf / co;
        arc_r    = buf % co;
        buf   = ans_r * var_r1;
        quot  = ans_r * var_q;
        quot += 62500L * (buf / va);
        buf2  = 250 * (buf % va) + ans_r * var_r2;
        quot += 250 * (buf2 / va);
        buf   = 250 * (buf2 % va) + ans_r * var_r3 + *(arc+i);
        *(ans+i) -= buf / va + quot;
        ans_r     = buf % va;
      }
      for (j = ex.max; i <= j; i++) {
        buf   = ans_r * var_r1;
        quot  = ans_r * var_q;
        quot += 62500L * (buf / va);
        buf2  = 250 * (buf % va) + ans_r * var_r2;
        quot += 250 * (buf2 / va);
        buf   = 250 * (buf2 % va) + ans_r * var_r3;
        *(ans+i) -= buf / va + quot;
        ans_r     = buf % va;
      }
    } else {             /* arctan(1/5)以外 */
      for (co = con, va = var, arc_r = 0, ans_r = 0,
           i = ex.head, j = ex.tail; i <= j; i++) {
        buf      = *(arc+i) + arc_r * con_r;
        *(arc+i) = buf / co + arc_r * con_q;
        arc_r    = buf % co;
        buf   = ans_r * var_r1;
        quot  = ans_r * var_q;
        quot += 62500L * (buf / va);
        buf2  = 250 * (buf % va) + ans_r * var_r2;
        quot += 250 * (buf2 / va);
        buf   = 250 * (buf2 % va) + ans_r * var_r3 + *(arc+i);
        *(ans+i) -= buf / va + quot;
        ans_r     = buf % va;
      }
      for (j = ex.max; i <= j; i++) {
        buf   = ans_r * var_r1;
        quot  = ans_r * var_q;
        quot += 62500L * (buf / va);
        buf2  = 250 * (buf % va) + ans_r * var_r2;
        quot += 250 * (buf2 / va);
        buf   = 250 * (buf2 % va) + ans_r * var_r3;
        *(ans+i) -= buf / va + quot;
        ans_r     = buf % va;
      }
    }
  }
}

void regular( int mode ) {    /* ans[]を0〜TEN8-1の範囲に収める */
  int i;
  long buf, carry;
  long *ans = ex.ans;
  buf   = 0;
  carry = 0;

  if (mode == OFF) {    /* 計算中は高速化の為、負の処理をしない */
    for (i = ex.max; i >= 0; i--) {
      buf = *(ans + i) + carry;
      *(ans + i) = buf % TEN8;
      carry      = buf / TEN8;
    }
  } else {              /* 負の処理も行う */
    for (i = ex.max; i >= 0; i--) {
      buf = *(ans + i) + carry;
      *(ans + i) = buf % TEN8;
      carry      = buf / TEN8;
      if ( *(ans + i) < 0) {    /* 負の% 演算の処理 */
        *(ans + i) += TEN8;
          carry--;
      }
    }
  }
}

void load( char *fname ) {
  FILE *fp;
  char buf[80];
  int m;

  if ( ( fp = fopen(fname, "r") ) == NULL ) {
    fputs( "file not found", stderr );
    usage( !DATA );
  }
  ex.load = ON; /* ロードフラグON */

  fgets( buf, 79, fp );  /* まずコメント行をダミーリード */
  fgets( buf, 79, fp );  /* 次に配列に変数群を読み込む */
  if ( sscanf( buf, "%d%ld%ld%ld%d%d%d%d%ld%ld%ld",
               &ex.ver,   /* Version  */
               &ex.sec,   /* 計算時間 */
               &ex.cent,  /* 〃1/100  */
               &ex.keta,  /* 計算桁数 */
               &ex.max,   /* 静的最終 */
               &ex.type,  /* 計算公式 */
               &ex.head,  /* 動的先頭 */
               &ex.tail,  /* 動的最終 */
               &ex.save,  /* 退避flag */
               &ex.con,   /* base^2   */
               &ex.var )  /* arc変数  */ 
      != 11                              /* 変数が足りない   */
      || (ex.max != (int)ex.keta/8 + 2)  /* 関係式がおかしい */
      || ! (ex.save && ex.con && ex.var) /* どれかが0 である */
      ) {                                /* どれでも成立すれば異常 */
      fclose( fp );
      fputs( "Not regular file", stderr );
      exit( 3 );
  }
  get_memory();

  for ( m = 0; m <= ex.max; m++ ) {
    fgets( buf, 79, fp );
    if ( sscanf( buf, "%ld%ld", ex.ans+m, ex.arc+m ) != 2 ) {
      fclose( fp );
      fputs( "Not regular file\n", stderr );
      exit( 3 );
    }
  }
  fclose( fp );
}

void save( void ) {
  FILE *fp;
  int  m;
  time_t cur_time; /* 現在時刻記録用 */
  char fname[13];  /* 8+1+3+NUL = 13 */
  char buf[80];

  time_end();
  sprintf( fname, "%ld.cal", ex.keta );
  if ( ( fp = fopen(fname, "w") ) == NULL ) {
    time( &cur_time );
    fprintf( stderr, "file can't be opened at %s", ctime( &cur_time ) );
    time_set();  /* セーブ失敗でも計算は続行 */
    return;
  }

  regular( ON );  /* 既にregular()しているのでこれはズルでない */
  fputs( "PI.C(Ver 2.0)(Normal mode) ここはコメント行（79文字以下）\n", fp );
  sprintf( buf, "%d %ld %ld %ld %d %d %d %d %ld %ld %ld\n",
           ex.ver,   /* Version  */
           ex.sec,   /* 計算時間 */
           ex.cent,  /* 〃1/100  */
           ex.keta,  /* 計算桁数 */
           ex.max,   /* 静的最終 */
           ex.type,  /* 計算公式 */
           ex.head,  /* 動的先頭 */
           ex.tail,  /* 動的最終 */
           ex.save,  /* 退避flag */
           ex.con,   /* base^2   */
           ex.var ); /* arc変数  */
  fputs( buf, fp );

  for ( m = 0; m <= ex.max; m++ ) {
    sprintf( buf, "%ld %ld\n", *(ex.ans+m), *(ex.arc+m) );
    fputs( buf, fp );
  }
  fclose( fp );
  time_set();  /* 計算開始時刻を再セット */
}

void time_set(void) {  /* 開始時間セット   */
  time( &ex.s_sec );   /* 計算開始時刻(秒) */
  ex.s_cent = clock(); /* 〃1/100 秒       */ 
}

void time_end( void ) {  /* 経過時間更新 */
  long    sec, cent;
  time_t  e_sec;
  clock_t e_cent;

  time( &e_sec );
  e_cent = clock();
  sec  = (long)difftime( e_sec, ex.s_sec );
  cent = ( e_cent - ex.s_cent + 86400 * (long)CLK_TCK ) 
          % ( 86400 * (long)CLK_TCK );

  sec /= 86400;
  sec *= 86400;                  /* 一日未満の秒は一旦切捨て、 */
  sec += (cent / (long)CLK_TCK); /*  cent より求めて再加算する */
  cent *= 100;            /* 1/100 秒単位の数にする為 */
  cent /= (long)CLK_TCK;  /* これでなっているハズ    */
  cent  = cent % 100;     /* 1 秒未満 */

  ex.sec  += sec;      /* 計算時間更新 */
  ex.cent += cent;     /* 〃           */
  if (ex.cent >= 100) {
    ex.cent -= 100;
    ex.sec++;
  }
  ex.s_sec  = e_sec;    /* 再開時刻更新(time_setと同じ) */
  ex.s_cent = e_cent;   /* 〃1/100 秒 */ 
}

void show_ans( void ) {
  int i = 0;    /* 配列count  */
  int l = 0;    /* 文字列位置 */
  char buf[60];
  char *buf_p = buf;

  fprintf( stderr, "計算時間は%ld時間", ex.sec/3600 );
  fprintf( stderr, "%0.2ld分%0.2ld.%0.2ld秒でした\n\n", 
                   (ex.sec/60)%60, ex.sec%60, ex.cent );
  printf("πの小数点以下%d(＋α）桁の値です。\n", ex.keta );

  printf( "00000001:%0.1ld.", *ex.ans );
  for ( i = 1; i <= ex.max; i++) {
    sprintf( buf_p, "%0.8ld\0", *(ex.ans + i) );
    buf_p += 8;
    if ( strlen(buf) >= 50 || i == ex.max ) {
      for ( l = 0; l < 50 && buf[l]; l++) {
        putchar( buf[l] );
        if ( l % 10 == 9 )
          putchar( ' ' );
      }
      if (l != 50)
        break;
      printf( "\n%0.8d:  ", (i*8 / 10)*10 + 1 );
      strcpy( buf, buf + 50 );  /* 50桁分詰める */
      buf_p -= 50;              /* 50桁分戻る   */
    }
  }
  printf( "\n計算時間は%ld時間", ex.sec/3600 );
  printf( "%0.2ld分%0.2ld.%0.2ld秒でした\n", 
           (ex.sec/60)%60, ex.sec%60, ex.cent );
}

void usage( int mode ) {
  if (mode == DATA) {
    puts("TECNICAL DATA････通常の使い方は  A>RUN386 PI で表示\n");
    puts("TOWNS 2F(初代) 1万桁  10万桁  100万桁");
    puts("ｳｪｲﾄ無調整     96秒   165分   333時間");
    puts("公式別時間     ﾏﾁﾝ    ｶﾞｳｽ    ｼｭﾃﾙﾏｰ  ﾗｻﾞﾌｫｰﾄﾞ");
    puts("               96秒   131秒   141秒   141秒");
    puts("\nセーブファイルについて");
    puts("  サイズ：  　桁数×(1.5〜2.1)倍");
    puts("  ファイル名：セーブ時：桁数.CAL（.CALは自動的に付く）");
    puts("              ロード時：拡張子も含め任意の名が可能");
    puts("\n計算再開時にセーブ時間の再変更が可能");
    puts("セーブ指定有りでは、各arctan()終了毎にもセーブする。");
    puts("\n計算進行表示は約１分毎だが途中多少の変動あり");
    puts("\nまた、第２オプションは複数記述が可能");
    exit(0);
  }
  puts("RUN386 PI 100000     ←第１オプション：桁数 or ファイル名");
  puts("          100000.CAL   拡張子必要：ファイルを読み込んで計算再開");
  puts("          (others)     無効オプション： Tecnical data  表示");
  puts("RUN386 PI 100000 G20 ←第２オプション：公式 or セーブ時間間隔(分)");
  puts("                       g:ｶﾞｳｽ m:ﾏﾁﾝ(ﾃﾞﾌｫﾙﾄ) r:ﾗｻﾞﾌｫｰﾄﾞ s:ｼｭﾃﾙﾏｰ");
  puts("                       例はガウス公式で20分おきにセーブする設定");
  puts("約１分毎の計算進行表示時、 Ctrl-C による中断が可能である");
  exit(0);
}