Re: bdiv vs redc

2014-04-11 Thread Niels Möller
ni...@lysator.liu.se (Niels Möller) writes: I've checked in my changes in a new repository, ssh://shell.gmplib.org//var/hg/gmp-proj/gmp-bdiv Patch also appended below. The code passes make check now. I just remember this old sub-project: Changing bdiv convention from U = Q D + R B^n

Re: 2-adic roots (Re: bdiv vs redc)

2012-07-31 Thread Niels Möller
Replying to myself yet again... For the power x^n, I currently use mpn_powlo. But the least significant half is known from the previous iteration, so wraparound would be desirable. To me it would make some sense with a pow function for (mod (2^k - 1)), i.e., using mpn_mulmod_bnm1 and

Re: 2-adic roots (Re: bdiv vs redc)

2012-07-30 Thread Niels Möller
ni...@lysator.liu.se (Niels Möller) writes: Correction: Each iteration gets from \ell to 2 \ell-2, and it needs 2 \ell-1 bits of precision for intermediate values. Yet another correction, after more work on the implementation. For numbers x = 1 (mod 8), there are four square roots mod 2^k. If

Re: 2-adic roots (Re: bdiv vs redc)

2012-07-30 Thread David Harvey
On Jul 31, 2012, at 7:42 AM, Niels Möller wrote: We currently have modular exponentation, powlo and regular powering with no reduction of any kind. I'm suggesting a pow_modbnm1. For euclidean square root, and for mpfr, it might also be useful with a pow_high, keeping only the n most

Re: bdiv vs redc

2012-07-17 Thread Torbjorn Granlund
ni...@lysator.liu.se (Niels Möller) writes: ni...@lysator.liu.se (Niels Möller) writes: Some other test cases fail, though: t-root, t-perfpow, t-divis, and t-cong. I'm not familiar with this code, I had a quick look at the root code but didn't see any obvious dependency on bdiv.

Re: bdiv vs redc

2012-07-17 Thread Torbjorn Granlund
ni...@lysator.liu.se (Niels Möller) writes: How do you define hensel square root with remainder? Given a and n, if there exists an x such that x^2 = a (mod B^n), that seems like the reasonable definition of the square root. But what if no such x exists; where should we put the remainder

Re: bdiv vs redc

2012-07-17 Thread Niels Möller
Torbjorn Granlund t...@gmplib.org writes: I don't think a remainder is meaningful here. Ok, so the return value should be a proper root mod B^k (B = 2^{GMP_NUMB_BITS, as usual), or an indication that no root exists. When a square root exists (and the input is non-zero), then there are two

Re: bdiv vs redc

2012-07-16 Thread Niels Möller
ni...@lysator.liu.se (Niels Möller) writes: Some other test cases fail, though: t-root, t-perfpow, t-divis, and t-cong. I'm not familiar with this code, I had a quick look at the root code but didn't see any obvious dependency on bdiv. Found it, it's the divisibility check in perfpow. With

Re: bdiv vs redc

2012-07-10 Thread Niels Möller
Torbjorn Granlund t...@gmplib.org writes: What are your plans with these bdiv changes, do you intend to debug things, or are you hoping that I and Marco debug it? ;-) I don't think I have a real plan. But the most practical way forward I see is to 1. Adapt bdiv_mu to the new convention by

Re: bdiv vs redc

2012-06-29 Thread Niels Möller
Torbjorn Granlund t...@gmplib.org writes: hi = mpn_addmul_1 (np, dp, dn, q) hi += cy cy = hi cy // can this be true? Unless something interesting can be inferred from the hypothesis previous iteration generated a carry, I'm pretty sure carry is possible here. An

Re: bdiv vs redc

2012-06-28 Thread Torbjorn Granlund
Torbjorn Granlund t...@gmplib.org writes: iterate nn-dn times q = np[0] * dinv hi = mpn_addmul_1 (np, dp, dn, q) hi += cy cy = hi cy // can this be true? hi += np[dn] cy2 += hi np[dn] np[dn] = hi np += 1; I wrote a possible start

bdiv vs redc

2012-06-27 Thread Niels Möller
Let us review the differences between bdiv and redc, and focus on the balanced case (since that's what redc uses). We have an odd modulo M of size n limbs, and a numerator U of size 2n limbs. bdiv computes Q = U M^-1 (mod B^n) R = B^{-n} (U - Q M) where Q is n limbs, and R is n limbs plus