Short & Sweet Math Challenge #21: Powers that be

11022016, 01:13 AM
(This post was last modified: 11022016 09:48 PM by Valentin Albillo.)
Post: #1




Short & Sweet Math Challenge #21: Powers that be
Short & Sweet Math Challenge #21: Powers that be Hi all, At the end of April 2008 and after 20 releases I concluded my "Short & Sweet Math Challenges" series which allowed you to dust off and show off both your favorite programmable HP calcs and your math/HPprogramming skills while introducing and discussing interesting and unusual mathematical topics. The series was well received and now, some 8 years later, I'm releasing "season 2", so to say, beginning with this brandnew S&SMC #21, in the tradition and style of the preceding ones. Iif this second installment gets a similarly good reception, I intend to eventually post another full 20strong batch. Now as for this present new challenge, first a short prelude. Many of you will surely have heard about Phi, the ubiquitous socalled Golden Ratio constant which can be promptly evaluated as (1+Sqrt(5))/2, namely: Phi = 1.61803+ Phi is further characterized by being the largest root in absolute value of the polynomial x^2x1 (i.e.: Phi^2Phi1 = 0), which happens to be the lowestdegree polynomial for which Phi is a root and thus it gets called the "minimal polynomial" of Phi. As it happens, Phi has a plethora of very interesting, sometimes unique properties, one of them being that its increasing powers are ever nearer to integer values, for instance: Phi^20 = 15126.99993+ Phi^50 = 28143753122.99999999996+ Phi^100 = 792[...]126.999999999999999999998+ Phi^200 = 627[...]126.999999999999999999999999999999999999999998+ Alas, while certainly interesting, this quasiinteger property of the powers of Phi is far from unique, there are some other constants that boast this very property, as for instance (let's call it "Hpi" for the time being): Hpi = 1.32471+ which happens to be the largest root in absolute value of the polynomial x^3x1 (i.e.: Hpi^3Hpi1 = 0), which is its minimal polynomial. Indeed, we do have: Hpi^20 = 276.992+ Hpi^50 = 1276942.001+ Hpi^100 = 163[...]001.9999995+ Hpi^200 = 265[...]250.000000000001+ Hpi^500 = 115[...]876.9999999999999999999999999999994+ demonstrating that Hpi's powers also get nearer and nearer to integer values, like Phi's do. That said, now go and get a programmable HP calc of your choice (say an HP10C or better, hardware/software based HPcalc emulators also welcomed) and use that calc's language of your choice (say RPN, FOCAL, RPL, BASIC, FORTH, etc.) to try and meet: The Challenge: Write a program to find and output the constants that have the aforementioned quasiintegerpowers property, subject to the requirements that your program must:  find ALL such constants in the range 1 < constant < 2, for minimal polynomials up to degree 8 having all their coefficients equal to either +1, 0, or 1, the leading coefficient in particular being always +1.  output both each constant AND the minimal polynomial for which it is a root, sorted by increasing numerical value, with a final tally of how many constants were found (hint: more than 30) and, optionally, timing. Of course, the faster your program runs the better, smaller program sizes being important but secondary. If your chosen HP calc runs too slowly you might consider, in order: . using a better algorithm and/or optimizing your solution for speed, . coding in a faster language (say FORTH or assembler, if available), . using a faster HP calc, . using an HPcalc emulator running on a faster platform, or . just be patient and have a meal (or two) while it merrily runs. That's all. Within a few days I'll post and comment my 12line (60statement, 585byte) solution for the HP71B, which runs like this (no spoilers): >RUN ... some lines of output ... 1.57367896839 x^8x^7x^6+x^21 1.59000537390 x^7x^5x^4x^3x^2x1 1.60134733379 x^7x^6x^4x^21 1.61803398875 x^2x1 ... more lines of output ... (number of constants found) I'll also comment on the underlying theory and algorithm used as well as trivia and problems I found. Now for the small print:  you can use any ROMs, libraries or binary files for your calc as long as they are readily available, preferably for free, as well as extra RAM if needed. For instance, for the HP41 and emulators you can use the Math ROM, Advantage ROM, card reader ROM, printer ROM and many others. My HP71B solutions typically may use the Math ROM, HPIL ROM, JPCROM, STRNGLEX binary and additional RAM modules.  you can use any hard/soft emulator for a given HP calc model but solutions given for nonHP calcs (say written in Excel, Visual Basic, C#, for vintage SHARP or TI machines, etc) won't be considered to have met the Challenge and in fact may spoil the fun for all.  ditto for using the web to search for solutions. You will ruin the fun for yourself and other people attempting the challenge so you'd better try and solve it by your own efforts first, not mercilessly scavenging the web. Enough said, now let's see your solutions ! (and please do not use CODE sections in your posts, they don't print well). Best regards. V. . All My Articles & other Materials here: Valentin Albillo's HP Collection 

11022016, 05:47 AM
Post: #2




RE: Short & Sweet Math Challenge #21: Powers that be
HP 50g, 181.5 bytes, 14 seconds. Only 6 constants, though.
Best regards, Gerson. %%HP: T(3)A(R)F(.); \<< 1. 6. FOR i i X SWAP ^ [ 1. 1. 1. ] X PEVAL * [ 1. 0. 1. ] X PEVAL + EVAL DUP 'X' ZEROS DUP TYPE IF THEN OBJ\> DROP NIP END NEXT \>> TEVAL 'X^3.+0.*X^2.+1.*X1.' 1.32471795724 'X^4.+1.*X^3.+0.*X^2.1.' 1.3802775691 'X^5.+1.*X^4.+1.*X^3.+X^2.1.' 1.44326879127 'X^6.+1.*X^5.+1.*X^4.+X^2.1.' 1.50159480354 'X^7.+1.*X^6.+1.*X^5.+X^2.1.' 1.54521564973 'X^8.+1.*X^7.+1.*X^6.+X^2.1.' 1.57367896839 s:14.0446 

11022016, 07:16 AM
Post: #3




RE: Short & Sweet Math Challenge #21: Powers that be
Welcome back Valentin and S&SMC!
Greetings, Massimo +×÷ ↔ left is right and right is wrong 

11022016, 01:48 PM
(This post was last modified: 11022016 01:55 PM by JF Garnier.)
Post: #4




RE: Short & Sweet Math Challenge #21: Powers that be
Hi Valentin,
I'm very happy to read you again in this forum ! One question: From the output examples here: (11022016 01:13 AM)Valentin Albillo Wrote:it seems that null polynomial coefficients are allowed, not only +1 or 1, Is it correct? JF 

11022016, 02:00 PM
Post: #5




RE: Short & Sweet Math Challenge #21: Powers that be
Hello Valentin,
Welcome back and thanks for another very interesting challenge. Some of the newer forum members may not be familiar with your past challenges. Back in 2005, I assembled the first eleven of these into a single PDF file. I never did update it. You are up to 21 now. I guess I need to spend a little time getting it updated to 21. Please see the following forum post for a link to the PDF file of the first 11: SSMC PDF Bill Smithville, NJ 

11022016, 10:08 PM
Post: #6




RE: Short & Sweet Math Challenge #21: Powers that be
.
Hi, Gerson ! I'm truly glad to get your always valuable contributions to one of my S&SMC's, as in the good old times. A few comments: (11022016 05:47 AM)Gerson W. Barbosa Wrote: HP 50g, 181.5 bytes, 14 seconds. Only 6 constants, though. Why only 6 ? Not being versed in RPL I don't fully understand your code but I also don't see any reference to the maximum degree for the polynomials, which should be 8. Quote:'X^3.+0.*X^2.+1.*X1.' 1.32471795724 Also, even within this limited range, your program is missing two constants, one between 1.50.. and 1.54.. and another between 1.54.. and 1.57.., perhaps due to the typo that JF pointed out (which I duly corrected) about the coefficients being 1, 0, or +1 (the '0' was missing in my OP). Thnaks for your continued interest and best regards. V. . All My Articles & other Materials here: Valentin Albillo's HP Collection 

11022016, 10:14 PM
Post: #7




RE: Short & Sweet Math Challenge #21: Powers that be
Hi, Massimo !: (11022016 07:16 AM)Massimo Gnerucci Wrote: Welcome back Valentin and S&SMC! Thank you very much for your kind welcome, much appreciated. Perhaps you'd consider contributing your very own solution to the present challenge ? It's just a warmer, there are many interesting ones ahead ! Best regards. V. All My Articles & other Materials here: Valentin Albillo's HP Collection 

11022016, 10:22 PM
Post: #8




RE: Short & Sweet Math Challenge #21: Powers that be
Hi, JeanFrançois !: Thanks a lot for your kind comment and continued interest in my posts, much appreciated. (11022016 01:48 PM)JF Garnier Wrote: One question: Yes, of course, my bad, thanks for pointing it out to me, I've already corrected it in my original post. I was going to write that the absolute value of the integer coefficients had to be up to and including 1 but decided instead to just enumerate them and the '0' was simply left out. Thanks and here's hoping for your own solution to this challenge, it would be an amazing way to start the "second season" ! ... 8D Best regards. V. All My Articles & other Materials here: Valentin Albillo's HP Collection 

11022016, 10:39 PM
Post: #9




RE: Short & Sweet Math Challenge #21: Powers that be
Hi, Bill ! Thank for your kind welcome and appreciation, I'm glad you like my S&SMC's (11022016 02:00 PM)Bill (Smithville NJ) Wrote: Back in 2005, I assembled the first eleven of these into a single PDF file. I never did update it. You are up to 21 now. I guess I need to spend a little time getting it updated to 21. It would be great if you'd eventually find the time to update your PDF file with the whole collection but you might consider waiting till the "second season" is completed instead of doing incremental updates. Whatever suits you. Quote:Please see the following forum post for a link to the PDF file of the first 11: SSMC PDF I already downloaded your PDF file back then and found it extremely useful to me because I was missing some of the earliest challenges and your PDF put them back in my hands in a most convenient format so I was very obliged to you for it. If you'd like to try your hand at providing your very own solution (or any comments) to the present challenge, it would be my pleasure to have a look at it. Best regards. V. All My Articles & other Materials here: Valentin Albillo's HP Collection 

11022016, 11:19 PM
(This post was last modified: 11022016 11:34 PM by Gerson W. Barbosa.)
Post: #10




RE: Short & Sweet Math Challenge #21: Powers that be
(11022016 10:08 PM)Valentin Albillo Wrote: . Hello, Valentin! Thanks for starting Season 2. I'm looking forward for the next episodes. I've always appreciated your insightful and wellthought S&SMC series, even though most of the time I was able to solve only the easier items. Regarding this particular problem, I remember Phi belongs to a special set of numbers named after a French mathematician which produce nearintegers when raised to high powers. I misspelled his name, but Google pointed me to the right reference wherein I found three polynomials that generate such numbers. The RPL program uses only the first generating polynomial. PEVAL creates a symbolic polynomial expression from a variable name and a coefficents vector. For instance, [ 1 1 1 ] 'X' PEVAL > '1+(1+X)*x' FACTOR > 'X²+X+1' PROOT might be a better alternative to ZEROS, but I couldn't find a builtin inverse to PEVAL so PROOT could be used. Anyway, I guess this is not the kind of solution your looking for, so I won't proceed with this approach any longer. (11022016 10:08 PM)Valentin Albillo Wrote: Also, even within this limited range, your program is missing two constants, one between 1.50.. and 1.54.. and another between 1.54.. and 1.57.., perhaps due to the typo that JF pointed out (which I duly corrected) about the coefficients being 1, 0, or +1 (the '0' was missing in my OP). Yes, I was aware of the missing constants. By using the third generating polynomial in the aforementioned reference, a few more can be found: %%HP: T(3)A(R)F(.); \<< 3. 8. FOR n X n ^ X n 1. + ^ 1.  X 2. ^ 1.  /  DUP 'X' ZEROS DUP SIZE GET NEXT \>> TEVAL 'X^3.(X^4.1.)/(X^2.1.)' 1.46557123188 'X^4.(X^5.1.)/(X^2.1.)' 1.53415774491 'X^5.(X^6.1.)/(X^2.1.)' 1.5701473122 'X^6.(X^7.1.)/(X^2.1.)' 1.5900053739 'X^7.(X^8.1.)/(X^2.1.)' 1.60134733379 'X^8.(X^9.1.)/(X^2.1.)' 1.60798272793 :s: 18.1289 Not in the required format, though. Also, the polynomial have yet to be simplified and checked whether the degrees are no greater than 8. Best regards, Gerson. Edited to fix a typo 

11062016, 05:21 PM
(This post was last modified: 11062016 05:22 PM by JF Garnier.)
Post: #11




RE: Short & Sweet Math Challenge #21: Powers that be
(11022016 10:22 PM)Valentin Albillo Wrote: ... hoping for your own solution to this challenge, it would be an amazing way to start the "second season" ! ... 8D Well, I tried first on Emu71, then switched to Free42 to have higher computing accuracy, but I'm afraid I wasn't able to build a proper solution with either tool. With Emu71, I was able to find the roots, but wasn't able to identify all the constants with the desired property due to the limited accuracy. And with Free42, I had the right accuracy, but had difficulties to find the roots of the polynomials. Waiting for your solution... JF 

11072016, 10:23 PM
(This post was last modified: 11072016 10:23 PM by Valentin Albillo.)
Post: #12




RE: Short & Sweet Math Challenge #21: Powers that be
Hi, JF: (11062016 05:21 PM)JF Garnier Wrote: Well, I tried first on Emu71, then switched to Free42 to have higher computing accuracy, but I'm afraid I wasn't able to build a proper solution with either tool. How many did you identify ? How do you know they aren't all ? For the particular limits of this challenge, i.e.: minimal polynomials up to degree 8 (or less) and with coefficients 1,0,+1, I found no accuracy problems at all (though perhaps there are and I just didn't find them ... 8D ) Quote:And with Free42, I had the right accuracy, but had difficulties to find the roots of the polynomials. How so ? Details ? Quote:Waiting for your solution... I'll post it within two days, give or take a day. It does take a significant amount of time to write down the solution post and regrettably it seems there wasn't much interest at all, no one but you and Gerson made any attempt at a solution or posted comments, let alone post actual code. Thanks for your interest. Best regards. V. All My Articles & other Materials here: Valentin Albillo's HP Collection 

11082016, 08:28 AM
(This post was last modified: 11082016 10:11 AM by JF Garnier.)
Post: #13




RE: Short & Sweet Math Challenge #21: Powers that be
(11072016 10:23 PM)Valentin Albillo Wrote: How many did you identify ? How do you know they aren't all ? I didn't attempt to count them, but I know that I didn't identify all constants with Emu71, because I missed at least one: the 1.3802775691 value that Gerson reported. When trying to check if the powers of this value get closer and closer to an integer, I got: I 1.3802775691^I 51 13749532.9553 52 18978171.9238 53 26195145.0089 54 36156571.0751 55 49906104.0305 56 68884275.9545 57 95079420.9637 58 131235992.039 59 181142096.07 60 250026372.026 61 345105792.99 62 476341785.031 63 657483881.103 64 907510253.132 65 1252616046.13 66 1728957831.16 67 2386441712.27 68 3293951965.42 69 4546568011.56 70 6275525842.74 ... The closest values are around power 60, then the lack of accuracy makes the fractional part no more significant. My criteria was that the value must be close to an integer with 0.01 tolerance for 3 successive power values, to have a good confidence. I could have relax my criteria, but what if other values had an even worst behaviour? With Free42, I had the right accuracy, with 35 digits. I first solved the equation x^4x^31=0 with the solver to have the root X with full accuracy, then the powers (frac part) are: FP(X^75) = 0.98147851 FP(X^76) = 0.98907003 FP(X^77) = 0.01158426 FP(X^78) = 0.01475045 FP(X^79) = 0.99622896 FP(X^80) = 0.98529899 FP(X^81) = 0.99688326 FP(X^82) = 0.01163371 FP(X^83) = 0.00786267 FP(X^84) = 0.99316166 FP(X^85) = 0.99004492 FP(X^86) = 0.00167862 FP(X^87) = 0.00954129 FP(X^88) = 0.00270295 FP(X^89) = 0.99274787 FP(X^90) = 0.99442649 FP(X^91) = 0.00396778 FP(X^92) = 0.00667073 FP(X^93) = 0.99941859 FP(X^94) = 0.99384508 FP(X^95) = 0.99781286 Here my criteria is met around power 85. But since there is no polynomial root solver on the HP42S (and I didn't want to write one), I had to use the equation solver, and I wasn't sure that the reported root was the correct, largest one. Regarding the participation to the challenge, maybe others met the same difficulties. Unless there is another, better way to solve the challenge :) JF 

11082016, 01:38 PM
(This post was last modified: 11082016 01:40 PM by Paul Dale.)
Post: #14




RE: Short & Sweet Math Challenge #21: Powers that be
I'd figured out an approach to this problem. I've not implemented it on any hardware.
It is essentially a brute force approach to the problem: For each length of polynomial (1 .. 8):


11082016, 02:01 PM
Post: #15




RE: Short & Sweet Math Challenge #21: Powers that be  
11082016, 02:03 PM
Post: #16




RE: Short & Sweet Math Challenge #21: Powers that be  
11082016, 03:15 PM
Post: #17




RE: Short & Sweet Math Challenge #21: Powers that be
Here are the results of my search on HP71/Emu71:
Order 2 1.61803398875 ok x^2x1 Order 3 1.83928675521 ok x^3x^2x1 1.46557123188 ok x^3x^21 1.32471795724 ok x^3x1 Order 4 1.92756197548 ok x^4x^3x^2x1 1.75487766625 ok x^4x^3x^21 1.61803398875 ok x^4x^3x1 1.46557123188 ok x^4x^2x1 Order 5 1.61803398875 ok x^5x^4x^3+x^2x1 1.32471795724 ok x^5x^41 1.32471795724 ok x^5x^2x1 Order 6 1.83928675521 ok x^6x^5x^4x^2x1 1.61803398875 ok x^6x^5x^4+x^2x1 1.46557123188 ok x^6x^5x^4+x^3x^2+x1 1.61803398875 ok x^6x^5x^3x1 1.46557123188 ok x^6x^5x^21 1.46557123188 ok x^6x^5+x^4x^3x^2x1 1.32471795724 ok x^6x^4x1 Order 7 1.83928675521 ok x^7x^6x^5x^4+x^3x^2x1 1.61803398875 ok x^7x^6x^5+x^2x1 ... 1.32471795724 ok x^8x^6x^2x1 1.32471795724 ok x^8x^5x^4x1 884 candidates (1<root<2) 59 candidates (quasiinteger powers) Quite disappointing since I got only ... 7 unique constants. Here is my HP71 working program: 10 !  SSMC21  20 OPTION BASE 0 @ DIM A(10) 30 C=0 @ C2=0 40 FOR D=2 TO 8 50 DISP "Order";D 60 DIM A(D) @ COMPLEX R(D1) 70 A(0)=1 80 A(D)=1 ! assumed... 90 K=3^(D1) ! numbers of coefficient combinaisons 100 FOR J=0 TO K1 110 L=J 120 ! build the coefficients 130 FOR I=D1 TO 1 STEP 1 140 A(I)=MOD(L,3)1 @ L=L DIV 3 150 NEXT I 160 ! find roots of polynomia A 170 MAT R=PROOT(A) 180 ! DISP "Polynomia";J 190 ! MAT DISP A 200 ! DISP "Roots" 210 ! MAT DISP R 220 X=REPT(R(D1)) 230 IF IMPT(R(D1))=0 AND X>1.000001 AND X<2 THEN GOSUB 300 240 ! PAUSE 250 NEXT J ! next polynomia of order D 260 NEXT D ! next order polynomiae 270 DISP C;"candidates (1<root<2)" 280 DISP C2;"candidates (quasiinteger powers)" 290 END 300 ! evaluate candidate 310 'T': 320 C=C+1 330 ! DISP "Candidate x=";X 340 ! MAT DISP A 350 F=0 ! flag candidate found 360 N=20 370 Y=X^N 380 IF ABS(FP(Y).5)>.49 THEN F=F+1 ELSE F=0 390 N=N+1 400 IF Y<1E10 AND N<80 AND F<3 THEN 370 ! no need to go beyong 1E10 or power 80 410 ! IF X=1.38027756910 THEN PAUSE 420 IF F<3 THEN 530 430 C2=C2+1 440 FIX 11 @ DISP X;"ok"; @ STD 450 DISP " x^";STR$(D); 460 FOR I=1 TO D1 470 IF A(I)=1 THEN DISP "+"; 480 IF A(I)=1 THEN DISP ""; 490 IF A(I)<>0 THEN DISP "x"; 500 IF A(I)<>0 AND DI<>1 THEN DISP "^";STR$(DI); 510 NEXT I 520 DISP STR$(A(D)) 530 ! PAUSE 540 RETURN 

11082016, 05:24 PM
Post: #18




RE: Short & Sweet Math Challenge #21: Powers that be
(11082016 03:29 PM)Mike (Stgt) Wrote: And for a(0) = +1, well... then the restriction of the leading coefficient to +1 is thwarted as we look for the maximum _value_ of the root in the range ]1..2[. No?Well, not necessarily. For instance, 1.32471795724 (one of the constants) is a solution of x^4x^3x^2+1. Anyway, in my code above I assumed a(0)=1: 80 A(D)=1 ! assumed... because all examples posted by Valentin and Gerson are so :) JF 

11082016, 06:52 PM
Post: #19




RE: Short & Sweet Math Challenge #21: Powers that be
(11082016 02:03 PM)JF Garnier Wrote:(11082016 01:38 PM)Paul Dale Wrote: If there is one real root > 1 and all of the remaining (possibly complex) roots have root < 1, then the polynomial is of the required form.It may be the element I was missing, but can you explain or justify this statement? The following is a test for 4th order polynomials, using Pauli's and Mike's conditions. I have assumed the solutions given by PROOT are ordered by magnitude, but I am not sure about that. But this implementation is probably wrong since it gives only two solutions (three when the program is modified for 3rd order polynomials). Hopefully no typing errors since the EMU71 version I have doesn't work on Windows 10 64bit and I don't know how to copy and paste the listing in DosBox.  3 DESTROY ALL 5 INTEGER A,B,C,D,E,N,T 7 A=1 @ E=1 10 OPTION BASE 1 15 INTEGER P(5) 20 COMPLEX R(4) 25 FOR B=1 TO 1 30 FOR C=1 TO 1 35 FOR D=1 TO 1 40 P(1)=A @ P(2)=B @ P(3)=C @ P(4)=D @ P(5)=E 45 MAT R=PROOT(P) 50 IF IMPT(R(4))=0 THEN GOSUB 1000 55 NEXT D 60 NEXT C 70 NEXT B 75 END 1000 IF REPT(R(4))>1 THEN GOSUB 2000 1005 RETURN 2000 T=0 2005 FOR N=1 TO 3 2010 IF ABS(R(N))<1 THEN T=T+1 2015 NEXT N 2020 IF T=3 THEN PRINT REPT(R(4));A;B;C;D;E 2025 RETURN >RUN 1.92756197548 1 1 1 1 1 1.3802775691 1 1 0 0 1  3rd order polynomial solutions: 1.83928675521 1 1 1 1 1.46557123188 1 1 0 1 1.32471795721 1 0 1 1  

11082016, 10:24 PM
Post: #20




RE: Short & Sweet Math Challenge #21: Powers that be
(11082016 01:38 PM)Paul Dale Wrote: I'd figured out an approach to this problem. I've not implemented it on any hardware. My approach is similar to Paul's. One subroutine sequentially cycles through all coefficients  1, 0, 1, starting with 1 1 1 1 1 1 1 1 1. The other routine uses the solver to check for roots in the proscribed range. Coded for the 15C, the routines appear to work independently, but not sure how they will behave together across the whole range. 

« Next Oldest  Next Newest »

User(s) browsing this thread: 1 Guest(s)