[MUG] isqrt command

[MUG] Re: isqrt command
Author: Maple User Group    Posted: Wed, 27 Nov 2002 10:33:13 -0500
>> From: Maple User Group "maple_gr" -=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=- Date: Mon, 25 Nov 2002 07:28:27 -0500 (EST) From: Carl Devore "devore" To: "maple-list" Subject: isqrt command On Fri, 22 Nov 2002, Adalberto Ayjara Dornelles Filho wrote: | I would know how "isqrt" command work internally. It is a built-in | function and no algorithm are available. Anybody can tell me that or | send me a reference? I cannot tell you what the builtin isqrt does exactly, but if you are merely interested in an algorithm rather than the precise implementation, I think the following is pretty good. It uses pure integer arithmetic: Isqrt:= proc(n::posint) local r,r0,`r^2+r`; r0:= 0; r:= Scale2(1,iquo(ilog2(n),2)); while CopySign(r-r0,1) > 1 do r0:= r; r:= iquo(r0^2+n, Scale2(r0,1)) od; `r^2+r`:= r*(r+1); while `r^2+r` < n do r:= r+1; `r^2+r`:= `r^2+r`+Scale2(r,1) od; r end proc: Of course, isqrt will be much faster because it is builtin. -=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=- From: Adalberto Ayjara Dornelles Filho "ayjara" To: "maple-list" Carl Devore "devore" Date: Tue, 26 Nov 2002 16:42:21 -0200 Subject: isqrt command Thank you Mr. Devore, It's seems a Newton-Raphson iteration with a good start guess r. I'm right? But how isqrt work exactly? There are some bit manipulation? Thanks, **************************************************** Adalberto Ayjara Dornelles Filho - "ayjara" Universidade de Caxias do Sul - www.ucs.br Departamento de Matematica e Estatmstica Rua Francisco Getulio Vargas, 1130 95070-970 Caxias do Sul RS Fone (54) 212-1133 R.2103 **************************************************** -=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=- From: Adalberto Ayjara Dornelles Filho "ayjara" To: "maple-list" Carl Devore "devore" Date: Tue, 26 Nov 2002 16:12:03 -0200 Subject: isqrt command bug? Hi, It's MY computer failure or a MAPLE bug?? > isqrt(67890); 261 > floor(sqrt(67890)); 260 Bye, **************************************************** Adalberto Ayjara Dornelles Filho - "ayjara" Universidade de Caxias do Sul - www.ucs.br Departamento de Matematica e Estatistica Rua Francisco Getulio Vargas, 1130 95070-970 Caxias do Sul RS Fone (54) 212-1133 R.2103 **************************************************** -=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=-=*=- Date: Wed, 27 Nov 2002 01:57:05 -0500 (EST) From: Carl Devore "devore" To: Adalberto Ayjara Dornelles Filho "ayjara" Subject: several messages On Tue, 26 Nov 2002, Adalberto Ayjara Dornelles Filho wrote: > It's seems a Newton-Raphson iteration with a good start guess r. I'm right? Yes, but for the case of square roots of integers I believe that the algorithm I gave has been known since ancient times. > But how isqrt work exactly? I don't know. Someone from WMI would have to answer that. My guess would be that it does the same thing as the algorithm I gave except that (1) it probably falls back on floating point when arguemt is word sized, (2) iquo's with a divisor of 2 are done as bit shifts, and (3) the main iteration is broken into two steps r:= iquo(r0 + iquo(n,r0), 2); which avoids the need to square. Doing it this way is a bit slower in high-level Maple code because there are two calls to iquo. > There are some bit manipulation? I coded it to show the bit manipulations. ilog2, CopySign, and Scale2 are low-level bit manipulation commands. > It's MY computer failure or a MAPLE bug?? > > isqrt(67890); > 261 > > floor(sqrt(67890)); > 260 I believe that isqrt(n) returns round(sqrt(n)). The documentation only claims that it returns either floor(sqrt(n)) or ceil(sqrt(n)). However, I've tested it for hundreds of thousands of integers, both word-sized and very large, and in every case it returned round(sqrt(n)). I think that it would be more useful if it always returned floor(sqrt(n)).

Previous by date: [MUG] Re: animation problems, Maple User Group
Next by date: [MUG] why the plot is not consistent with evaluated value?, Xiaoyan Li
Previous thread: [MUG] more about implicitdiff:complexity and recursiveness,  Charles James Leonardo Quarra Cappiello
Next thread: [MUG] why the plot is not consistent with evaluated value?, Xiaoyan Li