Subversion Repositories planix.SVN

Rev

Rev 22 | Show entire file | Ignore whitespace | Details | Blame | Last modification | View Log | RSS feed

Rev 22 Rev 26
Line 7... Line 7...
7 
// counts down rather than up.
 7 
// counts down rather than up.
8 
 
 8 
 
9 
void
 9 
void
10 
mpdiv(mpint *dividend, mpint *divisor, mpint *quotient, mpint *remainder)
 10 
mpdiv(mpint *dividend, mpint *divisor, mpint *quotient, mpint *remainder)
11 
{
 11 
{
12 
	int j, s, vn, sign;
 12 
	int j, s, vn, sign, qsign, rsign;
13 
	mpdigit qd, *up, *vp, *qp;
 13 
	mpdigit qd, *up, *vp, *qp;
14 
	mpint *u, *v, *t;
 14 
	mpint *u, *v, *t;
- 
 
 15 
 
- 
 
 16 
	assert(quotient != remainder);
- 
 
 17 
	assert(divisor->flags & MPnorm);
15 
 
 18 
 
16 
	// divide bv zero
 19 
	// divide bv zero
17 
	if(divisor->top == 0)
 20 
	if(divisor->top == 0)
18 
		abort();
 21 
		abort();
- 
 
 22 
 
- 
 
 23 
	// division by one or small powers of two
- 
 
 24 
	if(divisor->top == 1 && (divisor->p[0] & divisor->p[0]-1) == 0){
- 
 
 25 
		vlong r = 0;
- 
 
 26 
		if(dividend->top > 0)
- 
 
 27 
			r = (vlong)dividend->sign * (dividend->p[0] & divisor->p[0]-1);
- 
 
 28 
		if(quotient != nil){
- 
 
 29 
			sign = divisor->sign;
- 
 
 30 
			for(s = 0; ((divisor->p[0] >> s) & 1) == 0; s++)
- 
 
 31 
				;
- 
 
 32 
			mpright(dividend, s, quotient);
- 
 
 33 
			if(sign < 0)
- 
 
 34 
				quotient->sign ^= (-mpmagcmp(quotient, mpzero) >> 31) << 1;
- 
 
 35 
		}
- 
 
 36 
		if(remainder != nil){
- 
 
 37 
			remainder->flags |= dividend->flags & MPtimesafe;
- 
 
 38 
			vtomp(r, remainder);
- 
 
 39 
		}
- 
 
 40 
		return;
- 
 
 41 
	}
- 
 
 42 
	assert((dividend->flags & MPtimesafe) == 0);
19 
 
 43 
 
20 
	// quick check
 44 
	// quick check
21 
	if(mpmagcmp(dividend, divisor) < 0){
 45 
	if(mpmagcmp(dividend, divisor) < 0){
22 
		if(remainder != nil)
 46 
		if(remainder != nil)
23 
			mpassign(dividend, remainder);
 47 
			mpassign(dividend, remainder);
24 
		if(quotient != nil)
 48 
		if(quotient != nil)
25 
			mpassign(mpzero, quotient);
 49 
			mpassign(mpzero, quotient);
26 
		return;
 50 
		return;
27 
	}
 51 
	}
- 
 
 52
- 
 
 53 
	qsign = divisor->sign * dividend->sign;
- 
 
 54 
	rsign = dividend->sign;
28 
 
 55 
 
29 
	// D1: shift until divisor, v, has hi bit set (needed to make trial
 56 
	// D1: shift until divisor, v, has hi bit set (needed to make trial
30 
	//     divisor accurate)
 57 
	//     divisor accurate)
31 
	qd = divisor->p[divisor->top-1];
 58 
	qd = divisor->p[divisor->top-1];
32 
	for(s = 0; (qd & mpdighi) == 0; s++)
 59 
	for(s = 0; (qd & mpdighi) == 0; s++)
Line 93... Line 120...
93 
		// push top of u down one
 120 
		// push top of u down one
94 
		u->top--;
 121 
		u->top--;
95 
		*up-- = 0;
 122 
		*up-- = 0;
96 
	}
 123 
	}
97 
	if(qp != nil){
 124 
	if(qp != nil){
- 
 
 125 
		assert((quotient->flags & MPtimesafe) == 0);
98 
		mpnorm(quotient);
 126 
		mpnorm(quotient);
99 
		if(dividend->sign != divisor->sign)
 127 
		if(quotient->top != 0)
100 
			quotient->sign = -1;
 128 
			quotient->sign = qsign;
101 
	}
 129 
	}
102 
 
 130 
 
103 
	if(remainder != nil){
 131 
	if(remainder != nil){
- 
 
 132 
		assert((remainder->flags & MPtimesafe) == 0);
104 
		mpright(u, s, remainder);	// u is the remainder shifted
 133 
		mpright(u, s, remainder);	// u is the remainder shifted
- 
 
 134 
		if(remainder->top != 0)
105 
		remainder->sign = dividend->sign;
 135 
			remainder->sign = rsign;
106 
	}
 136 
	}
107 
 
 137 
 
108 
	mpfree(t);
 138 
	mpfree(t);
109 
	mpfree(u);
 139 
	mpfree(u);
110 
	if(v != divisor)
 140 
	if(v != divisor)