Search Apps Documentation Source Content File Folder Download Copy Actions Download

arithmetic.gno

2.59 Kb ยท 146 lines
  1package int256
  2
  3import "math/bits"
  4
  5func udivrem(quot, u []uint64, d *Int) (rem Int) {
  6	var dLen int
  7	for i := len(d) - 1; i >= 0; i-- {
  8		if d[i] != 0 {
  9			dLen = i + 1
 10			break
 11		}
 12	}
 13
 14	shift := uint(bits.LeadingZeros64(d[dLen-1]))
 15
 16	var dnStorage Int
 17	dn := dnStorage[:dLen]
 18	for i := dLen - 1; i > 0; i-- {
 19		dn[i] = (d[i] << shift) | (d[i-1] >> (64 - shift))
 20	}
 21	dn[0] = d[0] << shift
 22
 23	var uLen int
 24	for i := len(u) - 1; i >= 0; i-- {
 25		if u[i] != 0 {
 26			uLen = i + 1
 27			break
 28		}
 29	}
 30
 31	if uLen < dLen {
 32		copy(rem[:], u)
 33		return rem
 34	}
 35
 36	var unStorage [9]uint64
 37	un := unStorage[:uLen+1]
 38	un[uLen] = u[uLen-1] >> (64 - shift)
 39	for i := uLen - 1; i > 0; i-- {
 40		un[i] = (u[i] << shift) | (u[i-1] >> (64 - shift))
 41	}
 42	un[0] = u[0] << shift
 43
 44	if dLen == 1 {
 45		r := udivremBy1(quot, un, dn[0])
 46		rem.SetUint64(r >> shift)
 47		return rem
 48	}
 49
 50	udivremKnuth(quot, un, dn)
 51
 52	for i := 0; i < dLen-1; i++ {
 53		rem[i] = (un[i] >> shift) | (un[i+1] << (64 - shift))
 54	}
 55	rem[dLen-1] = un[dLen-1] >> shift
 56
 57	return rem
 58}
 59
 60func udivremBy1(quot, u []uint64, d uint64) (rem uint64) {
 61	reciprocal := reciprocal2by1(d)
 62	rem = u[len(u)-1]
 63	for j := len(u) - 2; j >= 0; j-- {
 64		quot[j], rem = udivrem2by1(rem, u[j], d, reciprocal)
 65	}
 66	return rem
 67}
 68
 69func udivremKnuth(quot, u, d []uint64) {
 70	dLen := len(d)
 71	dh := d[dLen-1]
 72	dl := d[dLen-2]
 73	reciprocal := reciprocal2by1(dh)
 74
 75	for j := len(u) - dLen - 1; j >= 0; j-- {
 76		u2 := u[j+dLen]
 77		u1 := u[j+dLen-1]
 78		u0 := u[j+dLen-2]
 79
 80		var qhat, rhat uint64
 81		if u2 >= dh {
 82			qhat = MAX_UINT64
 83		} else {
 84			qhat, rhat = udivrem2by1(u2, u1, dh, reciprocal)
 85			ph, pl := bits.Mul64(qhat, dl)
 86			if ph > rhat || (ph == rhat && pl > u0) {
 87				qhat--
 88			}
 89		}
 90
 91		borrow := subMulTo(u[j:], d, qhat)
 92		u[j+dLen] = u2 - borrow
 93		if u2 < borrow {
 94			qhat--
 95			u[j+dLen] += addTo(u[j:], d)
 96		}
 97
 98		quot[j] = qhat
 99	}
100}
101
102func reciprocal2by1(d uint64) uint64 {
103	reciprocal, _ := bits.Div64(^d, MAX_UINT64, d)
104	return reciprocal
105}
106
107func udivrem2by1(uh, ul, d, reciprocal uint64) (quot, rem uint64) {
108	qh, ql := bits.Mul64(reciprocal, uh)
109	ql, carry := bits.Add64(ql, ul, 0)
110	qh += uh + carry
111	qh++
112
113	r := ul - qh*d
114
115	if r > ql {
116		qh--
117		r += d
118	}
119
120	if r >= d {
121		qh++
122		r -= d
123	}
124
125	return qh, r
126}
127
128func subMulTo(x, y []uint64, multiplier uint64) uint64 {
129	var borrow uint64
130	for i := 0; i < len(y); i++ {
131		s, carry1 := bits.Sub64(x[i], borrow, 0)
132		ph, pl := bits.Mul64(y[i], multiplier)
133		t, carry2 := bits.Sub64(s, pl, 0)
134		x[i] = t
135		borrow = ph + carry1 + carry2
136	}
137	return borrow
138}
139
140func addTo(x, y []uint64) uint64 {
141	var carry uint64
142	for i := 0; i < len(y); i++ {
143		x[i], carry = bits.Add64(x[i], y[i], carry)
144	}
145	return carry
146}