Calkin-Wilf sequence: Difference between revisions
m (change appearance of a improper fraction.) |
m (→{{header|REXX}}: changed some comments.) |
||
Line 996: | Line 996: | ||
CW_frac: procedure; parse arg p '/' q .; say |
CW_frac: procedure; parse arg p '/' q .; say |
||
if q=='' then do; p= 83116; q= 51639; end |
if q=='' then do; p= 83116; q= 51639; end |
||
n= rle2dec( |
n= rle2dec( frac2cf(p q) ); @CWS= 'the Calkin─Wilf sequence' |
||
say 'for ' p"/"q', the element number for' @CWS "is: " commas(n)th(n) |
say 'for ' p"/"q', the element number for' @CWS "is: " commas(n)th(n) |
||
return |
return |
||
Line 1,013: | Line 1,013: | ||
return |
return |
||
/*──────────────────────────────────────────────────────────────────────────────────────*/ |
/*──────────────────────────────────────────────────────────────────────────────────────*/ |
||
frac2cf: procedure; parse arg p q; a0= p % q; cfrac= a0; a.0= 0; m= q |
|||
p= p - a0*q; |
p= p - a0*q; n= p; if p==0 then return cFrac |
||
do k=1 until n==0; |
do k=1 until n==0; a.k= m % n |
||
m= m - a.k * n; parse value n m with m n |
m= m - a.k * n; parse value n m with m n |
||
end /*k*/ |
end /*k*/ |
||
Line 1,021: | Line 1,021: | ||
if k//2 then do; a.k= a.k - 1; a.k= 1; k= k + 1; end |
if k//2 then do; a.k= a.k - 1; a.k= 1; k= k + 1; end |
||
a.0= k |
a.0= k |
||
do k=1 for a.0; |
do k=1 for a.0; cFrac= cFrac a.k; end /*k*/ |
||
return |
return cFrac |
||
/*──────────────────────────────────────────────────────────────────────────────────────*/ |
/*──────────────────────────────────────────────────────────────────────────────────────*/ |
||
rle2dec: procedure; parse arg rle; obin=; parse var rle f1 rle |
rle2dec: procedure; parse arg rle; obin=; parse var rle f1 rle |
||
Line 1,030: | Line 1,030: | ||
parse var rle f1 rle; obin= copies(1, f1)obin |
parse var rle f1 rle; obin= copies(1, f1)obin |
||
end /*until*/ |
end /*until*/ |
||
return x2d( b2x(obin) )</lang> |
return x2d( b2x(obin) ) /*RLE2DEC: Run Length Encoding ──► decimal.*/</lang> |
||
{{out|output|text= when using the default inputs:}} |
{{out|output|text= when using the default inputs:}} |
||
<pre> |
<pre> |
Revision as of 07:25, 17 April 2021
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.
- See also
- Wikipedia entry: Calkin-Wilf tree
- Continued fraction
- Continued fraction/Arithmetic/Construct from rational number
AppleScript
<lang applescript>-- Return the first n terms of the sequence. Tree generation. Faster for this purpose. on CalkinWilfSequence(n)
script o property sequence : Template: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</lang>
- Output:
<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 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"</lang>
C++
<lang cpp>#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";
}</lang>
- 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.
F#
The Function
<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 </lang>
The Tasks
- first 20
<lang fsharp> cW |> Seq.take 20 |> Seq.iter(fun(n,g)->printf "%d/%d " n g);printfn "" </lang>
- 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
<lang fsharp> printfn "%d" (1+Seq.findIndex(fun n->n=(83116,51639)) cW) </lang>
- Output:
123456789
Factor
<lang 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</lang>
- 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.
FreeBASIC
Uses the code from Greatest common divisor#FreeBASIC as an include.
<lang freebasic>#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."</lang>
- 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.
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. <lang go>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))
}</lang>
- 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
<lang 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)</lang>
- 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 <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 </lang>
Julia
<lang 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.")
</lang>
- 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.
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.
<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."</lang>
- 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.
Perl
<lang perl>use strict; use warnings; use feature 'say';
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);</lang>
- 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
<lang Phix>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 to_continued_fraction(sequence r)
integer {a,b} = r sequence res = {} while true do res &= floor(a/b) {a, b} = {b, remainder(a,b)} if a=1 then exit end if end while return res
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 return bits_to_int(b)
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 -- 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 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)})</lang>
- 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.
Python
<lang 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.")</lang>
- 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.
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.
<lang perl6>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: '/';
}</lang>
- 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. <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). */ 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) then call CW_frac /*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) 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; a0= p % q; cfrac= a0; a.0= 0; m= q
p= p - a0*q; n= p; if p==0 then return cFrac do k=1 until n==0; a.k= m % n m= m - a.k * n; parse value n m with m n end /*k*/ /*for inverse(c_w), K must be even*/ if k//2 then do; a.k= a.k - 1; a.k= 1; k= k + 1; end a.0= k do k=1 for a.0; cFrac= cFrac a.k; end /*k*/ return cFrac
/*──────────────────────────────────────────────────────────────────────────────────────*/ rle2dec: procedure; parse arg rle; obin=; parse var rle f1 rle
obin= copies(1, f1) /*F1 may be zero*/ do until rle==; parse var rle f0 rle /*zeros next. */ obin= copies(0, f0)obin; if rle== then leave parse var rle f1 rle; obin= copies(1, f1)obin end /*until*/ return x2d( b2x(obin) ) /*RLE2DEC: Run Length Encoding ──► decimal.*/</lang>
- 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
Rust
<lang 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));
}</lang>
- 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.
Wren
<lang ecmascript>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)</lang>
- 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.