符号なし除算(回復法)

[←]　[↑]　[→]

---

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

a / d = q + e

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

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

商q：　符号なし、整数部IQビット、小数部FQビット (0≦q≦2^IQ-2^-FQ)

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

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

　　　　　a/dの最小=2^-FA/(2^ID-2^-FD)＞2^-(FA+ID)だから小数部はFQ=FA+IDで十分

　　　　(2) 商の範囲がa/d≦2^IQ-2^-FQと限られている場合、整数部はIQ<IA+FDで十分

　　　　　例) 平均値Σx_i/Nはx_iのとりうる範囲内に限られる

誤差e：丸め -2^-FQ-1≦e＜2^-FQ-1

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

#include <ap_int.h>

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

 * division a/d = q + e

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

 *   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 LIMITED_RANGE_OF_Q

  const int IQ=IA+FD-1,FQ=FA+ID;       // IQ<IA+FD

#else

  const int IQ=IA+FD,FQ=FA+ID;  // IQ=IA+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 = (double) i * pow(2,-FA);

      d = (double) 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;

}

---

[←]　[↑]　[→]