Hofstadter Q sequence
The Hofstadter Q sequence is defined as:
It is defined like the Fibonacci sequence, but whereas the next term in the Fibonacci sequence is the sum of the previous two terms, in the Q sequence the previous two terms tell you how far to go back in the Q sequence to find the two numbers to sum to make the next term of the sequence.
- Task
- Confirm and display that the first ten terms of the sequence are: 1, 1, 2, 3, 3, 4, 5, 5, 6, and 6
- Confirm and display that the 1000'th term is: 502
- Optional extra credit
- Count and display how many times a member of the sequence is less than its preceeding term for terms up to and including the 100,000'th term.
Common Lisp
<lang lisp>(defparameter *mm* (make-hash-table :test #'equal))
- generic memoization macro
(defmacro defun-memoize (f (&rest args) &body body)
(defmacro hash () `(gethash (cons ',f (list ,@args)) *mm*))
(let ((h (gensym)))
`(defun ,f (,@args)
(let ((,h (hash)))
(if ,h ,h (setf (hash) (progn ,@body)))))))
- def q
(defun-memoize q (n)
(if (<= n 2) 1
(+ (q (- n (q (- n 1))))
(q (- n (q (- n 2)))))))
- test
(format t "First of Q: ~a~%Q(1000): ~a~%Bumps up to 100000: ~a~%" (loop for i from 1 to 10 collect (q i)) (q 1000) (loop with c = 0 with last-q = (q 1) for i from 2 to 100000 do (let ((next-q (q i))) (if (< next-q last-q) (incf c)) (setf last-q next-q)) finally (return c)))</lang>output<lang>First of Q: (1 1 2 3 3 4 5 5 6 6) Q(1000): 502 Bumps up to 100000: 49798</lang>
Although the above definition of q is more general, for this specific problem the following is faster:<lang lisp>(let ((cc (make-array 3 :element-type 'integer
:initial-element 1
:adjustable t
:fill-pointer 3)))
(defun q (n)
(when (>= n (length cc)) (loop for i from (length cc) below n do (q i)) (vector-push-extend (+ (aref cc (- n (aref cc (- n 1)))) (aref cc (- n (aref cc (- n 2))))) cc)) (aref cc n)))</lang>
D
<lang d>import std.stdio, std.algorithm, std.functional, std.range;
long Q(long n) {
assert(n > 0);
alias memoize!Q mQ;
if (n == 1 || n == 2)
return 1;
else
return mQ(n - mQ(n - 1)) + mQ(n - mQ(n - 2));
}
void main() {
writeln("Q(n) for n = [1..10] is: ", map!Q(iota(1, 11)));
writeln("Q(1000) = ", Q(1000));
writefln("Q(i) is less than Q(i-1) for i [2..100_000] %d times.",
count!((i){ return Q(i) < Q(i-1); })(iota(2, 100_001)));
}</lang> Output:
Q(n) for n = [1..10] is: [1, 1, 2, 3, 3, 4, 5, 5, 6, 6] Q(1000) = 502 Q(i) is less than Q(i-1) for i [2..100_000] 49798 times.
Icon and Unicon
<lang Icon>link printf
procedure main()
V := [1, 1, 2, 3, 3, 4, 5, 5, 6, 6] every i := 1 to *V do
if Q(i) ~= V[i] then stop("Assertion failure for position ",i)
printf("Q(1 to %d) - verified.\n",*V)
q := Q(n := 1000) v := 502 printf("Q[%d]=%d - %s.\n",n,v,if q = v then "verified" else "failed")
invcount := 0 every i := 2 to (n := 100000) do
if Q(i) < Q(i-1) then {
printf("Q(%d)=%d < Q(%d)=%d\n",i,Q(i),i-1,Q(i-1))
invcount +:= 1
}
printf("There were %d inversions in Q up to %d\n",invcount,n) end
procedure Q(n) #: Hofstader Q sequence static S initial S := [1,1]
if q := S[n] then return q else {
q := Q(n - Q(n - 1)) + Q(n - Q(n - 2))
if *S = n - 1 then {
put(S,q)
return q
}
else
runerr(500,n)
}
end</lang>
printf.icn provides formatting
Output:
Q(1 to 10) - verified. Q[1000]=502 - verified. Q(16)=9 < Q(15)=10 Q(25)=14 < Q(24)=16 Q(32)=17 < Q(31)=20 Q(36)=19 < Q(35)=21 ... Q(99996)=48252 < Q(99995)=50276 Q(99999)=48456 < Q(99998)=50901 Q(100000)=48157 < Q(99999)=48456 There were 49798 inversions in Q up to 100000
J
<lang j>Qs=:0 1 1 Q=: verb define
n=. >./,y
while. n>:#Qs do.
Qs=: Qs,+/((#Qs)-_2{.Qs){Qs
end.
y{Qs
)</lang>
Examples:
<lang j> Q 1+i.10 1 1 2 3 3 4 5 5 6 6
Q 1000
502
+/2>/\ Q 1+i.100000
49798</lang>
PicoLisp
<lang PicoLisp>(de q (N)
(cache '(NIL) (pack (char (hash N)) N)
(if (>= 2 N)
1
(+
(q (- N (q (dec N))))
(q (- N (q (- N 2)))) ) ) ) )</lang>
Test: <lang PicoLisp>: (mapcar q (range 1 10)) -> (1 1 2 3 3 4 5 5 6 6)
- (q 1000)
-> 502
- (let L (mapcar q (range 1 100000))
(cnt < (cdr L) L) )
-> 49798</lang>
Python
<lang python>def q(n):
if n < 1 or type(n) != int: raise ValueError("n must be an int >= 1")
try:
return q.seq[n]
except IndexError:
ans = q(n - q(n - 1)) + q(n - q(n - 2))
q.seq.append(ans)
return ans
q.seq = [None, 1, 1]
if __name__ == '__main__':
first10 = [q(i) for i in range(1,11)]
assert first10 == [1, 1, 2, 3, 3, 4, 5, 5, 6, 6], "Q() value error(s)"
print("Q(n) for n = [1..10] is:", ', '.join(str(i) for i in first10))
assert q(1000) == 502, "Q(1000) value error"
print("Q(1000) =", q(1000))</lang>
- Extra credit
If you try and initially compute larger values of n then you tend to hit the Python recursion limit.
The function q1 gets around this by calling function q to extend the Q series in increments below the recursion limit.
The following code is to be concatenated to the code above: <lang python>from sys import getrecursionlimit
def q1(n):
if n < 1 or type(n) != int: raise ValueError("n must be an int >= 1")
try:
return q.seq[n]
except IndexError:
len_q, rlimit = len(q.seq), getrecursionlimit()
if (n - len_q) > (rlimit // 5):
for i in range(len_q, n, rlimit // 5):
q(i)
ans = q(n - q(n - 1)) + q(n - q(n - 2))
q.seq.append(ans)
return ans
if __name__ == '__main__':
tmp = q1(100000)
print("Q(i+1) < Q(i) for i [1..100000] is true %i times." %
sum(k1 < k0 for k0, k1 in zip(q.seq[1:], q.seq[2:])))</lang>
- Combined output
Q(n) for n = [1..10] is: 1, 1, 2, 3, 3, 4, 5, 5, 6, 6 Q(1000) = 502 Q(i+1) < Q(i) for i [1..10000] is true 49798 times.
Alternative
<lang python>def q(n):
l = len(q.seq)
while l <= n:
q.seq.append(q.seq[l - q.seq[l - 1]] + q.seq[l - q.seq[l - 2]])
l += 1
return q.seq[n]
q.seq = [None, 1, 1]
print("Q(n) for n = [1..10] is:", [q(i) for i in range(1, 11)]) print("Q(1000) =", q(1000)) q(100000) print("Q(i+1) < Q(i) for i [1..100000] is true %i times." %
sum([q.seq[i] > q.seq[i + 1] for i in range(1, 100000)]))</lang>
Scheme
I wish there were a portable way to define-syntax, or to resize arrays, or to do formated output--anything to make the code less silly looking while still run under more than one interpreter.
<lang lisp>(define qc '#(0 1 1))
(define filled 3)
(define len 3)
- chicken scheme
- vector-resize!
- gambit
- vector-append
(define (extend-qc)
(let* ((new-len (* 2 len))
(new-qc (make-vector new-len)))
(let copy ((n 0))
(if (< n len)
(begin (vector-set! new-qc n (vector-ref qc n)) (copy (+ 1 n)))))
(set! len new-len) (set! qc new-qc)))
(define (q n)
(let loop ()
(if (>= filled len) (extend-qc))
(if (>= n filled)
(begin
(vector-set! qc filled (+ (q (- filled (q (- filled 1)))) (q (- filled (q (- filled 2)))))) (set! filled (+ 1 filled)) (loop))
(vector-ref qc n))))
(display "Q(1 .. 10): ") (let loop ((i 1))
;; (print) behave differently regarding newline across compilers (display (q i)) (display " ") (if (< i 10) (loop (+ 1 i)) (newline)))
(display "Q(1000): ") (display (q 1000)) (newline)
(display "bumps up to 100000: ") (display
(let loop ((s 0) (i 1))
(if (>= i 100000) s
(loop (+ s (if (> (q i) (q (+ 1 i))) 1 0)) (+ 1 i)))))
(newline)</lang>output<lang>Q(1 .. 10): 1 1 2 3 3 4 5 5 6 6 Q(1000): 502 bumps up to 100000: 49798</lang>
Tcl
<lang tcl>package require Tcl 8.5
- Index 0 is not used, but putting it in makes the code a bit shorter
set tcl::mathfunc::Qcache {Q:-> 1 1} proc tcl::mathfunc::Q {n} {
variable Qcache
if {$n >= [llength $Qcache]} {
lappend Qcache [expr {Q($n - Q($n-1)) + Q($n - Q($n-2))}]
} return [lindex $Qcache $n]
}
- Demonstration code
for {set i 1} {$i <= 10} {incr i} {
puts "Q($i) == [expr {Q($i)}]"
}
- This runs very close to recursion limit...
puts "Q(1000) == [expr Q(1000)]"
- This code is OK, because the calculations are done step by step
set q [expr Q(1)] for {set i 2} {$i <= 100000} {incr i} {
incr count [expr {$q > [set q [expr {Q($i)}]]}]
} puts "Q(i)<Q(i-1) for i \[2..100000\] is true $count times"</lang> Output:
Q(1) == 1 Q(2) == 1 Q(3) == 2 Q(4) == 3 Q(5) == 3 Q(6) == 4 Q(7) == 5 Q(8) == 5 Q(9) == 6 Q(10) == 6 Q(1000) == 502 Q(i)<Q(i-1) for i [2..100000] is true 49798 times