NAME
Math::NumSeq::SternDiatomic  Stern's diatomic sequence
SYNOPSIS
use Math::NumSeq::SternDiatomic;
my $seq = Math::NumSeq::SternDiatomic>new;
my ($i, $value) = $seq>next;
DESCRIPTION
This is Moritz Stern's diatomic sequence
0, 1, 1, 2, 1, 3, 2, 3, ...
starting i=0
It's constructed by successive levels with a recurrence
D(0) = 0
D(1) = 1
D(2*i) = D(i)
D(2*i+1) = D(i) + D(i+1)
So the sequence is extended by copying the previous level to the next spead out to even indices, and at the odd indices fill in the sum of adjacent terms,
0, i=0
1, i=1
1, 2, i=2 to 3
1, 3, 2, 3, i=4 to 7
1,4,3,5,2,5,3,4, i=8 to 15
For example the i=4to7 row is a copy of the preceding row values 1,2 with sums 1+2 and 2+1 interleaved.
For the new value at the end of each row the sum wraps around so as to take the last copied value and the first value of the next row, which is always 1. This means the last value in each row increments 1,2,3,4,5,etc.
Odd and Even
The sequence makes a repeating pattern even,odd,odd,
0, 1, 1, 2, 1, 3, 2, 3
E O O E O O E ...
This can be seen from the copying in the recurrence above. For example the i=8 to 15 row copying to i=16 to 31,
O . E . O . O . E . O . O . E . spread
O O E O O E O O sum adjacent
Adding adjacent terms odd+even and even+odd are both odd. Adding adjacent odd+odd gives even. So the pattern E,O,O in the original row when spread and added gives E,O,O again in the next row.
FUNCTIONS
See "FUNCTIONS" in Math::NumSeq for behaviour common to all sequence classes.
Random Access
$value = $seq>ith($i)

Return the
$i
'th value of the sequence. ($v0, $v1) = $seq>ith_pair($i)

Return two values
ith($i)
andith($i+1)
from the sequence. As described below ("Ith Pair") two values can be calculated with the same work as one. $bool = $seq>pred($value)

Return true if
$value
occurs in the sequence, which means simply integer$value>=0
.
FORMULAS
Next
The sequence is iterated using a method by Moshe Newman in
"Recounting the Rationals, Continued", answers to problem 10906 posed by Donald E. Knuth, C. P. Rupert, Alex Smith and Richard Stong, American Mathematical Monthly, volume 110, number 7, AugSep 2003, pages 642643, http://www.jstor.org/stable/3647762
Two successive sequence values are maintained and are advanced by a simple operation.
maintaining p = seq[i], q = seq[i+1]
initial p = seq[0] = 0
q = seq[1] = 1
p_next = seq[i+1] = q
q_next = seq[i+2] = p+q  2*(p mod q)
where the mod operation rounds towards zero
0 <= (p mod q) <= q1
The form by Newman uses a floor operation. That suits expressing the iteration in terms of a rational x[i]=p/q.
p_next 1
 = 
q_next 1 + 2*floor(p/q)  p/q
For separate p,q a little rearrangement gives it in terms of the remainder p mod q, as per formula a(n)=a(n2)+a(n1)2*(a(n2) mod a(n1)) by Mike Stay OEIS A002487, November 2006.
division p = q*floor(p/q) + rem where rem = (p mod q)
then
p_next/q_next = 1 / (1 + 2*floor(p/q)  p/q) per Newman
= q / (2*q*floor(p/q) + q  p)
= q / (2*(p  rem) + q  p)
= q / (p+q  2*rem) using p,q
In terms of the CalkinWilf tree this method works because the number of trailing right leg steps can be found by m=floor(p/q), then take a step across, then back down again by m many left leg steps. When at the rightmost edge of the tree the step across goes down by one extra left, thereby automatically wrapping around at the end of each row.
seek_to_i()
is implemented by calculating a pair of new p,q values with an ith_pair()
per below.
Ith Pair
For ith_pair()
the two sequence values at an arbitrary i,i+1 can be calculated from the bits of i,
p = 0
q = 1
for each bit of i from high to low
if bit=1 then p += q
if bit=0 then q += p
return p,q # are ith(i) and ith(i+1)
For example i=6 is binary "110" so
p,q

initial 0,1
high bit=1 so p+=q 1,1
next bit=1 so p+=q 2,1
low bit=0 so q+=p 2,3 is ith(6),ith(7)
This is the same as the CalkinWilf tree descent, per "CalkinWilf Tree" in Math::PlanePath::RationalsTree. Its X/Y fractions are successive Stern diatomic sequence values.
Ith Alone
If only a single ith() value is desired then a variation can be made on the "Ith Pair" above. Taking the bits of i from low to high (instead of high to low) gives p=ith(i), but q is not ith(i+1). Low zero bits can be ignored for this approach since initial p=0 means the steps q+=p for bit=0 do nothing.
SEE ALSO
Math::PlanePath::RationalsTree
HOME PAGE
http://user42.tuxfamily.org/mathnumseq/index.html
LICENSE
Copyright 2011, 2012, 2013, 2014 Kevin Ryde
MathNumSeq is free software; you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation; either version 3, or (at your option) any later version.
MathNumSeq is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details.
You should have received a copy of the GNU General Public License along with MathNumSeq. If not, see <http://www.gnu.org/licenses/>.