Opened 6 years ago
Closed 9 months ago
#19120 closed enhancement (fixed)
efficient algorithm to compute continued fraction of a sum of continued fractions
Reported by:  sstarosta  Owned by:  sstarosta 

Priority:  minor  Milestone:  sage9.3 
Component:  number theory  Keywords:  continued_fraction 
Cc:  vdelecroix, sstarosta  Merged in:  
Authors:  Miroslav Kovar  Reviewers:  Vincent Delecroix, Frédéric Chapoton 
Report Upstream:  N/A  Work issues:  
Branch:  8cafbb4 (Commits, GitHub, GitLab)  Commit:  8cafbb46fc9f1a5256e80f53ea37b6712d01652b 
Dependencies:  Stopgaps: 
Description (last modified by )
Implement a method which, taking continued fraction x as an input, computes the continued fraction of (a*x+b)/(c*x+d) using Gosper's algorithm. Call this method with proper arguments upon
sage: 3*x sage: x+1
(overlaps with the following TODO item of sage/rings/continued_fraction.py: Add Gosper’s algorithm to compute the continued fraction of (ax + b)/(cx + d) knowing the one of x (see Gosper (1972, http://www.inwap.com/pdp10/hbaker/hakmem/cf.html), Knuth (1998, TAOCP vol 2, Exercise 4.5.3.15), Fowler (1999). See also Liardet, P. and Stambul, P. “Algebraic Computation with Continued Fractions.” J. Number Th. 73, 92121, 1998.)
Change History (39)
comment:1 Changed 5 years ago by
 Branch set to u/mirgee/my_sage_1
comment:2 Changed 5 years ago by
 Commit set to 3fcc597e97e028f66e1fc1048338d6dd8f20a836
comment:3 Changed 5 years ago by
 Commit changed from 3fcc597e97e028f66e1fc1048338d6dd8f20a836 to 51f21d56f009aa01c8fd39a6409956ee935c0e9f
Branch pushed to git repo; I updated commit sha1. New commits:
51f21d5  Better check coefficients using Integer.

comment:4 Changed 5 years ago by
Hello,
Some remarks:
 You would better edit the description of this ticket to explain what you want to do.
 Your branch does not merge cleanly with the last beta. The reason is that you added a line in
lazy_list.pyx
:from sage.rings.all import Infinity
. It is always better to base your work on the latest available version of the source code. You can consult the developer manual.
 The commit message should reflect what is inside the commit. Not your mood at the time you did the commit.
 The file with gosper algorithm should not be in
sage/misc/
. Put it in the same directory ascontinued_fraction.py
. And you would better call itcontinued_fraction_gosper.py
.
 In
continued_fraction.py
you added the following in the comparison codeI understand why you did it but it is definitely not a solution. What you can do is to raise a@@ 530,6 +530,8 @@ class ContinuedFraction_base(SageObject): if a == ZZ_0 and b == ZZ_0 and i: # rational case return 0 i += 1 + if i == 300: + return 0
RuntimeError
after a certain limit.
 All methods must be documented as explained in details in the developer manual.
Vincent
comment:5 Changed 5 years ago by
And last but not the least: it would be very nice to have Gosper algorithm within Sage!!
comment:6 followup: ↓ 10 Changed 5 years ago by
Hi,
thanks for the great input.
As for 5.: We need a way to compare instances of _infinite, and even though this is not a valid test, if enough terms are the same, the two continued fractions could be regarded as equal for most practical purposes. Another option would be to compare values, as in _real, but the problem is that
 value() method is by documentation optional,
 in case of _infinite, this would delegate the computation of value to RLF, which itself may use comparison,
 this solution would effectively do the same thing  compute the value of the operands with certain precision and then compare the intervals.
What do you think?
comment:7 Changed 5 years ago by
 Description modified (diff)
comment:8 Changed 5 years ago by
 Description modified (diff)
comment:9 Changed 5 years ago by
 Commit changed from 51f21d56f009aa01c8fd39a6409956ee935c0e9f to fa1173a4a5d7a2570212b7dbaa62faf7a9669b37
Branch pushed to git repo; I updated commit sha1. New commits:
fa1173a  Moved and renamed the continued_fraction_gosper. Added some documentation.

comment:10 in reply to: ↑ 6 Changed 5 years ago by
Replying to mirgee:
Hi,
thanks for the great input.
As for 5.: We need a way to compare instances of _infinite, and even though this is not a valid test, if enough terms are the same, the two continued fractions could be regarded as equal for most practical purposes. Another option would be to compare values, as in _real, but the problem is that
 value() method is by documentation optional,
 in case of _infinite, this would delegate the computation of value to RLF, which itself may use comparison,
 this solution would effectively do the same thing  compute the value of the operands with certain precision and then compare the intervals.
What do you think?
Everything would be better than an silent approximate comparison. If equality can not be decided just raise an error (see Python exceptions). I think that in that case an ArithmeticError
would make sense.
If you want to check equality up to some precision, you should do
sage: abs(cf1  cf2) < 2**100
or
sage: R = RealIntervalField(256) sage: R(cf1).overlaps(R(cf2))
comment:11 Changed 5 years ago by
For the sake of this ticket I would not implement something like 1+cf
or 3*cf
which needs some work. Having a method apply_homography
would be enough for a first ticket. In order to include this first part in Sage you need to
 make it compatible with the last version of Sage source code (see 2. in 4). By the way why did you make this change in
lazy_list.pyx
?
 add doctests to each method or function
comment:12 Changed 5 years ago by
The change in lazy_list()
is unnecessary and can be deleted. Personally I don't think it is necessary to add doctests for the simple helper functions, but I will comply with the standards and add them, no problem. I will also look into the the comparison.
comment:13 Changed 5 years ago by
 Commit changed from fa1173a4a5d7a2570212b7dbaa62faf7a9669b37 to fea8eb997260a2a4d056fc26e4e72752a2e13702
comment:14 Changed 5 years ago by
 Commit changed from fea8eb997260a2a4d056fc26e4e72752a2e13702 to 92b8aa03f51c1d1abef29b5b4e2f657d267db9eb
Branch pushed to git repo; I updated commit sha1. New commits:
92b8aa0  Added more unit tests, replaced overwritten method quotients() with get preperiod length, so that all doctests pass. Added a few more doctests. Now all methods are accompanied with doctests.

comment:15 Changed 5 years ago by
 Commit changed from 92b8aa03f51c1d1abef29b5b4e2f657d267db9eb to e063f2fbdb0d43ee4d75b13d23c1f88eb1d125cb
Branch pushed to git repo; I updated commit sha1. New commits:
e063f2f  Using a set of tuples instead of list of dictionaries for state tracking, should improve performance. Using lazy_list instead of Word for instances of infinite continued fractions.

comment:16 Changed 2 years ago by
 Branch changed from u/mirgee/my_sage_1 to public/ticket/19120
 Commit changed from e063f2fbdb0d43ee4d75b13d23c1f88eb1d125cb to 3308c92d78c4ef22373f2cf7cccaf5d2273f6a79
rebase, squashed, py3compatible
New commits:
3308c92  Implement Gosper algorithm for homographic action on continued fractions

comment:17 Changed 2 years ago by
 Commit changed from 3308c92d78c4ef22373f2cf7cccaf5d2273f6a79 to 031088f7ca9b4a18530fb972729472f8e8d71a27
Branch pushed to git repo; I updated commit sha1. New commits:
031088f  fix for py2: also define next

comment:18 Changed 10 months ago by
 Keywords continued_fraction added; continued fractions removed
comment:19 Changed 10 months ago by
 Commit changed from 031088f7ca9b4a18530fb972729472f8e8d71a27 to aeb6da12b14ff152dd1ec02abdc0b8cfa35970af
Branch pushed to git repo; I updated commit sha1. New commits:
aeb6da1  Merge branch 'public/ticket/19120' in 9.2.b12

comment:20 Changed 10 months ago by
 Status changed from new to needs_review
comment:21 Changed 10 months ago by
 Milestone changed from sagepending to sage9.3
comment:22 Changed 10 months ago by
 Commit changed from aeb6da12b14ff152dd1ec02abdc0b8cfa35970af to fe2cd08a88e946b7dd9878d3c13d951fb159f376
Branch pushed to git repo; I updated commit sha1. New commits:
fe2cd08  fix one returns

comment:23 Changed 10 months ago by
 Status changed from needs_review to needs_work
not suitable yet. I will provide a commit.
comment:24 followup: ↓ 25 Changed 10 months ago by
Author? (name in the git log is John Doe <miroslav.kovar@concur.com>
)
comment:25 in reply to: ↑ 24 ; followup: ↓ 27 Changed 10 months ago by
Replying to vdelecroix:
Author? (name in the git log is
John Doe <miroslav.kovar@concur.com>
)
I am the author: Miroslav Kovar (miroslavkovar@…)
comment:26 Changed 10 months ago by
 Commit changed from fe2cd08a88e946b7dd9878d3c13d951fb159f376 to 81e069389c5944f1b26ecc7aac111bda7bc9b1c2
Branch pushed to git repo; I updated commit sha1. New commits:
81e0693  fixes, documentation, copyright, optimization

comment:27 in reply to: ↑ 25 Changed 10 months ago by
Replying to mirgee:
Replying to vdelecroix:
Author? (name in the git log is
John Doe <miroslav.kovar@concur.com>
)I am the author: Miroslav Kovar (miroslavkovar@…)
Thanks. I put your name and mail in the copyright header of the file.
comment:28 Changed 10 months ago by
 Reviewers set to Vincent Delecroix
 Status changed from needs_work to needs_review
comment:29 Changed 10 months ago by
patchbot plugins complain
comment:30 Changed 10 months ago by
Le reference [LS1998] a une allure bizarre, avec un P isolé
comment:31 Changed 10 months ago by
 Commit changed from 81e069389c5944f1b26ecc7aac111bda7bc9b1c2 to 3323f09f75c96787cd0d7496aa4b761cd32fc1a0
Branch pushed to git repo; I updated commit sha1. New commits:
3323f09  enoding, import, references

comment:32 Changed 10 months ago by
Merci. Pyflakes est encore pas content:
src/sage/rings/continued_fraction_gosper.py:37:1 'sage.rings.real_mpfr.RR' imported but unused
comment:33 Changed 10 months ago by
 Commit changed from 3323f09f75c96787cd0d7496aa4b761cd32fc1a0 to 7161eafc6a96ba597d3b6165257cac1070f56ac1
Branch pushed to git repo; I updated commit sha1. New commits:
7161eaf  remove unused import

comment:34 Changed 10 months ago by
 Status changed from needs_review to needs_work
There is something wrong in the way infinity is handled. In the situation where we apply x > 1 / (qx  p)
to x=p/q
Gosper iterator is empty. But currently continued_fraction([])
is wrongly 0
.
comment:35 Changed 10 months ago by
 Commit changed from 7161eafc6a96ba597d3b6165257cac1070f56ac1 to 8cafbb46fc9f1a5256e80f53ea37b6712d01652b
Branch pushed to git repo; I updated commit sha1. New commits:
8cafbb4  fix infinities in continued fraction and do not forward values

comment:37 Changed 10 months ago by
ok, this seems to be ready for a positive review ?
comment:38 Changed 10 months ago by
 Reviewers changed from Vincent Delecroix to Vincent Delecroix, Frédéric Chapoton
 Status changed from needs_review to positive_review
comment:39 Changed 9 months ago by
 Branch changed from public/ticket/19120 to 8cafbb46fc9f1a5256e80f53ea37b6712d01652b
 Resolution set to fixed
 Status changed from positive_review to closed
Branch pushed to git repo; I updated commit sha1. Last 10 new commits:
Solved the big bug with the negative value of infinite continued fraction! :)
Don't know what changed. Backing up.
Added a few lines to check for problem with singularity.
Latest version... Don't know what changed :)
Done some code refactoring, tried to make the code simpler and more readable. I hope I didn't break it.
Of course I broke something. Bug for infinite fractions fixed and code shortened.
Added a test and fixed associated bugs (hopefully).
Added a test and fixed associated bugs (hopefully).
Added check for arguments' validity.
Moved Gosper iterator.