Calkin-Wilf sequence: Difference between revisions
m (→{{header|Factor}}: remove template) |
|||
(72 intermediate revisions by 36 users not shown) | |||
Line 1: | Line 1: | ||
{{ |
{{task}} |
||
The '''Calkin-Wilf sequence''' contains every nonnegative rational number exactly once. |
|||
{{Clarified-review|Task description updated, all examples need updating!}} |
|||
It can be calculated recursively as follows: |
|||
{{math|a<sub>1</sub>}} = {{math|1}} |
<big> {{math|a<sub>1</sub>}} = {{math|1}} </big> |
||
<big> {{math|a<sub>n+1</sub>}} = {{math|1/(2⌊a<sub>n</sub>⌋+1-a<sub>n</sub>)}} for n > 1 </big> |
|||
{{math|a<sub>n+1</sub>}} = {{math|1/(2⌊a<sub>n</sub>⌋+1-a<sub>n</sub>)}} for n > 1 |
|||
;Task part 1: |
|||
* Show on this page terms 1 through 20 of the Calkin-Wilf sequence.<br> |
|||
* Show on this page terms 1 through 20 of the Calkin-Wilf sequence. |
|||
To avoid floating point error, you may want to use a rational number data type. |
|||
<br>To avoid floating point error, you may want to use a rational number data type. |
|||
It is also possible, given a nonnegative rational number, to determine where it appears in the sequence without calculating the sequence. The procedure is to get the continued fraction representation of the rational and use it as the run-length encoding of the binary representation of the term number, beginning from the end of the continued fraction. |
|||
It is also possible, given a non-negative rational number, to determine where it appears in the sequence without calculating the sequence. The procedure is to get the continued fraction representation of the rational and use it as the run-length encoding of the binary representation of the term number, beginning from the end of the continued fraction. |
|||
It only works if the number of terms in the continued fraction is odd- use either of the two equivalent representations to achieve this: |
It only works if the number of terms in the continued fraction is odd- use either of the two equivalent representations to achieve this: |
||
{{math|[a<sub>0</sub>; a<sub>1</sub>, a<sub>2</sub>, ..., a<sub>n</sub>]}} = {{math|[a<sub>0</sub>; a<sub>1</sub>, a<sub>2</sub> ,..., a<sub>n</sub>-1, 1]}} |
<big> {{math|[a<sub>0</sub>; a<sub>1</sub>, a<sub>2</sub>, ..., a<sub>n</sub>]}} = {{math|[a<sub>0</sub>; a<sub>1</sub>, a<sub>2</sub> ,..., a<sub>n</sub>-1, 1]}} </big> |
||
Thus, for example, the fraction 9/4 has odd continued fraction representation {{math|2; 3, 1}}, giving a binary representation of 100011, which means 9/4 appears as the 35th term of the sequence. |
|||
;Example: |
|||
* Find the position of the number {{math|83116/51639}} in the Calkin-Wilf sequence. |
|||
The fraction '''9/4''' has odd continued fraction representation <big> {{math|2; 3, 1}},</big> giving a binary representation of '''100011''', |
|||
<br>which means '''9/4''' appears as the '''35th''' term of the sequence. |
|||
;See also: |
|||
;Task part 2: |
|||
* Find the position of the number <big>'''<sup>83116</sup>'''<big>'''/'''</big>'''<sub>51639</sub>'''</big> in the Calkin-Wilf sequence. |
|||
;Related tasks: |
|||
:* [[Fusc sequence]]. |
|||
;See also: |
|||
* Wikipedia entry: [[wp:Calkin%E2%80%93Wilf_tree|Calkin-Wilf tree]] |
* Wikipedia entry: [[wp:Calkin%E2%80%93Wilf_tree|Calkin-Wilf tree]] |
||
* [[Continued fraction]] |
* [[Continued fraction]] |
||
* [[Continued fraction/Arithmetic/Construct from rational number]] |
* [[Continued fraction/Arithmetic/Construct from rational number]] |
||
<br><br> |
|||
=={{header|11l}}== |
|||
{{trans|Nim}} |
|||
<syntaxhighlight lang="11l">T CalkinWilf |
|||
n = 1 |
|||
d = 1 |
|||
F ()() |
|||
V r = (.n, .d) |
|||
.n = 2 * (.n I/ .d) * .d + .d - .n |
|||
swap(&.n, &.d) |
|||
R r |
|||
print(‘The first 20 terms of the Calkwin-Wilf sequence are:’) |
|||
V cw = CalkinWilf() |
|||
[String] seq |
|||
L 20 |
|||
V (n, d) = cw() |
|||
seq.append(I d == 1 {String(n)} E n‘/’d) |
|||
print(seq.join(‘, ’)) |
|||
cw = CalkinWilf() |
|||
V index = 1 |
|||
L cw() != (83116, 51639) |
|||
index++ |
|||
print("\nThe element 83116/51639 is at position "index‘ in the sequence.’)</syntaxhighlight> |
|||
{{out}} |
|||
<pre> |
|||
The first 20 terms of the Calkwin-Wilf sequence are: |
|||
1, 1/2, 2, 1/3, 3/2, 2/3, 3, 1/4, 4/3, 3/5, 5/2, 2/5, 5/3, 3/4, 4, 1/5, 5/4, 4/7, 7/3, 3/8 |
|||
The element 83116/51639 is at position 123456789 in the sequence. |
|||
</pre> |
|||
=={{header|ALGOL 68}}== |
|||
Uses code from the [[Arithmetic/Rational]] and [[Continued fraction/Arithmetic/Construct from rational number]] tasks. |
|||
<syntaxhighlight lang="algol68">BEGIN |
|||
# Show elements 1-20 of the Calkin-Wilf sequence as rational numbers # |
|||
# also show the position of a specific element in the seuence # |
|||
# Uses code from the Arithmetic/Rational # |
|||
# & Continued fraction/Arithmetic/Construct from rational number tasks # |
|||
# Code from the Arithmetic/Rational task # |
|||
# ============================================================== # |
|||
MODE FRAC = STRUCT( INT num #erator#, den #ominator#); |
|||
PROC gcd = (INT a, b) INT: # greatest common divisor # |
|||
(a = 0 | b |: b = 0 | a |: ABS a > ABS b | gcd(b, a MOD b) | gcd(a, b MOD a)); |
|||
PROC lcm = (INT a, b)INT: # least common multiple # |
|||
a OVER gcd(a, b) * b; |
|||
PRIO // = 9; # higher then the ** operator # |
|||
OP // = (INT num, den)FRAC: ( # initialise and normalise # |
|||
INT common = gcd(num, den); |
|||
IF den < 0 THEN |
|||
( -num OVER common, -den OVER common) |
|||
ELSE |
|||
( num OVER common, den OVER common) |
|||
FI |
|||
); |
|||
OP + = (FRAC a, b)FRAC: ( |
|||
INT common = lcm(den OF a, den OF b); |
|||
FRAC result := ( common OVER den OF a * num OF a + common OVER den OF b * num OF b, common ); |
|||
num OF result//den OF result |
|||
); |
|||
OP - = (FRAC a, b)FRAC: a + -b, |
|||
* = (FRAC a, b)FRAC: ( |
|||
INT num = num OF a * num OF b, |
|||
den = den OF a * den OF b; |
|||
INT common = gcd(num, den); |
|||
(num OVER common) // (den OVER common) |
|||
); |
|||
OP - = (FRAC frac)FRAC: (-num OF frac, den OF frac); |
|||
# ============================================================== # |
|||
# end code from the Arithmetic/Rational task # |
|||
# code from the Continued fraction/Arithmetic/Construct from rational number task # |
|||
# ================================================================================# |
|||
# returns the quotient of numerator over denominator and sets # |
|||
# numerator and denominator to the next values for # |
|||
# the continued fraction # |
|||
PROC r2cf = ( REF INT numerator, REF INT denominator )INT: |
|||
IF denominator = 0 |
|||
THEN 0 |
|||
ELSE INT quotient := numerator OVER denominator; |
|||
INT prev numerator = numerator; |
|||
numerator := denominator; |
|||
denominator := prev numerator MOD denominator; |
|||
quotient |
|||
FI # r2cf # ; |
|||
# ====================================================================================# |
|||
# end code from the Continued fraction/Arithmetic/Construct from rational number task # |
|||
# Additional FRACrelated operators # |
|||
OP * = ( INT a, FRAC b )FRAC: ( num OF b * a ) // den OF b; |
|||
OP / = ( FRAC a, b )FRAC: ( num OF a * den OF b ) // ( num OF b * den OF a ); |
|||
OP FLOOR = ( FRAC a )INT: num OF a OVER den OF a; |
|||
OP + = ( INT a, FRAC b )FRAC: ( a // 1 ) + b; |
|||
FRAC one = 1 // 1; |
|||
# returns the first n elements of the Calkin-Wilf sequence # |
|||
PROC calkin wilf = ( INT n )[]FRAC: |
|||
BEGIN |
|||
[ 1 : n ]FRAC q; |
|||
IF n > 0 THEN |
|||
q[ 1 ] := 1 // 1; |
|||
FOR i FROM 2 TO UPB q DO |
|||
q[ i ] := one / ( ( 2 * FLOOR q[ i - 1 ] ) + one - q[ i - 1 ] ) |
|||
OD |
|||
FI; |
|||
q |
|||
END # calkin wilf # ; |
|||
# returns the position of a FRAC in the Calkin-Wilf sequence by computing its # |
|||
# continued fraction representation and converting that to a bit string # |
|||
# - the position must fit in a 2-bit number # |
|||
PROC position in calkin wilf sequence = ( FRAC f )INT: |
|||
IF INT result := 0; |
|||
[ 1 : 32 ]INT cf; FOR i FROM LWB cf TO UPB cf DO cf[ i ] := 0 OD; |
|||
INT num := num OF f; |
|||
INT den := den OF f; |
|||
INT cf length := 0; |
|||
FOR i FROM LWB cf WHILE den /= 0 DO |
|||
cf[ cf length := i ] := r2cf( num, den ) |
|||
OD; |
|||
NOT ODD cf length |
|||
THEN # the continued fraction does not have an odd length # |
|||
-1 |
|||
ELSE # the continued fraction has an odd length so we can compute the seuence length # |
|||
# build the number by alternating d 1s and 0s where d is the digits of the # |
|||
# continued fraction, starting at the least significant # |
|||
INT digit := 1; |
|||
FOR d pos FROM cf length BY -1 TO 1 DO |
|||
FOR i TO cf[ d pos ] DO |
|||
result *:= 2 +:= digit |
|||
OD; |
|||
digit := IF digit = 0 THEN 1 ELSE 0 FI |
|||
OD; |
|||
result |
|||
FI # position in calkin wilf sequence # ; |
|||
BEGIN # task # |
|||
# get and show the first 20 Calkin-Wilf sequence numbers # |
|||
[]FRAC cw = calkin wilf( 20 ); |
|||
print( ( "The first 20 elements of the Calkin-Wilf sequence are:", newline, " " ) ); |
|||
FOR n FROM LWB cw TO UPB cw DO |
|||
FRAC sn = cw[ n ]; |
|||
print( ( " ", whole( num OF sn, 0 ), "/", whole( den OF sn, 0 ) ) ) |
|||
OD; |
|||
print( ( newline ) ); |
|||
# show the position of a specific element in the sequence # |
|||
print( ( "Position of 83116/51639 in the sequence: " |
|||
, whole( position in calkin wilf sequence( 83116//51639 ), 0 ) |
|||
) |
|||
) |
|||
END |
|||
END</syntaxhighlight> |
|||
{{out}} |
|||
<pre> |
|||
The first 20 elements of the Calkin-Wilf sequence are: |
|||
1/1 1/2 2/1 1/3 3/2 2/3 3/1 1/4 4/3 3/5 5/2 2/5 5/3 3/4 4/1 1/5 5/4 4/7 7/3 3/8 |
|||
Position of 83116/51639 in the sequence: 123456789 |
|||
</pre> |
|||
=={{header|AppleScript}}== |
=={{header|AppleScript}}== |
||
<syntaxhighlight lang="applescript">-- Return the first n terms of the sequence. Tree generation. Faster for this purpose. |
|||
<lang applescript>on CalkinWilfSequence(n) |
|||
on CalkinWilfSequence(n) |
|||
script o |
script o |
||
property |
property sequence : {{1, 1}} -- Initialised with the first term ({numerator, denominator}). |
||
end script |
end script |
||
-- Work through the growing sequence list, adding the two children of each term to the end and |
|||
set i to 0 |
|||
-- converting each term to text representing the vulgar fraction. Stop adding children halfway through. |
|||
set counter to 1 |
|||
set halfway to n div 2 |
|||
repeat with position from 1 to n |
|||
set { |
set {numerator, denominator} to item position of o's sequence |
||
if (position ≤ halfway) then |
|||
tell numerator + denominator |
|||
set end of o's |
set end of o's sequence to {numerator, it} |
||
if ((position < halfway) or (position * 2 < n)) then set end of o's sequence to {it, denominator} |
|||
set counter to counter + 1 |
|||
end tell |
|||
end if |
|||
set |
set item position of o's sequence to (numerator as text) & "/" & denominator |
||
end repeat |
|||
repeat with i from (i + 1) to counter |
|||
set {a, b} to item i of o's output |
|||
set item i of o's output to (a as text) & "/" & b |
|||
end repeat |
end repeat |
||
return o's |
return o's sequence |
||
end CalkinWilfSequence |
end CalkinWilfSequence |
||
-- Alternatively, return terms pos1 to pos2. Binary run-length encoding. Doesn't need to work from the beginning of the sequence. |
|||
on CalkinWilfSequencePos(numerator, denominator) |
|||
on CalkinWilfSequence2(pos1, pos2) |
|||
script o |
|||
property sequence : {} |
|||
end script |
|||
repeat with position from pos1 to pos2 |
|||
-- Build a continued fraction list from the binary run-length encoding of this position index. |
|||
-- There's no need to put the last value into the list as it's used immediately. |
|||
set continuedFraction to {} |
|||
set bitValue to 1 |
|||
set runLength to 0 |
|||
repeat until (position = 0) |
|||
if (position mod 2 = bitValue) then |
|||
set runLength to runLength + 1 |
|||
else |
|||
set end of continuedFraction to runLength |
|||
set bitValue to (bitValue + 1) mod 2 |
|||
set runLength to 1 |
|||
end if |
|||
set position to position div 2 |
|||
end repeat |
|||
-- Work out the numerator and denominator from the continued fraction and derive text representing the vulgar fraction. |
|||
set numerator to runLength |
|||
set denominator to 1 |
|||
repeat with i from (count continuedFraction) to 1 by -1 |
|||
tell numerator |
|||
set numerator to numerator * (item i of continuedFraction) + denominator |
|||
set denominator to it |
|||
end tell |
|||
end repeat |
|||
set end of o's sequence to (numerator as text) & "/" & denominator |
|||
end repeat |
|||
return o's sequence |
|||
end CalkinWilfSequence2 |
|||
-- Return the sequence position of the term with the given numerator and denominator. |
|||
on CalkinWilfSequencePosition(numerator, denominator) |
|||
-- Build a continued fraction list from the input. |
|||
set continuedFraction to {} |
set continuedFraction to {} |
||
repeat until (denominator is 0) |
repeat until (denominator is 0) |
||
Line 62: | Line 279: | ||
set {numerator, denominator} to {denominator, numerator mod denominator} |
set {numerator, denominator} to {denominator, numerator mod denominator} |
||
end repeat |
end repeat |
||
-- If it has an even number of entries, convert to the equivalent odd number. |
|||
if ((count continuedFraction) mod 2 is 0) then |
if ((count continuedFraction) mod 2 is 0) then |
||
set item |
set last item of continuedFraction to (last item of continuedFraction) - 1 |
||
set end of continuedFraction to 1 |
set end of continuedFraction to 1 |
||
end if |
end if |
||
-- "Binary run-length decode" the entries to get the position index. |
|||
set position to 0 |
set position to 0 |
||
set |
set bitValue to 1 |
||
repeat with i from (count continuedFraction) to 1 by -1 |
repeat with i from (count continuedFraction) to 1 by -1 |
||
repeat (item i of continuedFraction) times |
repeat (item i of continuedFraction) times |
||
set position to position * 2 + |
set position to position * 2 + bitValue |
||
end repeat |
end repeat |
||
set |
set bitValue to (bitValue + 1) mod 2 |
||
end repeat |
end repeat |
||
return position |
return position |
||
end CalkinWilfSequencePosition |
|||
end CalkinWilfSequencePos |
|||
-- Task code: |
-- Task code: |
||
local |
local sequenceResult1, sequenceResult2, positionResult, output, astid |
||
set |
set sequenceResult1 to CalkinWilfSequence(20) |
||
set |
set sequenceResult2 to CalkinWilfSequence2(1, 20) |
||
set positionResult to CalkinWilfSequencePosition(83116, 51639) |
|||
set astid to AppleScript's text item delimiters |
set astid to AppleScript's text item delimiters |
||
set AppleScript's text item delimiters to ", " |
set AppleScript's text item delimiters to ", " |
||
set output to "First twenty terms |
set output to "First twenty terms of sequence using tree generation:" & (linefeed & sequenceResult1) |
||
set output to output & (linefeed & "Ditto using binary run-length encoding:") & (linefeed & sequenceResult1) |
|||
set AppleScript's text item delimiters to astid |
set AppleScript's text item delimiters to astid |
||
set output to output & linefeed & "83116/51639 is term number " & positionResult |
set output to output & (linefeed & "83116/51639 is term number " & positionResult) |
||
return output</ |
return output</syntaxhighlight> |
||
{{output}} |
{{output}} |
||
< |
<syntaxhighlight lang="applescript">"First twenty terms of sequence using tree generation: |
||
1/1, 1/2, 2/1, 1/3, 3/2, 2/3, 3/1, 1/4, 4/3, 3/5, 5/2, 2/5, 5/3, 3/4, 4/1, 1/5, 5/4, 4/7, 7/3, 3/8 |
1/1, 1/2, 2/1, 1/3, 3/2, 2/3, 3/1, 1/4, 4/3, 3/5, 5/2, 2/5, 5/3, 3/4, 4/1, 1/5, 5/4, 4/7, 7/3, 3/8 |
||
Ditto using binary run-length encoding: |
|||
83116/51639 is term number 123456789"</lang> |
|||
1/1, 1/2, 2/1, 1/3, 3/2, 2/3, 3/1, 1/4, 4/3, 3/5, 5/2, 2/5, 5/3, 3/4, 4/1, 1/5, 5/4, 4/7, 7/3, 3/8 |
|||
83116/51639 is term number 123456789"</syntaxhighlight> |
|||
=={{header|Arturo}}== |
|||
<syntaxhighlight lang="rebol">n: new 1 |
|||
d: new 1 |
|||
calkinWilf: function [] .export:[n,d] [ |
|||
n: (d - n) + 2 * (n/d) * d |
|||
tmp: d |
|||
d: n |
|||
n: tmp |
|||
return @[n d] |
|||
] |
|||
first20: [[1 1]] ++ map 1..19 => calkinWilf |
|||
print "The first 20 terms of the Calkwin-Wilf sequence are:" |
|||
print map first20 'f -> ~"|f\0|/|f\1|" |
|||
n: new 1 |
|||
d: new 1 |
|||
indx: new 1 |
|||
target: [83116, 51639] |
|||
while ø [ |
|||
inc 'indx |
|||
if target = calkinWilf -> break |
|||
] |
|||
print "" |
|||
print ["The element" ~"|target\0|/|target\1|" "is at position" indx "in the sequence."]</syntaxhighlight> |
|||
{{out}} |
|||
<pre>The first 20 terms of the Calkwin-Wilf sequence are: |
|||
1/1 1/2 2/1 1/3 3/2 2/3 3/1 1/4 4/3 3/5 5/2 2/5 5/3 3/4 4/1 1/5 5/4 4/7 7/3 3/8 |
|||
The element 83116/51639 is at position 123456789 in the sequence.</pre> |
|||
=={{header|BQN}}== |
|||
BQN does not have rational number arithmetic yet, so it is manually implemented. |
|||
Part 2 runs in ~150 secs on [[CBQN]]. |
|||
<code>GCD</code> and <code>_while_</code> are idioms from [https://mlochbaum.github.io/bqncrate/ BQNcrate]. |
|||
<syntaxhighlight lang="bqn">GCD ← {m 𝕊⍟(0<m←𝕨|𝕩) 𝕨} |
|||
_while_ ← {𝔽⍟𝔾∘𝔽_𝕣_𝔾∘𝔽⍟𝔾𝕩} |
|||
Sim ← { # Simplify a fraction |
|||
x𝕊1: 𝕨‿1; |
|||
0𝕊y: 0‿𝕩; |
|||
⌊𝕨‿𝕩 ÷ 𝕨 GCD 𝕩 |
|||
} |
|||
Add ← { # Add two fractions |
|||
0‿b 𝕊 𝕩: 𝕩; |
|||
𝕨 𝕊 0‿y: 𝕨; |
|||
a‿b 𝕊 x‿y: |
|||
((a×y)+x×b) Sim b×y |
|||
} |
|||
Next ← {n‿d: ⌽(2×⌊÷´n‿d)‿1 Add (d-n)‿d} # Next term |
|||
Cal ← {Next⍟𝕩 1‿1} |
|||
•Show Cal 1+↕20 |
|||
•Show { |
|||
cnt‿fr: |
|||
⟨cnt+1,Next fr⟩ |
|||
} _while_ { |
|||
cnt‿fr: |
|||
fr ≢ 83116‿51639 |
|||
} ⟨1,1‿1⟩</syntaxhighlight> |
|||
<syntaxhighlight lang="bqn">⟨ ⟨ 1 2 ⟩ ⟨ 2 1 ⟩ ⟨ 1 3 ⟩ ⟨ 3 2 ⟩ ⟨ 2 3 ⟩ ⟨ 3 1 ⟩ ⟨ 1 4 ⟩ ⟨ 4 3 ⟩ ⟨ 3 5 ⟩ ⟨ 5 2 ⟩ ⟨ 2 5 ⟩ ⟨ 5 3 ⟩ ⟨ 3 4 ⟩ ⟨ 4 1 ⟩ ⟨ 1 5 ⟩ ⟨ 5 4 ⟩ ⟨ 4 7 ⟩ ⟨ 7 3 ⟩ ⟨ 3 8 ⟩ ⟨ 8 5 ⟩ ⟩ |
|||
⟨ 123456789 ⟨ 83116 51639 ⟩ ⟩</syntaxhighlight> |
|||
You can try Part 1 [https://mlochbaum.github.io/BQN/try.html#code=R0NEIOKGkCB7bSDwnZWK4o2fKDA8beKGkPCdlah88J2VqSkg8J2VqH0KX3doaWxlXyDihpAge/CdlL3ijZ/wnZS+4oiY8J2UvV/wnZWjX/CdlL7iiJjwnZS94o2f8J2UvvCdlal9ClNpbSDihpAgewogIHjwnZWKMTog8J2VqOKAvzE7CiAgMPCdlYp5OiAw4oC/8J2VqTsKICDijIrwnZWo4oC/8J2VqSDDtyDwnZWoIEdDRCDwnZWpCn0KQWRkIOKGkCB7CiAgMOKAv2Ig8J2ViiDwnZWpOiDwnZWpOwogIPCdlagg8J2ViiAw4oC/eTog8J2VqDsKICBh4oC/YiDwnZWKIHjigL95OgogICgoYcOXeSkreMOXYikgU2ltIGLDl3kKfQpOZXh0IOKGkCB7buKAv2Q6IOKMvSgyw5fijIrDt8K0buKAv2Qp4oC/MSBBZGQgKGQtbinigL9kfQpDYWwg4oaQIHtOZXh04o2f8J2VqSAx4oC/MX0KCuKAolNob3cgQ2FsIDEr4oaVMjA= here.] Second part can and will hang your browser, so it is best to try locally on [[CBQN]]. |
|||
=={{header|Bracmat}}== |
|||
{{trans|Python}} |
|||
<syntaxhighlight lang="bracmat">( 1:?a |
|||
& 0:?i |
|||
& whl |
|||
' ( 1+!i:<20:?i |
|||
& (2*div$(!a,1)+1+-1*!a)^-1:?a |
|||
& out$!a |
|||
) |
|||
& ( r2cf |
|||
= floor |
|||
. div$(!arg,1):?floor |
|||
& ( !floor:!arg |
|||
| !floor r2cf$((!arg+-1*!floor)^-1) |
|||
) |
|||
) |
|||
& ( get-term-num |
|||
= ans dig pwr |
|||
. (0,1,1):(?ans,?dig,?pwr) |
|||
& r2cf$!arg:?n |
|||
& map |
|||
$ ( ( |
|||
= |
|||
. whl |
|||
' ( !arg+-1:~<0:?arg |
|||
& !dig*!pwr+!ans:?ans |
|||
& 2*!pwr:?pwr |
|||
) |
|||
& 1+-1*!dig:?dig |
|||
) |
|||
. !n |
|||
) |
|||
& !ans |
|||
) |
|||
& out$(get-term-num$83116/51639) |
|||
);</syntaxhighlight> |
|||
{{out}} |
|||
<pre>1/2 |
|||
2 |
|||
1/3 |
|||
3/2 |
|||
2/3 |
|||
3 |
|||
1/4 |
|||
4/3 |
|||
3/5 |
|||
5/2 |
|||
2/5 |
|||
5/3 |
|||
3/4 |
|||
4 |
|||
1/5 |
|||
5/4 |
|||
4/7 |
|||
7/3 |
|||
3/8 |
|||
123456789</pre> |
|||
=={{header|C++}}== |
=={{header|C++}}== |
||
{{libheader|Boost}} |
{{libheader|Boost}} |
||
< |
<syntaxhighlight lang="cpp">#include <iostream> |
||
#include <vector> |
#include <vector> |
||
#include <boost/rational.hpp> |
#include <boost/rational.hpp> |
||
Line 149: | Line 499: | ||
rational r(83116, 51639); |
rational r(83116, 51639); |
||
std::cout << r << " is the " << term_number(r) << "th term of the sequence.\n"; |
std::cout << r << " is the " << term_number(r) << "th term of the sequence.\n"; |
||
}</ |
}</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
Line 175: | Line 525: | ||
20: 3/8 |
20: 3/8 |
||
83116/51639 is the 123456789th term of the sequence. |
83116/51639 is the 123456789th term of the sequence. |
||
</pre> |
|||
=={{header|EasyLang}}== |
|||
{{trans|Nim}} |
|||
<syntaxhighlight> |
|||
subr first |
|||
n = 1 ; d = 1 |
|||
. |
|||
proc next . . |
|||
n = 2 * (n div d) * d + d - n |
|||
swap n d |
|||
. |
|||
print "The first 20 terms of the Calkwin-Wilf sequence are:" |
|||
first |
|||
for i to 20 |
|||
write n & "/" & d & " " |
|||
next |
|||
. |
|||
print "" |
|||
# |
|||
first |
|||
i = 1 |
|||
while n <> 83116 or d <> 51639 |
|||
next |
|||
i += 1 |
|||
. |
|||
print "83116/51639 is at position " & i |
|||
</syntaxhighlight> |
|||
{{out}} |
|||
<pre> |
|||
The first 20 terms of the Calkwin-Wilf sequence are: |
|||
1/1 1/2 2/1 1/3 3/2 2/3 3/1 1/4 4/3 3/5 5/2 2/5 5/3 3/4 4/1 1/5 5/4 4/7 7/3 3/8 |
|||
83116/51639 is at position 123456789 |
|||
</pre> |
</pre> |
||
=={{header|EDSAC order code}}== |
|||
===Find first n terms=== |
|||
{{trans|Pascal}} |
|||
<syntaxhighlight lang="edsac"> |
|||
[For Rosetta Code. EDSAC program, Initial Orders 2. |
|||
Prints the first 20 terms of the Calkin-Wilf sequence. |
|||
Uses term{n} to calculate term{n + 1}.] |
|||
[Print subroutine for non-negative 17-bit integers. |
|||
Parameters: 0F = integer to be printed (not preserved) |
|||
1F = character for leading zero (preserved) |
|||
Workspace: 4F, 5F. Even address; 40 locations] |
|||
T 56 K [define load address] |
|||
GKA3FT34@A1FT35@S37@T36@AFT5FT4FH38#@V4DH30@ |
|||
R32FR16FYFE23@O35@A2FT36@T5FV4DYFL8FT4DA5FL1024F |
|||
UFA36@G16@OFTFT35@A36@G17@ZFPFPFP4FT1714FZ219D |
|||
[Main routine] |
|||
T 100 K [define load address] |
|||
G K [set up relative addressing via @ (theta)] |
|||
[Constants] |
|||
[0] P 10 F [maximum index = 20, edit ad lib.] |
|||
[1] P D [constant 1] |
|||
[Teleprinter characters] |
|||
[2] # F [set figures mode] |
|||
[3] C F [colon (in figures mode)] |
|||
[4] X F [slash (in figures mode)] |
|||
[5] ! F [space] |
|||
[6] @ F [carriage return] |
|||
[7] & F [line feed] |
|||
[8] K 4096 F [null] |
|||
[Variables] |
|||
[9] P F [index] |
|||
[10] P F [a (where term = a/b)] |
|||
[11] P F [b] |
|||
[Enter with acc = 0] |
|||
[12] O 2 @ [set teleprinter to figures] |
|||
A 1 @ [acc := 1] |
|||
U 9 @ [index := 1] |
|||
U 10 @ [a := 1] |
|||
T 11 @ [b := 1 (and clear acc)] |
|||
E 34 @ [jump to print first term] |
|||
[Loop back here if not yet printed enough terms] |
|||
[18] A @ [restore index after test] |
|||
A 1 @ [add 1] |
|||
T 9 @ [update index] |
|||
[Calculate next term. New b := a + b - 2(a mod b). |
|||
Code below calculates c := (a mod b) - b, then new b := a - b - 2*c] |
|||
A 10 @ [acc := a] |
|||
[22] S 11 @ [subtract b] |
|||
E 22 @ [if acc >= 0, subtract again] |
|||
T F [result c < 0, store in 0F] |
|||
A 10 @ [acc := a] |
|||
S 11 @ [subtract b] |
|||
S F [subtract c] |
|||
S F [subtract c] |
|||
T F [new b = a - b - 2*c; store in 0F] |
|||
A 11 @ [acc := old b] |
|||
T 10 @ [copy to a] |
|||
A F [acc := new b] |
|||
T 11 @ [copy to b] |
|||
[Print index and a/b. Assume acc = 0 here.] |
|||
[34] A 5 @ [space to replace leading 0's] |
|||
T 1 F [pass to print subroutine] |
|||
A 9 @ [acc := index] |
|||
T F [pass to print subroutine] |
|||
[38] A 38 @ [for return from subroutine] |
|||
G 56 F [call subroutine, clears acc] |
|||
O 3 @ [print colon] |
|||
O 5 @ [print space] |
|||
A 8 @ [null to replace leading 0's] |
|||
T 1 F [pass to print subroutine] |
|||
A10@ TF A46@ G56F O4@ [print a followed by slash] |
|||
A11@ TF A51@ G56F O6@ O7@ [print b followed by CR LF] |
|||
[Test whether enough terms have been printed] |
|||
A 9 @ [acc := index] |
|||
S @ [subtract maximum index] |
|||
G 18 @ [loop back if acc < 0] |
|||
[Exit] |
|||
O 8 @ [print null to flush teleprinter buffer] |
|||
Z F [stop] |
|||
E 12 Z [relative address of entry point] |
|||
P F [enter with acc = 0] |
|||
[end] |
|||
</syntaxhighlight> |
|||
{{out}} |
|||
<pre> |
|||
1: 1/1 |
|||
2: 1/2 |
|||
3: 2/1 |
|||
4: 1/3 |
|||
5: 3/2 |
|||
6: 2/3 |
|||
7: 3/1 |
|||
8: 1/4 |
|||
9: 4/3 |
|||
10: 3/5 |
|||
11: 5/2 |
|||
12: 2/5 |
|||
13: 5/3 |
|||
14: 3/4 |
|||
15: 4/1 |
|||
16: 1/5 |
|||
17: 5/4 |
|||
18: 4/7 |
|||
19: 7/3 |
|||
20: 3/8 |
|||
</pre> |
|||
===Find index of a given term=== |
|||
{{trans|Pascal}} |
|||
<syntaxhighlight lang="edsac"> |
|||
[For Rosetta Code. EDSAC program, Initial Orders 2.] |
|||
[Finds the index of a given rational in the Calkin-Wilf series.] |
|||
[Library subroutine R2: input of positive integers. |
|||
Runs during input of the program, and is then overwritten. |
|||
Allows integers to be written in decimal, rather than as "pseudo-orders". |
|||
See Wilkes, Williams & Gill, 1951 edn, pp. 96-97, 148.] |
|||
T 54 K [to access integers via C parameter] |
|||
P 110 F [where to load integers] |
|||
GKT20FVDL8FA40DUDTFI40FA40FS39FG@S2FG23FA5@T5@E4@E13Z |
|||
T #C [tell R2 where to load integers] |
|||
[F after each integer except the last, and # after the last.] |
|||
83116F51639# |
|||
[Modified library subroutine P7. |
|||
Prints signed integer up to 10 digits, left-justified. |
|||
Input: Number to be printed is at 0D. |
|||
54 locations. Load at even address. Workspace 4D.] |
|||
T 56 K |
|||
GKA3FT42@A49@T31@ADE10@T31@A48@T31@SDTDH44#@NDYFLDT4DS43@ |
|||
TFH17@S17@A43@G23@UFS43@T1FV4DAFG50@SFLDUFXFOFFFSFL4FT4DA49@ |
|||
T31@A1FA43@G20@XFP1024FP610D@524D!FO46@O26@XFSFL8FT4DE39@ |
|||
[Main routine.] |
|||
T 120 K [define load address (must be even)] |
|||
G K [set up relative addressing via @ (theta)] |
|||
[Put 35-bit values first, to ensure each is at an even address] |
|||
[Variables] |
|||
[0] P F P F [a] |
|||
[2] P F P F [b] |
|||
[4] P F P F [power of 2] |
|||
[6] P F P F [calculated index] |
|||
[Constants] |
|||
T8#Z PF T8Z [clears sandwich digit between 8 and 9] |
|||
[8] P D P F [35-bit constant 1] |
|||
[Teleprinter characters] |
|||
[10] # F [set figures mode] |
|||
[11] X F [slash (in figures mode)] |
|||
[12] K 2048 F [set letters mode] |
|||
[13] I F [letter I] |
|||
[14] R F [letter R] |
|||
[15] ! F [space] |
|||
[16] @ F [carriage return] |
|||
[17] & F [line feed] |
|||
[18] K 4096 F [null char] |
|||
[Enter with acc = 0] |
|||
[19] A #C [acc := initial a] |
|||
T #@ [copy to variable] |
|||
A 2#C [acc := initial b] |
|||
T 2#@ [copy to variable] |
|||
[23] A 8#@ [acc := 1] |
|||
[24] T 4#@ [initialize power of 2] |
|||
T 6#@ [initialize index to 0] |
|||
[Loop] |
|||
[26] A #@ [acc := a] |
|||
[27] S 2#@ [subtract b] |
|||
[28] E 33 @ [jump if a >= b] |
|||
[Here if a < b] |
|||
T D [store a - b in 0D] |
|||
S D [negate] |
|||
T 2#@ [b := b - a] |
|||
E 40 @ [join common code] |
|||
[Here if a >= b] |
|||
[33] S 8#@ [acc = a - b; test for a = b] |
|||
G 45 @ [jump out of loop if so] |
|||
A 8#@ [restore a - b] |
|||
T #@ [a := a - b] |
|||
A 6#@ [acc := index] |
|||
A 4#@ [inc index by power of 2] |
|||
T 6#@ |
|||
[Code common to both cases] |
|||
[40] A 4#@ [acc := power of 2] |
|||
L D [shift left] |
|||
G 76 @ |
|||
T 4#@ [update power of 2] |
|||
E 26 @ [loop back] |
|||
[Exit from loop.] |
|||
[45] T D [dump acc to clear it] |
|||
A 6#@ [acc := index] |
|||
A 4#@ [add power of 2 ] |
|||
T 6#@ [store final value of index] |
|||
[Finished calcualting index, now do printing] |
|||
O 10 @ [set teleprinter to figures] |
|||
A #C [acc := initial a] |
|||
T D [to 0D for printing] |
|||
[52] A 52 @ [for return from subroutine] |
|||
G 56 F [call print subroutine, clears acc] |
|||
O 11 @ [print slash] |
|||
A 2#C [print initial b similarly] |
|||
T D |
|||
[57] A 57 @ |
|||
G 56 F |
|||
O 12 @ [set teleprinter to letters and print ' IS AT '] |
|||
O15@ O13@ O27@ O15@ O23@ O24@ O15@ |
|||
O 10 @ [set teleprinter to figures] |
|||
A 6#@ [acc := calculated index] |
|||
T D [send to print subroutine] |
|||
[70] A 70 @ |
|||
G 56 F |
|||
[72] O16@ O17@ [print CR, LF] |
|||
O 18 @ [print null to flush teleprinter buffer] |
|||
Z F [stop] |
|||
[Here if power of 2 goes negative (accumulator overflow)] |
|||
[76] O 12 @ [set teleprinter to letters] |
|||
O28@ O14@ O14@ O76@ O14@ [print'ERROR'] |
|||
G 72 @ [jump to common exit] |
|||
E 19 Z [relative address of entry point] |
|||
P F [enter with acc = 0] |
|||
</syntaxhighlight> |
|||
{{out}} |
|||
<pre> |
|||
83116/51639 IS AT 123456789 |
|||
</pre> |
|||
=={{header|F_Sharp|F#}}== |
|||
===The Function=== |
|||
<syntaxhighlight lang="fsharp"> |
|||
// Calkin Wilf Sequence. Nigel Galloway: January 9th., 2021 |
|||
let cW=Seq.unfold(fun(n)->Some(n,seq{for n,g in n do yield (n,n+g); yield (n+g,g)}))(seq[(1,1)])|>Seq.concat |
|||
</syntaxhighlight> |
|||
===The Tasks=== |
|||
; first 20 |
|||
<syntaxhighlight lang="fsharp"> |
|||
cW |> Seq.take 20 |> Seq.iter(fun(n,g)->printf "%d/%d " n g);printfn "" |
|||
</syntaxhighlight> |
|||
{{out}} |
|||
<pre> |
|||
1/1 1/2 2/1 1/3 3/2 2/3 3/1 1/4 4/3 3/5 5/2 2/5 5/3 3/4 4/1 1/5 5/4 4/7 7/3 3/8 |
|||
</pre> |
|||
; Indexof 83116/51639 |
|||
<syntaxhighlight lang="fsharp"> |
|||
printfn "%d" (1+Seq.findIndex(fun n->n=(83116,51639)) cW) |
|||
</syntaxhighlight> |
|||
{{out}} |
|||
<pre> |
|||
123456789 |
|||
</pre> |
|||
=={{header|Factor}}== |
=={{header|Factor}}== |
||
{{works with|Factor|0.99 2020-08-14}} |
{{works with|Factor|0.99 2020-08-14}} |
||
< |
<syntaxhighlight lang="factor">USING: formatting io kernel lists lists.lazy math |
||
math.continued-fractions math.functions math.parser prettyprint |
math.continued-fractions math.functions math.parser prettyprint |
||
sequences strings vectors ; |
sequences strings vectors ; |
||
Line 199: | Line 832: | ||
20 calkin-wilf ltake [ pprint bl ] leach nl nl |
20 calkin-wilf ltake [ pprint bl ] leach nl nl |
||
83116/51639 cw-index "83116/51639 is at index %d.\n" printf</ |
83116/51639 cw-index "83116/51639 is at index %d.\n" printf</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
<pre> |
<pre> |
||
Line 207: | Line 840: | ||
83116/51639 is at index 123456789. |
83116/51639 is at index 123456789. |
||
</pre> |
</pre> |
||
=={{header|Forth}}== |
|||
{{works with|gforth|0.7.3}} |
|||
<syntaxhighlight lang="forth">\ Calkin-Wilf sequence |
|||
: frac. swap . ." / " . ; |
|||
: cw-next ( num den -- num den ) 2dup / over * 2* over + rot - ; |
|||
: cw-seq ( n -- ) |
|||
1 1 rot |
|||
0 do |
|||
cr 2dup frac. cw-next |
|||
loop 2drop ; |
|||
variable index |
|||
variable bit-state |
|||
variable bit-position |
|||
: r2cf-next ( num1 den1 -- num2 den2 u ) swap over >r s>d r> sm/rem ; |
|||
: n2bitlength ( n -- ) |
|||
bit-state @ if |
|||
1 swap lshift 1- bit-position @ lshift index +! |
|||
else drop then ; |
|||
: index-init true bit-state ! 0 bit-position ! 0 index ! ; |
|||
: index-build ( n -- ) |
|||
dup n2bitlength bit-position +! bit-state @ invert bit-state ! ; |
|||
: index-finish ( n 0 -- ) 2drop -1 bit-position +! 1 index-build ; |
|||
: cw-index ( num den -- ) |
|||
index-init |
|||
begin r2cf-next index-build dup 0<> while repeat |
|||
index-finish ; |
|||
: cw-demo |
|||
20 cw-seq |
|||
cr 83116 51639 2dup frac. cw-index index @ . ; |
|||
cw-demo</syntaxhighlight> |
|||
{{out}} |
|||
<pre>1 / 1 |
|||
1 / 2 |
|||
2 / 1 |
|||
1 / 3 |
|||
3 / 2 |
|||
2 / 3 |
|||
3 / 1 |
|||
1 / 4 |
|||
4 / 3 |
|||
3 / 5 |
|||
5 / 2 |
|||
2 / 5 |
|||
5 / 3 |
|||
3 / 4 |
|||
4 / 1 |
|||
1 / 5 |
|||
5 / 4 |
|||
4 / 7 |
|||
7 / 3 |
|||
3 / 8 |
|||
83116 / 51639 123456789 ok</pre> |
|||
=={{header|FreeBASIC}}== |
=={{header|FreeBASIC}}== |
||
Uses the code from [[Greatest common divisor#FreeBASIC]] as an include. |
Uses the code from [[Greatest common divisor#FreeBASIC]] as an include. |
||
< |
<syntaxhighlight lang="freebasic">#include "gcd.bas" |
||
type rational |
type rational |
||
Line 319: | Line 1,012: | ||
q.num = 83116 |
q.num = 83116 |
||
q.den = 51639 |
q.den = 51639 |
||
print disp_rational(q)+" is the "+str(frac_to_int(q))+"th term."</ |
print disp_rational(q)+" is the "+str(frac_to_int(q))+"th term."</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
Line 344: | Line 1,037: | ||
20 3/8 |
20 3/8 |
||
83116/51639 is the 123456789th term.</pre> |
83116/51639 is the 123456789th term.</pre> |
||
=={{header|Fōrmulæ}}== |
|||
{{FormulaeEntry|page=https://formulae.org/?script=examples/Calkin-Wilf_correspondence}} |
|||
=={{header|Go}}== |
=={{header|Go}}== |
||
{{trans|Wren}} |
{{trans|Wren}} |
||
Go just has arbitrary precision rational numbers which we use here whilst assuming the numbers needed for this task can be represented exactly by the 64 bit built-in types. |
Go just has arbitrary precision rational numbers which we use here whilst assuming the numbers needed for this task can be represented exactly by the 64 bit built-in types. |
||
< |
<syntaxhighlight lang="go">package main |
||
import ( |
import ( |
||
Line 438: | Line 1,134: | ||
tn := getTermNumber(cf) |
tn := getTermNumber(cf) |
||
fmt.Printf("%s is the %sth term of the sequence.\n", r.RatString(), commatize(tn)) |
fmt.Printf("%s is the %sth term of the sequence.\n", r.RatString(), commatize(tn)) |
||
}</ |
}</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
Line 466: | Line 1,162: | ||
83116/51639 is the 123,456,789th term of the sequence. |
83116/51639 is the 123,456,789th term of the sequence. |
||
</pre> |
</pre> |
||
=={{header|Haskell}}== |
=={{header|Haskell}}== |
||
<syntaxhighlight lang="haskell">import Control.Monad (forM_) |
|||
{{needs-review|Haskell|Now starts a1 = 1}} |
|||
<lang haskell>import Control.Monad (forM_) |
|||
import Data.Bool (bool) |
import Data.Bool (bool) |
||
import Data.List.NonEmpty (NonEmpty, fromList, toList, unfoldr) |
import Data.List.NonEmpty (NonEmpty, fromList, toList, unfoldr) |
||
import Text.Printf (printf) |
import Text.Printf (printf) |
||
-- The infinite Calkin-Wilf sequence. |
-- The infinite Calkin-Wilf sequence, a(n), starting with a(1) = 1. |
||
calkinWilfs :: [Rational] |
calkinWilfs :: [Rational] |
||
calkinWilfs = iterate ( |
calkinWilfs = iterate (recip . succ . ((-) =<< (2 *) . fromIntegral . floor)) 1 |
||
-- The index into the Calkin-Wilf sequence of a given rational number, starting |
-- The index into the Calkin-Wilf sequence of a given rational number, starting |
||
-- with |
-- with 1 at index 1. |
||
calkinWilfIdx :: Rational -> Integer |
calkinWilfIdx :: Rational -> Integer |
||
calkinWilfIdx = rld . cfo |
calkinWilfIdx = rld . cfo |
||
Line 492: | Line 1,186: | ||
cf :: Rational -> NonEmpty Int |
cf :: Rational -> NonEmpty Int |
||
cf = unfoldr step |
cf = unfoldr step |
||
where |
|||
where step r = case properFraction r of |
|||
step r = |
|||
(n, 1) -> (1+n, Nothing) |
|||
case properFraction r of |
|||
(n, 0) -> ( n, Nothing) |
|||
(n, 1) -> (succ n, Nothing) |
|||
(n, 0) -> (n, Nothing) |
|||
(n, f) -> (n, Just (recip f)) |
|||
-- Ensure a continued fraction has an odd length. |
-- Ensure a continued fraction has an odd length. |
||
oddLen :: NonEmpty Int -> NonEmpty Int |
oddLen :: NonEmpty Int -> NonEmpty Int |
||
oddLen = fromList . go . toList |
oddLen = fromList . go . toList |
||
where |
where |
||
go [x, y] = [x, pred y, 1] |
|||
go (x:y:zs) = x : y : go zs |
|||
go xs = xs |
|||
-- Run-length decode a continued fraction. |
-- Run-length decode a continued fraction. |
||
rld :: NonEmpty Int -> Integer |
rld :: NonEmpty Int -> Integer |
||
rld = snd . foldr step (True, 0) |
rld = snd . foldr step (True, 0) |
||
where |
|||
where step i (b, n) = let p = 2^i in (not b, n*p + bool 0 (p-1) b) |
|||
step i (b, n) = |
|||
let p = 2 ^ i |
|||
in (not b, n * p + bool 0 (pred p) b) |
|||
main :: IO () |
main :: IO () |
||
main = do |
main = do |
||
forM_ (take |
forM_ (take 20 $ zip [1 :: Int ..] calkinWilfs) $ |
||
printf "%2d %s\n" i (show r) |
\(i, r) -> printf "%2d %s\n" i (show r) |
||
let r = 83116 / 51639 |
let r = 83116 / 51639 |
||
printf |
|||
printf "\n%s is at index %d of the Calkin-Wilf sequence.\n" |
|||
"\n%s is at index %d of the Calkin-Wilf sequence.\n" |
|||
(show r) (calkinWilfIdx r)</lang> |
|||
(show r) |
|||
(calkinWilfIdx r)</syntaxhighlight> |
|||
{{out}} |
{{out}} |
||
<pre> |
<pre> |
||
0 0 % 1 |
|||
1 1 % 1 |
1 1 % 1 |
||
2 1 % 2 |
2 1 % 2 |
||
Line 542: | Line 1,242: | ||
83116 % 51639 is at index 123456789 of the Calkin-Wilf sequence. |
83116 % 51639 is at index 123456789 of the Calkin-Wilf sequence. |
||
</pre> |
|||
=={{header|J}}== |
|||
<pre> |
|||
cw_next_term^:(<20)1x |
|||
1 1r2 2 1r3 3r2 2r3 3 1r4 4r3 3r5 5r2 2r5 5r3 3r4 4 1r5 5r4 4r7 7r3 3r8 |
|||
(,. index_cw_term&>) 3r4 53r37 83116r51639 |
|||
3r4 14 |
|||
53r37 1081 |
|||
83116r51639 123456789 |
|||
</pre> |
|||
given definitions |
|||
<syntaxhighlight lang="j"> |
|||
cw_next_term=: [: % +:@<. + -. |
|||
ccf =: compute_continued_fraction=: 3 :0 |
|||
if. 0 -: y do. |
|||
, 0 |
|||
else. |
|||
result=. i. 0 |
|||
remainder=. % y |
|||
whilst. remainder do. |
|||
remainder=. % remainder |
|||
integer_part=. <. remainder |
|||
remainder=. remainder - integer_part |
|||
result=. result , integer_part |
|||
end. |
|||
end. |
|||
) |
|||
molcf =: make_odd_length_continued_fraction=: (}: , 1 ,~ <:@{:)^:(0 -: 2 | #) |
|||
NB. base 2 @ reverse @ the cf's representation copies of 1 0 1 0 ... |
|||
index_cw_term=: #.@|.@(# 1 0 $~ #)@molcf@ccf |
|||
</syntaxhighlight> |
|||
Note that <code>ccf</code> could be expressed more concisely: |
|||
<syntaxhighlight lang=J>ccf=: _1 {"1 |.@(0 1 #: %@{.)^:(0~:{.)^:a:</syntaxhighlight> |
|||
=={{header|Java}}== |
|||
<syntaxhighlight lang="java"> |
|||
import java.util.ArrayDeque; |
|||
import java.util.Deque; |
|||
public final class CalkinWilfSequence { |
|||
public static void main(String[] aArgs) { |
|||
Rational term = Rational.ONE; |
|||
System.out.println("First 20 terms of the Calkin-Wilf sequence are:"); |
|||
for ( int i = 1; i <= 20; i++ ) { |
|||
System.out.println(String.format("%2d", i) + ": " + term); |
|||
term = nextCalkinWilf(term); |
|||
} |
|||
System.out.println(); |
|||
Rational rational = new Rational(83_116, 51_639); |
|||
System.out.println(" " + rational + " is the " + termIndex(rational) + "th term of the sequence."); |
|||
} |
|||
private static Rational nextCalkinWilf(Rational aTerm) { |
|||
Rational divisor = Rational.TWO.multiply(aTerm.floor()).add(Rational.ONE).subtract(aTerm); |
|||
return Rational.ONE.divide(divisor); |
|||
} |
|||
private static long termIndex(Rational aRational) { |
|||
long result = 0; |
|||
long binaryDigit = 1; |
|||
long power = 0; |
|||
for ( long term : continuedFraction(aRational) ) { |
|||
for ( long i = 0; i < term; power++, i++ ) { |
|||
result |= ( binaryDigit << power ); |
|||
} |
|||
binaryDigit = ( binaryDigit == 0 ) ? 1 : 0; |
|||
} |
|||
return result; |
|||
} |
|||
private static Deque<Long> continuedFraction(Rational aRational) { |
|||
long numerator = aRational.numerator(); |
|||
long denominator = aRational.denominator(); |
|||
Deque<Long> result = new ArrayDeque<Long>(); |
|||
while ( numerator != 1 ) { |
|||
result.addLast(numerator / denominator); |
|||
long copyNumerator = numerator; |
|||
numerator = denominator; |
|||
denominator = copyNumerator % denominator; |
|||
} |
|||
if ( ! result.isEmpty() && result.size() % 2 == 0 ) { |
|||
final long back = result.removeLast(); |
|||
result.addLast(back - 1); |
|||
result.addLast(1L); |
|||
} |
|||
return result; |
|||
} |
|||
} |
|||
final class Rational { |
|||
public Rational(long aNumerator, long aDenominator) { |
|||
if ( aDenominator == 0 ) { |
|||
throw new ArithmeticException("Denominator cannot be zero"); |
|||
} |
|||
if ( aNumerator == 0 ) { |
|||
aDenominator = 1; |
|||
} |
|||
if ( aDenominator < 0 ) { |
|||
numer = -aNumerator; |
|||
denom = -aDenominator; |
|||
} else { |
|||
numer = aNumerator; |
|||
denom = aDenominator; |
|||
} |
|||
final long gcd = gcd(numer, denom); |
|||
numer = numer / gcd; |
|||
denom = denom / gcd; |
|||
} |
|||
public Rational add(Rational aOther) { |
|||
return new Rational(numer * aOther.denom + aOther.numer * denom, denom * aOther.denom); |
|||
} |
|||
public Rational subtract(Rational aOther) { |
|||
return new Rational(numer * aOther.denom - aOther.numer * denom, denom * aOther.denom); |
|||
} |
|||
public Rational multiply(Rational aOther) { |
|||
return new Rational(numer * aOther.numer, denom * aOther.denom); |
|||
} |
|||
public Rational divide(Rational aOther) { |
|||
return new Rational(numer * aOther.denom, denom * aOther.numer); |
|||
} |
|||
public Rational floor() { |
|||
return new Rational(numer / denom, 1); |
|||
} |
|||
public long numerator() { |
|||
return numer; |
|||
} |
|||
public long denominator() { |
|||
return denom; |
|||
} |
|||
@Override |
|||
public String toString() { |
|||
return numer + "/" + denom; |
|||
} |
|||
public static final Rational ONE = new Rational(1, 1); |
|||
public static final Rational TWO = new Rational(2, 1); |
|||
private long gcd(long aOne, long aTwo) { |
|||
if ( aTwo == 0 ) { |
|||
return aOne; |
|||
} |
|||
return gcd(aTwo, aOne % aTwo); |
|||
} |
|||
private long numer; |
|||
private long denom; |
|||
} |
|||
</syntaxhighlight> |
|||
{{ out }} |
|||
<pre> |
|||
First 20 terms of the Calkin-Wilf sequence are: |
|||
1: 1/1 |
|||
2: 1/2 |
|||
3: 2/1 |
|||
4: 1/3 |
|||
5: 3/2 |
|||
6: 2/3 |
|||
7: 3/1 |
|||
8: 1/4 |
|||
9: 4/3 |
|||
10: 3/5 |
|||
11: 5/2 |
|||
12: 2/5 |
|||
13: 5/3 |
|||
14: 3/4 |
|||
15: 4/1 |
|||
16: 1/5 |
|||
17: 5/4 |
|||
18: 4/7 |
|||
19: 7/3 |
|||
20: 3/8 |
|||
83116/51639 is the 123456789th term of the sequence. |
|||
</pre> |
|||
=={{header|jq}}== |
|||
'''Adapted from [[#Wren|Wren]]''' |
|||
{{works with|jq}} |
|||
'''Also works with gojq, the Go implementation of jq, and with fq''' |
|||
See [[Arithmetic/Rational#jq]] for the Rational module included by the `include` directive. |
|||
In this module, rationals are represented by JSON objects of the form {n, d}, where .n and .d are |
|||
the numerator and denominator respectively. r(n;d) is the constructor function, |
|||
and r(n;d) is pretty-printed as `n // d`. |
|||
<syntaxhighlight lang=jq> |
|||
include "rational"; # see [[Arithmetic/Rational#jq]] |
|||
### Generic Utilities |
|||
# counting from 0 |
|||
def enumerate(s): foreach s as $x (-1; .+1; [., $x]); |
|||
# input string is converted from "base" to an integer, within limits |
|||
# of the underlying arithmetic operations, and without error-checking: |
|||
def to_i(base): |
|||
explode |
|||
| reverse |
|||
| map(if . > 96 then . - 87 else . - 48 end) # "a" ~ 97 => 10 ~ 87 |
|||
| reduce .[] as $c |
|||
# state: [power, ans] |
|||
([1,0]; (.[0] * base) as $b | [$b, .[1] + (.[0] * $c)]) |
|||
| .[1]; |
|||
### The Calkin-Wilf Sequence |
|||
# Emit an array of $n terms |
|||
def calkinWilf($n): |
|||
reduce range(1;$n) as $i ( [r(1;1)]; |
|||
radd(1; rminus( rmult(2; (.[$i-1]|rfloor)); .[$i-1])) as $t |
|||
| .[$i] = rdiv(r(1;1) ; $t)) ; |
|||
# input: a Rational |
|||
def toContinued: |
|||
{ a: .n, |
|||
b: .d, |
|||
res: [] } |
|||
| until( .break; |
|||
.res += [.a / .b | floor] |
|||
| (.a % .b) as $t |
|||
| .a = .b |
|||
| .b = $t |
|||
| .break = (.a == 1) ) |
|||
| if .res|length % 2 == 0 |
|||
then # ensure always odd |
|||
.res[-1] += -1 |
|||
| .res += [1] |
|||
else . |
|||
end |
|||
| .res; |
|||
# input: an array representing a continued fraction |
|||
def getTermNumber: |
|||
reduce .[] as $n ( {b: "", d: "1"}; |
|||
.b = (.d * $n) + .b |
|||
| .d = (if .d == "1" then "0" else "1" end)) |
|||
| .b | to_i(2) ; |
|||
# input: a Rational in the Calkin-Wilf sequence |
|||
def getTermNumber: |
|||
reduce .[] as $n ( {b: "", d: "1"}; |
|||
.b = (.d * $n) + .b |
|||
| .d = (if .d == "1" then "0" else "1" end)) |
|||
| .b | to_i(2) ; |
|||
def task(r): |
|||
"The first 20 terms of the Calkin-Wilf sequence are:", |
|||
(enumerate(calkinWilf(20)[]) | "\(1+.[0]): \(.[1]|rpp)" ), |
|||
"", |
|||
"\(r|rpp) is term # \(r|toContinued|getTermNumber) of the sequence."; |
|||
task( r(83116; 51639) ) |
|||
</syntaxhighlight> |
|||
'''Invocation''': jq -nrf calkin-wilf-sequence.jq |
|||
{{output}} |
|||
<pre> |
|||
The first 20 terms of the Calkin-Wilf sequence are: |
|||
1: 1 // 1 |
|||
2: 1 // 2 |
|||
3: 2 // 1 |
|||
4: 1 // 3 |
|||
5: 3 // 2 |
|||
6: 2 // 3 |
|||
7: 3 // 1 |
|||
8: 1 // 4 |
|||
9: 4 // 3 |
|||
10: 3 // 5 |
|||
11: 5 // 2 |
|||
12: 2 // 5 |
|||
13: 5 // 3 |
|||
14: 3 // 4 |
|||
15: 4 // 1 |
|||
16: 1 // 5 |
|||
17: 5 // 4 |
|||
18: 4 // 7 |
|||
19: 7 // 3 |
|||
20: 3 // 8 |
|||
83116 // 51639 is term # 123456789 of the sequence. |
|||
</pre> |
</pre> |
||
=={{header|Julia}}== |
=={{header|Julia}}== |
||
{{needs-review|Julia|Now starts a1 = 1}} |
|||
{{trans|Wren}} |
{{trans|Wren}} |
||
< |
<syntaxhighlight lang="julia">function calkin_wilf(n) |
||
cw = zeros(Rational, n + 1) |
cw = zeros(Rational, n + 1) |
||
for i in 2:n + 1 |
for i in 2:n + 1 |
||
Line 553: | Line 1,554: | ||
cw[i] = 1 // t |
cw[i] = 1 // t |
||
end |
end |
||
return cw |
return cw[2:end] |
||
end |
end |
||
Line 577: | Line 1,578: | ||
const cw = calkin_wilf(20) |
const cw = calkin_wilf(20) |
||
println("The first |
println("The first 20 terms of the Calkin-Wilf sequence are: $cw") |
||
const r = 83116 // 51639 |
const r = 83116 // 51639 |
||
Line 583: | Line 1,584: | ||
const tn = term_number(cf) |
const tn = term_number(cf) |
||
println("$r is the $tn-th term of the sequence.") |
println("$r is the $tn-th term of the sequence.") |
||
</ |
</syntaxhighlight>{{out}} |
||
<pre> |
<pre> |
||
The first |
The first 20 terms of the Calkin-Wilf sequence are: Rational[1//1, 1//2, 2//1, 1//3, 3//2, 2//3, 3//1, 1//4, 4//3, 3//5, 5//2, 2//5, 5//3, 3//4, 4//1, 1//5, 5//4, 4//7, 7//3, 3//8] |
||
83116//51639 is the 123456789-th term of the sequence. |
83116//51639 is the 123456789-th term of the sequence. |
||
</pre> |
</pre> |
||
=={{header| |
=={{header|Little Man Computer}}== |
||
Runs in a home-made simulator, which is mostly compatible with Peter Higginson's online simulator. Only, for better control of the output format, I've added an instruction OTX (extended output). To run the code in PH's simulator, replace OTX and its parameter with OUT and no parameter. |
|||
<lang Phix>function calkin_wilf(integer len) |
|||
===Find first n terms=== |
|||
sequence cw = repeat(0,len) |
|||
{{trans|Pascal}} |
|||
integer n=0, d=1 |
|||
<syntaxhighlight lang="little man computer"> |
|||
for i=1 to len do |
|||
// Little Man Computer, for Rosetta Code. |
|||
{n,d} = {d,(floor(n/d)*2+1)*d-n} |
|||
// Displays terms of Calkin-Wilf sequence up to the given index. |
|||
cw[i] = {n,d} |
|||
// The chosen algorithm calculates the i-th term directly from i |
|||
end for |
|||
// (i.e. not using any previous terms). |
|||
return cw |
|||
input INP // get number of terms from user |
|||
end function |
|||
BRZ exit // exit if 0 |
|||
STA max_i // store maximum index |
|||
LDA c1 // index := 1 |
|||
next_i STA i |
|||
// Write index followed by '->' |
|||
OTX 3 // non-standard: minimum width 3, no new line |
|||
LDA asc_hy |
|||
OTC |
|||
LDA asc_gt |
|||
OTC |
|||
// Find greatest power of 2 not exceeding i, |
|||
// and count the number of binary digits in i. |
|||
LDA c1 |
|||
STA pwr2 |
|||
loop2 STA nrDigits |
|||
LDA i |
|||
SUB pwr2 |
|||
SUB pwr2 |
|||
BRP double |
|||
BRA part2 // jump out if next power of 2 would exceed i |
|||
double LDA pwr2 |
|||
ADD pwr2 |
|||
STA pwr2 |
|||
LDA nrDigits |
|||
ADD c1 |
|||
BRA loop2 |
|||
// The nth term a/b is calculated from the binary digits of i. |
|||
// The leading 1 is not used. |
|||
part2 LDA c1 |
|||
STA a // a := 1 |
|||
STA b // b := 1 |
|||
LDA i |
|||
SUB pwr2 |
|||
STA diff |
|||
// Pre-decrement count, since leading 1 is not used |
|||
dec_ct LDA nrDigits // count down the number of digits |
|||
SUB c1 |
|||
BRZ output // if all digits done, output the result |
|||
STA nrDigits |
|||
// We now want to compare diff with pwr2/2. |
|||
// Since division is awkward in LMC, we compare 2*diff with pwr2. |
|||
LDA diff // diff := 2*diff |
|||
ADD diff |
|||
STA diff |
|||
SUB pwr2 // is diff >= pwr2 ? |
|||
BRP digit_1 // binary digit is 1 if yes, 0 if no |
|||
// If binary digit is 0 then set b := a + b |
|||
LDA a |
|||
ADD b |
|||
STA b |
|||
BRA dec_ct |
|||
// If binary digit is 1 then update diff and set a := a + b |
|||
digit_1 STA diff |
|||
LDA a |
|||
ADD b |
|||
STA a |
|||
BRA dec_ct |
|||
// Now have nth term a/b. Write it to the output. |
|||
output LDA a // write a |
|||
OTX 1 // non-standard: minimum width 1; no new line |
|||
LDA asc_sl // write slash |
|||
OTC |
|||
LDA b // write b |
|||
OTX 11 // non-standard: minimum width 1; add new line |
|||
LDA i // have we done maximum i yet? |
|||
SUB max_i |
|||
BRZ exit // if yes, exit |
|||
LDA i // if no, increment i and loop back |
|||
ADD c1 |
|||
BRA next_i |
|||
exit HLT |
|||
// Constants |
|||
c1 DAT 1 |
|||
asc_hy DAT 45 |
|||
asc_gt DAT 62 |
|||
asc_sl DAT 47 |
|||
// Variables |
|||
i DAT |
|||
max_i DAT |
|||
pwr2 DAT |
|||
nrDigits DAT |
|||
diff DAT |
|||
a DAT |
|||
b DAT |
|||
// end |
|||
</syntaxhighlight> |
|||
{{out}} |
|||
<pre> |
|||
1->1/1 |
|||
2->1/2 |
|||
3->2/1 |
|||
4->1/3 |
|||
5->3/2 |
|||
6->2/3 |
|||
7->3/1 |
|||
8->1/4 |
|||
9->4/3 |
|||
10->3/5 |
|||
11->5/2 |
|||
12->2/5 |
|||
13->5/3 |
|||
14->3/4 |
|||
15->4/1 |
|||
16->1/5 |
|||
17->5/4 |
|||
18->4/7 |
|||
19->7/3 |
|||
20->3/8 |
|||
</pre> |
|||
===Find index of a given term=== |
|||
{{trans|Pascal}} |
|||
The numbers in part 2 of the task are too large for LMC, so the demo program just confirms the example, that 9/4 is the 35th term. |
|||
<syntaxhighlight lang="little man computer"> |
|||
// Little Man Computer, for Rosetta Code. |
|||
// Calkin-Wilf sequence: displays index of term entered by user. |
|||
INP // get numerator from user |
|||
BRZ exit // exit if 0 |
|||
STA num |
|||
STA a // initialize a := numerator |
|||
INP // get denominator from user |
|||
BRZ exit // exit if 0 |
|||
STA den |
|||
STA b // initialize b := denominator |
|||
LDA c0 // initialize index := 0 |
|||
STA index |
|||
LDA c1 // initialize power of 2 := 1 |
|||
STA pwr2 |
|||
// Build binary digits of the index |
|||
loop LDA a // is a = b yet? |
|||
SUB b |
|||
BRZ break // if yes, break out of loop |
|||
BRP a_gt_b // jump if a > b |
|||
// If a < b then b := b - a, binary digit is 0 |
|||
LDA b |
|||
SUB a |
|||
STA b |
|||
BRA double |
|||
// If a > b then a := a - b, binary digit is 1 |
|||
a_gt_b STA a |
|||
LDA index |
|||
ADD pwr2 |
|||
STA index |
|||
// In either case, on to next power of 2 |
|||
double LDA pwr2 |
|||
ADD pwr2 |
|||
STA pwr2 |
|||
BRA loop |
|||
// Out of loop, add leading binary digit 1 |
|||
break LDA index |
|||
ADD pwr2 |
|||
STA index |
|||
// Output the result |
|||
LDA num |
|||
OTX 1 // non-standard: minimum width = 1, no new line |
|||
LDA asc_sl |
|||
OTC |
|||
LDA den |
|||
OTX 1 |
|||
LDA asc_lt // write '<-' after fraction |
|||
OTC |
|||
LDA asc_hy |
|||
OTC |
|||
LDA index |
|||
OTX 11 // non-standard: minimum width = 1, add new line |
|||
exit HLT |
|||
// Constants |
|||
c0 DAT 0 |
|||
c1 DAT 1 |
|||
asc_sl DAT 47 |
|||
asc_lt DAT 60 |
|||
asc_hy DAT 45 |
|||
// Variables |
|||
num DAT |
|||
den DAT |
|||
a DAT |
|||
b DAT |
|||
pwr2 DAT |
|||
index DAT |
|||
// end |
|||
</syntaxhighlight> |
|||
{{out}} |
|||
<pre> |
|||
9/4<-35 |
|||
</pre> |
|||
=={{header|Mathematica}} / {{header|Wolfram Language}}== |
|||
<syntaxhighlight lang="mathematica">ClearAll[a] |
|||
a[1] = 1; |
|||
a[n_?(GreaterThan[1])] := a[n] = 1/(2 Floor[a[n - 1]] + 1 - a[n - 1]) |
|||
a /@ Range[20] |
|||
ClearAll[a] |
|||
function to_continued_fraction(sequence r) |
|||
a = 1; |
|||
integer {a,b} = r |
|||
n = 1; |
|||
sequence res = {} |
|||
Dynamic[n] |
|||
while true do |
|||
done = False; |
|||
res &= floor(a/b) |
|||
While[! done, |
|||
{a, b} = {b, remainder(a,b)} |
|||
a = 1/(2 Floor[a] + 1 - a); |
|||
if a=1 then exit end if |
|||
n++; |
|||
end while |
|||
If[a == 83116/51639, |
|||
return res |
|||
Print[n]; |
|||
end function |
|||
Break[]; |
|||
] |
|||
function get_term_number(sequence cf) |
|||
]</syntaxhighlight> |
|||
sequence b = {} |
|||
{{out}} |
|||
integer d = 1 |
|||
<pre>{1, 1/2, 2, 1/3, 3/2, 2/3, 3, 1/4, 4/3, 3/5, 5/2, 2/5, 5/3, 3/4, 4, 1/5, 5/4, 4/7, 7/3, 3/8} |
|||
for i=1 to length(cf) do |
|||
123456789</pre> |
|||
b &= repeat(d,cf[i]) |
|||
d = 1-d |
|||
end for |
|||
return bits_to_int(b) |
|||
end function |
|||
=={{header|Maxima}}== |
|||
-- additional verification methods (2 of) |
|||
<syntaxhighlight lang="maxima"> |
|||
function i_to_cf(integer i) |
|||
/* The function fusc is related to Calkin-Wilf sequence */ |
|||
-- sequence b = trim_tail(int_to_bits(i,32),0)&2, |
|||
fusc(n):=block( |
|||
sequence b = int_to_bits(i)&2, |
|||
[k:n,a:1,b:0], |
|||
cf = iff(b[1]=0?{0}:{}) |
|||
while k>0 do (if evenp(k) then (k:k/2,a:a+b) else (k:(k-1)/2,b:a+b)), |
|||
while length(b)>1 do |
|||
b)$ |
|||
for j=2 to length(b) do |
|||
if b[j]!=b[1] then |
|||
cf &= j-1 |
|||
b = b[j..$] |
|||
exit |
|||
end if |
|||
end for |
|||
end while |
|||
-- replace even length with odd length equivalent: |
|||
if remainder(length(cf),2)=0 then |
|||
cf[$] -= 1 |
|||
cf &= 1 |
|||
end if |
|||
return cf |
|||
end function |
|||
function |
/* Calkin-Wilf function using fusc */ |
||
calkin_wilf(n):=fusc(n)/fusc(n+1)$ |
|||
integer n=0, d=1 |
|||
for i=length(cf) to 2 by -1 do |
|||
{n,d} = {d,n+d*cf[i]} |
|||
end for |
|||
return {n+cf[1]*d,d} |
|||
end function |
|||
/* Function that given a nonnegative rational returns its position in the Calkin-Wilf sequence */ |
|||
function prettyr(sequence r) |
|||
cf_bin(fracti):=block( |
|||
integer {n,d} = r |
|||
cf_list:cf(fracti), |
|||
return iff(d=1?sprintf("%d",n):sprintf("%d/%d",{n,d})) |
|||
cf_len:length(cf_list), |
|||
end function |
|||
if oddp(cf_len) then cf_list:reverse(cf_list) else cf_list:reverse(append(at(cf_list,[cf_list[cf_len]=cf_list[cf_len]-1]),[1])), |
|||
makelist(lambda([x],if oddp(x) then makelist(1,j,1,cf_list[x]) else makelist(0,j,1,cf_list[x]))(i),i,1,length(cf_list)), /* decoding part begins here */ |
|||
apply(append,%%), |
|||
apply("+",makelist(2^i,i,0,length(%%)-1)*reverse(%%)))$ |
|||
/* Test cases */ |
|||
sequence cw = calkin_wilf(20) |
|||
/* 20 first terms of the sequence */ |
|||
makelist(calkin_wilf(i),i,1,20); |
|||
for i=1 to 20 do |
|||
string s = prettyr(cw[i]), |
|||
/* Position of 83116/51639 in Calkin-Wilf sequence */ |
|||
r = prettyr(cf2r(i_to_cf(i))) |
|||
83116/51639$ |
|||
integer t = get_term_number(to_continued_fraction(cw[i])) |
|||
cf_bin(%); |
|||
printf(1,"%2d: %-4s [==> %2d: %-3s]\n", {i, s, t, r}) |
|||
</syntaxhighlight> |
|||
end for |
|||
{{out}} |
|||
printf(1,"\n") |
|||
<pre> |
|||
sequence r = {83116,51639} |
|||
[1,1/2,2,1/3,3/2,2/3,3,1/4,4/3,3/5,5/2,2/5,5/3,3/4,4,1/5,5/4,4/7,7/3,3/8] |
|||
sequence cf = to_continued_fraction(r) |
|||
integer tn = get_term_number(cf) |
|||
123456789 |
|||
printf(1,"%d/%d is the %,d%s term of the sequence.\n", r&{tn,ord(tn)})</lang> |
|||
</pre> |
|||
=={{header|Nim}}== |
|||
We ignored the standard module “rationals” which is slow and have rather defined a fraction as a tuple of two 32 bits unsigned integers (slightly faster than 64 bits signed integers and sufficient for this task). Moreover, we didn’t do operations on fractions and computed directly the numerator and denominator values at each step. The fractions built this way are irreducible (which avoids to compute a GCD which is a slow operation). |
|||
With these optimizations, the program runs in less than 1.3 s on our laptop. |
|||
<syntaxhighlight lang="nim">type Fraction = tuple[num, den: uint32] |
|||
iterator calkinWilf(): Fraction = |
|||
## Yield the successive values of the sequence. |
|||
var n, d = 1u32 |
|||
yield (n, d) |
|||
while true: |
|||
n = 2 * (n div d) * d + d - n |
|||
swap n, d |
|||
yield (n, d) |
|||
proc `$`(fract: Fraction): string = |
|||
## Return the representation of a fraction. |
|||
$fract.num & '/' & $fract.den |
|||
func `==`(a, b: Fraction): bool {.inline.} = |
|||
## Compare two fractions. Slightly faster than comparison of tuples. |
|||
a.num == b.num and a.den == b.den |
|||
when isMainModule: |
|||
echo "The first 20 terms of the Calkwin-Wilf sequence are:" |
|||
var count = 0 |
|||
for an in calkinWilf(): |
|||
inc count |
|||
stdout.write $an & ' ' |
|||
if count == 20: break |
|||
stdout.write '\n' |
|||
const Target: Fraction = (83116u32, 51639u32) |
|||
var index = 0 |
|||
for an in calkinWilf(): |
|||
inc index |
|||
if an == Target: break |
|||
echo "\nThe element ", $Target, " is at position ", $index, " in the sequence."</syntaxhighlight> |
|||
{{out}} |
|||
<pre>The first 20 terms of the Calkwin-Wilf sequence are: |
|||
1/1 1/2 2/1 1/3 3/2 2/3 3/1 1/4 4/3 3/5 5/2 2/5 5/3 3/4 4/1 1/5 5/4 4/7 7/3 3/8 |
|||
The element 83116/51639 is at position 123456789 in the sequence.</pre> |
|||
=={{header|PARI/GP}}== |
|||
{{trans|Mathematica_/_Wolfram_Language}} |
|||
<syntaxhighlight lang="PARI/GP"> |
|||
\\ This function assumes the existence of a global variable 'an' for 'a[n]' |
|||
a(n) = if(n==1, 1, 1 / (2 * floor(an[n-1]) + 1 - an[n-1])); |
|||
\\ We will use a vector to hold the values and compute them iteratively to avoid stack overflow |
|||
an = vector(20); |
|||
an[1] = 1; |
|||
for(i=2, 20, an[i] = a(i)); |
|||
\\ Now we print the vector |
|||
print(an); |
|||
\\ Initialize variables for the while loop |
|||
a = 1; |
|||
n = 1; |
|||
\\ Loop until the condition is met |
|||
while(a != 83116/51639,{ |
|||
a = 1/(2 * floor(a) + 1 - a); |
|||
if(n>=123456789,print(n)); |
|||
n++; |
|||
}); |
|||
\\ Output the number of iterations needed to reach 83116/51639 |
|||
print(n); |
|||
</syntaxhighlight> |
|||
{{out}} |
|||
<pre> |
|||
[1, 1/2, 2, 1/3, 3/2, 2/3, 3, 1/4, 4/3, 3/5, 5/2, 2/5, 5/3, 3/4, 4, 1/5, 5/4, 4/7, 7/3, 3/8] |
|||
123456789 |
|||
</pre> |
|||
=={{header|Pascal}}== |
|||
These programs were written in Free Pascal, using the Lazarus IDE and the Free Pascal compiler version 3.2.0. They are based on the Wikipedia article "Calkin-Wilf tree", rather than the algorithms in the task description. |
|||
<syntaxhighlight lang="pascal"> |
|||
program CWTerms; |
|||
{------------------------------------------------------------------------------- |
|||
FreePascal command-line program. |
|||
Calculates the Calkin-Wilf sequence up to the specified maximum index, |
|||
where the first term 1/1 has index 1. |
|||
Command line format is: CWTerms <max_index> |
|||
The program demonstrates 3 algorithms for calculating the sequence: |
|||
(1) Calculate term[2n] and term[2n + 1] from term[n] |
|||
(2) Calculate term[n + 1] from term[n] |
|||
(3) Calculate term[n] directly from n, without using other terms |
|||
Algorithm 1 is called first, and stores the terms in an array. |
|||
Then the program calls Algorithms 2 and 3, and checks that they agree |
|||
with Algorithm 1. |
|||
-------------------------------------------------------------------------------} |
|||
uses SysUtils; |
|||
type TRational = record |
|||
Num, Den : integer; |
|||
end; |
|||
var |
|||
terms : array of TRational; |
|||
max_index, k : integer; |
|||
// Routine to calculate array of terms up the the maiximum index |
|||
procedure CalcTerms_algo_1(); |
|||
var |
|||
j, k : integer; |
|||
begin |
|||
SetLength( terms, max_index + 1); |
|||
j := 1; // index to earlier term, from which current term is calculated |
|||
k := 1; // index to current term |
|||
terms[1].Num := 1; |
|||
terms[1].Den := 1; |
|||
while (k < max_index) do begin |
|||
inc(k); |
|||
if (k and 1) = 0 then begin // or could write "if not Odd(k)" |
|||
terms[k].Num := terms[j].Num; |
|||
terms[k].Den := terms[j].Num + terms[j].Den; |
|||
end |
|||
else begin |
|||
terms[k].Num := terms[j].Num + terms[j].Den; |
|||
terms[k].Den := terms[j].Den; |
|||
inc(j); |
|||
end; |
|||
end; |
|||
end; |
|||
// Method to get each term from the preceding term. |
|||
// a/b --> b/(a + b - 2(a mod b)); |
|||
function CheckTerms_algo_2() : boolean; |
|||
var |
|||
index, a, b, temp : integer; |
|||
begin |
|||
result := true; |
|||
index := 1; |
|||
a := 1; |
|||
b := 1; |
|||
while (index <= max_index) do begin |
|||
if (a <> terms[index].Num) or (b <> terms[index].Den) then |
|||
result := false; |
|||
temp := a + b - 2*(a mod b); |
|||
a := b; |
|||
b := temp; |
|||
inc( index) |
|||
end; |
|||
end; |
|||
// Mathod to calcualte each term from its index, without using other terms. |
|||
function CheckTerms_algo_3() : boolean; |
|||
var |
|||
index, a, b, pwr2, idiv2 : integer; |
|||
begin |
|||
result := true; |
|||
for index := 1 to max_index do begin |
|||
idiv2 := index div 2; |
|||
pwr2 := 1; |
|||
while (pwr2 <= idiv2) do pwr2 := pwr2 shl 1; |
|||
a := 1; |
|||
b := 1; |
|||
while (pwr2 > 1) do begin |
|||
pwr2 := pwr2 shr 1; |
|||
if (pwr2 and index) = 0 then |
|||
inc( b, a) |
|||
else |
|||
inc( a, b); |
|||
end; |
|||
if (a <> terms[index].Num) or (b <> terms[index].Den) then |
|||
result := false; |
|||
end; |
|||
end; |
|||
begin |
|||
// Read and validate maximum index |
|||
max_index := SysUtils.StrToIntDef( paramStr(1), -1); // -1 if not an integer |
|||
if (max_index <= 0) then begin |
|||
WriteLn( 'Maximum index must be a positive integer'); |
|||
exit; |
|||
end; |
|||
// Calculate terms by algo 1, then check that algos 2 and 3 agree. |
|||
CalcTerms_algo_1(); |
|||
if not CheckTerms_algo_2() then begin |
|||
WriteLn( 'Algorithm 2 failed'); |
|||
exit; |
|||
end; |
|||
if not CheckTerms_algo_3() then begin |
|||
WriteLn( 'Algorithm 3 failed'); |
|||
exit; |
|||
end; |
|||
// Display the terms |
|||
for k := 1 to max_index do |
|||
with terms[k] do |
|||
WriteLn( SysUtils.Format( '%8d: %d/%d', [k, Num, Den])); |
|||
end. |
|||
</syntaxhighlight> |
|||
{{out}} |
|||
<pre> |
|||
1: 1/1 |
|||
2: 1/2 |
|||
3: 2/1 |
|||
4: 1/3 |
|||
5: 3/2 |
|||
6: 2/3 |
|||
7: 3/1 |
|||
8: 1/4 |
|||
9: 4/3 |
|||
10: 3/5 |
|||
11: 5/2 |
|||
12: 2/5 |
|||
13: 5/3 |
|||
14: 3/4 |
|||
15: 4/1 |
|||
16: 1/5 |
|||
17: 5/4 |
|||
18: 4/7 |
|||
19: 7/3 |
|||
20: 3/8 |
|||
</pre> |
|||
<syntaxhighlight lang="pascal"> |
|||
program CWIndex; |
|||
{------------------------------------------------------------------------------- |
|||
FreePascal command-line program. |
|||
Calculates index of a rational number in the Calkin-Wilf sequence, |
|||
where the first term 1/1 has index 1. |
|||
Command line format is |
|||
CWIndex <numerator> <denominator> |
|||
e.g. for the Rosetta Code example |
|||
CWIndex 83116 51639 |
|||
-------------------------------------------------------------------------------} |
|||
uses SysUtils; |
|||
var |
|||
num, den : integer; |
|||
a, b : integer; |
|||
pwr2, index : qword; // 64-bit unsiged |
|||
begin |
|||
// Read and validate input. |
|||
num := SysUtils.StrToIntDef( paramStr(1), -1); // return -1 if not an integer |
|||
den := SysUtils.StrToIntDef( paramStr(2), -1); |
|||
if (num <= 0) or (den <= 0) then begin |
|||
WriteLn( 'Numerator and denominator must be positive integers'); |
|||
exit; |
|||
end; |
|||
// Input OK, calculate and display index of num/den |
|||
// The index may overflow 64 bits, so turn on overflow detection |
|||
{$Q+} |
|||
a := num; |
|||
b := den; |
|||
pwr2 := 1; |
|||
index := 0; |
|||
try |
|||
while (a <> b) do begin |
|||
if (a < b) then |
|||
dec( b, a) |
|||
else begin |
|||
dec( a, b); |
|||
inc( index, pwr2); |
|||
end; |
|||
pwr2 := 2*pwr2; |
|||
end; |
|||
inc( index, pwr2); |
|||
WriteLn( SysUtils.Format( 'Index of %d/%d is %u', [num, den, index])); |
|||
except |
|||
WriteLn( 'Index is too large for 64 bits'); |
|||
end; |
|||
end. |
|||
</syntaxhighlight> |
|||
{{out}} |
|||
<pre> |
|||
Index of 83116/51639 is 123456789 |
|||
</pre> |
|||
=={{header|Perl}}== |
|||
{{trans|Raku}} |
|||
{{libheader|ntheory}} |
|||
<syntaxhighlight lang="perl">use strict; |
|||
use warnings; |
|||
use feature qw(say state); |
|||
use ntheory 'fromdigits'; |
|||
use List::Lazy 'lazy_list'; |
|||
use Math::AnyNum ':overload'; |
|||
my $calkin_wilf = lazy_list { state @cw = 1; push @cw, 1 / ( (2 * int $cw[0]) + 1 - $cw[0] ); shift @cw }; |
|||
sub r2cf { |
|||
my($num, $den) = @_; |
|||
my($n, @cf); |
|||
my $f = sub { return unless $den; |
|||
my $q = int($num/$den); |
|||
($num, $den) = ($den, $num - $q*$den); |
|||
$q; |
|||
}; |
|||
push @cf, $n while defined($n = $f->()); |
|||
reverse @cf; |
|||
} |
|||
sub r2cw { |
|||
my($num, $den) = @_; |
|||
my $bits; |
|||
my @f = r2cf($num, $den); |
|||
$bits .= ($_%2 ? 0 : 1) x $f[$_] for 0..$#f; |
|||
fromdigits($bits, 2); |
|||
} |
|||
say 'First twenty terms of the Calkin-Wilf sequence:'; |
|||
printf "%s ", $calkin_wilf->next() for 1..20; |
|||
say "\n\n83116/51639 is at index: " . r2cw(83116,51639);</syntaxhighlight> |
|||
{{out}} |
|||
<pre>First twenty terms of the Calkin-Wilf sequence: |
|||
1 1/2 2 1/3 3/2 2/3 3 1/4 4/3 3/5 5/2 2/5 5/3 3/4 4 1/5 5/4 4/7 7/3 3/8 |
|||
83116/51639 is at index: 123456789</pre> |
|||
=={{header|Phix}}== |
|||
<!--<syntaxhighlight lang="phix">(phixonline)--> |
|||
<span style="color: #008080;">with</span> <span style="color: #008080;">javascript_semantics</span> |
|||
<span style="color: #7060A8;">requires</span><span style="color: #0000FF;">(</span><span style="color: #008000;">"1.0.0"</span><span style="color: #0000FF;">)</span> <span style="color: #000080;font-style:italic;">-- (new even() builtin)</span> |
|||
<span style="color: #008080;">function</span> <span style="color: #000000;">calkin_wilf</span><span style="color: #0000FF;">(</span><span style="color: #004080;">integer</span> <span style="color: #000000;">len</span><span style="color: #0000FF;">)</span> |
|||
<span style="color: #004080;">sequence</span> <span style="color: #000000;">cw</span> <span style="color: #0000FF;">=</span> <span style="color: #7060A8;">repeat</span><span style="color: #0000FF;">(</span><span style="color: #000000;">0</span><span style="color: #0000FF;">,</span><span style="color: #000000;">len</span><span style="color: #0000FF;">)</span> |
|||
<span style="color: #004080;">integer</span> <span style="color: #000000;">n</span><span style="color: #0000FF;">=</span><span style="color: #000000;">0</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">d</span><span style="color: #0000FF;">=</span><span style="color: #000000;">1</span> |
|||
<span style="color: #008080;">for</span> <span style="color: #000000;">i</span><span style="color: #0000FF;">=</span><span style="color: #000000;">1</span> <span style="color: #008080;">to</span> <span style="color: #000000;">len</span> <span style="color: #008080;">do</span> |
|||
<span style="color: #0000FF;">{</span><span style="color: #000000;">n</span><span style="color: #0000FF;">,</span><span style="color: #000000;">d</span><span style="color: #0000FF;">}</span> <span style="color: #0000FF;">=</span> <span style="color: #0000FF;">{</span><span style="color: #000000;">d</span><span style="color: #0000FF;">,(</span><span style="color: #7060A8;">floor</span><span style="color: #0000FF;">(</span><span style="color: #000000;">n</span><span style="color: #0000FF;">/</span><span style="color: #000000;">d</span><span style="color: #0000FF;">)*</span><span style="color: #000000;">2</span><span style="color: #0000FF;">+</span><span style="color: #000000;">1</span><span style="color: #0000FF;">)*</span><span style="color: #000000;">d</span><span style="color: #0000FF;">-</span><span style="color: #000000;">n</span><span style="color: #0000FF;">}</span> |
|||
<span style="color: #000000;">cw</span><span style="color: #0000FF;">[</span><span style="color: #000000;">i</span><span style="color: #0000FF;">]</span> <span style="color: #0000FF;">=</span> <span style="color: #0000FF;">{</span><span style="color: #000000;">n</span><span style="color: #0000FF;">,</span><span style="color: #000000;">d</span><span style="color: #0000FF;">}</span> |
|||
<span style="color: #008080;">end</span> <span style="color: #008080;">for</span> |
|||
<span style="color: #008080;">return</span> <span style="color: #000000;">cw</span> |
|||
<span style="color: #008080;">end</span> <span style="color: #008080;">function</span> |
|||
<span style="color: #008080;">function</span> <span style="color: #000000;">odd_length</span><span style="color: #0000FF;">(</span><span style="color: #004080;">sequence</span> <span style="color: #000000;">cf</span><span style="color: #0000FF;">)</span> |
|||
<span style="color: #000080;font-style:italic;">-- replace even length continued fraction with odd length equivalent |
|||
-- if remainder(length(cf),2)=0 then</span> |
|||
<span style="color: #008080;">if</span> <span style="color: #7060A8;">even</span><span style="color: #0000FF;">(</span><span style="color: #7060A8;">length</span><span style="color: #0000FF;">(</span><span style="color: #000000;">cf</span><span style="color: #0000FF;">))</span> <span style="color: #008080;">then</span> |
|||
<span style="color: #000000;">cf</span><span style="color: #0000FF;">[$]</span> <span style="color: #0000FF;">-=</span> <span style="color: #000000;">1</span> |
|||
<span style="color: #000000;">cf</span> <span style="color: #0000FF;">&=</span> <span style="color: #000000;">1</span> |
|||
<span style="color: #008080;">end</span> <span style="color: #008080;">if</span> |
|||
<span style="color: #008080;">return</span> <span style="color: #000000;">cf</span> |
|||
<span style="color: #008080;">end</span> <span style="color: #008080;">function</span> |
|||
<span style="color: #008080;">function</span> <span style="color: #000000;">to_continued_fraction</span><span style="color: #0000FF;">(</span><span style="color: #004080;">sequence</span> <span style="color: #000000;">r</span><span style="color: #0000FF;">)</span> |
|||
<span style="color: #004080;">integer</span> <span style="color: #0000FF;">{</span><span style="color: #000000;">a</span><span style="color: #0000FF;">,</span><span style="color: #000000;">b</span><span style="color: #0000FF;">}</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">r</span> |
|||
<span style="color: #004080;">sequence</span> <span style="color: #000000;">cf</span> <span style="color: #0000FF;">=</span> <span style="color: #0000FF;">{}</span> |
|||
<span style="color: #008080;">while</span> <span style="color: #004600;">true</span> <span style="color: #008080;">do</span> |
|||
<span style="color: #000000;">cf</span> <span style="color: #0000FF;">&=</span> <span style="color: #7060A8;">floor</span><span style="color: #0000FF;">(</span><span style="color: #000000;">a</span><span style="color: #0000FF;">/</span><span style="color: #000000;">b</span><span style="color: #0000FF;">)</span> |
|||
<span style="color: #0000FF;">{</span><span style="color: #000000;">a</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">b</span><span style="color: #0000FF;">}</span> <span style="color: #0000FF;">=</span> <span style="color: #0000FF;">{</span><span style="color: #000000;">b</span><span style="color: #0000FF;">,</span> <span style="color: #7060A8;">remainder</span><span style="color: #0000FF;">(</span><span style="color: #000000;">a</span><span style="color: #0000FF;">,</span><span style="color: #000000;">b</span><span style="color: #0000FF;">)}</span> |
|||
<span style="color: #008080;">if</span> <span style="color: #000000;">a</span><span style="color: #0000FF;">=</span><span style="color: #000000;">1</span> <span style="color: #008080;">then</span> <span style="color: #008080;">exit</span> <span style="color: #008080;">end</span> <span style="color: #008080;">if</span> |
|||
<span style="color: #008080;">end</span> <span style="color: #008080;">while</span> |
|||
<span style="color: #000000;">cf</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">odd_length</span><span style="color: #0000FF;">(</span><span style="color: #000000;">cf</span><span style="color: #0000FF;">)</span> |
|||
<span style="color: #008080;">return</span> <span style="color: #000000;">cf</span> |
|||
<span style="color: #008080;">end</span> <span style="color: #008080;">function</span> |
|||
<span style="color: #008080;">function</span> <span style="color: #000000;">get_term_number</span><span style="color: #0000FF;">(</span><span style="color: #004080;">sequence</span> <span style="color: #000000;">cf</span><span style="color: #0000FF;">)</span> |
|||
<span style="color: #004080;">sequence</span> <span style="color: #000000;">b</span> <span style="color: #0000FF;">=</span> <span style="color: #0000FF;">{}</span> |
|||
<span style="color: #004080;">integer</span> <span style="color: #000000;">d</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">1</span> |
|||
<span style="color: #008080;">for</span> <span style="color: #000000;">i</span><span style="color: #0000FF;">=</span><span style="color: #000000;">1</span> <span style="color: #008080;">to</span> <span style="color: #7060A8;">length</span><span style="color: #0000FF;">(</span><span style="color: #000000;">cf</span><span style="color: #0000FF;">)</span> <span style="color: #008080;">do</span> |
|||
<span style="color: #000000;">b</span> <span style="color: #0000FF;">&=</span> <span style="color: #7060A8;">repeat</span><span style="color: #0000FF;">(</span><span style="color: #000000;">d</span><span style="color: #0000FF;">,</span><span style="color: #000000;">cf</span><span style="color: #0000FF;">[</span><span style="color: #000000;">i</span><span style="color: #0000FF;">])</span> |
|||
<span style="color: #000000;">d</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">1</span><span style="color: #0000FF;">-</span><span style="color: #000000;">d</span> |
|||
<span style="color: #008080;">end</span> <span style="color: #008080;">for</span> |
|||
<span style="color: #004080;">integer</span> <span style="color: #000000;">t</span> <span style="color: #0000FF;">=</span> <span style="color: #7060A8;">bits_to_int</span><span style="color: #0000FF;">(</span><span style="color: #000000;">b</span><span style="color: #0000FF;">)</span> |
|||
<span style="color: #008080;">return</span> <span style="color: #000000;">t</span> |
|||
<span style="color: #008080;">end</span> <span style="color: #008080;">function</span> |
|||
<span style="color: #000080;font-style:italic;">-- additional verification methods (2 of)</span> |
|||
<span style="color: #008080;">function</span> <span style="color: #000000;">i_to_cf</span><span style="color: #0000FF;">(</span><span style="color: #004080;">integer</span> <span style="color: #000000;">i</span><span style="color: #0000FF;">)</span> |
|||
<span style="color: #000080;font-style:italic;">-- sequence b = trim_tail(int_to_bits(i,32),0)&2,</span> |
|||
<span style="color: #004080;">sequence</span> <span style="color: #000000;">b</span> <span style="color: #0000FF;">=</span> <span style="color: #7060A8;">int_to_bits</span><span style="color: #0000FF;">(</span><span style="color: #000000;">i</span><span style="color: #0000FF;">)&</span><span style="color: #000000;">2</span><span style="color: #0000FF;">,</span> |
|||
<span style="color: #000000;">cf</span> <span style="color: #0000FF;">=</span> <span style="color: #008080;">iff</span><span style="color: #0000FF;">(</span><span style="color: #000000;">b</span><span style="color: #0000FF;">[</span><span style="color: #000000;">1</span><span style="color: #0000FF;">]=</span><span style="color: #000000;">0</span><span style="color: #0000FF;">?{</span><span style="color: #000000;">0</span><span style="color: #0000FF;">}:{})</span> |
|||
<span style="color: #008080;">while</span> <span style="color: #7060A8;">length</span><span style="color: #0000FF;">(</span><span style="color: #000000;">b</span><span style="color: #0000FF;">)></span><span style="color: #000000;">1</span> <span style="color: #008080;">do</span> |
|||
<span style="color: #008080;">for</span> <span style="color: #000000;">j</span><span style="color: #0000FF;">=</span><span style="color: #000000;">2</span> <span style="color: #008080;">to</span> <span style="color: #7060A8;">length</span><span style="color: #0000FF;">(</span><span style="color: #000000;">b</span><span style="color: #0000FF;">)</span> <span style="color: #008080;">do</span> |
|||
<span style="color: #008080;">if</span> <span style="color: #000000;">b</span><span style="color: #0000FF;">[</span><span style="color: #000000;">j</span><span style="color: #0000FF;">]!=</span><span style="color: #000000;">b</span><span style="color: #0000FF;">[</span><span style="color: #000000;">1</span><span style="color: #0000FF;">]</span> <span style="color: #008080;">then</span> |
|||
<span style="color: #000000;">cf</span> <span style="color: #0000FF;">&=</span> <span style="color: #000000;">j</span><span style="color: #0000FF;">-</span><span style="color: #000000;">1</span> |
|||
<span style="color: #000000;">b</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">b</span><span style="color: #0000FF;">[</span><span style="color: #000000;">j</span><span style="color: #0000FF;">..$]</span> |
|||
<span style="color: #008080;">exit</span> |
|||
<span style="color: #008080;">end</span> <span style="color: #008080;">if</span> |
|||
<span style="color: #008080;">end</span> <span style="color: #008080;">for</span> |
|||
<span style="color: #008080;">end</span> <span style="color: #008080;">while</span> |
|||
<span style="color: #000000;">cf</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">odd_length</span><span style="color: #0000FF;">(</span><span style="color: #000000;">cf</span><span style="color: #0000FF;">)</span> |
|||
<span style="color: #008080;">return</span> <span style="color: #000000;">cf</span> |
|||
<span style="color: #008080;">end</span> <span style="color: #008080;">function</span> |
|||
<span style="color: #008080;">function</span> <span style="color: #000000;">cf2r</span><span style="color: #0000FF;">(</span><span style="color: #004080;">sequence</span> <span style="color: #000000;">cf</span><span style="color: #0000FF;">)</span> |
|||
<span style="color: #004080;">integer</span> <span style="color: #000000;">n</span><span style="color: #0000FF;">=</span><span style="color: #000000;">0</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">d</span><span style="color: #0000FF;">=</span><span style="color: #000000;">1</span> |
|||
<span style="color: #008080;">for</span> <span style="color: #000000;">i</span><span style="color: #0000FF;">=</span><span style="color: #7060A8;">length</span><span style="color: #0000FF;">(</span><span style="color: #000000;">cf</span><span style="color: #0000FF;">)</span> <span style="color: #008080;">to</span> <span style="color: #000000;">2</span> <span style="color: #008080;">by</span> <span style="color: #0000FF;">-</span><span style="color: #000000;">1</span> <span style="color: #008080;">do</span> |
|||
<span style="color: #0000FF;">{</span><span style="color: #000000;">n</span><span style="color: #0000FF;">,</span><span style="color: #000000;">d</span><span style="color: #0000FF;">}</span> <span style="color: #0000FF;">=</span> <span style="color: #0000FF;">{</span><span style="color: #000000;">d</span><span style="color: #0000FF;">,</span><span style="color: #000000;">n</span><span style="color: #0000FF;">+</span><span style="color: #000000;">d</span><span style="color: #0000FF;">*</span><span style="color: #000000;">cf</span><span style="color: #0000FF;">[</span><span style="color: #000000;">i</span><span style="color: #0000FF;">]}</span> |
|||
<span style="color: #008080;">end</span> <span style="color: #008080;">for</span> |
|||
<span style="color: #008080;">return</span> <span style="color: #0000FF;">{</span><span style="color: #000000;">n</span><span style="color: #0000FF;">+</span><span style="color: #000000;">cf</span><span style="color: #0000FF;">[</span><span style="color: #000000;">1</span><span style="color: #0000FF;">]*</span><span style="color: #000000;">d</span><span style="color: #0000FF;">,</span><span style="color: #000000;">d</span><span style="color: #0000FF;">}</span> |
|||
<span style="color: #008080;">end</span> <span style="color: #008080;">function</span> |
|||
<span style="color: #008080;">function</span> <span style="color: #000000;">prettyr</span><span style="color: #0000FF;">(</span><span style="color: #004080;">sequence</span> <span style="color: #000000;">r</span><span style="color: #0000FF;">)</span> |
|||
<span style="color: #004080;">integer</span> <span style="color: #0000FF;">{</span><span style="color: #000000;">n</span><span style="color: #0000FF;">,</span><span style="color: #000000;">d</span><span style="color: #0000FF;">}</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">r</span> |
|||
<span style="color: #008080;">return</span> <span style="color: #008080;">iff</span><span style="color: #0000FF;">(</span><span style="color: #000000;">d</span><span style="color: #0000FF;">=</span><span style="color: #000000;">1</span><span style="color: #0000FF;">?</span><span style="color: #7060A8;">sprintf</span><span style="color: #0000FF;">(</span><span style="color: #008000;">"%d"</span><span style="color: #0000FF;">,</span><span style="color: #000000;">n</span><span style="color: #0000FF;">):</span><span style="color: #7060A8;">sprintf</span><span style="color: #0000FF;">(</span><span style="color: #008000;">"%d/%d"</span><span style="color: #0000FF;">,{</span><span style="color: #000000;">n</span><span style="color: #0000FF;">,</span><span style="color: #000000;">d</span><span style="color: #0000FF;">}))</span> |
|||
<span style="color: #008080;">end</span> <span style="color: #008080;">function</span> |
|||
<span style="color: #004080;">sequence</span> <span style="color: #000000;">cw</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">calkin_wilf</span><span style="color: #0000FF;">(</span><span style="color: #000000;">20</span><span style="color: #0000FF;">)</span> |
|||
<span style="color: #7060A8;">printf</span><span style="color: #0000FF;">(</span><span style="color: #000000;">1</span><span style="color: #0000FF;">,</span><span style="color: #008000;">"The first 20 terms of the Calkin-Wilf sequence are:\n"</span><span style="color: #0000FF;">)</span> |
|||
<span style="color: #008080;">for</span> <span style="color: #000000;">i</span><span style="color: #0000FF;">=</span><span style="color: #000000;">1</span> <span style="color: #008080;">to</span> <span style="color: #000000;">20</span> <span style="color: #008080;">do</span> |
|||
<span style="color: #004080;">string</span> <span style="color: #000000;">s</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">prettyr</span><span style="color: #0000FF;">(</span><span style="color: #000000;">cw</span><span style="color: #0000FF;">[</span><span style="color: #000000;">i</span><span style="color: #0000FF;">]),</span> |
|||
<span style="color: #000000;">r</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">prettyr</span><span style="color: #0000FF;">(</span><span style="color: #000000;">cf2r</span><span style="color: #0000FF;">(</span><span style="color: #000000;">i_to_cf</span><span style="color: #0000FF;">(</span><span style="color: #000000;">i</span><span style="color: #0000FF;">)))</span> |
|||
<span style="color: #004080;">integer</span> <span style="color: #000000;">t</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">get_term_number</span><span style="color: #0000FF;">(</span><span style="color: #000000;">to_continued_fraction</span><span style="color: #0000FF;">(</span><span style="color: #000000;">cw</span><span style="color: #0000FF;">[</span><span style="color: #000000;">i</span><span style="color: #0000FF;">]))</span> |
|||
<span style="color: #7060A8;">printf</span><span style="color: #0000FF;">(</span><span style="color: #000000;">1</span><span style="color: #0000FF;">,</span><span style="color: #008000;">"%2d: %-4s [==> %2d: %-3s]\n"</span><span style="color: #0000FF;">,</span> <span style="color: #0000FF;">{</span><span style="color: #000000;">i</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">s</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">t</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">r</span><span style="color: #0000FF;">})</span> |
|||
<span style="color: #008080;">end</span> <span style="color: #008080;">for</span> |
|||
<span style="color: #7060A8;">printf</span><span style="color: #0000FF;">(</span><span style="color: #000000;">1</span><span style="color: #0000FF;">,</span><span style="color: #008000;">"\n"</span><span style="color: #0000FF;">)</span> |
|||
<span style="color: #004080;">sequence</span> <span style="color: #000000;">r</span> <span style="color: #0000FF;">=</span> <span style="color: #0000FF;">{</span><span style="color: #000000;">83116</span><span style="color: #0000FF;">,</span><span style="color: #000000;">51639</span><span style="color: #0000FF;">}</span> |
|||
<span style="color: #004080;">sequence</span> <span style="color: #000000;">cf</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">to_continued_fraction</span><span style="color: #0000FF;">(</span><span style="color: #000000;">r</span><span style="color: #0000FF;">)</span> |
|||
<span style="color: #004080;">integer</span> <span style="color: #000000;">tn</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">get_term_number</span><span style="color: #0000FF;">(</span><span style="color: #000000;">cf</span><span style="color: #0000FF;">)</span> |
|||
<span style="color: #7060A8;">printf</span><span style="color: #0000FF;">(</span><span style="color: #000000;">1</span><span style="color: #0000FF;">,</span><span style="color: #008000;">"%d/%d is the %,d%s term of the sequence.\n"</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">r</span><span style="color: #0000FF;">&{</span><span style="color: #000000;">tn</span><span style="color: #0000FF;">,</span><span style="color: #7060A8;">ord</span><span style="color: #0000FF;">(</span><span style="color: #000000;">tn</span><span style="color: #0000FF;">)})</span> |
|||
<!--</syntaxhighlight>--> |
|||
{{out}} |
{{out}} |
||
<pre> |
<pre> |
||
Line 695: | Line 2,285: | ||
83116/51639 is the 123,456,789th term of the sequence. |
83116/51639 is the 123,456,789th term of the sequence. |
||
</pre> |
</pre> |
||
=={{header|Prolog}}== |
|||
<syntaxhighlight lang="prolog"> |
|||
% John Devou: 26-Nov-2021 |
|||
% g(N,X):- consecutively generate in X the first N elements of the Calkin-Wilf sequence |
|||
g(N,[A/B|_]-_,A/B):- N > 0. |
|||
g(N,[A/B|Ls]-[A/C,C/B|Ys],X):- N > 1, M is N-1, C is A+B, g(M,Ls-Ys,X). |
|||
g(N,X):- g(N,[1/1|Ls]-Ls,X). |
|||
% t(A/B,X):- generate in X the index of A/B in the Calkin-Wilf sequence |
|||
t(A/1,S,C,X):- X is C*(2**(A-1+S)-S). |
|||
t(A/B,S,C,X):- B > 1, divmod(A,B,M,N), T is 1-S, D is C*2**M, t(B/N,T,D,Y), X is Y + S*C*(2**M-1). |
|||
t(A/B,X):- t(A/B,1,1,X), !. |
|||
</syntaxhighlight> |
|||
{{out}} |
|||
<pre> |
|||
?- findall(X, g(20,X), L), write(L). |
|||
[1/1,1/2,2/1,1/3,3/2,2/3,3/1,1/4,4/3,3/5,5/2,2/5,5/3,3/4,4/1,1/5,5/4,4/7,7/3,3/8] |
|||
L = [1/1, 1/2, 2/1, 1/3, 3/2, 2/3, 3/1, 1/4, ... / ...|...]. |
|||
?- t(83116/51639,X). |
|||
X = 123456789. |
|||
</pre> |
|||
=={{header|Python}}== |
=={{header|Python}}== |
||
< |
<syntaxhighlight lang="python">from fractions import Fraction |
||
from math import floor |
from math import floor |
||
from itertools import islice, groupby |
from itertools import islice, groupby |
||
Line 727: | Line 2,341: | ||
print('TERMS 1..20: ', ', '.join(str(x) for x in islice(cw(), 20))) |
print('TERMS 1..20: ', ', '.join(str(x) for x in islice(cw(), 20))) |
||
x = Fraction(83116, 51639) |
x = Fraction(83116, 51639) |
||
print(f"\n{x} is the {get_term_num(x):_}'th term.")</ |
print(f"\n{x} is the {get_term_num(x):_}'th term.")</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
Line 733: | Line 2,347: | ||
83116/51639 is the 123_456_789'th term.</pre> |
83116/51639 is the 123_456_789'th term.</pre> |
||
=={{header|Quackery}}== |
|||
<code>cf</code> is defined at [[Continued fraction/Arithmetic/Construct from rational number#Quackery]]. |
|||
=={{header|Raku}}== |
|||
{{needs-review|Raku|Now starts a1 = 1}} |
|||
Technically, the Calkin-Wilf sequence should begin with 1, but start with 0 as that is what the task specifies. |
|||
<syntaxhighlight lang="quackery"> [ $ "bigrat.qky" loadfile ] now! |
|||
Conveniently, having the bogus first term shifts the indices up by one, making the ordinal position and index match. |
|||
[ ' [ [ 1 1 ] ] |
|||
swap 1 - times |
|||
[ dup -1 peek do |
|||
2dup proper 2drop |
|||
2 * n->v |
|||
2swap -v 1 n->v v+ v+ |
|||
1/v join nested join ] ] is calkin-wilf ( n --> [ ) |
|||
[ 1 & ] is odd ( n --> b ) |
|||
[ dup size odd not if |
|||
[ -1 split do |
|||
1 - join |
|||
1 join ] ] is oddcf ( [ --> [ ) |
|||
[ 0 swap |
|||
reverse witheach |
|||
[ i odd iff |
|||
<< done |
|||
dup dip << |
|||
bit 1 - | ] ] is rl->n ( [ --> n ) |
|||
[ cf oddcf rl->n ] is cw-term ( n/d --> n ) |
|||
20 calkin-wilf |
|||
witheach |
|||
[ do vulgar$ echo$ sp ] |
|||
cr cr |
|||
83116 51639 cw-term echo</syntaxhighlight> |
|||
{{out}} |
|||
<pre>1/1 1/2 2/1 1/3 3/2 2/3 3/1 1/4 4/3 3/5 5/2 2/5 5/3 3/4 4/1 1/5 5/4 4/7 7/3 3/8 |
|||
123456789</pre> |
|||
=={{header|Raku}}== |
|||
In Raku, arrays are indexed from 0. The Calkin-Wilf sequence does not have a term defined at 0. |
|||
This implementation includes a bogus undefined value at position 0, having the bogus first term shifts the indices up by one, making the ordinal position and index match. Useful due to how reversibility function works. |
|||
Only show the first twenty terms that are ''actually'' in the sequence. |
|||
<lang |
<syntaxhighlight lang="raku" line>my @calkin-wilf = Any, 1, {1 / (.Int × 2 + 1 - $_)} … *; |
||
# Rational to Calkin-Wilf index |
# Rational to Calkin-Wilf index |
||
Line 773: | Line 2,424: | ||
return $num.numerator if $num.denominator == 1; |
return $num.numerator if $num.denominator == 1; |
||
$num.nude.join: '/'; |
$num.nude.join: '/'; |
||
}</ |
}</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
<pre>First twenty terms of the Calkin-Wilf sequence: 1, 1/2, 2, 1/3, 3/2, 2/3, 3, 1/4, 4/3, 3/5, 5/2, 2/5, 5/3, 3/4, 4, 1/5, 5/4, 4/7, 7/3, 3/8 |
<pre>First twenty terms of the Calkin-Wilf sequence: 1, 1/2, 2, 1/3, 3/2, 2/3, 3, 1/4, 4/3, 3/5, 5/2, 2/5, 5/3, 3/4, 4, 1/5, 5/4, 4/7, 7/3, 3/8 |
||
Line 782: | Line 2,433: | ||
83116/51639 is at index: 123456789</pre> |
83116/51639 is at index: 123456789</pre> |
||
=={{header|REXX}}== |
|||
The meat of this REXX program was provided by Paul Kislanko. |
|||
<syntaxhighlight lang="rexx">/*REXX pgm finds the Nth value of the Calkin─Wilf sequence (which will be a fraction),*/ |
|||
/*────────────────────── or finds which sequence number contains a specified fraction). */ |
|||
numeric digits 2000 /*be able to handle ginormic integers. */ |
|||
parse arg LO HI te . /*obtain optional arguments from the CL*/ |
|||
if LO=='' | LO=="," then LO= 1 /*Not specified? Then use the default.*/ |
|||
if HI=='' | HI=="," then HI= 20 /* " " " " " " */ |
|||
if te=='' | te=="," then te= '/' /* " " " " " " */ |
|||
if datatype(LO, 'W') then call CW_terms /*Is LO numeric? Then show some terms.*/ |
|||
if pos('/', te)>0 then call CW_frac te /*Does TE have a / ? Then find term #*/ |
|||
exit 0 |
|||
/*──────────────────────────────────────────────────────────────────────────────────────*/ |
|||
commas: parse arg ?; do jc=length(?)-3 to 1 by -3; ?=insert(',', ?, jc); end; return ? |
|||
th: parse arg th; return word('th st nd rd', 1+(th//10) *(th//100%10\==1) *(th//10<4)) |
|||
/*──────────────────────────────────────────────────────────────────────────────────────*/ |
|||
CW_frac: procedure; parse arg p '/' q .; say |
|||
if q=='' then do; p= 83116; q= 51639; end |
|||
n= rle2dec( frac2cf(p q) ); @CWS= 'the Calkin─Wilf sequence' |
|||
say 'for ' p"/"q', the element number for' @CWS "is: " commas(n)th(n) |
|||
if length(n)<10 then return |
|||
say; say 'The above number has ' commas(length(n)) " decimal digits." |
|||
return |
|||
/*──────────────────────────────────────────────────────────────────────────────────────*/ |
|||
CW_term: procedure; parse arg z; dd= 1; nn= 0 |
|||
do z |
|||
parse value dd dd*(2*(nn%dd)+1)-nn with nn dd |
|||
end /*z*/ |
|||
return nn'/'dd |
|||
/*──────────────────────────────────────────────────────────────────────────────────────*/ |
|||
CW_terms: $=; if LO\==0 then do j=LO to HI; $= $ CW_term(j)',' |
|||
end /*j*/ |
|||
if $=='' then return |
|||
say 'Calkin─Wilf sequence terms for ' commas(LO) " ──► " commas(HI) ' are:' |
|||
say strip( strip($), 'T', ",") |
|||
return |
|||
/*──────────────────────────────────────────────────────────────────────────────────────*/ |
|||
frac2cf: procedure; parse arg p q; if q=='' then return p; cf= p % q; m= q |
|||
p= p - cf*q; n= p; if p==0 then return cf |
|||
do k=1 until n==0; @.k= m % n |
|||
m= m - @.k * n; parse value n m with m n /*swap N M*/ |
|||
end /*k*/ |
|||
/*for inverse Calkin─Wilf, K must be even.*/ |
|||
if k//2 then do; @.k= @.k - 1; k= k + 1; @.k= 1; end |
|||
do k=1 for k; cf= cf @.k; end /*k*/ |
|||
return cf |
|||
/*──────────────────────────────────────────────────────────────────────────────────────*/ |
|||
rle2dec: procedure; parse arg f1 rle; obin= copies(1, f1) |
|||
do until rle==''; parse var rle f0 f1 rle |
|||
obin= copies(1, f1)copies(0, f0)obin |
|||
end /*until*/ |
|||
return x2d( b2x(obin) ) /*RLE2DEC: Run Length Encoding ──► decimal*/</syntaxhighlight> |
|||
{{out|output|text= when using the default inputs:}} |
|||
<pre> |
|||
Calkin─Wilf sequence terms for 1 ──► 20 are: |
|||
1/1, 1/2, 2/1, 1/3, 3/2, 2/3, 3/1, 1/4, 4/3, 3/5, 5/2, 2/5, 5/3, 3/4, 4/1, 1/5, 5/4, 4/7, 7/3, 3/8 |
|||
for 83116/51639, the element number for the Calkin─Wilf sequence is: 123,456,789th |
|||
</pre> |
|||
=={{header|RPL}}== |
|||
{{works with|HP|49g}} |
|||
≪ { } SWAP |
|||
'''WHILE''' DUP '''REPEAT ''' |
|||
ROT OVER IDIV2 |
|||
4 ROLL ROT + SWAP |
|||
'''END''' ROT DROP2 |
|||
≫ '<span style="color:blue">CONTFRAC</span>' STO |
|||
≪ {1} |
|||
'''WHILE''' DUP2 SIZE > '''REPEAT''' |
|||
DUP DUP SIZE GET |
|||
DUP IP R→I 2 * 1 + SWAP - INV EVAL + |
|||
'''END''' NIP |
|||
≫ ≫ '<span style="color:blue">CWILF</span>' STO |
|||
≪ |
|||
<span style="color:blue">CONTFRAC</span> DUP SIZE |
|||
'''IF''' DUP MOD '''THEN''' DROP '''ELSE''' |
|||
DUP2 GET 1 - PUT 1 + '''END''' |
|||
1 → frac pow2 |
|||
≪ 0 |
|||
1 frac SIZE '''FOR''' j |
|||
frac j GET |
|||
'''WHILE''' DUP '''REPEAT''' |
|||
'''IF''' j 2 MOD '''THEN''' SWAP pow2 + SWAP '''END''' |
|||
2 'pow2' STO* |
|||
1 - |
|||
'''END''' DROP |
|||
'''NEXT''' |
|||
≫ ≫ '<span style="color:blue">CWPOS</span>' STO |
|||
20 <span style="color:blue">CWILF</span> |
|||
83116 51639 <span style="color:blue">CWPOS</span> |
|||
{{out}} |
|||
<pre> |
|||
2: {1 '1/2' 2 '1/3' '3/2' '2/3' 3 '1/4' '4/3' '3/5' '5/2' '2/5' '5/3' '3/4' 4 '1/5' '5/4' '4/7' '7/3' '3/8'} |
|||
1: 123456789 |
|||
</pre> |
|||
=={{header|Ruby}}== |
|||
{{trans|Python}} |
|||
<syntaxhighlight lang="ruby">cw = Enumerator.new do |y| |
|||
y << a = 1.to_r |
|||
loop { y << a = 1/(2*a.floor + 1 - a) } |
|||
end |
|||
def term_num(rat) |
|||
num, den, res, pwr, dig = rat.numerator, rat.denominator, 0, 0, 1 |
|||
while den > 0 |
|||
num, (digit, den) = den, num.divmod(den) |
|||
digit.times do |
|||
res |= dig << pwr |
|||
pwr += 1 |
|||
end |
|||
dig ^= 1 |
|||
end |
|||
res |
|||
end |
|||
puts cw.take(20).join(", ") |
|||
puts term_num (83116/51639r) |
|||
</syntaxhighlight> |
|||
<pre>1/1, 1/2, 2/1, 1/3, 3/2, 2/3, 3/1, 1/4, 4/3, 3/5, 5/2, 2/5, 5/3, 3/4, 4/1, 1/5, 5/4, 4/7, 7/3, 3/8 |
|||
123456789 |
|||
</pre> |
|||
=={{header|Rust}}== |
=={{header|Rust}}== |
||
< |
<syntaxhighlight lang="rust">// [dependencies] |
||
// num = "0.3" |
// num = "0.3" |
||
Line 837: | Line 2,613: | ||
let r = Rational::new(83116, 51639); |
let r = Rational::new(83116, 51639); |
||
println!("{} is the {}th term of the sequence.", r, term_number(&r)); |
println!("{} is the {}th term of the sequence.", r, term_number(&r)); |
||
}</ |
}</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
Line 863: | Line 2,639: | ||
20: 3/8 |
20: 3/8 |
||
83116/51639 is the 123456789th term of the sequence. |
83116/51639 is the 123456789th term of the sequence. |
||
</pre> |
|||
=={{header|Scheme}}== |
|||
{{works with|Chez Scheme}} |
|||
'''Continued Fraction support''' |
|||
<syntaxhighlight lang="scheme">; Create a terminating Continued Fraction generator for the given rational number. |
|||
; Returns one term per call; returns #f when no more terms remaining. |
|||
(define make-continued-fraction-gen |
|||
(lambda (rat) |
|||
(let ((num (numerator rat)) (den (denominator rat))) |
|||
(lambda () |
|||
(if (= den 0) |
|||
#f |
|||
(let ((ret (quotient num den)) |
|||
(rem (modulo num den))) |
|||
(set! num den) |
|||
(set! den rem) |
|||
ret)))))) |
|||
; Return the continued fraction representation of a rational number as a list of terms. |
|||
(define rat->cf-list |
|||
(lambda (rat) |
|||
(let ((cf (make-continued-fraction-gen rat)) |
|||
(lst '())) |
|||
(let loop ((term (cf))) |
|||
(when term |
|||
(set! lst (append lst (list term))) |
|||
(loop (cf)))) |
|||
lst))) |
|||
; Enforce the length of the given continued fraction list to be odd. |
|||
; Changes the list in situ (if needed), and returns its possibly changed value. |
|||
(define continued-fraction-list-enforce-odd-length! |
|||
(lambda (cf) |
|||
(when (even? (length cf)) |
|||
(let ((cf-last-cons (list-tail cf (1- (length cf))))) |
|||
(set-car! cf-last-cons (1- (car cf-last-cons))) |
|||
(set-cdr! cf-last-cons (cons 1 '())))) |
|||
cf))</syntaxhighlight> |
|||
'''Calkin-Wilf sequence''' |
|||
<syntaxhighlight lang="scheme">; Create a Calkin-Wilf sequence generator. |
|||
(define make-calkin-wilf-gen |
|||
(lambda () |
|||
(let ((an 1)) |
|||
(lambda () |
|||
(let ((ret an)) |
|||
(set! an (/ 1 (+ (* 2 (floor an)) 1 (- an)))) |
|||
ret))))) |
|||
; Return the position in the Calkin-Wilf sequence of the given rational number. |
|||
(define calkin-wilf-position |
|||
(lambda (rat) |
|||
; Run-length encodes binary value. Assumes first run is 1's. Args: initial value, |
|||
; starting place value (a power of 2), and list of run lengths (list must be odd length). |
|||
(define encode-list-of-runs |
|||
(lambda (value placeval lstruns) |
|||
; Encode a single run in a binary value. Args: initial value, bit value (0 or 1), |
|||
; starting place value (a power of 2), number of places (bits) to encode. |
|||
; Returns multiple values: the encoded value, and the new place value. |
|||
(define encode-run |
|||
(lambda (value bitval placeval places) |
|||
(if (= places 1) |
|||
(values (+ value (* bitval placeval)) (* 2 placeval)) |
|||
(encode-run (+ value (* bitval placeval)) bitval (* 2 placeval) (1- places))))) |
|||
; Loop through the list of runs two at a time. If list of length 1, do a final |
|||
; '1'-bit encode and return the value. Otherwise, do a '1'-bit then '0'-bit encode, |
|||
; and recurse to do the next two runs. |
|||
(let-values (((value-1 placeval-1) (encode-run value 1 placeval (car lstruns)))) |
|||
(if (= 1 (length lstruns)) |
|||
value-1 |
|||
(let-values (((value-2 placeval-2) (encode-run value-1 0 placeval-1 (cadr lstruns)))) |
|||
(encode-list-of-runs value-2 placeval-2 (cddr lstruns))))))) |
|||
; Return the run-length binary encoding from the odd-length Calkin-Wilf sequence of the |
|||
; given rational number. This is equal to the number's position in the sequence. |
|||
(encode-list-of-runs 0 1 (continued-fraction-list-enforce-odd-length! (rat->cf-list rat)))))</syntaxhighlight> |
|||
'''The Task''' |
|||
<syntaxhighlight lang="scheme">(let ((count 20) |
|||
(cw (make-calkin-wilf-gen))) |
|||
(printf "~%First ~a terms of the Calkin-Wilf sequence:~%" count) |
|||
(do ((num 1 (1+ num))) |
|||
((> num count)) |
|||
(printf "~2d : ~a~%" num (cw)))) |
|||
(printf "~%Positions in Calkin-Wilf sequence of given numbers:~%") |
|||
(let ((num 9/4)) |
|||
(printf "~a @ ~a~%" num (calkin-wilf-position num))) |
|||
(let ((num 83116/51639)) |
|||
(printf "~a @ ~a~%" num (calkin-wilf-position num)))</syntaxhighlight> |
|||
{{out}} |
|||
<pre> |
|||
First 20 terms of the Calkin-Wilf sequence: |
|||
1 : 1 |
|||
2 : 1/2 |
|||
3 : 2 |
|||
4 : 1/3 |
|||
5 : 3/2 |
|||
6 : 2/3 |
|||
7 : 3 |
|||
8 : 1/4 |
|||
9 : 4/3 |
|||
10 : 3/5 |
|||
11 : 5/2 |
|||
12 : 2/5 |
|||
13 : 5/3 |
|||
14 : 3/4 |
|||
15 : 4 |
|||
16 : 1/5 |
|||
17 : 5/4 |
|||
18 : 4/7 |
|||
19 : 7/3 |
|||
20 : 3/8 |
|||
Positions in Calkin-Wilf sequence of given numbers: |
|||
9/4 @ 35 |
|||
83116/51639 @ 123456789 |
|||
</pre> |
|||
=={{header|Sidef}}== |
|||
<syntaxhighlight lang="ruby">func calkin_wilf(n) is cached { |
|||
return 1 if (n == 1) |
|||
1/(2*floor(__FUNC__(n-1)) + 1 - __FUNC__(n-1)) |
|||
} |
|||
func r2cw(r) { |
|||
var cfrac = r.as_cfrac |
|||
cfrac.len.is_odd || return nil |
|||
Num(cfrac.flip.map_kv {|k,v| (k.is_odd ? '0' : '1') * v }.join, 2) |
|||
} |
|||
with (20) {|n| |
|||
say "First #{n} terms of the Calkin-Wilf sequence:" |
|||
say calkin_wilf.map(1..n) |
|||
} |
|||
with (83116/51639) {|r| |
|||
say ("\n#{r.as_rat} is at index: ", r2cw(r)) |
|||
}</syntaxhighlight> |
|||
{{out}} |
|||
<pre> |
|||
First 20 terms of the Calkin-Wilf sequence: |
|||
[1, 1/2, 2, 1/3, 3/2, 2/3, 3, 1/4, 4/3, 3/5, 5/2, 2/5, 5/3, 3/4, 4, 1/5, 5/4, 4/7, 7/3, 3/8] |
|||
83116/51639 is at index: 123456789 |
|||
</pre> |
|||
=={{header|V (Vlang)}}== |
|||
{{trans|Go}}s. |
|||
<syntaxhighlight lang="v (vlang)">import math.fractions |
|||
import math |
|||
import strconv |
|||
fn calkin_wilf(n int) []fractions.Fraction { |
|||
mut cw := []fractions.Fraction{len: n+1} |
|||
cw[0] = fractions.fraction(1, 1) |
|||
one := fractions.fraction(1, 1) |
|||
two := fractions.fraction(2, 1) |
|||
for i in 1..n { |
|||
mut t := cw[i-1] |
|||
mut f := t.f64() |
|||
f = math.floor(f) |
|||
t = fractions.approximate(f) |
|||
t*=two |
|||
t-= cw[i-1] |
|||
t+=one |
|||
t=t.reciprocal() |
|||
cw[i] = t |
|||
} |
|||
return cw |
|||
} |
|||
fn to_continued(r fractions.Fraction) []int { |
|||
idx := r.str().index('/') or {0} |
|||
mut a := r.str()[..idx].i64() |
|||
mut b := r.str()[idx+1..].i64() |
|||
mut res := []int{} |
|||
for { |
|||
res << int(a/b) |
|||
t := a % b |
|||
a, b = b, t |
|||
if a == 1 { |
|||
break |
|||
} |
|||
} |
|||
le := res.len |
|||
if le%2 == 0 { // ensure always odd |
|||
res[le-1]-- |
|||
res << 1 |
|||
} |
|||
return res |
|||
} |
|||
fn get_term_number(cf []int) ?int { |
|||
mut b := "" |
|||
mut d := "1" |
|||
for n in cf { |
|||
b = d.repeat(n)+b |
|||
if d == "1" { |
|||
d = "0" |
|||
} else { |
|||
d = "1" |
|||
} |
|||
} |
|||
i := strconv.parse_int(b, 2, 64)? |
|||
return int(i) |
|||
} |
|||
fn commatize(n int) string { |
|||
mut s := "$n" |
|||
if n < 0 { |
|||
s = s[1..] |
|||
} |
|||
le := s.len |
|||
for i := le - 3; i >= 1; i -= 3 { |
|||
s = s[0..i] + "," + s[i..] |
|||
} |
|||
if n >= 0 { |
|||
return s |
|||
} |
|||
return "-" + s |
|||
} |
|||
fn main() { |
|||
cw := calkin_wilf(20) |
|||
println("The first 20 terms of the Calkin-Wilf sequnence are:") |
|||
for i := 1; i <= 20; i++ { |
|||
println("${i:2}: ${cw[i-1]}") |
|||
} |
|||
println('') |
|||
r := fractions.fraction(83116, 51639) |
|||
cf := to_continued(r) |
|||
tn := get_term_number(cf) or {0} |
|||
println("$r is the ${commatize(tn)}th term of the sequence.") |
|||
}</syntaxhighlight> |
|||
{{out}} |
|||
<pre> |
|||
The first 20 terms of the Calkin-Wilf sequnence are: |
|||
1: 1/1 |
|||
2: 1/2 |
|||
3: 2/1 |
|||
4: 1/3 |
|||
5: 3/2 |
|||
6: 2/3 |
|||
7: 3/1 |
|||
8: 1/4 |
|||
9: 4/3 |
|||
10: 3/5 |
|||
11: 5/2 |
|||
12: 2/5 |
|||
13: 5/3 |
|||
14: 3/4 |
|||
15: 4/1 |
|||
16: 1/5 |
|||
17: 5/4 |
|||
18: 4/7 |
|||
19: 7/3 |
|||
20: 3/8 |
|||
83116/51639 is the 123,456,789th term of the sequence. |
|||
</pre> |
</pre> |
||
Line 868: | Line 2,902: | ||
{{libheader|Wren-rat}} |
{{libheader|Wren-rat}} |
||
{{libheader|Wren-fmt}} |
{{libheader|Wren-fmt}} |
||
< |
<syntaxhighlight lang="wren">import "./rat" for Rat |
||
import "/fmt" for Fmt, Conv |
import "./fmt" for Fmt, Conv |
||
var calkinWilf = Fn.new { |n| |
var calkinWilf = Fn.new { |n| |
||
Line 917: | Line 2,951: | ||
var cf = toContinued.call(r) |
var cf = toContinued.call(r) |
||
var tn = getTermNumber.call(cf) |
var tn = getTermNumber.call(cf) |
||
Fmt.print("$s is the $,r term of the sequence.", r, tn)</ |
Fmt.print("$s is the $,r term of the sequence.", r, tn)</syntaxhighlight> |
||
{{out}} |
{{out}} |
Latest revision as of 10:48, 3 April 2024
You are encouraged to solve this task according to the task description, using any language you may know.
The Calkin-Wilf sequence contains every nonnegative rational number exactly once.
It can be calculated recursively as follows:
a1 = 1 an+1 = 1/(2⌊an⌋+1-an) for n > 1
- Task part 1
- Show on this page terms 1 through 20 of the Calkin-Wilf sequence.
To avoid floating point error, you may want to use a rational number data type.
It is also possible, given a non-negative rational number, to determine where it appears in the sequence without calculating the sequence. The procedure is to get the continued fraction representation of the rational and use it as the run-length encoding of the binary representation of the term number, beginning from the end of the continued fraction.
It only works if the number of terms in the continued fraction is odd- use either of the two equivalent representations to achieve this:
[a0; a1, a2, ..., an] = [a0; a1, a2 ,..., an-1, 1]
- Example
The fraction 9/4 has odd continued fraction representation 2; 3, 1, giving a binary representation of 100011,
which means 9/4 appears as the 35th term of the sequence.
- Task part 2
- Find the position of the number 83116/51639 in the Calkin-Wilf sequence.
- Related tasks
- See also
- Wikipedia entry: Calkin-Wilf tree
- Continued fraction
- Continued fraction/Arithmetic/Construct from rational number
11l
T CalkinWilf
n = 1
d = 1
F ()()
V r = (.n, .d)
.n = 2 * (.n I/ .d) * .d + .d - .n
swap(&.n, &.d)
R r
print(‘The first 20 terms of the Calkwin-Wilf sequence are:’)
V cw = CalkinWilf()
[String] seq
L 20
V (n, d) = cw()
seq.append(I d == 1 {String(n)} E n‘/’d)
print(seq.join(‘, ’))
cw = CalkinWilf()
V index = 1
L cw() != (83116, 51639)
index++
print("\nThe element 83116/51639 is at position "index‘ in the sequence.’)
- Output:
The first 20 terms of the Calkwin-Wilf sequence are: 1, 1/2, 2, 1/3, 3/2, 2/3, 3, 1/4, 4/3, 3/5, 5/2, 2/5, 5/3, 3/4, 4, 1/5, 5/4, 4/7, 7/3, 3/8 The element 83116/51639 is at position 123456789 in the sequence.
ALGOL 68
Uses code from the Arithmetic/Rational and Continued fraction/Arithmetic/Construct from rational number tasks.
BEGIN
# Show elements 1-20 of the Calkin-Wilf sequence as rational numbers #
# also show the position of a specific element in the seuence #
# Uses code from the Arithmetic/Rational #
# & Continued fraction/Arithmetic/Construct from rational number tasks #
# Code from the Arithmetic/Rational task #
# ============================================================== #
MODE FRAC = STRUCT( INT num #erator#, den #ominator#);
PROC gcd = (INT a, b) INT: # greatest common divisor #
(a = 0 | b |: b = 0 | a |: ABS a > ABS b | gcd(b, a MOD b) | gcd(a, b MOD a));
PROC lcm = (INT a, b)INT: # least common multiple #
a OVER gcd(a, b) * b;
PRIO // = 9; # higher then the ** operator #
OP // = (INT num, den)FRAC: ( # initialise and normalise #
INT common = gcd(num, den);
IF den < 0 THEN
( -num OVER common, -den OVER common)
ELSE
( num OVER common, den OVER common)
FI
);
OP + = (FRAC a, b)FRAC: (
INT common = lcm(den OF a, den OF b);
FRAC result := ( common OVER den OF a * num OF a + common OVER den OF b * num OF b, common );
num OF result//den OF result
);
OP - = (FRAC a, b)FRAC: a + -b,
* = (FRAC a, b)FRAC: (
INT num = num OF a * num OF b,
den = den OF a * den OF b;
INT common = gcd(num, den);
(num OVER common) // (den OVER common)
);
OP - = (FRAC frac)FRAC: (-num OF frac, den OF frac);
# ============================================================== #
# end code from the Arithmetic/Rational task #
# code from the Continued fraction/Arithmetic/Construct from rational number task #
# ================================================================================#
# returns the quotient of numerator over denominator and sets #
# numerator and denominator to the next values for #
# the continued fraction #
PROC r2cf = ( REF INT numerator, REF INT denominator )INT:
IF denominator = 0
THEN 0
ELSE INT quotient := numerator OVER denominator;
INT prev numerator = numerator;
numerator := denominator;
denominator := prev numerator MOD denominator;
quotient
FI # r2cf # ;
# ====================================================================================#
# end code from the Continued fraction/Arithmetic/Construct from rational number task #
# Additional FRACrelated operators #
OP * = ( INT a, FRAC b )FRAC: ( num OF b * a ) // den OF b;
OP / = ( FRAC a, b )FRAC: ( num OF a * den OF b ) // ( num OF b * den OF a );
OP FLOOR = ( FRAC a )INT: num OF a OVER den OF a;
OP + = ( INT a, FRAC b )FRAC: ( a // 1 ) + b;
FRAC one = 1 // 1;
# returns the first n elements of the Calkin-Wilf sequence #
PROC calkin wilf = ( INT n )[]FRAC:
BEGIN
[ 1 : n ]FRAC q;
IF n > 0 THEN
q[ 1 ] := 1 // 1;
FOR i FROM 2 TO UPB q DO
q[ i ] := one / ( ( 2 * FLOOR q[ i - 1 ] ) + one - q[ i - 1 ] )
OD
FI;
q
END # calkin wilf # ;
# returns the position of a FRAC in the Calkin-Wilf sequence by computing its #
# continued fraction representation and converting that to a bit string #
# - the position must fit in a 2-bit number #
PROC position in calkin wilf sequence = ( FRAC f )INT:
IF INT result := 0;
[ 1 : 32 ]INT cf; FOR i FROM LWB cf TO UPB cf DO cf[ i ] := 0 OD;
INT num := num OF f;
INT den := den OF f;
INT cf length := 0;
FOR i FROM LWB cf WHILE den /= 0 DO
cf[ cf length := i ] := r2cf( num, den )
OD;
NOT ODD cf length
THEN # the continued fraction does not have an odd length #
-1
ELSE # the continued fraction has an odd length so we can compute the seuence length #
# build the number by alternating d 1s and 0s where d is the digits of the #
# continued fraction, starting at the least significant #
INT digit := 1;
FOR d pos FROM cf length BY -1 TO 1 DO
FOR i TO cf[ d pos ] DO
result *:= 2 +:= digit
OD;
digit := IF digit = 0 THEN 1 ELSE 0 FI
OD;
result
FI # position in calkin wilf sequence # ;
BEGIN # task #
# get and show the first 20 Calkin-Wilf sequence numbers #
[]FRAC cw = calkin wilf( 20 );
print( ( "The first 20 elements of the Calkin-Wilf sequence are:", newline, " " ) );
FOR n FROM LWB cw TO UPB cw DO
FRAC sn = cw[ n ];
print( ( " ", whole( num OF sn, 0 ), "/", whole( den OF sn, 0 ) ) )
OD;
print( ( newline ) );
# show the position of a specific element in the sequence #
print( ( "Position of 83116/51639 in the sequence: "
, whole( position in calkin wilf sequence( 83116//51639 ), 0 )
)
)
END
END
- Output:
The first 20 elements of the Calkin-Wilf sequence are: 1/1 1/2 2/1 1/3 3/2 2/3 3/1 1/4 4/3 3/5 5/2 2/5 5/3 3/4 4/1 1/5 5/4 4/7 7/3 3/8 Position of 83116/51639 in the sequence: 123456789
AppleScript
-- Return the first n terms of the sequence. Tree generation. Faster for this purpose.
on CalkinWilfSequence(n)
script o
property sequence : {{1, 1}} -- Initialised with the first term ({numerator, denominator}).
end script
-- Work through the growing sequence list, adding the two children of each term to the end and
-- converting each term to text representing the vulgar fraction. Stop adding children halfway through.
set halfway to n div 2
repeat with position from 1 to n
set {numerator, denominator} to item position of o's sequence
if (position ≤ halfway) then
tell numerator + denominator
set end of o's sequence to {numerator, it}
if ((position < halfway) or (position * 2 < n)) then set end of o's sequence to {it, denominator}
end tell
end if
set item position of o's sequence to (numerator as text) & "/" & denominator
end repeat
return o's sequence
end CalkinWilfSequence
-- Alternatively, return terms pos1 to pos2. Binary run-length encoding. Doesn't need to work from the beginning of the sequence.
on CalkinWilfSequence2(pos1, pos2)
script o
property sequence : {}
end script
repeat with position from pos1 to pos2
-- Build a continued fraction list from the binary run-length encoding of this position index.
-- There's no need to put the last value into the list as it's used immediately.
set continuedFraction to {}
set bitValue to 1
set runLength to 0
repeat until (position = 0)
if (position mod 2 = bitValue) then
set runLength to runLength + 1
else
set end of continuedFraction to runLength
set bitValue to (bitValue + 1) mod 2
set runLength to 1
end if
set position to position div 2
end repeat
-- Work out the numerator and denominator from the continued fraction and derive text representing the vulgar fraction.
set numerator to runLength
set denominator to 1
repeat with i from (count continuedFraction) to 1 by -1
tell numerator
set numerator to numerator * (item i of continuedFraction) + denominator
set denominator to it
end tell
end repeat
set end of o's sequence to (numerator as text) & "/" & denominator
end repeat
return o's sequence
end CalkinWilfSequence2
-- Return the sequence position of the term with the given numerator and denominator.
on CalkinWilfSequencePosition(numerator, denominator)
-- Build a continued fraction list from the input.
set continuedFraction to {}
repeat until (denominator is 0)
set end of continuedFraction to numerator div denominator
set {numerator, denominator} to {denominator, numerator mod denominator}
end repeat
-- If it has an even number of entries, convert to the equivalent odd number.
if ((count continuedFraction) mod 2 is 0) then
set last item of continuedFraction to (last item of continuedFraction) - 1
set end of continuedFraction to 1
end if
-- "Binary run-length decode" the entries to get the position index.
set position to 0
set bitValue to 1
repeat with i from (count continuedFraction) to 1 by -1
repeat (item i of continuedFraction) times
set position to position * 2 + bitValue
end repeat
set bitValue to (bitValue + 1) mod 2
end repeat
return position
end CalkinWilfSequencePosition
-- Task code:
local sequenceResult1, sequenceResult2, positionResult, output, astid
set sequenceResult1 to CalkinWilfSequence(20)
set sequenceResult2 to CalkinWilfSequence2(1, 20)
set positionResult to CalkinWilfSequencePosition(83116, 51639)
set astid to AppleScript's text item delimiters
set AppleScript's text item delimiters to ", "
set output to "First twenty terms of sequence using tree generation:" & (linefeed & sequenceResult1)
set output to output & (linefeed & "Ditto using binary run-length encoding:") & (linefeed & sequenceResult1)
set AppleScript's text item delimiters to astid
set output to output & (linefeed & "83116/51639 is term number " & positionResult)
return output
- Output:
"First twenty terms of sequence using tree generation:
1/1, 1/2, 2/1, 1/3, 3/2, 2/3, 3/1, 1/4, 4/3, 3/5, 5/2, 2/5, 5/3, 3/4, 4/1, 1/5, 5/4, 4/7, 7/3, 3/8
Ditto using binary run-length encoding:
1/1, 1/2, 2/1, 1/3, 3/2, 2/3, 3/1, 1/4, 4/3, 3/5, 5/2, 2/5, 5/3, 3/4, 4/1, 1/5, 5/4, 4/7, 7/3, 3/8
83116/51639 is term number 123456789"
Arturo
n: new 1
d: new 1
calkinWilf: function [] .export:[n,d] [
n: (d - n) + 2 * (n/d) * d
tmp: d
d: n
n: tmp
return @[n d]
]
first20: [[1 1]] ++ map 1..19 => calkinWilf
print "The first 20 terms of the Calkwin-Wilf sequence are:"
print map first20 'f -> ~"|f\0|/|f\1|"
n: new 1
d: new 1
indx: new 1
target: [83116, 51639]
while ø [
inc 'indx
if target = calkinWilf -> break
]
print ""
print ["The element" ~"|target\0|/|target\1|" "is at position" indx "in the sequence."]
- Output:
The first 20 terms of the Calkwin-Wilf sequence are: 1/1 1/2 2/1 1/3 3/2 2/3 3/1 1/4 4/3 3/5 5/2 2/5 5/3 3/4 4/1 1/5 5/4 4/7 7/3 3/8 The element 83116/51639 is at position 123456789 in the sequence.
BQN
BQN does not have rational number arithmetic yet, so it is manually implemented.
Part 2 runs in ~150 secs on CBQN.
GCD
and _while_
are idioms from BQNcrate.
GCD ← {m 𝕊⍟(0<m←𝕨|𝕩) 𝕨}
_while_ ← {𝔽⍟𝔾∘𝔽_𝕣_𝔾∘𝔽⍟𝔾𝕩}
Sim ← { # Simplify a fraction
x𝕊1: 𝕨‿1;
0𝕊y: 0‿𝕩;
⌊𝕨‿𝕩 ÷ 𝕨 GCD 𝕩
}
Add ← { # Add two fractions
0‿b 𝕊 𝕩: 𝕩;
𝕨 𝕊 0‿y: 𝕨;
a‿b 𝕊 x‿y:
((a×y)+x×b) Sim b×y
}
Next ← {n‿d: ⌽(2×⌊÷´n‿d)‿1 Add (d-n)‿d} # Next term
Cal ← {Next⍟𝕩 1‿1}
•Show Cal 1+↕20
•Show {
cnt‿fr:
⟨cnt+1,Next fr⟩
} _while_ {
cnt‿fr:
fr ≢ 83116‿51639
} ⟨1,1‿1⟩
⟨ ⟨ 1 2 ⟩ ⟨ 2 1 ⟩ ⟨ 1 3 ⟩ ⟨ 3 2 ⟩ ⟨ 2 3 ⟩ ⟨ 3 1 ⟩ ⟨ 1 4 ⟩ ⟨ 4 3 ⟩ ⟨ 3 5 ⟩ ⟨ 5 2 ⟩ ⟨ 2 5 ⟩ ⟨ 5 3 ⟩ ⟨ 3 4 ⟩ ⟨ 4 1 ⟩ ⟨ 1 5 ⟩ ⟨ 5 4 ⟩ ⟨ 4 7 ⟩ ⟨ 7 3 ⟩ ⟨ 3 8 ⟩ ⟨ 8 5 ⟩ ⟩
⟨ 123456789 ⟨ 83116 51639 ⟩ ⟩
You can try Part 1 here. Second part can and will hang your browser, so it is best to try locally on CBQN.
Bracmat
( 1:?a
& 0:?i
& whl
' ( 1+!i:<20:?i
& (2*div$(!a,1)+1+-1*!a)^-1:?a
& out$!a
)
& ( r2cf
= floor
. div$(!arg,1):?floor
& ( !floor:!arg
| !floor r2cf$((!arg+-1*!floor)^-1)
)
)
& ( get-term-num
= ans dig pwr
. (0,1,1):(?ans,?dig,?pwr)
& r2cf$!arg:?n
& map
$ ( (
=
. whl
' ( !arg+-1:~<0:?arg
& !dig*!pwr+!ans:?ans
& 2*!pwr:?pwr
)
& 1+-1*!dig:?dig
)
. !n
)
& !ans
)
& out$(get-term-num$83116/51639)
);
- Output:
1/2 2 1/3 3/2 2/3 3 1/4 4/3 3/5 5/2 2/5 5/3 3/4 4 1/5 5/4 4/7 7/3 3/8 123456789
C++
#include <iostream>
#include <vector>
#include <boost/rational.hpp>
using rational = boost::rational<unsigned long>;
unsigned long floor(const rational& r) {
return r.numerator()/r.denominator();
}
rational calkin_wilf_next(const rational& term) {
return 1UL/(2UL * floor(term) + 1UL - term);
}
std::vector<unsigned long> continued_fraction(const rational& r) {
unsigned long a = r.numerator();
unsigned long b = r.denominator();
std::vector<unsigned long> result;
do {
result.push_back(a/b);
unsigned long c = a;
a = b;
b = c % b;
} while (a != 1);
if (result.size() > 0 && result.size() % 2 == 0) {
--result.back();
result.push_back(1);
}
return result;
}
unsigned long term_number(const rational& r) {
unsigned long result = 0;
unsigned long d = 1;
unsigned long p = 0;
for (unsigned long n : continued_fraction(r)) {
for (unsigned long i = 0; i < n; ++i, ++p)
result |= (d << p);
d = !d;
}
return result;
}
int main() {
rational term = 1;
std::cout << "First 20 terms of the Calkin-Wilf sequence are:\n";
for (int i = 1; i <= 20; ++i) {
std::cout << std::setw(2) << i << ": " << term << '\n';
term = calkin_wilf_next(term);
}
rational r(83116, 51639);
std::cout << r << " is the " << term_number(r) << "th term of the sequence.\n";
}
- Output:
First 20 terms of the Calkin-Wilf sequence are: 1: 1/1 2: 1/2 3: 2/1 4: 1/3 5: 3/2 6: 2/3 7: 3/1 8: 1/4 9: 4/3 10: 3/5 11: 5/2 12: 2/5 13: 5/3 14: 3/4 15: 4/1 16: 1/5 17: 5/4 18: 4/7 19: 7/3 20: 3/8 83116/51639 is the 123456789th term of the sequence.
EasyLang
subr first
n = 1 ; d = 1
.
proc next . .
n = 2 * (n div d) * d + d - n
swap n d
.
print "The first 20 terms of the Calkwin-Wilf sequence are:"
first
for i to 20
write n & "/" & d & " "
next
.
print ""
#
first
i = 1
while n <> 83116 or d <> 51639
next
i += 1
.
print "83116/51639 is at position " & i
- Output:
The first 20 terms of the Calkwin-Wilf sequence are: 1/1 1/2 2/1 1/3 3/2 2/3 3/1 1/4 4/3 3/5 5/2 2/5 5/3 3/4 4/1 1/5 5/4 4/7 7/3 3/8 83116/51639 is at position 123456789
EDSAC order code
Find first n terms
[For Rosetta Code. EDSAC program, Initial Orders 2.
Prints the first 20 terms of the Calkin-Wilf sequence.
Uses term{n} to calculate term{n + 1}.]
[Print subroutine for non-negative 17-bit integers.
Parameters: 0F = integer to be printed (not preserved)
1F = character for leading zero (preserved)
Workspace: 4F, 5F. Even address; 40 locations]
T 56 K [define load address]
GKA3FT34@A1FT35@S37@T36@AFT5FT4FH38#@V4DH30@
R32FR16FYFE23@O35@A2FT36@T5FV4DYFL8FT4DA5FL1024F
UFA36@G16@OFTFT35@A36@G17@ZFPFPFP4FT1714FZ219D
[Main routine]
T 100 K [define load address]
G K [set up relative addressing via @ (theta)]
[Constants]
[0] P 10 F [maximum index = 20, edit ad lib.]
[1] P D [constant 1]
[Teleprinter characters]
[2] # F [set figures mode]
[3] C F [colon (in figures mode)]
[4] X F [slash (in figures mode)]
[5] ! F [space]
[6] @ F [carriage return]
[7] & F [line feed]
[8] K 4096 F [null]
[Variables]
[9] P F [index]
[10] P F [a (where term = a/b)]
[11] P F [b]
[Enter with acc = 0]
[12] O 2 @ [set teleprinter to figures]
A 1 @ [acc := 1]
U 9 @ [index := 1]
U 10 @ [a := 1]
T 11 @ [b := 1 (and clear acc)]
E 34 @ [jump to print first term]
[Loop back here if not yet printed enough terms]
[18] A @ [restore index after test]
A 1 @ [add 1]
T 9 @ [update index]
[Calculate next term. New b := a + b - 2(a mod b).
Code below calculates c := (a mod b) - b, then new b := a - b - 2*c]
A 10 @ [acc := a]
[22] S 11 @ [subtract b]
E 22 @ [if acc >= 0, subtract again]
T F [result c < 0, store in 0F]
A 10 @ [acc := a]
S 11 @ [subtract b]
S F [subtract c]
S F [subtract c]
T F [new b = a - b - 2*c; store in 0F]
A 11 @ [acc := old b]
T 10 @ [copy to a]
A F [acc := new b]
T 11 @ [copy to b]
[Print index and a/b. Assume acc = 0 here.]
[34] A 5 @ [space to replace leading 0's]
T 1 F [pass to print subroutine]
A 9 @ [acc := index]
T F [pass to print subroutine]
[38] A 38 @ [for return from subroutine]
G 56 F [call subroutine, clears acc]
O 3 @ [print colon]
O 5 @ [print space]
A 8 @ [null to replace leading 0's]
T 1 F [pass to print subroutine]
A10@ TF A46@ G56F O4@ [print a followed by slash]
A11@ TF A51@ G56F O6@ O7@ [print b followed by CR LF]
[Test whether enough terms have been printed]
A 9 @ [acc := index]
S @ [subtract maximum index]
G 18 @ [loop back if acc < 0]
[Exit]
O 8 @ [print null to flush teleprinter buffer]
Z F [stop]
E 12 Z [relative address of entry point]
P F [enter with acc = 0]
[end]
- Output:
1: 1/1 2: 1/2 3: 2/1 4: 1/3 5: 3/2 6: 2/3 7: 3/1 8: 1/4 9: 4/3 10: 3/5 11: 5/2 12: 2/5 13: 5/3 14: 3/4 15: 4/1 16: 1/5 17: 5/4 18: 4/7 19: 7/3 20: 3/8
Find index of a given term
[For Rosetta Code. EDSAC program, Initial Orders 2.]
[Finds the index of a given rational in the Calkin-Wilf series.]
[Library subroutine R2: input of positive integers.
Runs during input of the program, and is then overwritten.
Allows integers to be written in decimal, rather than as "pseudo-orders".
See Wilkes, Williams & Gill, 1951 edn, pp. 96-97, 148.]
T 54 K [to access integers via C parameter]
P 110 F [where to load integers]
GKT20FVDL8FA40DUDTFI40FA40FS39FG@S2FG23FA5@T5@E4@E13Z
T #C [tell R2 where to load integers]
[F after each integer except the last, and # after the last.]
83116F51639#
[Modified library subroutine P7.
Prints signed integer up to 10 digits, left-justified.
Input: Number to be printed is at 0D.
54 locations. Load at even address. Workspace 4D.]
T 56 K
GKA3FT42@A49@T31@ADE10@T31@A48@T31@SDTDH44#@NDYFLDT4DS43@
TFH17@S17@A43@G23@UFS43@T1FV4DAFG50@SFLDUFXFOFFFSFL4FT4DA49@
T31@A1FA43@G20@XFP1024FP610D@524D!FO46@O26@XFSFL8FT4DE39@
[Main routine.]
T 120 K [define load address (must be even)]
G K [set up relative addressing via @ (theta)]
[Put 35-bit values first, to ensure each is at an even address]
[Variables]
[0] P F P F [a]
[2] P F P F [b]
[4] P F P F [power of 2]
[6] P F P F [calculated index]
[Constants]
T8#Z PF T8Z [clears sandwich digit between 8 and 9]
[8] P D P F [35-bit constant 1]
[Teleprinter characters]
[10] # F [set figures mode]
[11] X F [slash (in figures mode)]
[12] K 2048 F [set letters mode]
[13] I F [letter I]
[14] R F [letter R]
[15] ! F [space]
[16] @ F [carriage return]
[17] & F [line feed]
[18] K 4096 F [null char]
[Enter with acc = 0]
[19] A #C [acc := initial a]
T #@ [copy to variable]
A 2#C [acc := initial b]
T 2#@ [copy to variable]
[23] A 8#@ [acc := 1]
[24] T 4#@ [initialize power of 2]
T 6#@ [initialize index to 0]
[Loop]
[26] A #@ [acc := a]
[27] S 2#@ [subtract b]
[28] E 33 @ [jump if a >= b]
[Here if a < b]
T D [store a - b in 0D]
S D [negate]
T 2#@ [b := b - a]
E 40 @ [join common code]
[Here if a >= b]
[33] S 8#@ [acc = a - b; test for a = b]
G 45 @ [jump out of loop if so]
A 8#@ [restore a - b]
T #@ [a := a - b]
A 6#@ [acc := index]
A 4#@ [inc index by power of 2]
T 6#@
[Code common to both cases]
[40] A 4#@ [acc := power of 2]
L D [shift left]
G 76 @
T 4#@ [update power of 2]
E 26 @ [loop back]
[Exit from loop.]
[45] T D [dump acc to clear it]
A 6#@ [acc := index]
A 4#@ [add power of 2 ]
T 6#@ [store final value of index]
[Finished calcualting index, now do printing]
O 10 @ [set teleprinter to figures]
A #C [acc := initial a]
T D [to 0D for printing]
[52] A 52 @ [for return from subroutine]
G 56 F [call print subroutine, clears acc]
O 11 @ [print slash]
A 2#C [print initial b similarly]
T D
[57] A 57 @
G 56 F
O 12 @ [set teleprinter to letters and print ' IS AT ']
O15@ O13@ O27@ O15@ O23@ O24@ O15@
O 10 @ [set teleprinter to figures]
A 6#@ [acc := calculated index]
T D [send to print subroutine]
[70] A 70 @
G 56 F
[72] O16@ O17@ [print CR, LF]
O 18 @ [print null to flush teleprinter buffer]
Z F [stop]
[Here if power of 2 goes negative (accumulator overflow)]
[76] O 12 @ [set teleprinter to letters]
O28@ O14@ O14@ O76@ O14@ [print'ERROR']
G 72 @ [jump to common exit]
E 19 Z [relative address of entry point]
P F [enter with acc = 0]
- Output:
83116/51639 IS AT 123456789
F#
The Function
// Calkin Wilf Sequence. Nigel Galloway: January 9th., 2021
let cW=Seq.unfold(fun(n)->Some(n,seq{for n,g in n do yield (n,n+g); yield (n+g,g)}))(seq[(1,1)])|>Seq.concat
The Tasks
- first 20
cW |> Seq.take 20 |> Seq.iter(fun(n,g)->printf "%d/%d " n g);printfn ""
- Output:
1/1 1/2 2/1 1/3 3/2 2/3 3/1 1/4 4/3 3/5 5/2 2/5 5/3 3/4 4/1 1/5 5/4 4/7 7/3 3/8
- Indexof 83116/51639
printfn "%d" (1+Seq.findIndex(fun n->n=(83116,51639)) cW)
- Output:
123456789
Factor
USING: formatting io kernel lists lists.lazy math
math.continued-fractions math.functions math.parser prettyprint
sequences strings vectors ;
: next-cw ( x -- y ) [ floor dup + ] [ 1 swap - + recip ] bi ;
: calkin-wilf ( -- list ) 1 [ next-cw ] lfrom-by ;
: >continued-fraction ( x -- seq )
1vector [ dup last integer? ] [ dup next-approx ] until
dup length even? [ unclip-last 1 - suffix! 1 suffix! ] when ;
: cw-index ( x -- n )
>continued-fraction <reversed>
[ even? CHAR: 1 CHAR: 0 ? <string> ] map-index concat bin> ;
! Task
"First 20 terms of the Calkin-Wilf sequence:" print
20 calkin-wilf ltake [ pprint bl ] leach nl nl
83116/51639 cw-index "83116/51639 is at index %d.\n" printf
- Output:
First 20 terms of the Calkin-Wilf sequence: 1 1/2 2 1/3 1+1/2 2/3 3 1/4 1+1/3 3/5 2+1/2 2/5 1+2/3 3/4 4 1/5 1+1/4 4/7 2+1/3 3/8 83116/51639 is at index 123456789.
Forth
\ Calkin-Wilf sequence
: frac. swap . ." / " . ;
: cw-next ( num den -- num den ) 2dup / over * 2* over + rot - ;
: cw-seq ( n -- )
1 1 rot
0 do
cr 2dup frac. cw-next
loop 2drop ;
variable index
variable bit-state
variable bit-position
: r2cf-next ( num1 den1 -- num2 den2 u ) swap over >r s>d r> sm/rem ;
: n2bitlength ( n -- )
bit-state @ if
1 swap lshift 1- bit-position @ lshift index +!
else drop then ;
: index-init true bit-state ! 0 bit-position ! 0 index ! ;
: index-build ( n -- )
dup n2bitlength bit-position +! bit-state @ invert bit-state ! ;
: index-finish ( n 0 -- ) 2drop -1 bit-position +! 1 index-build ;
: cw-index ( num den -- )
index-init
begin r2cf-next index-build dup 0<> while repeat
index-finish ;
: cw-demo
20 cw-seq
cr 83116 51639 2dup frac. cw-index index @ . ;
cw-demo
- Output:
1 / 1 1 / 2 2 / 1 1 / 3 3 / 2 2 / 3 3 / 1 1 / 4 4 / 3 3 / 5 5 / 2 2 / 5 5 / 3 3 / 4 4 / 1 1 / 5 5 / 4 4 / 7 7 / 3 3 / 8 83116 / 51639 123456789 ok
FreeBASIC
Uses the code from Greatest common divisor#FreeBASIC as an include.
#include "gcd.bas"
type rational
num as integer
den as integer
end type
dim shared as rational ONE, TWO
ONE.num = 1 : ONE.den = 1
TWO.num = 2 : TWO.den = 1
function simplify( byval a as rational ) as rational
dim as uinteger g = gcd( a.num, a.den )
a.num /= g : a.den /= g
if a.den < 0 then
a.den = -a.den
a.num = -a.num
end if
return a
end function
operator + ( a as rational, b as rational ) as rational
dim as rational ret
ret.num = a.num * b.den + b.num*a.den
ret.den = a.den * b.den
return simplify(ret)
end operator
operator - ( a as rational, b as rational ) as rational
dim as rational ret
ret.num = a.num * b.den - b.num*a.den
ret.den = a.den * b.den
return simplify(ret)
end operator
operator * ( a as rational, b as rational ) as rational
dim as rational ret
ret.num = a.num * b.num
ret.den = a.den * b.den
return simplify(ret)
end operator
operator / ( a as rational, b as rational ) as rational
dim as rational ret
ret.num = a.num * b.den
ret.den = a.den * b.num
return simplify(ret)
end operator
function floor( a as rational ) as rational
dim as rational ret
ret.den = 1
ret.num = a.num \ a.den
return ret
end function
function cw_nextterm( q as rational ) as rational
dim as rational ret = (TWO*floor(q))
ret = ret + ONE : ret = ret - q
return ONE / ret
end function
function frac_to_int( byval a as rational ) as uinteger
redim as uinteger cfrac(-1)
dim as integer lt = -1, ones = 1, ret = 0
do
lt += 1
redim preserve as uinteger cfrac(0 to lt)
cfrac(lt) = floor(a).num
a = a - floor(a) : a = ONE / a
loop until a.num = 0 or a.den = 0
if lt mod 2 = 1 and cfrac(lt) = 1 then
lt -= 1
cfrac(lt)+=1
redim preserve as uinteger cfrac(0 to lt)
end if
if lt mod 2 = 1 and cfrac(lt) > 1 then
cfrac(lt) -= 1
lt += 1
redim preserve as uinteger cfrac(0 to lt)
cfrac(lt) = 1
end if
for i as integer = lt to 0 step -1
for j as integer = 1 to cfrac(i)
ret *= 2
if ones = 1 then ret += 1
next j
ones = 1 - ones
next i
return ret
end function
function disp_rational( a as rational ) as string
if a.den = 1 or a.num= 0 then return str(a.num)
return str(a.num)+"/"+str(a.den)
end function
dim as rational q
q.num = 1
q.den = 1
for i as integer = 1 to 20
print i, disp_rational(q)
q = cw_nextterm(q)
next i
q.num = 83116
q.den = 51639
print disp_rational(q)+" is the "+str(frac_to_int(q))+"th term."
- Output:
1 1 2 1/2 3 2 4 1/3 5 3/2 6 2/3 7 3 8 1/4 9 4/3 10 3/5 11 5/2 12 2/5 13 5/3 14 3/4 15 4 16 1/5 17 5/4 18 4/7 19 7/3 20 3/8 83116/51639 is the 123456789th term.
Fōrmulæ
Fōrmulæ programs are not textual, visualization/edition of programs is done showing/manipulating structures but not text. Moreover, there can be multiple visual representations of the same program. Even though it is possible to have textual representation —i.e. XML, JSON— they are intended for storage and transfer purposes more than visualization and edition.
Programs in Fōrmulæ are created/edited online in its website.
In this page you can see and run the program(s) related to this task and their results. You can also change either the programs or the parameters they are called with, for experimentation, but remember that these programs were created with the main purpose of showing a clear solution of the task, and they generally lack any kind of validation.
Go
Go just has arbitrary precision rational numbers which we use here whilst assuming the numbers needed for this task can be represented exactly by the 64 bit built-in types.
package main
import (
"fmt"
"math"
"math/big"
"strconv"
"strings"
)
func calkinWilf(n int) []*big.Rat {
cw := make([]*big.Rat, n+1)
cw[0] = big.NewRat(1, 1)
one := big.NewRat(1, 1)
two := big.NewRat(2, 1)
for i := 1; i < n; i++ {
t := new(big.Rat).Set(cw[i-1])
f, _ := t.Float64()
f = math.Floor(f)
t.SetFloat64(f)
t.Mul(t, two)
t.Sub(t, cw[i-1])
t.Add(t, one)
t.Inv(t)
cw[i] = new(big.Rat).Set(t)
}
return cw
}
func toContinued(r *big.Rat) []int {
a := r.Num().Int64()
b := r.Denom().Int64()
var res []int
for {
res = append(res, int(a/b))
t := a % b
a, b = b, t
if a == 1 {
break
}
}
le := len(res)
if le%2 == 0 { // ensure always odd
res[le-1]--
res = append(res, 1)
}
return res
}
func getTermNumber(cf []int) int {
b := ""
d := "1"
for _, n := range cf {
b = strings.Repeat(d, n) + b
if d == "1" {
d = "0"
} else {
d = "1"
}
}
i, _ := strconv.ParseInt(b, 2, 64)
return int(i)
}
func commatize(n int) string {
s := fmt.Sprintf("%d", n)
if n < 0 {
s = s[1:]
}
le := len(s)
for i := le - 3; i >= 1; i -= 3 {
s = s[0:i] + "," + s[i:]
}
if n >= 0 {
return s
}
return "-" + s
}
func main() {
cw := calkinWilf(20)
fmt.Println("The first 20 terms of the Calkin-Wilf sequnence are:")
for i := 1; i <= 20; i++ {
fmt.Printf("%2d: %s\n", i, cw[i-1].RatString())
}
fmt.Println()
r := big.NewRat(83116, 51639)
cf := toContinued(r)
tn := getTermNumber(cf)
fmt.Printf("%s is the %sth term of the sequence.\n", r.RatString(), commatize(tn))
}
- Output:
The first 20 terms of the Calkin-Wilf sequnence are: 1: 1 2: 1/2 3: 2 4: 1/3 5: 3/2 6: 2/3 7: 3 8: 1/4 9: 4/3 10: 3/5 11: 5/2 12: 2/5 13: 5/3 14: 3/4 15: 4 16: 1/5 17: 5/4 18: 4/7 19: 7/3 20: 3/8 83116/51639 is the 123,456,789th term of the sequence.
Haskell
import Control.Monad (forM_)
import Data.Bool (bool)
import Data.List.NonEmpty (NonEmpty, fromList, toList, unfoldr)
import Text.Printf (printf)
-- The infinite Calkin-Wilf sequence, a(n), starting with a(1) = 1.
calkinWilfs :: [Rational]
calkinWilfs = iterate (recip . succ . ((-) =<< (2 *) . fromIntegral . floor)) 1
-- The index into the Calkin-Wilf sequence of a given rational number, starting
-- with 1 at index 1.
calkinWilfIdx :: Rational -> Integer
calkinWilfIdx = rld . cfo
-- A continued fraction representation of a given rational number, guaranteed
-- to have an odd length.
cfo :: Rational -> NonEmpty Int
cfo = oddLen . cf
-- The canonical (i.e. shortest) continued fraction representation of a given
-- rational number.
cf :: Rational -> NonEmpty Int
cf = unfoldr step
where
step r =
case properFraction r of
(n, 1) -> (succ n, Nothing)
(n, 0) -> (n, Nothing)
(n, f) -> (n, Just (recip f))
-- Ensure a continued fraction has an odd length.
oddLen :: NonEmpty Int -> NonEmpty Int
oddLen = fromList . go . toList
where
go [x, y] = [x, pred y, 1]
go (x:y:zs) = x : y : go zs
go xs = xs
-- Run-length decode a continued fraction.
rld :: NonEmpty Int -> Integer
rld = snd . foldr step (True, 0)
where
step i (b, n) =
let p = 2 ^ i
in (not b, n * p + bool 0 (pred p) b)
main :: IO ()
main = do
forM_ (take 20 $ zip [1 :: Int ..] calkinWilfs) $
\(i, r) -> printf "%2d %s\n" i (show r)
let r = 83116 / 51639
printf
"\n%s is at index %d of the Calkin-Wilf sequence.\n"
(show r)
(calkinWilfIdx r)
- Output:
1 1 % 1 2 1 % 2 3 2 % 1 4 1 % 3 5 3 % 2 6 2 % 3 7 3 % 1 8 1 % 4 9 4 % 3 10 3 % 5 11 5 % 2 12 2 % 5 13 5 % 3 14 3 % 4 15 4 % 1 16 1 % 5 17 5 % 4 18 4 % 7 19 7 % 3 20 3 % 8 83116 % 51639 is at index 123456789 of the Calkin-Wilf sequence.
J
cw_next_term^:(<20)1x 1 1r2 2 1r3 3r2 2r3 3 1r4 4r3 3r5 5r2 2r5 5r3 3r4 4 1r5 5r4 4r7 7r3 3r8 (,. index_cw_term&>) 3r4 53r37 83116r51639 3r4 14 53r37 1081 83116r51639 123456789
given definitions
cw_next_term=: [: % +:@<. + -.
ccf =: compute_continued_fraction=: 3 :0
if. 0 -: y do.
, 0
else.
result=. i. 0
remainder=. % y
whilst. remainder do.
remainder=. % remainder
integer_part=. <. remainder
remainder=. remainder - integer_part
result=. result , integer_part
end.
end.
)
molcf =: make_odd_length_continued_fraction=: (}: , 1 ,~ <:@{:)^:(0 -: 2 | #)
NB. base 2 @ reverse @ the cf's representation copies of 1 0 1 0 ...
index_cw_term=: #.@|.@(# 1 0 $~ #)@molcf@ccf
Note that ccf
could be expressed more concisely:
ccf=: _1 {"1 |.@(0 1 #: %@{.)^:(0~:{.)^:a:
Java
import java.util.ArrayDeque;
import java.util.Deque;
public final class CalkinWilfSequence {
public static void main(String[] aArgs) {
Rational term = Rational.ONE;
System.out.println("First 20 terms of the Calkin-Wilf sequence are:");
for ( int i = 1; i <= 20; i++ ) {
System.out.println(String.format("%2d", i) + ": " + term);
term = nextCalkinWilf(term);
}
System.out.println();
Rational rational = new Rational(83_116, 51_639);
System.out.println(" " + rational + " is the " + termIndex(rational) + "th term of the sequence.");
}
private static Rational nextCalkinWilf(Rational aTerm) {
Rational divisor = Rational.TWO.multiply(aTerm.floor()).add(Rational.ONE).subtract(aTerm);
return Rational.ONE.divide(divisor);
}
private static long termIndex(Rational aRational) {
long result = 0;
long binaryDigit = 1;
long power = 0;
for ( long term : continuedFraction(aRational) ) {
for ( long i = 0; i < term; power++, i++ ) {
result |= ( binaryDigit << power );
}
binaryDigit = ( binaryDigit == 0 ) ? 1 : 0;
}
return result;
}
private static Deque<Long> continuedFraction(Rational aRational) {
long numerator = aRational.numerator();
long denominator = aRational.denominator();
Deque<Long> result = new ArrayDeque<Long>();
while ( numerator != 1 ) {
result.addLast(numerator / denominator);
long copyNumerator = numerator;
numerator = denominator;
denominator = copyNumerator % denominator;
}
if ( ! result.isEmpty() && result.size() % 2 == 0 ) {
final long back = result.removeLast();
result.addLast(back - 1);
result.addLast(1L);
}
return result;
}
}
final class Rational {
public Rational(long aNumerator, long aDenominator) {
if ( aDenominator == 0 ) {
throw new ArithmeticException("Denominator cannot be zero");
}
if ( aNumerator == 0 ) {
aDenominator = 1;
}
if ( aDenominator < 0 ) {
numer = -aNumerator;
denom = -aDenominator;
} else {
numer = aNumerator;
denom = aDenominator;
}
final long gcd = gcd(numer, denom);
numer = numer / gcd;
denom = denom / gcd;
}
public Rational add(Rational aOther) {
return new Rational(numer * aOther.denom + aOther.numer * denom, denom * aOther.denom);
}
public Rational subtract(Rational aOther) {
return new Rational(numer * aOther.denom - aOther.numer * denom, denom * aOther.denom);
}
public Rational multiply(Rational aOther) {
return new Rational(numer * aOther.numer, denom * aOther.denom);
}
public Rational divide(Rational aOther) {
return new Rational(numer * aOther.denom, denom * aOther.numer);
}
public Rational floor() {
return new Rational(numer / denom, 1);
}
public long numerator() {
return numer;
}
public long denominator() {
return denom;
}
@Override
public String toString() {
return numer + "/" + denom;
}
public static final Rational ONE = new Rational(1, 1);
public static final Rational TWO = new Rational(2, 1);
private long gcd(long aOne, long aTwo) {
if ( aTwo == 0 ) {
return aOne;
}
return gcd(aTwo, aOne % aTwo);
}
private long numer;
private long denom;
}
- Output:
First 20 terms of the Calkin-Wilf sequence are: 1: 1/1 2: 1/2 3: 2/1 4: 1/3 5: 3/2 6: 2/3 7: 3/1 8: 1/4 9: 4/3 10: 3/5 11: 5/2 12: 2/5 13: 5/3 14: 3/4 15: 4/1 16: 1/5 17: 5/4 18: 4/7 19: 7/3 20: 3/8 83116/51639 is the 123456789th term of the sequence.
jq
Adapted from Wren
Also works with gojq, the Go implementation of jq, and with fq
See Arithmetic/Rational#jq for the Rational module included by the `include` directive. In this module, rationals are represented by JSON objects of the form {n, d}, where .n and .d are the numerator and denominator respectively. r(n;d) is the constructor function, and r(n;d) is pretty-printed as `n // d`.
include "rational"; # see [[Arithmetic/Rational#jq]]
### Generic Utilities
# counting from 0
def enumerate(s): foreach s as $x (-1; .+1; [., $x]);
# input string is converted from "base" to an integer, within limits
# of the underlying arithmetic operations, and without error-checking:
def to_i(base):
explode
| reverse
| map(if . > 96 then . - 87 else . - 48 end) # "a" ~ 97 => 10 ~ 87
| reduce .[] as $c
# state: [power, ans]
([1,0]; (.[0] * base) as $b | [$b, .[1] + (.[0] * $c)])
| .[1];
### The Calkin-Wilf Sequence
# Emit an array of $n terms
def calkinWilf($n):
reduce range(1;$n) as $i ( [r(1;1)];
radd(1; rminus( rmult(2; (.[$i-1]|rfloor)); .[$i-1])) as $t
| .[$i] = rdiv(r(1;1) ; $t)) ;
# input: a Rational
def toContinued:
{ a: .n,
b: .d,
res: [] }
| until( .break;
.res += [.a / .b | floor]
| (.a % .b) as $t
| .a = .b
| .b = $t
| .break = (.a == 1) )
| if .res|length % 2 == 0
then # ensure always odd
.res[-1] += -1
| .res += [1]
else .
end
| .res;
# input: an array representing a continued fraction
def getTermNumber:
reduce .[] as $n ( {b: "", d: "1"};
.b = (.d * $n) + .b
| .d = (if .d == "1" then "0" else "1" end))
| .b | to_i(2) ;
# input: a Rational in the Calkin-Wilf sequence
def getTermNumber:
reduce .[] as $n ( {b: "", d: "1"};
.b = (.d * $n) + .b
| .d = (if .d == "1" then "0" else "1" end))
| .b | to_i(2) ;
def task(r):
"The first 20 terms of the Calkin-Wilf sequence are:",
(enumerate(calkinWilf(20)[]) | "\(1+.[0]): \(.[1]|rpp)" ),
"",
"\(r|rpp) is term # \(r|toContinued|getTermNumber) of the sequence.";
task( r(83116; 51639) )
Invocation: jq -nrf calkin-wilf-sequence.jq
- Output:
The first 20 terms of the Calkin-Wilf sequence are: 1: 1 // 1 2: 1 // 2 3: 2 // 1 4: 1 // 3 5: 3 // 2 6: 2 // 3 7: 3 // 1 8: 1 // 4 9: 4 // 3 10: 3 // 5 11: 5 // 2 12: 2 // 5 13: 5 // 3 14: 3 // 4 15: 4 // 1 16: 1 // 5 17: 5 // 4 18: 4 // 7 19: 7 // 3 20: 3 // 8 83116 // 51639 is term # 123456789 of the sequence.
Julia
function calkin_wilf(n)
cw = zeros(Rational, n + 1)
for i in 2:n + 1
t = Int(floor(cw[i - 1])) * 2 - cw[i - 1] + 1
cw[i] = 1 // t
end
return cw[2:end]
end
function continued(r::Rational)
a, b = r.num, r.den
res = []
while true
push!(res, Int(floor(a / b)))
a, b = b, a % b
a == 1 && break
end
return res
end
function term_number(cf)
b, d = "", "1"
for n in cf
b = d^n * b
d = (d == "1") ? "0" : "1"
end
return parse(Int, b, base=2)
end
const cw = calkin_wilf(20)
println("The first 20 terms of the Calkin-Wilf sequence are: $cw")
const r = 83116 // 51639
const cf = continued(r)
const tn = term_number(cf)
println("$r is the $tn-th term of the sequence.")
- Output:
The first 20 terms of the Calkin-Wilf sequence are: Rational[1//1, 1//2, 2//1, 1//3, 3//2, 2//3, 3//1, 1//4, 4//3, 3//5, 5//2, 2//5, 5//3, 3//4, 4//1, 1//5, 5//4, 4//7, 7//3, 3//8] 83116//51639 is the 123456789-th term of the sequence.
Little Man Computer
Runs in a home-made simulator, which is mostly compatible with Peter Higginson's online simulator. Only, for better control of the output format, I've added an instruction OTX (extended output). To run the code in PH's simulator, replace OTX and its parameter with OUT and no parameter.
Find first n terms
// Little Man Computer, for Rosetta Code.
// Displays terms of Calkin-Wilf sequence up to the given index.
// The chosen algorithm calculates the i-th term directly from i
// (i.e. not using any previous terms).
input INP // get number of terms from user
BRZ exit // exit if 0
STA max_i // store maximum index
LDA c1 // index := 1
next_i STA i
// Write index followed by '->'
OTX 3 // non-standard: minimum width 3, no new line
LDA asc_hy
OTC
LDA asc_gt
OTC
// Find greatest power of 2 not exceeding i,
// and count the number of binary digits in i.
LDA c1
STA pwr2
loop2 STA nrDigits
LDA i
SUB pwr2
SUB pwr2
BRP double
BRA part2 // jump out if next power of 2 would exceed i
double LDA pwr2
ADD pwr2
STA pwr2
LDA nrDigits
ADD c1
BRA loop2
// The nth term a/b is calculated from the binary digits of i.
// The leading 1 is not used.
part2 LDA c1
STA a // a := 1
STA b // b := 1
LDA i
SUB pwr2
STA diff
// Pre-decrement count, since leading 1 is not used
dec_ct LDA nrDigits // count down the number of digits
SUB c1
BRZ output // if all digits done, output the result
STA nrDigits
// We now want to compare diff with pwr2/2.
// Since division is awkward in LMC, we compare 2*diff with pwr2.
LDA diff // diff := 2*diff
ADD diff
STA diff
SUB pwr2 // is diff >= pwr2 ?
BRP digit_1 // binary digit is 1 if yes, 0 if no
// If binary digit is 0 then set b := a + b
LDA a
ADD b
STA b
BRA dec_ct
// If binary digit is 1 then update diff and set a := a + b
digit_1 STA diff
LDA a
ADD b
STA a
BRA dec_ct
// Now have nth term a/b. Write it to the output.
output LDA a // write a
OTX 1 // non-standard: minimum width 1; no new line
LDA asc_sl // write slash
OTC
LDA b // write b
OTX 11 // non-standard: minimum width 1; add new line
LDA i // have we done maximum i yet?
SUB max_i
BRZ exit // if yes, exit
LDA i // if no, increment i and loop back
ADD c1
BRA next_i
exit HLT
// Constants
c1 DAT 1
asc_hy DAT 45
asc_gt DAT 62
asc_sl DAT 47
// Variables
i DAT
max_i DAT
pwr2 DAT
nrDigits DAT
diff DAT
a DAT
b DAT
// end
- Output:
1->1/1 2->1/2 3->2/1 4->1/3 5->3/2 6->2/3 7->3/1 8->1/4 9->4/3 10->3/5 11->5/2 12->2/5 13->5/3 14->3/4 15->4/1 16->1/5 17->5/4 18->4/7 19->7/3 20->3/8
Find index of a given term
The numbers in part 2 of the task are too large for LMC, so the demo program just confirms the example, that 9/4 is the 35th term.
// Little Man Computer, for Rosetta Code.
// Calkin-Wilf sequence: displays index of term entered by user.
INP // get numerator from user
BRZ exit // exit if 0
STA num
STA a // initialize a := numerator
INP // get denominator from user
BRZ exit // exit if 0
STA den
STA b // initialize b := denominator
LDA c0 // initialize index := 0
STA index
LDA c1 // initialize power of 2 := 1
STA pwr2
// Build binary digits of the index
loop LDA a // is a = b yet?
SUB b
BRZ break // if yes, break out of loop
BRP a_gt_b // jump if a > b
// If a < b then b := b - a, binary digit is 0
LDA b
SUB a
STA b
BRA double
// If a > b then a := a - b, binary digit is 1
a_gt_b STA a
LDA index
ADD pwr2
STA index
// In either case, on to next power of 2
double LDA pwr2
ADD pwr2
STA pwr2
BRA loop
// Out of loop, add leading binary digit 1
break LDA index
ADD pwr2
STA index
// Output the result
LDA num
OTX 1 // non-standard: minimum width = 1, no new line
LDA asc_sl
OTC
LDA den
OTX 1
LDA asc_lt // write '<-' after fraction
OTC
LDA asc_hy
OTC
LDA index
OTX 11 // non-standard: minimum width = 1, add new line
exit HLT
// Constants
c0 DAT 0
c1 DAT 1
asc_sl DAT 47
asc_lt DAT 60
asc_hy DAT 45
// Variables
num DAT
den DAT
a DAT
b DAT
pwr2 DAT
index DAT
// end
- Output:
9/4<-35
Mathematica / Wolfram Language
ClearAll[a]
a[1] = 1;
a[n_?(GreaterThan[1])] := a[n] = 1/(2 Floor[a[n - 1]] + 1 - a[n - 1])
a /@ Range[20]
ClearAll[a]
a = 1;
n = 1;
Dynamic[n]
done = False;
While[! done,
a = 1/(2 Floor[a] + 1 - a);
n++;
If[a == 83116/51639,
Print[n];
Break[];
]
]
- Output:
{1, 1/2, 2, 1/3, 3/2, 2/3, 3, 1/4, 4/3, 3/5, 5/2, 2/5, 5/3, 3/4, 4, 1/5, 5/4, 4/7, 7/3, 3/8} 123456789
Maxima
/* The function fusc is related to Calkin-Wilf sequence */
fusc(n):=block(
[k:n,a:1,b:0],
while k>0 do (if evenp(k) then (k:k/2,a:a+b) else (k:(k-1)/2,b:a+b)),
b)$
/* Calkin-Wilf function using fusc */
calkin_wilf(n):=fusc(n)/fusc(n+1)$
/* Function that given a nonnegative rational returns its position in the Calkin-Wilf sequence */
cf_bin(fracti):=block(
cf_list:cf(fracti),
cf_len:length(cf_list),
if oddp(cf_len) then cf_list:reverse(cf_list) else cf_list:reverse(append(at(cf_list,[cf_list[cf_len]=cf_list[cf_len]-1]),[1])),
makelist(lambda([x],if oddp(x) then makelist(1,j,1,cf_list[x]) else makelist(0,j,1,cf_list[x]))(i),i,1,length(cf_list)), /* decoding part begins here */
apply(append,%%),
apply("+",makelist(2^i,i,0,length(%%)-1)*reverse(%%)))$
/* Test cases */
/* 20 first terms of the sequence */
makelist(calkin_wilf(i),i,1,20);
/* Position of 83116/51639 in Calkin-Wilf sequence */
83116/51639$
cf_bin(%);
- Output:
[1,1/2,2,1/3,3/2,2/3,3,1/4,4/3,3/5,5/2,2/5,5/3,3/4,4,1/5,5/4,4/7,7/3,3/8] 123456789
Nim
We ignored the standard module “rationals” which is slow and have rather defined a fraction as a tuple of two 32 bits unsigned integers (slightly faster than 64 bits signed integers and sufficient for this task). Moreover, we didn’t do operations on fractions and computed directly the numerator and denominator values at each step. The fractions built this way are irreducible (which avoids to compute a GCD which is a slow operation). With these optimizations, the program runs in less than 1.3 s on our laptop.
type Fraction = tuple[num, den: uint32]
iterator calkinWilf(): Fraction =
## Yield the successive values of the sequence.
var n, d = 1u32
yield (n, d)
while true:
n = 2 * (n div d) * d + d - n
swap n, d
yield (n, d)
proc `$`(fract: Fraction): string =
## Return the representation of a fraction.
$fract.num & '/' & $fract.den
func `==`(a, b: Fraction): bool {.inline.} =
## Compare two fractions. Slightly faster than comparison of tuples.
a.num == b.num and a.den == b.den
when isMainModule:
echo "The first 20 terms of the Calkwin-Wilf sequence are:"
var count = 0
for an in calkinWilf():
inc count
stdout.write $an & ' '
if count == 20: break
stdout.write '\n'
const Target: Fraction = (83116u32, 51639u32)
var index = 0
for an in calkinWilf():
inc index
if an == Target: break
echo "\nThe element ", $Target, " is at position ", $index, " in the sequence."
- Output:
The first 20 terms of the Calkwin-Wilf sequence are: 1/1 1/2 2/1 1/3 3/2 2/3 3/1 1/4 4/3 3/5 5/2 2/5 5/3 3/4 4/1 1/5 5/4 4/7 7/3 3/8 The element 83116/51639 is at position 123456789 in the sequence.
PARI/GP
\\ This function assumes the existence of a global variable 'an' for 'a[n]'
a(n) = if(n==1, 1, 1 / (2 * floor(an[n-1]) + 1 - an[n-1]));
\\ We will use a vector to hold the values and compute them iteratively to avoid stack overflow
an = vector(20);
an[1] = 1;
for(i=2, 20, an[i] = a(i));
\\ Now we print the vector
print(an);
\\ Initialize variables for the while loop
a = 1;
n = 1;
\\ Loop until the condition is met
while(a != 83116/51639,{
a = 1/(2 * floor(a) + 1 - a);
if(n>=123456789,print(n));
n++;
});
\\ Output the number of iterations needed to reach 83116/51639
print(n);
- Output:
[1, 1/2, 2, 1/3, 3/2, 2/3, 3, 1/4, 4/3, 3/5, 5/2, 2/5, 5/3, 3/4, 4, 1/5, 5/4, 4/7, 7/3, 3/8] 123456789
Pascal
These programs were written in Free Pascal, using the Lazarus IDE and the Free Pascal compiler version 3.2.0. They are based on the Wikipedia article "Calkin-Wilf tree", rather than the algorithms in the task description.
program CWTerms;
{-------------------------------------------------------------------------------
FreePascal command-line program.
Calculates the Calkin-Wilf sequence up to the specified maximum index,
where the first term 1/1 has index 1.
Command line format is: CWTerms <max_index>
The program demonstrates 3 algorithms for calculating the sequence:
(1) Calculate term[2n] and term[2n + 1] from term[n]
(2) Calculate term[n + 1] from term[n]
(3) Calculate term[n] directly from n, without using other terms
Algorithm 1 is called first, and stores the terms in an array.
Then the program calls Algorithms 2 and 3, and checks that they agree
with Algorithm 1.
-------------------------------------------------------------------------------}
uses SysUtils;
type TRational = record
Num, Den : integer;
end;
var
terms : array of TRational;
max_index, k : integer;
// Routine to calculate array of terms up the the maiximum index
procedure CalcTerms_algo_1();
var
j, k : integer;
begin
SetLength( terms, max_index + 1);
j := 1; // index to earlier term, from which current term is calculated
k := 1; // index to current term
terms[1].Num := 1;
terms[1].Den := 1;
while (k < max_index) do begin
inc(k);
if (k and 1) = 0 then begin // or could write "if not Odd(k)"
terms[k].Num := terms[j].Num;
terms[k].Den := terms[j].Num + terms[j].Den;
end
else begin
terms[k].Num := terms[j].Num + terms[j].Den;
terms[k].Den := terms[j].Den;
inc(j);
end;
end;
end;
// Method to get each term from the preceding term.
// a/b --> b/(a + b - 2(a mod b));
function CheckTerms_algo_2() : boolean;
var
index, a, b, temp : integer;
begin
result := true;
index := 1;
a := 1;
b := 1;
while (index <= max_index) do begin
if (a <> terms[index].Num) or (b <> terms[index].Den) then
result := false;
temp := a + b - 2*(a mod b);
a := b;
b := temp;
inc( index)
end;
end;
// Mathod to calcualte each term from its index, without using other terms.
function CheckTerms_algo_3() : boolean;
var
index, a, b, pwr2, idiv2 : integer;
begin
result := true;
for index := 1 to max_index do begin
idiv2 := index div 2;
pwr2 := 1;
while (pwr2 <= idiv2) do pwr2 := pwr2 shl 1;
a := 1;
b := 1;
while (pwr2 > 1) do begin
pwr2 := pwr2 shr 1;
if (pwr2 and index) = 0 then
inc( b, a)
else
inc( a, b);
end;
if (a <> terms[index].Num) or (b <> terms[index].Den) then
result := false;
end;
end;
begin
// Read and validate maximum index
max_index := SysUtils.StrToIntDef( paramStr(1), -1); // -1 if not an integer
if (max_index <= 0) then begin
WriteLn( 'Maximum index must be a positive integer');
exit;
end;
// Calculate terms by algo 1, then check that algos 2 and 3 agree.
CalcTerms_algo_1();
if not CheckTerms_algo_2() then begin
WriteLn( 'Algorithm 2 failed');
exit;
end;
if not CheckTerms_algo_3() then begin
WriteLn( 'Algorithm 3 failed');
exit;
end;
// Display the terms
for k := 1 to max_index do
with terms[k] do
WriteLn( SysUtils.Format( '%8d: %d/%d', [k, Num, Den]));
end.
- Output:
1: 1/1 2: 1/2 3: 2/1 4: 1/3 5: 3/2 6: 2/3 7: 3/1 8: 1/4 9: 4/3 10: 3/5 11: 5/2 12: 2/5 13: 5/3 14: 3/4 15: 4/1 16: 1/5 17: 5/4 18: 4/7 19: 7/3 20: 3/8
program CWIndex;
{-------------------------------------------------------------------------------
FreePascal command-line program.
Calculates index of a rational number in the Calkin-Wilf sequence,
where the first term 1/1 has index 1.
Command line format is
CWIndex <numerator> <denominator>
e.g. for the Rosetta Code example
CWIndex 83116 51639
-------------------------------------------------------------------------------}
uses SysUtils;
var
num, den : integer;
a, b : integer;
pwr2, index : qword; // 64-bit unsiged
begin
// Read and validate input.
num := SysUtils.StrToIntDef( paramStr(1), -1); // return -1 if not an integer
den := SysUtils.StrToIntDef( paramStr(2), -1);
if (num <= 0) or (den <= 0) then begin
WriteLn( 'Numerator and denominator must be positive integers');
exit;
end;
// Input OK, calculate and display index of num/den
// The index may overflow 64 bits, so turn on overflow detection
{$Q+}
a := num;
b := den;
pwr2 := 1;
index := 0;
try
while (a <> b) do begin
if (a < b) then
dec( b, a)
else begin
dec( a, b);
inc( index, pwr2);
end;
pwr2 := 2*pwr2;
end;
inc( index, pwr2);
WriteLn( SysUtils.Format( 'Index of %d/%d is %u', [num, den, index]));
except
WriteLn( 'Index is too large for 64 bits');
end;
end.
- Output:
Index of 83116/51639 is 123456789
Perl
use strict;
use warnings;
use feature qw(say state);
use ntheory 'fromdigits';
use List::Lazy 'lazy_list';
use Math::AnyNum ':overload';
my $calkin_wilf = lazy_list { state @cw = 1; push @cw, 1 / ( (2 * int $cw[0]) + 1 - $cw[0] ); shift @cw };
sub r2cf {
my($num, $den) = @_;
my($n, @cf);
my $f = sub { return unless $den;
my $q = int($num/$den);
($num, $den) = ($den, $num - $q*$den);
$q;
};
push @cf, $n while defined($n = $f->());
reverse @cf;
}
sub r2cw {
my($num, $den) = @_;
my $bits;
my @f = r2cf($num, $den);
$bits .= ($_%2 ? 0 : 1) x $f[$_] for 0..$#f;
fromdigits($bits, 2);
}
say 'First twenty terms of the Calkin-Wilf sequence:';
printf "%s ", $calkin_wilf->next() for 1..20;
say "\n\n83116/51639 is at index: " . r2cw(83116,51639);
- Output:
First twenty terms of the Calkin-Wilf sequence: 1 1/2 2 1/3 3/2 2/3 3 1/4 4/3 3/5 5/2 2/5 5/3 3/4 4 1/5 5/4 4/7 7/3 3/8 83116/51639 is at index: 123456789
Phix
with javascript_semantics requires("1.0.0") -- (new even() builtin) function calkin_wilf(integer len) sequence cw = repeat(0,len) integer n=0, d=1 for i=1 to len do {n,d} = {d,(floor(n/d)*2+1)*d-n} cw[i] = {n,d} end for return cw end function function odd_length(sequence cf) -- replace even length continued fraction with odd length equivalent -- if remainder(length(cf),2)=0 then if even(length(cf)) then cf[$] -= 1 cf &= 1 end if return cf end function function to_continued_fraction(sequence r) integer {a,b} = r sequence cf = {} while true do cf &= floor(a/b) {a, b} = {b, remainder(a,b)} if a=1 then exit end if end while cf = odd_length(cf) return cf end function function get_term_number(sequence cf) sequence b = {} integer d = 1 for i=1 to length(cf) do b &= repeat(d,cf[i]) d = 1-d end for integer t = bits_to_int(b) return t end function -- additional verification methods (2 of) function i_to_cf(integer i) -- sequence b = trim_tail(int_to_bits(i,32),0)&2, sequence b = int_to_bits(i)&2, cf = iff(b[1]=0?{0}:{}) while length(b)>1 do for j=2 to length(b) do if b[j]!=b[1] then cf &= j-1 b = b[j..$] exit end if end for end while cf = odd_length(cf) return cf end function function cf2r(sequence cf) integer n=0, d=1 for i=length(cf) to 2 by -1 do {n,d} = {d,n+d*cf[i]} end for return {n+cf[1]*d,d} end function function prettyr(sequence r) integer {n,d} = r return iff(d=1?sprintf("%d",n):sprintf("%d/%d",{n,d})) end function sequence cw = calkin_wilf(20) printf(1,"The first 20 terms of the Calkin-Wilf sequence are:\n") for i=1 to 20 do string s = prettyr(cw[i]), r = prettyr(cf2r(i_to_cf(i))) integer t = get_term_number(to_continued_fraction(cw[i])) printf(1,"%2d: %-4s [==> %2d: %-3s]\n", {i, s, t, r}) end for printf(1,"\n") sequence r = {83116,51639} sequence cf = to_continued_fraction(r) integer tn = get_term_number(cf) printf(1,"%d/%d is the %,d%s term of the sequence.\n", r&{tn,ord(tn)})
- Output:
The first 20 terms of the Calkin-Wilf sequence are: 1: 1 [==> 1: 1 ] 2: 1/2 [==> 2: 1/2] 3: 2 [==> 3: 2 ] 4: 1/3 [==> 4: 1/3] 5: 3/2 [==> 5: 3/2] 6: 2/3 [==> 6: 2/3] 7: 3 [==> 7: 3 ] 8: 1/4 [==> 8: 1/4] 9: 4/3 [==> 9: 4/3] 10: 3/5 [==> 10: 3/5] 11: 5/2 [==> 11: 5/2] 12: 2/5 [==> 12: 2/5] 13: 5/3 [==> 13: 5/3] 14: 3/4 [==> 14: 3/4] 15: 4 [==> 15: 4 ] 16: 1/5 [==> 16: 1/5] 17: 5/4 [==> 17: 5/4] 18: 4/7 [==> 18: 4/7] 19: 7/3 [==> 19: 7/3] 20: 3/8 [==> 20: 3/8] 83116/51639 is the 123,456,789th term of the sequence.
Prolog
% John Devou: 26-Nov-2021
% g(N,X):- consecutively generate in X the first N elements of the Calkin-Wilf sequence
g(N,[A/B|_]-_,A/B):- N > 0.
g(N,[A/B|Ls]-[A/C,C/B|Ys],X):- N > 1, M is N-1, C is A+B, g(M,Ls-Ys,X).
g(N,X):- g(N,[1/1|Ls]-Ls,X).
% t(A/B,X):- generate in X the index of A/B in the Calkin-Wilf sequence
t(A/1,S,C,X):- X is C*(2**(A-1+S)-S).
t(A/B,S,C,X):- B > 1, divmod(A,B,M,N), T is 1-S, D is C*2**M, t(B/N,T,D,Y), X is Y + S*C*(2**M-1).
t(A/B,X):- t(A/B,1,1,X), !.
- Output:
?- findall(X, g(20,X), L), write(L). [1/1,1/2,2/1,1/3,3/2,2/3,3/1,1/4,4/3,3/5,5/2,2/5,5/3,3/4,4/1,1/5,5/4,4/7,7/3,3/8] L = [1/1, 1/2, 2/1, 1/3, 3/2, 2/3, 3/1, 1/4, ... / ...|...]. ?- t(83116/51639,X). X = 123456789.
Python
from fractions import Fraction
from math import floor
from itertools import islice, groupby
def cw():
a = Fraction(1)
while True:
yield a
a = 1 / (2 * floor(a) + 1 - a)
def r2cf(rational):
num, den = rational.numerator, rational.denominator
while den:
num, (digit, den) = den, divmod(num, den)
yield digit
def get_term_num(rational):
ans, dig, pwr = 0, 1, 0
for n in r2cf(rational):
for _ in range(n):
ans |= dig << pwr
pwr += 1
dig ^= 1
return ans
if __name__ == '__main__':
print('TERMS 1..20: ', ', '.join(str(x) for x in islice(cw(), 20)))
x = Fraction(83116, 51639)
print(f"\n{x} is the {get_term_num(x):_}'th term.")
- Output:
TERMS 1..20: 1, 1/2, 2, 1/3, 3/2, 2/3, 3, 1/4, 4/3, 3/5, 5/2, 2/5, 5/3, 3/4, 4, 1/5, 5/4, 4/7, 7/3, 3/8 83116/51639 is the 123_456_789'th term.
Quackery
cf
is defined at Continued fraction/Arithmetic/Construct from rational number#Quackery.
[ $ "bigrat.qky" loadfile ] now!
[ ' [ [ 1 1 ] ]
swap 1 - times
[ dup -1 peek do
2dup proper 2drop
2 * n->v
2swap -v 1 n->v v+ v+
1/v join nested join ] ] is calkin-wilf ( n --> [ )
[ 1 & ] is odd ( n --> b )
[ dup size odd not if
[ -1 split do
1 - join
1 join ] ] is oddcf ( [ --> [ )
[ 0 swap
reverse witheach
[ i odd iff
<< done
dup dip <<
bit 1 - | ] ] is rl->n ( [ --> n )
[ cf oddcf rl->n ] is cw-term ( n/d --> n )
20 calkin-wilf
witheach
[ do vulgar$ echo$ sp ]
cr cr
83116 51639 cw-term echo
- Output:
1/1 1/2 2/1 1/3 3/2 2/3 3/1 1/4 4/3 3/5 5/2 2/5 5/3 3/4 4/1 1/5 5/4 4/7 7/3 3/8 123456789
Raku
In Raku, arrays are indexed from 0. The Calkin-Wilf sequence does not have a term defined at 0.
This implementation includes a bogus undefined value at position 0, having the bogus first term shifts the indices up by one, making the ordinal position and index match. Useful due to how reversibility function works.
my @calkin-wilf = Any, 1, {1 / (.Int × 2 + 1 - $_)} … *;
# Rational to Calkin-Wilf index
sub r2cw (Rat $rat) { :2( join '', flat (flat (1,0) xx *) Zxx reverse r2cf $rat ) }
# The task
say "First twenty terms of the Calkin-Wilf sequence: ",
@calkin-wilf[1..20]».&prat.join: ', ';
say "\n99991st through 100000th: ",
(my @tests = @calkin-wilf[99_991 .. 100_000])».&prat.join: ', ';
say "\nCheck reversibility: ", @tests».Rat».&r2cw.join: ', ';
say "\n83116/51639 is at index: ", r2cw 83116/51639;
# Helper subs
sub r2cf (Rat $rat is copy) { # Rational to continued fraction
gather loop {
$rat -= take $rat.floor;
last if !$rat;
$rat = 1 / $rat;
}
}
sub prat ($num) { # pretty Rat
return $num unless $num ~~ Rat|FatRat;
return $num.numerator if $num.denominator == 1;
$num.nude.join: '/';
}
- Output:
First twenty terms of the Calkin-Wilf sequence: 1, 1/2, 2, 1/3, 3/2, 2/3, 3, 1/4, 4/3, 3/5, 5/2, 2/5, 5/3, 3/4, 4, 1/5, 5/4, 4/7, 7/3, 3/8 99991st through 100000th: 1085/303, 303/1036, 1036/733, 733/1163, 1163/430, 430/987, 987/557, 557/684, 684/127, 127/713 Check reversibility: 99991, 99992, 99993, 99994, 99995, 99996, 99997, 99998, 99999, 100000 83116/51639 is at index: 123456789
REXX
The meat of this REXX program was provided by Paul Kislanko.
/*REXX pgm finds the Nth value of the Calkin─Wilf sequence (which will be a fraction),*/
/*────────────────────── or finds which sequence number contains a specified fraction). */
numeric digits 2000 /*be able to handle ginormic integers. */
parse arg LO HI te . /*obtain optional arguments from the CL*/
if LO=='' | LO=="," then LO= 1 /*Not specified? Then use the default.*/
if HI=='' | HI=="," then HI= 20 /* " " " " " " */
if te=='' | te=="," then te= '/' /* " " " " " " */
if datatype(LO, 'W') then call CW_terms /*Is LO numeric? Then show some terms.*/
if pos('/', te)>0 then call CW_frac te /*Does TE have a / ? Then find term #*/
exit 0
/*──────────────────────────────────────────────────────────────────────────────────────*/
commas: parse arg ?; do jc=length(?)-3 to 1 by -3; ?=insert(',', ?, jc); end; return ?
th: parse arg th; return word('th st nd rd', 1+(th//10) *(th//100%10\==1) *(th//10<4))
/*──────────────────────────────────────────────────────────────────────────────────────*/
CW_frac: procedure; parse arg p '/' q .; say
if q=='' then do; p= 83116; q= 51639; end
n= rle2dec( frac2cf(p q) ); @CWS= 'the Calkin─Wilf sequence'
say 'for ' p"/"q', the element number for' @CWS "is: " commas(n)th(n)
if length(n)<10 then return
say; say 'The above number has ' commas(length(n)) " decimal digits."
return
/*──────────────────────────────────────────────────────────────────────────────────────*/
CW_term: procedure; parse arg z; dd= 1; nn= 0
do z
parse value dd dd*(2*(nn%dd)+1)-nn with nn dd
end /*z*/
return nn'/'dd
/*──────────────────────────────────────────────────────────────────────────────────────*/
CW_terms: $=; if LO\==0 then do j=LO to HI; $= $ CW_term(j)','
end /*j*/
if $=='' then return
say 'Calkin─Wilf sequence terms for ' commas(LO) " ──► " commas(HI) ' are:'
say strip( strip($), 'T', ",")
return
/*──────────────────────────────────────────────────────────────────────────────────────*/
frac2cf: procedure; parse arg p q; if q=='' then return p; cf= p % q; m= q
p= p - cf*q; n= p; if p==0 then return cf
do k=1 until n==0; @.k= m % n
m= m - @.k * n; parse value n m with m n /*swap N M*/
end /*k*/
/*for inverse Calkin─Wilf, K must be even.*/
if k//2 then do; @.k= @.k - 1; k= k + 1; @.k= 1; end
do k=1 for k; cf= cf @.k; end /*k*/
return cf
/*──────────────────────────────────────────────────────────────────────────────────────*/
rle2dec: procedure; parse arg f1 rle; obin= copies(1, f1)
do until rle==''; parse var rle f0 f1 rle
obin= copies(1, f1)copies(0, f0)obin
end /*until*/
return x2d( b2x(obin) ) /*RLE2DEC: Run Length Encoding ──► decimal*/
- output when using the default inputs:
Calkin─Wilf sequence terms for 1 ──► 20 are: 1/1, 1/2, 2/1, 1/3, 3/2, 2/3, 3/1, 1/4, 4/3, 3/5, 5/2, 2/5, 5/3, 3/4, 4/1, 1/5, 5/4, 4/7, 7/3, 3/8 for 83116/51639, the element number for the Calkin─Wilf sequence is: 123,456,789th
RPL
≪ { } SWAP WHILE DUP REPEAT ROT OVER IDIV2 4 ROLL ROT + SWAP END ROT DROP2 ≫ 'CONTFRAC' STO ≪ {1} WHILE DUP2 SIZE > REPEAT DUP DUP SIZE GET DUP IP R→I 2 * 1 + SWAP - INV EVAL + END NIP ≫ ≫ 'CWILF' STO ≪ CONTFRAC DUP SIZE IF DUP MOD THEN DROP ELSE DUP2 GET 1 - PUT 1 + END 1 → frac pow2 ≪ 0 1 frac SIZE FOR j frac j GET WHILE DUP REPEAT IF j 2 MOD THEN SWAP pow2 + SWAP END 2 'pow2' STO* 1 - END DROP NEXT ≫ ≫ 'CWPOS' STO
20 CWILF 83116 51639 CWPOS
- Output:
2: {1 '1/2' 2 '1/3' '3/2' '2/3' 3 '1/4' '4/3' '3/5' '5/2' '2/5' '5/3' '3/4' 4 '1/5' '5/4' '4/7' '7/3' '3/8'} 1: 123456789
Ruby
cw = Enumerator.new do |y|
y << a = 1.to_r
loop { y << a = 1/(2*a.floor + 1 - a) }
end
def term_num(rat)
num, den, res, pwr, dig = rat.numerator, rat.denominator, 0, 0, 1
while den > 0
num, (digit, den) = den, num.divmod(den)
digit.times do
res |= dig << pwr
pwr += 1
end
dig ^= 1
end
res
end
puts cw.take(20).join(", ")
puts term_num (83116/51639r)
1/1, 1/2, 2/1, 1/3, 3/2, 2/3, 3/1, 1/4, 4/3, 3/5, 5/2, 2/5, 5/3, 3/4, 4/1, 1/5, 5/4, 4/7, 7/3, 3/8 123456789
Rust
// [dependencies]
// num = "0.3"
use num::rational::Rational;
fn calkin_wilf_next(term: &Rational) -> Rational {
Rational::from_integer(1) / (Rational::from_integer(2) * term.floor() + 1 - term)
}
fn continued_fraction(r: &Rational) -> Vec<isize> {
let mut a = *r.numer();
let mut b = *r.denom();
let mut result = Vec::new();
loop {
let (q, r) = num::integer::div_rem(a, b);
result.push(q);
a = b;
b = r;
if a == 1 {
break;
}
}
let len = result.len();
if len != 0 && len % 2 == 0 {
result[len - 1] -= 1;
result.push(1);
}
result
}
fn term_number(r: &Rational) -> usize {
let mut result: usize = 0;
let mut d: usize = 1;
let mut p: usize = 0;
for n in continued_fraction(r) {
for _ in 0..n {
result |= d << p;
p += 1;
}
d ^= 1;
}
result
}
fn main() {
println!("First 20 terms of the Calkin-Wilf sequence are:");
let mut term = Rational::from_integer(1);
for i in 1..=20 {
println!("{:2}: {}", i, term);
term = calkin_wilf_next(&term);
}
let r = Rational::new(83116, 51639);
println!("{} is the {}th term of the sequence.", r, term_number(&r));
}
- Output:
First 20 terms of the Calkin-Wilf sequence are: 1: 1 2: 1/2 3: 2 4: 1/3 5: 3/2 6: 2/3 7: 3 8: 1/4 9: 4/3 10: 3/5 11: 5/2 12: 2/5 13: 5/3 14: 3/4 15: 4 16: 1/5 17: 5/4 18: 4/7 19: 7/3 20: 3/8 83116/51639 is the 123456789th term of the sequence.
Scheme
Continued Fraction support
; Create a terminating Continued Fraction generator for the given rational number.
; Returns one term per call; returns #f when no more terms remaining.
(define make-continued-fraction-gen
(lambda (rat)
(let ((num (numerator rat)) (den (denominator rat)))
(lambda ()
(if (= den 0)
#f
(let ((ret (quotient num den))
(rem (modulo num den)))
(set! num den)
(set! den rem)
ret))))))
; Return the continued fraction representation of a rational number as a list of terms.
(define rat->cf-list
(lambda (rat)
(let ((cf (make-continued-fraction-gen rat))
(lst '()))
(let loop ((term (cf)))
(when term
(set! lst (append lst (list term)))
(loop (cf))))
lst)))
; Enforce the length of the given continued fraction list to be odd.
; Changes the list in situ (if needed), and returns its possibly changed value.
(define continued-fraction-list-enforce-odd-length!
(lambda (cf)
(when (even? (length cf))
(let ((cf-last-cons (list-tail cf (1- (length cf)))))
(set-car! cf-last-cons (1- (car cf-last-cons)))
(set-cdr! cf-last-cons (cons 1 '()))))
cf))
Calkin-Wilf sequence
; Create a Calkin-Wilf sequence generator.
(define make-calkin-wilf-gen
(lambda ()
(let ((an 1))
(lambda ()
(let ((ret an))
(set! an (/ 1 (+ (* 2 (floor an)) 1 (- an))))
ret)))))
; Return the position in the Calkin-Wilf sequence of the given rational number.
(define calkin-wilf-position
(lambda (rat)
; Run-length encodes binary value. Assumes first run is 1's. Args: initial value,
; starting place value (a power of 2), and list of run lengths (list must be odd length).
(define encode-list-of-runs
(lambda (value placeval lstruns)
; Encode a single run in a binary value. Args: initial value, bit value (0 or 1),
; starting place value (a power of 2), number of places (bits) to encode.
; Returns multiple values: the encoded value, and the new place value.
(define encode-run
(lambda (value bitval placeval places)
(if (= places 1)
(values (+ value (* bitval placeval)) (* 2 placeval))
(encode-run (+ value (* bitval placeval)) bitval (* 2 placeval) (1- places)))))
; Loop through the list of runs two at a time. If list of length 1, do a final
; '1'-bit encode and return the value. Otherwise, do a '1'-bit then '0'-bit encode,
; and recurse to do the next two runs.
(let-values (((value-1 placeval-1) (encode-run value 1 placeval (car lstruns))))
(if (= 1 (length lstruns))
value-1
(let-values (((value-2 placeval-2) (encode-run value-1 0 placeval-1 (cadr lstruns))))
(encode-list-of-runs value-2 placeval-2 (cddr lstruns)))))))
; Return the run-length binary encoding from the odd-length Calkin-Wilf sequence of the
; given rational number. This is equal to the number's position in the sequence.
(encode-list-of-runs 0 1 (continued-fraction-list-enforce-odd-length! (rat->cf-list rat)))))
The Task
(let ((count 20)
(cw (make-calkin-wilf-gen)))
(printf "~%First ~a terms of the Calkin-Wilf sequence:~%" count)
(do ((num 1 (1+ num)))
((> num count))
(printf "~2d : ~a~%" num (cw))))
(printf "~%Positions in Calkin-Wilf sequence of given numbers:~%")
(let ((num 9/4))
(printf "~a @ ~a~%" num (calkin-wilf-position num)))
(let ((num 83116/51639))
(printf "~a @ ~a~%" num (calkin-wilf-position num)))
- Output:
First 20 terms of the Calkin-Wilf sequence: 1 : 1 2 : 1/2 3 : 2 4 : 1/3 5 : 3/2 6 : 2/3 7 : 3 8 : 1/4 9 : 4/3 10 : 3/5 11 : 5/2 12 : 2/5 13 : 5/3 14 : 3/4 15 : 4 16 : 1/5 17 : 5/4 18 : 4/7 19 : 7/3 20 : 3/8 Positions in Calkin-Wilf sequence of given numbers: 9/4 @ 35 83116/51639 @ 123456789
Sidef
func calkin_wilf(n) is cached {
return 1 if (n == 1)
1/(2*floor(__FUNC__(n-1)) + 1 - __FUNC__(n-1))
}
func r2cw(r) {
var cfrac = r.as_cfrac
cfrac.len.is_odd || return nil
Num(cfrac.flip.map_kv {|k,v| (k.is_odd ? '0' : '1') * v }.join, 2)
}
with (20) {|n|
say "First #{n} terms of the Calkin-Wilf sequence:"
say calkin_wilf.map(1..n)
}
with (83116/51639) {|r|
say ("\n#{r.as_rat} is at index: ", r2cw(r))
}
- Output:
First 20 terms of the Calkin-Wilf sequence: [1, 1/2, 2, 1/3, 3/2, 2/3, 3, 1/4, 4/3, 3/5, 5/2, 2/5, 5/3, 3/4, 4, 1/5, 5/4, 4/7, 7/3, 3/8] 83116/51639 is at index: 123456789
V (Vlang)
s.
import math.fractions
import math
import strconv
fn calkin_wilf(n int) []fractions.Fraction {
mut cw := []fractions.Fraction{len: n+1}
cw[0] = fractions.fraction(1, 1)
one := fractions.fraction(1, 1)
two := fractions.fraction(2, 1)
for i in 1..n {
mut t := cw[i-1]
mut f := t.f64()
f = math.floor(f)
t = fractions.approximate(f)
t*=two
t-= cw[i-1]
t+=one
t=t.reciprocal()
cw[i] = t
}
return cw
}
fn to_continued(r fractions.Fraction) []int {
idx := r.str().index('/') or {0}
mut a := r.str()[..idx].i64()
mut b := r.str()[idx+1..].i64()
mut res := []int{}
for {
res << int(a/b)
t := a % b
a, b = b, t
if a == 1 {
break
}
}
le := res.len
if le%2 == 0 { // ensure always odd
res[le-1]--
res << 1
}
return res
}
fn get_term_number(cf []int) ?int {
mut b := ""
mut d := "1"
for n in cf {
b = d.repeat(n)+b
if d == "1" {
d = "0"
} else {
d = "1"
}
}
i := strconv.parse_int(b, 2, 64)?
return int(i)
}
fn commatize(n int) string {
mut s := "$n"
if n < 0 {
s = s[1..]
}
le := s.len
for i := le - 3; i >= 1; i -= 3 {
s = s[0..i] + "," + s[i..]
}
if n >= 0 {
return s
}
return "-" + s
}
fn main() {
cw := calkin_wilf(20)
println("The first 20 terms of the Calkin-Wilf sequnence are:")
for i := 1; i <= 20; i++ {
println("${i:2}: ${cw[i-1]}")
}
println('')
r := fractions.fraction(83116, 51639)
cf := to_continued(r)
tn := get_term_number(cf) or {0}
println("$r is the ${commatize(tn)}th term of the sequence.")
}
- Output:
The first 20 terms of the Calkin-Wilf sequnence are: 1: 1/1 2: 1/2 3: 2/1 4: 1/3 5: 3/2 6: 2/3 7: 3/1 8: 1/4 9: 4/3 10: 3/5 11: 5/2 12: 2/5 13: 5/3 14: 3/4 15: 4/1 16: 1/5 17: 5/4 18: 4/7 19: 7/3 20: 3/8 83116/51639 is the 123,456,789th term of the sequence.
Wren
import "./rat" for Rat
import "./fmt" for Fmt, Conv
var calkinWilf = Fn.new { |n|
var cw = List.filled(n, null)
cw[0] = Rat.one
for (i in 1...n) {
var t = cw[i-1].floor * 2 - cw[i-1] + 1
cw[i] = Rat.one / t
}
return cw
}
var toContinued = Fn.new { |r|
var a = r.num
var b = r.den
var res = []
while (true) {
res.add((a/b).floor)
var t = a % b
a = b
b = t
if (a == 1) break
}
if (res.count%2 == 0) { // ensure always odd
res[-1] = res[-1] - 1
res.add(1)
}
return res
}
var getTermNumber = Fn.new { |cf|
var b = ""
var d = "1"
for (n in cf) {
b = (d * n) + b
d = (d == "1") ? "0" : "1"
}
return Conv.atoi(b, 2)
}
var cw = calkinWilf.call(20)
System.print("The first 20 terms of the Calkin-Wilf sequence are:")
Rat.showAsInt = true
for (i in 1..20) Fmt.print("$2d: $s", i, cw[i-1])
System.print()
var r = Rat.new(83116, 51639)
var cf = toContinued.call(r)
var tn = getTermNumber.call(cf)
Fmt.print("$s is the $,r term of the sequence.", r, tn)
- Output:
The first 20 terms of the Calkin-Wilf sequence are: 1: 1 2: 1/2 3: 2 4: 1/3 5: 3/2 6: 2/3 7: 3 8: 1/4 9: 4/3 10: 3/5 11: 5/2 12: 2/5 13: 5/3 14: 3/4 15: 4 16: 1/5 17: 5/4 18: 4/7 19: 7/3 20: 3/8 83116/51639 is the 123,456,789th term of the sequence.
- Programming Tasks
- Solutions by Programming Task
- 11l
- ALGOL 68
- AppleScript
- Arturo
- BQN
- Bracmat
- C++
- Boost
- EasyLang
- EDSAC order code
- F Sharp
- Factor
- Forth
- FreeBASIC
- Fōrmulæ
- Go
- Haskell
- J
- Java
- Jq
- Julia
- Little Man Computer
- Mathematica
- Wolfram Language
- Maxima
- Nim
- PARI/GP
- Pascal
- Perl
- Ntheory
- Phix
- Prolog
- Python
- Quackery
- Raku
- REXX
- RPL
- Ruby
- Rust
- Scheme
- Sidef
- V (Vlang)
- Wren
- Wren-rat
- Wren-fmt