[Hello Robert -- are you in Bristol yet?] I'm puzzled now. My comment on inefficiency and Steffen's were based on the code > for i from 0 <= i <= n/2: > if (i*i) % n == self.ivalue: > return self._new_c(i) >
but now when I do a.sqrt?? I do not see that code. Steffen, are you using 2.8.12? John On 09/11/2007, Steffen <[EMAIL PROTECTED]> wrote: > > Thx, I am wondering why I did not try the command a.sqrt?? on my own. > > However, it seems as the implemented algorithm is not the most > efficient one. My result from a.sqrt?? from the latest release: > > def sqrt(self, extend=True, all=False): > cdef int_fast32_t i, n = self.__modulus.int32 > if n > 100: > moduli = self._parent.factored_order() > # Unless the modulus is tiny, test to see if we're in the > really > # easy case of n prime, n = 3 mod 4. > if n > 100 and n % 4 == 3 and len(moduli) == 1 and moduli[0] > [1] == 1: > if jacobi_int(self.ivalue, self.__modulus.int32) == 1: > # it's a non-zero square, sqrt(a) = a^(p+1)/4 > i = mod_pow_int(self.ivalue, (self.__modulus.int32+1)/ > 4, n) > if i > n/2: > i = n-i > if all: > return [self._new_c(i), self._new_c(n-i)] > else: > return self._new_c(i) > elif self.ivalue == 0: > return self > elif not extend: > raise ValueError, "self must be a square" > # Now we use a heuristic to guess whether or not it will > # be faster to just brute-force search for squares in a c > loop... > # TODO: more tuning? > elif n <= 100 or n / (1 << len(moduli)) < 5000: > if all: > return [self._new_c(i) for i from 0 <= i < n if (i*i) > % n == self.ivalue] > else: > for i from 0 <= i <= n/2: > if (i*i) % n == self.ivalue: > return self._new_c(i) > if not extend: > raise ValueError, "self must be a square" > # Either it failed but extend was True, or the generic > algorithm is better > return IntegerMod_abstract.sqrt(self, extend=extend, all=all) > > The easier %4 == 3 case seems to be implemented efficiently, but the > %4 == 1 not. The algo from Tonelli and Shanks might be a good solution > here. Any thoughts on other/better algorithm? > > Steffen > > > > > On 9 Nov., 18:24, "John Cremona" <[EMAIL PROTECTED]> wrote: > > Use the ?? operator to see the algorithm: > > > > sage: a=GF(next_prime(10^6)).random_element()^2; > > sage: a.sqrt?? > > > > Type: builtin_function_or_method > > Base Class: <type 'builtin_function_or_method'> > > String Form: <built-in method sqrt of > > sage.rings.integer_mod.IntegerMod_int64 object at 0x99055a4> > > Namespace: Interactive > > Source: > > def sqrt(self, extend=True, all=False): > > r""" > > Returns square root or square roots of self modulo n. > > > > INPUT: > > extend -- bool (default: True); if True, return a square > > root in an extension ring, if necessary. Otherwise, > > raise a ValueError if the square is not in the base ring. > > all -- bool (default: False); if True, return *all* square > > roots of self, instead of just one. > > > > ALGORITHM: Calculates the square roots mod $p$ for each of the > > primes $p$ dividing the order of the ring, then lifts them > > p-adically and uses the CRT to find a square root mod $n$. > > > > See also \code{square_root_mod_prime_power} and > > \code{square_root_mod_prime} (in this module) for more > > algorithmic details. > > > > EXAMPLES: > > sage: mod(-1, 17).sqrt() > > 4 > > sage: mod(5, 389).sqrt() > > 86 > > sage: mod(7, 18).sqrt() > > 5 > > sage: a = mod(14, 5^60).sqrt() > > sage: a*a > > 14 > > sage: mod(15, 389).sqrt(extend=False) > > Traceback (most recent call last): > > ... > > ValueError: self must be a square > > sage: Mod(1/9, next_prime(2^40)).sqrt()^(-2) > > 9 > > sage: Mod(1/25, next_prime(2^90)).sqrt()^(-2) > > 25 > > > > sage: a = Mod(3,5); a > > 3 > > sage: x = Mod(-1, 360) > > sage: x.sqrt(extend=False) > > Traceback (most recent call last): > > ... > > ValueError: self must be a square > > sage: y = x.sqrt(); y > > sqrt359 > > sage: y.parent() > > Univariate Quotient Polynomial Ring in sqrt359 over Ring > > of integers modulo 360 with modulus x^2 + 1 > > sage: y^2 > > 359 > > > > We compute all square roots in several cases: > > sage: R = Integers(5*2^3*3^2); R > > Ring of integers modulo 360 > > sage: R(40).sqrt(all=True) > > [20, 160, 200, 340] > > sage: [x for x in R if x^2 == 40] # Brute force verification > > [20, 160, 200, 340] > > sage: R(1).sqrt(all=True) > > [1, 19, 71, 89, 91, 109, 161, 179, 181, 199, 251, 269, > > 271, 289, 341, 359] > > sage: R(0).sqrt(all=True) > > [0, 60, 120, 180, 240, 300] > > > > sage: R = Integers(5*13^3*37); R > > Ring of integers modulo 406445 > > sage: v = R(-1).sqrt(all=True); v > > [78853, 111808, 160142, 193097, 213348, 246303, 294637, 327592] > > sage: [x^2 for x in v] > > [406444, 406444, 406444, 406444, 406444, 406444, 406444, 406444] > > sage: v = R(169).sqrt(all=True); min(v), -max(v), len(v) > > (13, 13, 104) > > sage: all([x^2==169 for x in v]) > > True > > > > Modulo a power of 2: > > sage: R = Integers(2^7); R > > Ring of integers modulo 128 > > sage: a = R(17) > > sage: a.sqrt() > > 23 > > sage: a.sqrt(all=True) > > [23, 41, 87, 105] > > sage: [x for x in R if x^2==17] > > [23, 41, 87, 105] > > > > """ > > if self.is_one(): > > if all: > > return list(self.parent().square_roots_of_one()) > > else: > > return self > > > > if not self.is_square_c(): > > if extend: > > > > and so on! > > > > John > > > > On 09/11/2007, Steffen <[EMAIL PROTECTED]> wrote: > > > > > > > > > > > > > Hi, I need to find square roots in GF(prime). I did it like this: > > > > > y = sqrt(GF(prime)(ySquare)) > > > > > So far I am not quite happy with the calculation periods and I would > > > like to know which algorithm is used. In my case is prime % 4 == 1, > > > which is the hardest case to find the square root mod p. I am sage > > > newbie and probably its again only a command that I cant find. I tried > > > sqrt? and similar things, but only got the information that sqrt is a > > > symbolic function. Could somebody tell me which algo is implemented or > > > better how to find the implemented algo. > > > > > Cheers, Steffen > > > > -- > > John Cremona > > > > > -- John Cremona --~--~---------~--~----~------------~-------~--~----~ To post to this group, send email to sage-devel@googlegroups.com To unsubscribe from this group, send email to [EMAIL PROTECTED] For more options, visit this group at http://groups.google.com/group/sage-devel URLs: http://sage.scipy.org/sage/ and http://modular.math.washington.edu/sage/ -~----------~----~----~----~------~----~------~--~---