符号なし除算(回復法)

[] [] []


分子と分母が符号なし整数の場合の除算は符号なし除算(回復法)に書いたとおりで、これを固定小数点に直したものは次のとおり。余りは切捨でなく丸め(四捨五入)とする。

a / d = q + e

分子a:符号なし、整数部IAビット、小数部FAビット (0a2IA-2-FA)

分母d:符号なし、整数部IDビット、小数部FDビット (2-FDd2ID-2-FD)

q: 符号なし、整数部IQビット、小数部FQビット (0q2IQ-2-FQ)

    (1) a,dが任意の値をとりうる場合

     a/d(2IA-1-2-FA)/2-FD2IA+FDだから整数部はIQ=IA+FDだけ必要

     小数部はFQFA-FD

    (2) a/d2IQ-2-FQとわかっている場合はビット幅IQ<IA+FDで十分

     例) N個の値xi の平均値はxiの値の範囲内に必ず収まる

誤差e:丸めの場合 -2-FQ-1≦e2-FQ-1

VitisHLSで合成可能なCコードは

#include <ap_int.h>

/******************************************************

 * division a/d = q + e

 *   param range: IQ<=IA+FD FQ>=FA-FD

 *   input range: 0<d,0<=a<=d*(2**(IQ)-2**(-FQ))

 *   error range: -2**(-FQ-1)<=e<2**(-FQ-1) (rounding)

 *****************************************************/

template <int IA,int FA,int ID,int FD,int IQ,int FQ>

ap_ufixed<IQ+FQ,IQ> div_ufixed // (o) quotient

  (ap_ufixed<IA+FA,IA> a // (i) dividend

  ,ap_ufixed<ID+FD,ID> d // (i) divisor

  ) {

  ap_uint<IA+FD+FQ+1> a_tmp

    = (ap_uint<IA+FD+FQ+1>) a(IA+FA-1,0)<<(FQ-FA+FD+1);

  ap_uint<ID+FD+1> r;

  ap_uint<ID+FD+1> t;

  ap_ufixed<IQ+FQ,IQ> q;

  for(int i=IQ+FQ; i>=0; i--) {

    if(i==IQ+FQ) {

      r = a_tmp(IA+FD+FQ,i); // sign extension

    } else {

      r = (q[i]? t(ID+FD-1,0):r(ID+FD-1,0),a_tmp[i]);

    }

    t = r - d(ID+FD-1,0);

    if(i==0) {

      q(IQ+FQ-1,0) = q(IQ+FQ-1,0) + !t[ID+FD];

    } else {

      q[i-1] = !t[ID+FD];

    }

  }

  return q;

}

このコードは次のテストベンチで検証してある。

int main(int argc, char *argv[]) {

  const int IA=3,FA=1;

  const int ID=2,FD=2;

#ifdef FULL_IQ

  const int IQ=5,FQ=3;  // IQ=IA+FD FQ>=FA-FD

#else

  const int IQ=4,FQ=3;  // IQ<IA+FD FQ>=FA-FD

#endif

  ap_ufixed<IA+FA,IA> a;

  ap_ufixed<ID+FD,ID> d;

  ap_ufixed<IQ+FQ,IQ> q;

  double a_val,d_val,q_val,e_val;

  printf("-- IA=%-d FA=%-d ID=%-d FD=%-d IQ=%-d FQ=%-d\n",IA,FA,ID,FD,IQ,FQ);

  /* a,dを全範囲振って商と誤差を表示する */

  for(int j=1; j<pow(2,ID+FD); j++) {

    for(int i=0; i<pow(2,IA+FA); i++) {

      a = (float) i * pow(2,-FA);

      d = (float) j * pow(2,-FD);

      a_val = a.to_double();

      d_val = d.to_double();

      printf("-- a=%+f d=%+f",a_val,d_val);

      /* IQ<IA+FDのとき大きいaはスキップする */

      if(a_val<=d_val*(pow(2,IQ)-pow(2,-FQ))) {

        q = div_ufixed<IA,FA,ID,FD,IQ,FQ>(a,d);

        q_val = q.to_double();

        e_val = a_val/d_val - q_val;

        printf(" q=%f e=%+f",q_val,e_val);

        if(!(-pow(2,-FQ-1)<=e_val && e_val<pow(2,-FQ-1))) printf(" [FAIL]");

      }

      printf("\n");

    }

  }

  return 0;

}


[] [] []