Linux 0.12 OS. math - div.c

来源:互联网 发布:caffe的安装和配置 编辑:程序博客网 时间:2024/06/05 07:10
static void shift_left(int * c){__asm__ __volatile__("movl (%0),%%eax ; addl %%eax,(%0)\n\t""movl 4(%0),%%eax ; adcl %%eax,4(%0)\n\t""movl 8(%0),%%eax ; adcl %%eax,8(%0)\n\t""movl 12(%0),%%eax ; adcl %%eax,12(%0)"::"r" ((long) c));}static void shift_right(int * c){__asm__("shrl $1,12(%0) ; rcrl $1,8(%0) ; rcrl $1,4(%0) ; rcrl $1,(%0)"::"r" ((long) c));}static int try_sub(int * a, int * b){char ok;/*b = b -a, ok = ~CF.sbb会减去借位标志*/__asm__ __volatile__("movl (%1),%%eax ; subl %%eax,(%2)\n\t""movl 4(%1),%%eax ; sbbl %%eax,4(%2)\n\t""movl 8(%1),%%eax ; sbbl %%eax,8(%2)\n\t""movl 12(%1),%%eax ; sbbl %%eax,12(%2)\n\t""setae %%al":"=a" (ok):"c" ((long) a),"d" ((long) b)); /* setae: ~CF */return ok;}static void div64(int * a, int * b, int * c){int tmp[4];int i;unsigned int mask = 0;/*二进制的减法: c = a/b.c[1]=c[0]=0,因为最多只循环了64位.我们可以用4位来模拟:a = 1011 0011b = 1001 0100c =                 1 0011 01011001 0100 | 1011 0011           a-           1001 0100           ... b>>0            0001 1111 000       a-           0100 1010 0         ... b>>1-           0010 0101 00        ... b>>2-           0001 0010 100       ... b>>3            0000 1100 1000      a-           0000 1001 0100      ... b>>4            0000 0011 0100 00   a-           0000 0100 1010 0    ... b>>5-           0000 0010 0101 00   ... b>>6            0000 0000 1111 0000 a-           0000 0001 0010 100  ... b>>7-           0000 0000 1001 0100 ... b>>8            0000 0000 0101 1100 a当然,还可以继续知道c被填满为止,不过这里只取前8位*/c += 4; /* 商 */for (i = 0 ; i<64 ; i++) {if (!(mask >>= 1)) { /* mask与c来得到最终的商 */c--;mask = 0x80000000; /* 商用c数组来存储,每个元素表示32位,因此这里表示32位满了,需要动用下一个元素 */}tmp[0] = a[0]; tmp[1] = a[1];tmp[2] = a[2]; tmp[3] = a[3];if (try_sub(b,tmp)) { /* 如果a-b还够减的话,商或上当前mask位置的1,用余数更新a */*c |= mask;a[0] = tmp[0]; a[1] = tmp[1];a[2] = tmp[2]; a[3] = tmp[3];}shift_right(b); /* 除数移位 */}}void fdiv(const temp_real * src1, const temp_real * src2, temp_real * result){int i,sign;int a[4],b[4],tmp[4] = {0,0,0,0};sign = (src1->exponent ^ src2->exponent) & 0x8000; /* 异或取符号位 */if (!(src2->a || src2->b)) { /* 除0异常 */set_ZE(); /* 8.5.3 Divide-by-zero exception (#Z) */return;}i = (src1->exponent & 0x7fff) - (src2->exponent & 0x7fff) + 16383; /* 计算出指数部分 */if (i<0) { /* 很小趋近于0 */set_UE(); /* 8.5.5 Numeric underflow exception (#U) */result->exponent = sign;result->a = result->b = 0;return;}a[0] = a[1] = 0; /* 被除数放在a数组中 */a[2] = src1->a;a[3] = src1->b;b[0] = b[1] = 0; /* 除数放在b数组中 */b[2] = src2->a;b[3] = src2->b;while (b[3] >= 0) { /* 除数去掉前导0 */i++;  /* 除数指数--,结果是++ */shift_left(b);}div64(a,b,tmp); /* tmp = a/b */if (tmp[0] || tmp[1] || tmp[2] || tmp[3]) {while (i && tmp[3] >= 0) { /* 规范化 */i--;shift_left(tmp);}if (tmp[3] >= 0) /* 非规范化异常 */set_DE(); /* 8.5.2 Denormal operand exception (#D) */} else /* 结果为0 */i = 0;if (i>0x7fff) { /* 很大异常 */set_OE();return;}if (tmp[0] || tmp[1]) /* 不能整除 */set_PE(); /* 8.5.6 Inexact-result (Precision) exception (#P) */result->exponent = i | sign;result->a = tmp[2];result->b = tmp[3];}
原创粉丝点击