mersenneforum.org  

Go Back   mersenneforum.org > Math Stuff > Computer Science & Computational Number Theory

Reply
 
Thread Tools
Old 2009-09-28, 20:51   #1
bsquared
 
bsquared's Avatar
 
"Ben"
Feb 2007

5·727 Posts
Default integer nth root code

I'm going back to look at my integer nth root code per this, and I'm running into a problem.

I'm trying variations like so:
nroot(3^180,150)

The actual integer root is 3, and my code's first guess is 4. So far so good. I then compute the next term in the recurrence:

x_{k+1} = \frac{1}{n}\left[{x_{k} (n-1) + \frac{A}{x_k^{n-1}}}\right]

and since 4^149 > 3^180, the fraction in the bracket is zero and the next term in the sequence is 3. Great! But we can't stop there because this new term is not the same as the old one (hasn't converged).

So now with x_k = 3, we get a huge term for the fraction in the bracket, the next term in the sequence is now something like 4 trillion, and the sequence never converges.

What do I do in situations like this, when even a correct guess of the integer nth root produces a wildly inaccurate next guess in the recurrence?
bsquared is offline   Reply With Quote
Old 2009-09-29, 04:46   #2
wblipp
 
wblipp's Avatar
 
"William"
May 2003
New Haven

22·593 Posts
Default

This isn't modular arithemetic - you are looking for the rounded off (or perhaps down) value of the real root - right?

Section 9 of Numerical Recipes says you should always bracket a root and never let your iteration method get outside of the best known bounds at any stage. I've had very good experience with the Van Wijngaarden-Dekker-Brent method of Section 9.3 in a wide variety of numerical root finding problems.

If I understand the problem correctly, from a bracket of (3,4) you should move into an endgame of picking 3 or 4, not iterating further.

William

Last fiddled with by wblipp on 2009-09-29 at 04:47
wblipp is offline   Reply With Quote
Old 2009-09-29, 05:06   #3
paleseptember
 
paleseptember's Avatar
 
Jun 2008
Wollongong, .au

3×61 Posts
Default

Quote:
Originally Posted by wblipp View Post
This isn't modular arithemetic - you are looking for the rounded off (or perhaps down) value of the real root - right?
Yes, the integer root is the floor function of the nth root.

Sorry, my programming skills are minimal, so I can't actually assist with the mechanics of the computation.
paleseptember is offline   Reply With Quote
Old 2009-09-29, 12:46   #4
bsquared
 
bsquared's Avatar
 
"Ben"
Feb 2007

5×727 Posts
Default

Quote:
Originally Posted by wblipp View Post
This isn't modular arithemetic - you are looking for the rounded off (or perhaps down) value of the real root - right?
Yes, the floor of the real root, as paleseptember says.

Quote:
Originally Posted by wblipp View Post
Section 9 of Numerical Recipes says you should always bracket a root and never let your iteration method get outside of the best known bounds at any stage. I've had very good experience with the Van Wijngaarden-Dekker-Brent method of Section 9.3 in a wide variety of numerical root finding problems.

If I understand the problem correctly, from a bracket of (3,4) you should move into an endgame of picking 3 or 4, not iterating further.

William
I'll get out my NR and take a look, thanks. Basically, this seems to jive with what I was thinking... I know that the answer should be a certain number of bits (bits(input) / n, where n is the root, to be exact), so if the iteration goes outside that range after being inside that range, stop iterating and decide between the last couple values. I've already got code in place to detect oscillations by keeping a history of the last few iterations, so this shouldn't be hard to implement.

Last fiddled with by bsquared on 2009-09-29 at 12:46 Reason: value to range...
bsquared is offline   Reply With Quote
Old 2009-09-29, 13:52   #5
jasonp
Tribal Bullet
 
jasonp's Avatar
 
Oct 2004

DD916 Posts
Default

In the specific case of needing the closest integer to the root, Crandall and Pomerance show that you can stop the iteration when the sequence of root estimates stops decreasing. That leads to a skewed loop where the convergence check is in the middle, after A/x_k^(n-1) has been computed, and x_(k+1) only gets computed if the convergence test fails.
jasonp is offline   Reply With Quote
Old 2009-09-29, 19:39   #6
Gammatester
 
Gammatester's Avatar
 
Mar 2009

2·19 Posts
Default

As Jason said, you can stop the iteration if the sequence stops decreasing. But you have to choose a starting value greater than the real root. Normalyy you start with x0 = 2^ceil(bitsize(A)/n).

If n is large and the guesses are not very close to the root, convergence is only linear not quadratic. You should check the first few x_i for progress, otherwise it is better to use bisection steps until the iterates come close to the final result.

Wolfgang
Gammatester is offline   Reply With Quote
Old 2009-09-29, 19:54   #7
bsquared
 
bsquared's Avatar
 
"Ben"
Feb 2007

5×727 Posts
Default

Quote:
Originally Posted by Gammatester View Post
As Jason said, you can stop the iteration if the sequence stops decreasing. But you have to choose a starting value greater than the real root. Normalyy you start with x0 = 2^ceil(bitsize(A)/n).
Yep. This is implemented now and seems to work fine.

Quote:
Originally Posted by Gammatester View Post
If n is large and the guesses are not very close to the root, convergence is only linear not quadratic. You should check the first few x_i for progress, otherwise it is better to use bisection steps until the iterates come close to the final result.

Wolfgang
The Van Wijngaarden-Dekker-Brent method that wblipp suggested includes allowances for that, I think. Implementing that method in YAFU is an exercise for another day. Stuff like this is very interesting and a fun diversion but I don't believe many people use YAFU for everyday bigint calculations for which they need the fastest and most bulletproof methods. I'm still concentrating on factoring improvements for now. That said, by all means continue to find and report these things as they come up, it's very much appreciated!
bsquared is offline   Reply With Quote
Reply

Thread Tools


Similar Threads
Thread Thread Starter Forum Replies Last Post
k*b^n+/-c where b is an integer greater than 2 and c is an integer from 1 to b-1 jasong Miscellaneous Math 5 2016-04-24 03:40
Square root of 3 Damian Math 3 2010-01-01 01:56
Reference code for all-integer convolutions jasonp Software 17 2009-01-29 02:25
Cube root mfgoode Homework Help 10 2007-10-05 04:12
invalid root junky NFSNET Discussion 0 2004-03-26 02:14

All times are UTC. The time now is 21:27.


Mon Aug 15 21:27:55 UTC 2022 up 39 days, 16:15, 1 user, load averages: 1.21, 1.37, 1.37

Powered by vBulletin® Version 3.8.11
Copyright ©2000 - 2022, Jelsoft Enterprises Ltd.

This forum has received and complied with 0 (zero) government requests for information.

Permission is granted to copy, distribute and/or modify this document under the terms of the GNU Free Documentation License, Version 1.2 or any later version published by the Free Software Foundation.
A copy of the license is included in the FAQ.

≠ ± ∓ ÷ × · − √ ‰ ⊗ ⊕ ⊖ ⊘ ⊙ ≤ ≥ ≦ ≧ ≨ ≩ ≺ ≻ ≼ ≽ ⊏ ⊐ ⊑ ⊒ ² ³ °
∠ ∟ ° ≅ ~ ‖ ⟂ ⫛
≡ ≜ ≈ ∝ ∞ ≪ ≫ ⌊⌋ ⌈⌉ ∘ ∏ ∐ ∑ ∧ ∨ ∩ ∪ ⨀ ⊕ ⊗ 𝖕 𝖖 𝖗 ⊲ ⊳
∅ ∖ ∁ ↦ ↣ ∩ ∪ ⊆ ⊂ ⊄ ⊊ ⊇ ⊃ ⊅ ⊋ ⊖ ∈ ∉ ∋ ∌ ℕ ℤ ℚ ℝ ℂ ℵ ℶ ℷ ℸ 𝓟
¬ ∨ ∧ ⊕ → ← ⇒ ⇐ ⇔ ∀ ∃ ∄ ∴ ∵ ⊤ ⊥ ⊢ ⊨ ⫤ ⊣ … ⋯ ⋮ ⋰ ⋱
∫ ∬ ∭ ∮ ∯ ∰ ∇ ∆ δ ∂ ℱ ℒ ℓ
𝛢𝛼 𝛣𝛽 𝛤𝛾 𝛥𝛿 𝛦𝜀𝜖 𝛧𝜁 𝛨𝜂 𝛩𝜃𝜗 𝛪𝜄 𝛫𝜅 𝛬𝜆 𝛭𝜇 𝛮𝜈 𝛯𝜉 𝛰𝜊 𝛱𝜋 𝛲𝜌 𝛴𝜎𝜍 𝛵𝜏 𝛶𝜐 𝛷𝜙𝜑 𝛸𝜒 𝛹𝜓 𝛺𝜔