Sturmian word: Difference between revisions

m
(Added Wren)
m (→‎{{header|jq}}: shorten)
 
(12 intermediate revisions by 7 users not shown)
Line 1:
{{draft task}}
 
A [[wp:Sturmian_word|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.
Line 15:
(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 <math>\frac{\sqrt{5}-1}{2}</math>, that is, <math>a=5, b=1, m = -1, n = 2 </math>, which would just output the [[Fibonacci word]].
 
Stretch goal: calculate the Sturmian word for other kinds of definable real numbers, such as cubic roots.
Line 23:
First calculate the [[continued fraction convergents]] to <math>\sqrt a</math>. Let <math>\frac mn </math> be a convergent to <math>\sqrt a</math>, such that <math>n \geq k</math>, 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 <math>\frac{0}{k}, \frac{1}{k}, \dots</math>, the sequence would stride right across the gap <math>(m/n, 2x - m/n)</math>. Thus, we can take the largest <math>l</math> such that <math>l/k \leq m/n</math>, and we would know for sure that <math>l = floor(k\sqrt a)</math>.
 
In summary, <math display="block">floor(k\sqrt a) = floor(mk/n)</math> where <math>m/n</math> is the first continued fraction approximant to <math>\sqrt a</math> with a denominator <math>n \geq k</math>
 
 
=={{header|ALGOL 68}}==
where <math>m/n</math> is the first continued fraction approximant to <math>\sqrt a</math> with a denominator <math>n \geq k</math>
{{Trans|Python}}
<syntaxhighlight lang="algol68">
BEGIN # Sturmian word - based on the Python sample #
 
PROC sturmian word = ( INT m, n )STRING:
== {{header|Phix}} ==
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
</syntaxhighlight>
{{out}}
<pre>
01001010010010100101001001010010
</pre>
 
=={{header|BASIC}}==
==={{header|BASIC256}}===
{{trans|FreeBASIC}}
<syntaxhighlight lang="vbnet">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</syntaxhighlight>
{{out}}
<pre>01001010010010100101001001010010</pre>
 
==={{header|FreeBASIC}}===
{{trans|Python}}
<syntaxhighlight lang="vbnet">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</syntaxhighlight>
{{out}}
<pre>01001010010010100101001001010010</pre>
 
==={{header|Yabasic}}===
{{trans|FreeBASIC}}
<syntaxhighlight lang="vbnet">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</syntaxhighlight>
{{out}}
<pre>01001010010010100101001001010010</pre>
 
=={{header|jq}}==
'''Adapted from [[#Wren|Wren]]'''
 
'''Works with jq, the C implementation of jq'''
 
'''Works with gojq, the Go implementation of jq'''
<syntaxhighlight lang="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
</syntaxhighlight>
{{output}}
<pre>
01001010010010100101001001010010 from rational number 13/21
01001010010010100101001001010010 from quadratic number (√5 - 1)/2 (k = 8)
</pre>
 
=={{header|Julia}}==
{{trans|Phix}}
<syntaxhighlight lang="julia">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)")
</syntaxhighlight>{{out}}
<pre>
01001010010010100101001001010010 <== 13/21
01001010010010100101001001010010 <== 1/phi (8th convergent golden ratio)
</pre>
 
== {{header|Phix}} ==
Rational part translated from Python, quadratic part a modified copy of the [[Continued fraction convergents]] task.
<!--(phixonline)-->
<syntaxhighlight lang="Phix">
Line 63 ⟶ 397:
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))
</syntaxhighlight>
{{out}}
<pre>
" 01001010010010100101001001010010" <== 13/21
01001010010010100101001001010010 <== 1/phi (8th convergent)
</pre>
 
Line 112 ⟶ 466:
# 01001010010010100101001001010010
</syntaxhighlight>
 
=={{header|Raku}}==
{{trans|Julia}}
<syntaxhighlight lang="raku" line># 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)";</syntaxhighlight>
{{out}} Same as Julia example.
You may [https://ato.pxeger.com/run?1=ZVRNb5tAEJWqnvwrJhYRrA2YdVTVKiXxtWcfLTfCGAKyvcACtiLL-SO55ND-qf6azCyfSQ_Yq9k3b2beG3j9I_199fb2tyoja_Hv65ei2kJRVvKY-MI6p3JnaEcTNMHgMgIAGeKd-IwQiDgyu5S-KAzd0cG7B53rJv2os6MzSCIEwT1yuUBUx2cwNBkWmLvHJ5PhiYEHOqZxExyXMOc4OYQI28OEKsAtZiOFUzdTk2hBJSUm9qhdcqIqDYRqwItXVwDPaxIeHqgtuLnBJh29A9egGtMF99OpOl8HEhCtO7qORqRYlGxbKRqh1HQr8Yij4J8azCFBHF7XilIJc9vGcS4KqGCPnKlBVuJlmPuh7EqoqjQ4VkV4V_s7c1W4NQfvPvrE70yYc0QV_jOMLx3wCj9RFn43m_OxsgkpbRyrKKXhmD2hHcS-LBiEeR9zawGCKNgbmo_gLe0CbQz5OhBjmZmwzBn88GCNpHxjwpp83nTboCkXi1yWyMTIyy1MlaGz1k5S7Tf6PLT_HKe4I2iZtH-J0u1vDC2jJnKlfgObwDJbTyy-QWZ1mmMf_V3e3eXqrmVbZnZWFTExsi6WN7G8j6kR-Az3GqyGlv23OeumBbOtt-kcJRktXKzqUCKREvWbeh0sfOYmLDrzPr2Bg8S1QyMNA3zDGotnWZyAsShjCFJxCuVTKEp4Sg-7UID0yyRlY7f-FDRfhPbL8A4 Attempt This Online!]
 
=={{header|Wren}}==
{{trans|Python}}
{{libheader|Wren-assert}}
The 'rational number' function is a translation of the Python entry.
For rational numbers m / n:
 
The 'quadratic number' function is a modified version of the one used in the [[Continued fraction convergents]] task.
<syntaxhighlight lang="wren">import "./assert" for Assert
 
var SturmianWordSturmianWordRat = Fn.new { |m, n|
if (m > n) return SturmianWordSturmianWordRat.call(n, m).reduce("") { |acc, c|
return acc + (c == "0" ? "1" : "0")
}
Line 132 ⟶ 532:
}
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])
}
 
Line 146 ⟶ 562:
 
var fib = fibWord.call(10)
var sturmian = SturmianWordSturmianWordRat.call(13, 21)
Assert.equal(fib[0...sturmian.count], sturmian)
System.print("%(sturmian)< from rational number 13/syntaxhighlight>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)")</syntaxhighlight>
 
{{out}}
<pre>
01001010010010100101001001010010 from rational number 13/21
01001010010010100101001001010010 from quadratic number (√5 - 1)/2 (k = 8)
</pre>
2,442

edits