Sturmian word: Difference between revisions

From Rosetta Code
Content added Content deleted
m (→‎{{header|jq}}: shorten)
 
Line 255: Line 255:
def fibWord($n):
def fibWord($n):
{ sn1: "0", sn: "01" }
{ sn1: "0", sn: "01" }
| reduce range (2; 1+$n) as $i (.;
| reduce range (2; 1+$n) as $i (.; {sn: (.sn+.sn1), sn1: .sn})
.sn as $tmp
| .sn += .sn1
| .sn1 = $tmp)
| .sn;
| .sn;


Latest revision as of 07:36, 18 April 2024

Sturmian word is a draft programming task. It is not yet considered ready to be promoted as a complete task, for reasons that should be found in its talk page.

A Sturmian word is a binary sequence, finite or infinite, that makes up the cutting sequence for a positive real number x, as shown in the picture.

Example Sturmian word when x = 0.618..., the golden ratio.

The Sturmian word can be computed thus as an algorithm:

  • If , then it is the inverse of the Sturmian word for . So we have reduced to the case of .
  • Iterate over
  • If is an integer, then the program terminates. Else, if , then the program outputs 0, else, it outputs 10.

The problem:

  • Given a positive rational number , specified by two positive integers , output its entire Sturmian word.
  • Given a quadratic real number , specified by integers , where is not a perfect square, output the first letters of Sturmian words when given a positive number .

(If the programming language can represent infinite data structures, then that works too.)

A simple check is to do this for the inverse golden ratio , that is, , which would just output the Fibonacci word.

Stretch goal: calculate the Sturmian word for other kinds of definable real numbers, such as cubic roots.

The key difficulty is accurately calculating for large . Floating point arithmetic would lose precision. One can either do this simply by directly searching for some integer such that , or by more trickly methods, such as the continued fraction approach.

First calculate the continued fraction convergents to . Let be a convergent to , such that , then since the convergent sequence is the best rational approximant for denominators up to that point, we know for sure that, if we write out , the sequence would stride right across the gap . Thus, we can take the largest such that , and we would know for sure that .

In summary, where is the first continued fraction approximant to with a denominator


ALGOL 68

Translation of: Python
BEGIN # Sturmian word - based on the Python sample                    #

    PROC sturmian word = ( INT m, n )STRING:
         IF m > n THEN
            STRING sturmian := sturmian word( n, m );
            FOR s pos FROM LWB sturmian TO UPB sturmian DO
                CHAR c = sturmian[ s pos ];
                sturmian[ s pos ] := IF c = "0" THEN "1" ELIF c = "1" THEN "0" ELSE c FI
            OD;
            sturmian
         ELSE
            STRING sturmian := "";
            INT current floor := 0;
            FOR k WHILE ( k * m ) MOD n /= 0 DO
                INT previous floor = current floor;
                current floor := ( k * m ) OVER n;
                IF previous floor = current floor THEN
                    sturmian +:= "0"
                ELSE
                    sturmian +:= "10"
                FI
            OD;
            sturmian
         FI # sturmian word # ;

    # Checking that it works on the finite Fibonacci word:            #

    PROC fibonacci word = ( INT n )STRING:
         BEGIN
            STRING sn1 := "0", sn := "01";
            FOR i FROM 2 TO n DO
                STRING tmp = sn;
                sn +:= sn1;
                sn1 := tmp
            OD;
            sn
         END # fibonacci word # ;

    STRING f word   = fibonacci word( 10 );
    STRING sturmian = sturmian word( 13, 21 );
    IF f word[ : UPB sturmian[ AT 1 ] AT 1 ] = sturmian THEN
        print( ( sturmian ))
    ELSE
        print( ( "*** unexpected result: (", sturmian, ") (", f word, ")" ) )
    FI

END
Output:
01001010010010100101001001010010

BASIC

BASIC256

Translation of: FreeBASIC
fib = fibWord(10)
sturmian = sturmian_word(13, 21)
if left(fib, length(sturmian)) = sturmian then print sturmian
end

function invert(cadena)
	for i = 1 to length(cadena)
		b = mid(cadena, i, 1)
		if b = "0" then
			inverted += "1"
		else
			if b = "1" then
				inverted += "0"
			end if
		end if
	next
	return inverted
end function

function sturmian_word(m, n)
	sturmian = ""
	if m > n then return invert(sturmian_word(n, m))

	k = 1
	while True
		current_floor = int(k * m / n)
		previous_floor = int((k - 1) * m / n)
		if k * m mod n = 0 then exit while
		if previous_floor = current_floor then
			sturmian += "0"
		else
			sturmian += "10"
		end if
		k += 1
	end while
	return sturmian
end function

function fibWord(n)
	Sn_1 = "0"
	Sn = "01"

	for i = 2 to n
		tmp = Sn
		Sn += Sn_1
		Sn_1 = tmp
	next
	return Sn
end function
Output:
01001010010010100101001001010010

FreeBASIC

Translation of: Python
Function invert(cadena As String) As String
    Dim As String inverted, b
    For i As Integer = 1 To Len(cadena)
        b = Mid(cadena, i, 1)
        inverted += Iif(b = "0", "1", "0")
    Next
    Return inverted
End Function

Function sturmian_word(m As Integer, n As Integer) As String
    Dim sturmian As String
    If m > n Then Return invert(sturmian_word(n, m))
    
    Dim As Integer k = 1
    Do
        Dim As Integer current_floor = Int(k * m / n)
        Dim As Integer previous_floor = Int((k - 1) * m / n)
        If k * m Mod n = 0 Then Exit Do
        sturmian += Iif(previous_floor = current_floor, "0", "10")
        k += 1
    Loop
    Return sturmian
End Function

Function fibWord(n As Integer) As String
    Dim As String Sn_1 = "0"
    Dim As String Sn = "01"
    
    For i As Integer = 2 To n
        Dim As String tmp = Sn
        Sn += Sn_1
        Sn_1 = tmp
    Next
    Return Sn
End Function

Dim As String fib = fibWord(10)
Dim As String sturmian = sturmian_word(13, 21)
If Left(fib, Len(sturmian)) = sturmian Then Print sturmian

Sleep
Output:
01001010010010100101001001010010

Yabasic

Translation of: FreeBASIC
fib$ = fibword$(10)
sturmian$ = sturmian_word$(13,21)
if left$(fib$,len(sturmian$)) = sturmian$  print sturmian$
end
sub invert$(cadena$)
    for i = 1 to len(cadena$)
        b$ = mid$(cadena$,i,1)
        if b$ = "0" then
            inverted$ = inverted$ +"1"
        elsif b$ = "1" then
            inverted$ = inverted$ +"0"
        fi
    next i
    return inverted$
end sub
sub sturmian_word$(m,n)
    sturmian$ = ""
    if m > n  return invert$(sturmian_word(n,m))
    k = 1
    while true
        current_floor = int(k*m/n)
        previous_floor = int((k-1)*m/n)
        if mod(k*m, n) = 0  break
        if previous_floor = current_floor then
            sturmian$ = sturmian$ +"0"
        else
            sturmian$ = sturmian$ +"10"
        fi
        k = k +1
    end while
    return sturmian$
end sub
sub fibword$(n)
    sn_1$ = "0"
    sn$ = "01"
    for i = 2 to n
        tmp$ = sn$
        sn$ = sn$ + sn_1$
        sn_1$ = tmp$
    next i
    return sn$
end sub
Output:
01001010010010100101001001010010

jq

Adapted from Wren

Works with jq, the C implementation of jq

Works with gojq, the Go implementation of jq 

## Generic functions
def assertEq($a; $b):
  if $a == $b then .
  else . as $in
  | "assertion violation: \($a) != \($b)" | debug
  | $in
  end;

# Emit a stream of characters
def chars: explode[] | [.] | implode;

# Input: a string of "0"s and "1"s
def flip:
  {"1": "0", "0": "1"} as $flip
  | reduce chars as $c (""; . + $flip[$c]);

def fibWord($n):
  { sn1: "0", sn: "01" }
  | reduce range (2; 1+$n) as $i (.; {sn: (.sn+.sn1), sn1: .sn})
  | .sn;

### Sturmian words

def SturmianWordRat($m; $n):
  if $m > $n
  then SturmianWordRat($n; $m) | flip
  else {sturmian: "", k: 1}
  | until (.k * $m % $n == 0;
      ( (.k * $m / $n)|floor) as $curr
      | (((.k - 1) * $m / $n)|floor) as $prev
      | .sturmian += (if $prev == $curr then "0" else "10" end)
      | .k += 1 )
  | .sturmian
  end;

def SturmianWordQuad($a; $b; $m; $n; $k):
  def fraction: . - trunc;
    { p: [0, 1],
      q: [1, 0],
     rem: ((($a|sqrt) * $b + $m) / $n)
    }
    | reduce range(1; 1+$k) as $i (.;
          (.rem|trunc) as $whole
          | (.rem|fraction) as $frac
          | ($whole * .p[-1] + .p[-2]) as $pn
          | ($whole * .q[-1] + .q[-2]) as $qn
          | .p += [$pn]
          | .q += [$qn]
          | .rem = (1 / $frac) )
    | SturmianWordRat(.p[-1]; .q[-1]) ;

### The tasks

def tasks:
  fibWord(10) as $fib
  | SturmianWordRat(13; 21) as $sturmian
  | SturmianWordQuad(5; 1; -1; 2; 8) as $sturmian2

  | assertEq($fib[0:$sturmian|length]; $sturmian)
  | assertEq($sturmian; $sturmian2)
  
  | "\($sturmian) from rational number 13/21",
    "\($sturmian2) from quadratic number (√5 - 1)/2 (k = 8)";

tasks
Output:
01001010010010100101001001010010 from rational number 13/21
01001010010010100101001001010010 from quadratic number (√5 - 1)/2 (k = 8)

Julia

Translation of: Phix
function sturmian_word(m, n)
    m > n && return replace(sturmian_word(n, m), '0' => '1', '1' => '0')
    res = ""
    k, prev = 1, 0
    while rem(k * m, n) > 0
        curr = (k * m) ÷ n
        res *= prev == curr ? "0" : "10"
        prev = curr
        k += 1
    end
    return res
end

function fibWord(n)
    Sn_1, Sn, tmp = "0", "01", ""
    for _ in 2:n 
        tmp = Sn
        Sn *= Sn_1
        Sn_1 = tmp
    end
    return Sn
end

const fib = fibWord(7)
const sturmian = sturmian_word(13, 21)
@assert fib[begin:length(sturmian)] == sturmian
println(" $sturmian <== 13/21")

""" return the kth convergent """
function cfck(a, b, m, n, k)
    p = [0, 1]
    q = [1, 0]
    r = (sqrt(a) * b + m) / n
    for _ in 1:k
        whole = Int(trunc(r))
        pn = whole * p[end] + p[end-1]
        qn = whole * q[end] + q[end-1]
        push!(p, pn)
        push!(q, qn)
        r = 1/(r - whole)
    end
    return [p[end], q[end]]
end

println(" $(sturmian_word(cfck(5, 1, -1, 2, 8)...)) <== 1/phi (8th convergent golden ratio)")
Output:
 01001010010010100101001001010010 <== 13/21
 01001010010010100101001001010010 <== 1/phi (8th convergent golden ratio)

Phix

Rational part translated from Python, quadratic part a modified copy of the Continued fraction convergents task. 

with javascript_semantics
function sturmian_word(integer m, n)
    if m > n then
        return sq_sub('0'+'1',sturmian_word(n,m))
    end if
    string res = ""
    integer k = 1, prev = 0
    while remainder(k*m,n) do
        integer curr = floor(k*m/n)
        res &= iff(prev=curr?"0":"10")
        prev = curr
        k += 1
    end while
    return res
end function

function fibWord(integer n)
    string Sn_1 = "0",
             Sn = "01",
            tmp = ""
    for i=2 to n do
        tmp = Sn
        Sn &= Sn_1
        Sn_1 = tmp
    end for
    return Sn
end function

string fib = fibWord(7),
  sturmian = sturmian_word(13, 21)
assert(fib[1..length(sturmian)] == sturmian)
printf(1," %s <== 13/21\n",sturmian)

function cfck(integer a, b, m, n, k)
    -- return the kth convergent
    sequence p = {0, 1},
             q = {1, 0}
    atom rem = (sqrt(a)*b + m) / n
    for i=1 to k do
        integer whole = trunc(rem),
                pn = whole * p[-1] + p[-2],
                qn = whole * q[-1] + q[-2]
        p &= pn
        q &= qn
        rem = 1/(rem-whole)
    end for
    return {p[$],q[$]}
end function

integer {m,n} = cfck(5, 1, -1, 2, 8) -- (inverse golden ratio)
printf(1," %s <== 1/phi (8th convergent)\n",sturmian_word(m,n))
Output:
 01001010010010100101001001010010 <== 13/21
 01001010010010100101001001010010 <== 1/phi (8th convergent)

Python

For rational numbers:

def sturmian_word(m, n):
    sturmian = ""
    def invert(string):
      return ''.join(list(map(lambda b: {"0":"1", "1":"0"}[b], string)))
    if m > n:
      return invert(sturmian_word(n, m))

    k = 1
    while True:
        current_floor = int(k * m / n)
        previous_floor = int((k - 1) * m / n)
        if k * m % n == 0:
            break
        if previous_floor == current_floor:
            sturmian += "0"
        else:
            sturmian += "10"
        k += 1
    return sturmian

Checking that it works on the finite Fibonacci word:

def fibWord(n):
    Sn_1 = "0"
    Sn = "01"
    tmp = ""
    for i in range(2, n + 1):
        tmp = Sn
        Sn += Sn_1
        Sn_1 = tmp
    return Sn

fib = fibWord(10)
sturmian = sturmian_word(13, 21)
assert fib[:len(sturmian)] == sturmian
print(sturmian)

# Output:
# 01001010010010100101001001010010

Raku

Translation of: Julia
# 20240215 Raku programming solution

sub chenfoxlyndonfactorization(Str $s) {
 sub sturmian-word($m, $n) {
   return sturmian-word($n, $m).trans('0' => '1', '1' => '0') if $m > $n; 
   my ($res, $k, $prev) = '', 1, 0;
   while ($k * $m) % $n > 0 {
      my $curr = ($k * $m) div $n;
      $res ~= $prev == $curr ?? '0' !! '10';
      $prev = $curr;
      $k++;
   }
   return $res;
}

sub fib-word($n) {
   my ($Sn_1, $Sn) = '0', '01';
   for 2..$n { ($Sn, $Sn_1) = ($Sn~$Sn_1, $Sn) }
   return $Sn;
}

my $fib = fib-word(7);
my $sturmian = sturmian-word(13, 21);
say "{$sturmian} <== 13/21" if $fib.substr(0, $sturmian.chars) eq $sturmian;

sub cfck($a, $b, $m, $n, $k) {
   my (@p, @q) := [0, 1], [1, 0]; 
   my $r = (sqrt($a) * $b + $m) / $n;
   for ^$k {
      my $whole = $r.Int;
      my ($pn, $qn) = $whole * @p[*-1] + @p[*-2], $whole * @q[*-1] + @q[*-2];
      @p.push($pn);
      @q.push($qn);
      $r = 1/($r - $whole);
   }
   return [@p[*-1], @q[*-1]];
}

my $cfck-result = cfck(5, 1, -1, 2, 8);
say "{sturmian-word($cfck-result[0], $cfck-result[1])} <== 1/phi (8th convergent golden ratio)";
Output:

Same as Julia example.

You may Attempt This Online!

Wren

Library: Wren-assert

The 'rational number' function is a translation of the Python entry.

The 'quadratic number' function is a modified version of the one used in the Continued fraction convergents task. 

import "./assert" for Assert

var SturmianWordRat = Fn.new { |m, n|
    if (m > n) return SturmianWordRat.call(n, m).reduce("") { |acc, c| 
        return acc + (c == "0" ? "1" : "0") 
    }
    var sturmian = ""
    var k = 1
    while (k * m % n != 0) {
        var currFloor = (k * m / n).floor
        var prevFloor = ((k - 1) * m / n).floor
        sturmian = sturmian + (prevFloor == currFloor ? "0" : "10")
        k  = k + 1
    }
    return sturmian
}

var SturmianWordQuad = Fn.new { |a, b, m, n, k|
    var p = [0, 1]
    var q = [1, 0]
    var rem = (a.sqrt * b + m) / n
    for (i in 1..k) {
        var whole = rem.truncate
        var frac  = rem.fraction
        var pn = whole * p[-1] + p[-2]
        var qn = whole * q[-1] + q[-2]
        p.add(pn)
        q.add(qn)
        rem = 1 / frac
    }
    return SturmianWordRat.call(p[-1], q[-1])
}

var fibWord = Fn.new { |n|
    var sn1 = "0"
    var sn  = "01"
    for (i in 2..n) {
        var tmp = sn
        sn = sn + sn1
        sn1 = tmp
    }
    return sn
}

var fib = fibWord.call(10)
var sturmian = SturmianWordRat.call(13, 21)
Assert.equal(fib[0...sturmian.count], sturmian)
System.print("%(sturmian) from rational number 13/21")
var sturmian2 = SturmianWordQuad.call(5, 1, -1, 2, 8)
Assert.equal(sturmian, sturmian2)
System.print("%(sturmian2) from quadratic number (√5 - 1)/2 (k = 8)")
Output:
01001010010010100101001001010010 from rational number 13/21
01001010010010100101001001010010 from quadratic number (√5 - 1)/2 (k = 8)
Retrieved from "https://rosettacode.org/wiki/Sturmian_word?oldid=363774"