▼
Posts: 1,755
Threads: 112
Joined: Jan 2005
Hi all:
A few days ago, Csaba posted:
"what was the program, what used by Feigenbaum originally(!), to calculate this number on his HP65? [...] The greatest problem, how to isolate the bifurcation points on the diagram..."
Actually, it's quite simple. Have a look at this small 4line, 174byte adhoc program I've written for the HP71B (which can be very easily ported to most any other HP calc, Feigenbaum's own HP65 included):
10 DESTROY ALL @ DELAY 0,0 @ F=3.2 @ M=0 @ N=1
20 FOR K=2 TO 9 @ P=N+(NM)/F @ X=FNROOT(P,P,FNB(2^K,FVAR)) @ F=(NM)/(XN)
30 PRINT "K=";K;"F=";F @ M=N @ N=X @ NEXT K @ END
40 DEF FNB(K,A) @ B=0 @ FOR I=1 TO K @ B=AB*B @ NEXT I @ FNB=B @ END DEF
When run, it will perform 9 iterations of a procedure designed to calculate an increasingly refined value of Feigenbaum's constant, starting from a rough initial approximation (3.2). Just press [RUN] and you'll get
these results (K = #iteration, F = estimated value for
Feigenbaum's constant):
K F

2 3.21851142201
3 4.38567759894
4 4.60094927695
5 4.65513048088
6 4.66611150124
7 4.66855134603
8 4.66905077160
9 4.66921465563
After some 200 seconds, the last result at the end of the 9th iteration is F = 4.66921+, correct to 6 digits within 1 ulp (err = 0.00001+). For comparison purposes, the value of Feigenbaum's constant is 4.66920 16091029906718+.
Within the constraints of HP71B's 12digit precision and speed, this is more or less the best result you can get. A further 10th iteration would take much, much longer for no real improvement, as the running time grows exponentially with every iteration (there are a number of easy ways of speeding my 4liner, of course, but I prefer to keep it at its simplest for didactic purposes, if nothing else).
By the time it reaches the 9th iteration, it is finding a root of a 512thdegree (!) polynomial, which takes a lot of time and incurs in growing cumulative rounding errors when evaluating it a number of times as part of the search for its root. A 10th iteration would mean solving (i.e.: evaluating repeatedly) a 1,024thdegree (!!) polynomial, and cumulative rounding errors mean 6 correct digits is the best you can do in 12digit arithmetic.
However, if you have access to an arbitrary precision package either for a PC (Maple, Mathematica, etc) or even your HP71B, you can easily get more correct digits using this algorithm, albeit at a heavy time penalty.
Best regards from V.
Edited: 3 June 2004, 6:17 a.m.
▼
Posts: 275
Threads: 41
Joined: Mar 2010
Hi Valentin,
Btw I search actually as much as possible information on the HP8x series (processor, memory controller, ...).
I have an HP85A, an HP86B some 91{3/5}3 HD and a 9121 floppy. To use the 91xx I need an EMS rom, but can't find one. So I think about disassembling the rom, patching to allowed relocation and reassemble sources to allow them to be loaded as binary. I think it's feasible but I need some more information on the 8x series. Actually I have a python programm who can deassemble an HP8x binary in a quite usefull form. The assembler will be easier to write.
For example after searching the web I still have some questions:
 What does SAD and PAD opcodes exactly ? (octal 232 and 237)
 Is there an opcode for octal 326 and 336 ?
 What are the vectors for interrupts in the capricorn processor ?
 How the rom and ram paging works in an HP86 (preciselly)?
 I partially understood the structure of a binary program, but I miss informations for the rom structure ...
I get some informations from HP journals, from HP75C programs and from assemblers globals (HP85 and HP86)
Thanks in advance if you could help me.
Olivier
▼
Posts: 1,755
Threads: 112
Joined: Jan 2005
Hi, Olivier:
Olivier wrote:
"Btw I search actually as much as possible information on the HP8x series (processor, memory controller, ...)."
I did read your request the first time you posted it.
Unfortunately, it's been more than 1518 years since I
last had anything to do with an HP86 and all materials I once owned (manuals, ROMs, listings, binaries, etc) are now irretrievably missing. And, regrettably, my memory is of no
help for such accurate details after so long a time lapse.
Anyway, should I find any materials of information, by sheer luck, I'll let you know promptly.
Best of lucks in your plausible search and
Best regards from Valentin Albillo
Posts: 1,755
Threads: 112
Joined: Jan 2005
Hi,
Just to prove the point, here's a conversion of my 4liner HP71B program to the HP15C !
Despite the obvious (and large) differences in language level (lowlevel RPN vs. highlevel BASIC), numerical precision (10digit vs. 12digit) and above all, speed (1530x slower), the HP15C is perfectly up to the task, delivering Feigenbaum's constant correct to 4 digits in 18 min., which fares well compared to the 71B results (6 digits in 3.3 min.).
Also, the HP15C program is much more compact, at 60+ bytes (47 steps) versus 71B's 174 bytes (4 multistatement lines), i.e. just 1/3rd the size (RPN strikes again!). Listing follows with a few pertinent comments:
01 LBL A 14 LBL 0 25 RCL 1 36 GTO 0
02 MATRIX 1 15 2 26 RCL 4 37 RTN
03 CLX 16 STOx 5 27 RCL 0 38 LBL 1
04 STO 1 17 RCL 0 28 / 39 RCL 5
05 7 18 RCL 1 29 STO 2 40 STO 6
06 STO 3 19 RCL/ 2 30 R/S 41 CLX
07 2 20 RCL+ 0 31 RCL 0 42 LBL 2
08 STO 5 21 ENTER 32 STO 1 43 X^2
09 3.2 22 SOLVE 1 33 RCL 4 44 
12 STO 2 23 STO 4 34 STO 0 45 DSE 6
13 FIX 6 24 RCL 0 35 DSE 3 46 GTO 2
47 RND
To use, simply press RUN with no inputs. The program will compute and display the increasingly accurate value of Feigenbaum's constant for each iteration, stopping for you to see it. Press R/S to go on with the next iteration. After completing all 7 iterations, the program will stop displaying the last computed superstable a value, a = 1.401146+, which is the SOLVEd root of the corresponding 256thdegree polynomial.
Results:
Iteration F value Iteration time

1 3.218510 17 sec.
2 4.385687 50
3 4.600945 65
4 4.655086 83
5 4.666256 162
6 4.667830 234
7 4.669803 477
Last result is 4.669+, accurate to 4 digits within 1 ulp.
Notes:
 As the HP15C uses 10digit precision instead of HP71B's 12digit, we can't get a result accurate to more than 45 significant digits, so just 7 iterations will suffice. Besides, we don't need to find the highdegree polynomials' roots to full 10digit precision, thus steps 13 FIX 6 and 47 RND provide early termination for 22 SOLVE 1, saving much time at no cost in final accuracy.
 The highdegree polynomial that is being SOLVEd is computed at steps 3847. Notice the very tight loop 4246 which is once again a superb practical demonstration of classic RPN's stack capabilities at its best: no matter how high the degree of the polynomial, the repeated subtraction at step 44 never exhausts the provision of a values, thanks to toplevel replication, where T is replicated continuously as the stack drops and drops. If literally coded in RPL this loop would collapse the stack in no time !
 The ever bigger polynomials' roots are found using the powerful SOLVE instruction at step 22. In a running program, SOLVE behaves like a test: if a root is found, the next step is executed, else it is skipped. Normally this must be taken into account, but for this particular program, the root always exists and SOLVE always finds it for the number of iterations considered, so there's no need to code the possibility of it failing.
 Last but not least, notice how every advanced feature of the HP15C instruction set helps to reduce step count and increase speed, from the matrix instruction at step 02, to recall arithmetic everywhere, to being able to use any numbered register for loop indexing purposes.
In conclusion, I hope you'll agree that the mere fact that this kind of computation is feasible in the HP15C in reasonable times and with such a small, simple program, is but another proof of the outstanding quality of this wonderful machine. As "Moulin Rouge!"'s Zidler himself would say on behalf of the HP15C: "A magnificent, opulent, tremendous, stupendous, gargantuan success"
Best regards from V.
Edited: 4 June 2004, 7:55 a.m.
Posts: 113
Threads: 20
Joined: Sep 2013
I ported Valentin's Feigenbaum program to the HP49G+, and it runs a little faster (as it should, compared to the 71B), i.e., 9 iterations in 1 min 35 sec.
However, it seems to be less accurate (more prone to rounding error).
After 7 iterations, F= 4.66854950126,
after 8 iterations, F= 4.6690698451,
after 9 iterations, F= 4.66915720392. (only accurate to 4 digits)!
Even after 10 iterations, F=4.66903301492. The CAS mode (approximate or exact, Real or Complex) doesn't matter.
Strange, I thought the 49g+ was also using 12 digit precision internally, so the results should be identitical to those obtained on the 71B.
Here is the 49G+ program:
FEIG:
<< CLLCD TIME 'T' STO 3.2 'F' STO 0. 'M' STO 1. 'N' STO
2. 9. FOR k
N M  F / N + 'P' STO
<< 2. k ^ 'X' FNB EVAL >>
'X' P ROOT 'X' STO N M  X N  / 'F' STO
"K = " k + 1. DISP
"F = " F + 2. DISP
"Time = " TIME T HMS 100. * + 3. DISP
N 'M' STO X 'N' STO
NEXT 7. FREEZE
>>
FNB:
<< > K A
<< 0 1 K START DUPDUP * NEG A + NIP NEXT >>
>>
Any suggestions on how to improve the accuracy?
Posts: 97
Threads: 25
Joined: Jan 1970
Thank you!
I think I changing my mind about programming...
Csaba
