P-Adic square roots: Difference between revisions

Content added Content deleted
(→‎{{header|Wren}}: Updated in line with changes to FB example.)
(max. p-adic precision with standard data types)
Line 47: Line 47:




const emx = 64
const emx = 48
'exponent maximum
'exponent maximum


const amx = 6000
const amx = 700000
'tentative argument maximum
'tentative argument maximum


'------------------------------------------------
'------------------------------------------------
const Mxd = cdbl(2)^53 - 1
'max. float64 integer

const Pmax = 32749
const Pmax = 32749
'max. prime < 2^15
'max. prime < 2^15
Line 93: Line 96:
function padic.sqrt (byref g as ratio, byval sw as integer) as integer
function padic.sqrt (byref g as ratio, byval sw as integer) as integer
dim as longint a = g.a, b = g.b
dim as longint a = g.a, b = g.b
dim as longint f, q, x, pk
dim as longint q, x, pk
dim as long f1, r, s, t
dim as long f1, r, s, t
dim i as integer
dim i as integer, f as double
sqrt = 0
sqrt = 0


Line 152: Line 155:
'find root for small p
'find root for small p
for r = 1 to p1
for r = 1 to p1
f = b * r * r - a
q = b * r * r - a
if f mod p = 0 then exit for
if q mod p = 0 then exit for
next r
next r


Line 181: Line 184:
'evaluate f(x)
'evaluate f(x)
#macro evalf(x)
#macro evalf(x)
f = b * x * x - a
f = b * x * cdbl(x / pk)
q = f \ pk
f -= cdbl(a / pk)
'overflow
'overflow
if f - q * pk then exit for
if f > Mxd then exit for
q = clngint(f)
#endmacro
#endmacro


Line 229: Line 233:
function padic.crat (byval sw as integer) as ratio
function padic.crat (byval sw as integer) as ratio
dim as integer i, j, t = min(v, 0)
dim as integer i, j, t = min(v, 0)
dim r as ratio, f as double
dim as longint s, pk, pm
dim as long x, y, pk, p1
dim as long q, x, y
dim as longint q, s
dim as double f, h
dim r as ratio


'weighted digit sum
'weighted digit sum
s = 0: pk = 1
s = 0: pk = 1
for i = t to k - 1 + v
for i = t to k - 1 + v
p1 = pk: pk *= p
pm = pk: pk *= p


if pk \ p1 - p then
if pk \ pm - p then
'overflow
'overflow
pk = p1: exit for
pk = pm: exit for
end if
end if


s += d(i) * p1 '(mod pk)
s += d(i) * pm '(mod pk)
next i
next i


Line 249: Line 254:
dim as longint m(1) = {pk, s}
dim as longint m(1) = {pk, s}
dim as longint n(1) = {0, 1}
dim as longint n(1) = {0, 1}
⚫
'norm(v)^2
h = cdbl(s) * s + 1
i = 0: j = 1
i = 0: j = 1
⚫
s = s * s + 1 ' norm(v)^2


'Lagrange's algorithm
'Lagrange's algorithm
do
do
f = m(i) * cdbl(m(j) / s)
f = m(i) * (m(j) / h)
f += n(i) * cdbl(n(j) / s)
f += n(i) * (n(j) / h)


'Euclidean step
'Euclidean step
Line 262: Line 268:
n(i) -= q * n(j)
n(i) -= q * n(j)


q = s
f = h
s = m(i) * m(i) + n(i) * n(i)
h = cdbl(m(i)) * m(i)
h += cdbl(n(i)) * n(i)
'compare norms
'compare norms
if s < q then
if h < f then
'interchange vectors
'interchange vectors
swap i, j
swap i, j
Line 375: Line 382:
data 10496,497, 2,19
data 10496,497, 2,19


data 3141,5926, 3,15
data -577215,664901, 3,23
data 2718,281, 3,13
data 15403,26685, 3,18


data -1,1, 5,8
data -1,1, 5,8
Line 383: Line 390:


data 2,1, 7,8
data 2,1, 7,8
data -2645,28518, 7,9
data 11696,621467, 7,11
data 3029,4821, 7,9
data -27764,11521, 7,11
data 379,449, 7,8
data -27584,12953, 7,11


data 717,8, 11,7
data -166420,135131, 11,11
data 1414,213, 41,5
data 14142,135623, 5,15
data -255,256, 257,3
data -255,256, 257,3


Line 423: Line 430:
loop
loop


end
system
</lang>
</lang>

{{out|Examples}}
{{out|Examples}}
<pre>
<pre>
Line 478: Line 486:




3141/5926 + O(3^15)
-577215/664901 + O(3^23)
lift: 3406966 mod 3^15
lift: 44765673211 mod 3^23
sqrt +/-
sqrt +/-
... 2 0 1 0 2 0 0 2 1 1 0 2 2 1 0
... 1 0 2 1 1 1 2 2 1 0 1 1 1 2 1 1 0 1 0 2 0 1 0
... 0 2 1 2 0 2 2 0 1 1 2 0 0 2 0
... 1 2 0 1 1 1 0 0 1 2 1 1 1 0 1 1 2 1 2 0 2 2 0
sqrt^2
sqrt^2
2 1 1 1 2 2 0 0 0 2 0 1 1 0 0
0 1 2 1 1 0 0 2 1 1 0 0 1 1 0 2 0 1 1 0 1 0 0
-577215/664901
3141/5926




2718/281 + O(3^13)
15403/26685 + O(3^18)
lift: 693355 mod 3^13
lift: 238961782 mod 3^18
sqrt +/-
sqrt +/-
... 0 2 2 0 2 0 0 0 2 2 1 1 0
... 1 2 1 1 2 2 1 2 2 1 1 1 2 2 1 1 0. 1
... 2 0 0 2 0 2 2 2 0 0 1 2 0
... 1 0 1 1 0 0 1 0 0 1 1 1 0 0 1 1 2. 2
sqrt^2
sqrt^2
2 2 0 0 1 2 2 0 2 2 1 0 0
1 2 1 2 0 1 2 2 0 1 1 0 1 2 2 2. 0 1
15403/26685
2718/281




Line 538: Line 546:




-2645/28518 + O(7^9)
11696/621467 + O(7^11)
lift: 39082527 mod 7^9
lift: 967215488 mod 7^11
sqrt +/-
sqrt +/-
... 6 5 3 1 2 4 1 4. 1
... 3 2 6 5 3 1 2 4 1 4. 1
... 0 1 3 5 4 2 5 2. 6
... 3 4 0 1 3 5 4 2 5 2. 6
sqrt^2
sqrt^2
1 1 3 3 4 4 5. 1 1
3 3 1 1 3 3 4 4 5. 1 1
11696/621467
-2645/28518




3029/4821 + O(7^9)
-27764/11521 + O(7^11)
lift: 31104424 mod 7^9
lift: 453209984 mod 7^11
sqrt +/-
sqrt +/-
... 5 2 5 2 4 5 3 1 1
... 1 4 1 4 2 1 3 5 6 2 3
... 1 4 1 4 2 1 3 5 6
... 5 2 5 2 4 5 3 1 0 4 4
sqrt^2
sqrt^2
6 5 6 4 1 3 0 2 1
5 1 0 2 5 5 5 3 6 6 2
-27764/11521
3029/4821




379/449 + O(7^8)
-27584/12953 + O(7^11)
lift: 555087 mod 7^8
lift: 1239987302 mod 7^11
sqrt +/-
sqrt +/-
... 0 4 5 0 1 2 2 1
... 4 2 5 0 4 5 0 1 2 2 1
... 6 2 1 6 5 4 4 6
... 2 4 1 6 2 1 6 5 4 4 6
sqrt^2
sqrt^2
5 3 1 2 4 1 4 1
3 2 6 5 3 1 2 4 1 4 1
-27584/12953
379/449




717/8 + O(11^7)
-166420/135131 + O(11^11)
lift: 5280093 mod 11^7
lift: 34219617965 mod 11^11
sqrt +/-
sqrt +/-
... 2 10 8 7 0 1 5
... 1 3 5 7 0 0 9 10 5 0 5
... 8 0 2 3 10 9 6
... 9 7 5 3 10 10 1 0 5 10 6
sqrt^2
sqrt^2
1 4 1 4 2 1 3
1 4 1 4 2 1 3 5 6 2 3
-166420/135131
717/8




1414/213 + O(41^5)
14142/135623 + O(5^15)
lift: 25564902 mod 41^5
lift: 236397452 mod 5^15
sqrt +/-
sqrt +/-
... 9 1 38 6 8
... 0 0 0 4 4 1 0 0 4 2 0 4 3 0 2
... 31 39 2 34 33
... 4 4 4 0 0 3 4 4 0 2 4 0 1 4 3
sqrt^2
sqrt^2
12 36 31 15 23
4 2 1 3 3 4 3 4 3 4 2 3 2 0 4
14142/135623
1414/213




Retrieved from "https://rosettacode.org/wiki/P-Adic_square_roots"