/* From Don Reble, Oct 16 2017 */ /* %I A293160 %S 1 2 3 5 7 13 20 31 48 78 118 191 300 465 734 1175 1850 2926 4597 7296 11552 18278 28863 45832 72356 114742 %N Number of distinct terms in the n'th row of Stern's A049456. %I A135510 %S 2 3 4 6 6 14 20 28 38 54 90 150 216 350 506 876 1230 2034 3160 4470 7764 12190 18816 29952 43800 73968 112602 182210 %N Least positive number missing from row n of Stern's diatomic array (see A049456 or A002487). */ #include <iostream> #include <map> typedef unsigned long int uint32; /* Successive numerators yield Stern's diatomic sequence, A002487. Use 1/N to start row #N; it has length 2^(N-1). */ static void nextCalkinWilf(uint32& num, uint32& den) { uint32 aa = num / den; uint32 newden = (aa + aa + 1) * den - num; num = den; den = newden; } // Set representation: a sequence of [a,b) pairs. A std::map maps a -> b. typedef std::map<uint32,uint32> intset; typedef intset::value_type span; typedef intset::iterator setiter; static intset TheSet; // merges adjacent pairs when possible; // returns whether an insertion happened static bool insertTS(uint32 val) { span sp(val,val+1); if (TheSet.empty()) { TheSet.insert(sp); return true; } setiter post = TheSet.upper_bound(val); setiter prev = post; // decrement later if legal if (post != TheSet.end()) { if (post->first == sp.second) { sp.second = post->second; // will replace *post if (post == TheSet.begin()) { TheSet.erase(post); TheSet.insert(TheSet.begin(), sp); } else { // maybe merge --prev; if (prev->second == val) { prev->second = sp.second; } else { TheSet.insert(post, sp); } TheSet.erase(post); } return true; } } if (post == TheSet.begin()) { TheSet.insert(post, sp); return true; } --prev; if (prev->second > val) { return false; } if (prev->second == val) { prev->second = sp.second; } else { TheSet.insert(post, sp); } return true; } // return the number of elements in TheSet static uint32 countTS() { uint32 qu = 0; for (setiter it = TheSet.begin(); it != TheSet.end(); ++it) { qu += it->second - it->first; } return qu; } // prints N, A293160(number of distinct values), // A135510(first missing value), A000045(max value) // A000045 is just a double-check on program operation. static void printStuff(int nn) { setiter it = TheSet.end(); --it; std::cout << nn << ' ' << countTS() << ' ' << (TheSet.begin()->second) << ' ' << (it->second - 1) << std::endl; } int main() { static const int MaxN = 40; // Fib(47) < 2^32, so we can go that far. // (If ulongint has 32 bits and we're patient.) insertTS(1); printStuff(1); insertTS(2); printStuff(2); uint32 halfway = 1; for (int nn = 3; nn <= MaxN; nn += 1) { uint32 num = 1; // start of row N uint32 den = nn; /* Even-indexed values are old; second half is reverse of first; so go halfway, insert half of the time. */ for (uint32 kk = 0; kk < halfway; kk += 1) { nextCalkinWilf(num, den); insertTS(num); nextCalkinWilf(num, den); } printStuff(nn); halfway *= 2; } return 0; }