Short & Sweet Math Challenge #21: Powers that be

[Valentin Albillo]

 
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/HP-programming 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 brand-new 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 20-strong batch.

Now as for this present new challenge, first a short prelude. Many of you will surely have heard about Phi, the ubiquitous so-called 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^2-x-1 (i.e.: Phi^2-Phi-1 = 0), which happens to be the lowest-degree 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 quasi-integer 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^3-x-1 (i.e.: Hpi^3-Hpi-1 = 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 HP-10C or better, hardware/software based HP-calc 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 quasi-integer-powers 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 HP-calc 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 12-line (60-statement, 585-byte) solution for the HP-71B, which runs like this (no spoilers):

>RUN
       ... some lines of output ...

       1.57367896839    x^8-x^7-x^6+x^2-1
       1.59000537390    x^7-x^5-x^4-x^3-x^2-x-1
       1.60134733379    x^7-x^6-x^4-x^2-1
       1.61803398875    x^2-x-1

       ... 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 HP-71B solutions typically may use the Math ROM, HP-IL 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 non-HP 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.
.

[Gerson W. Barbosa]

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.*X-1.' 
                   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 

[Massimo Gnerucci]

Welcome back Valentin and S&SMC!

[J-F Garnier]

Hi Valentin,

I'm very happy to read you again in this forum !

One question:
From the output examples here:

(11-02-2016 01:13 AM)Valentin Albillo Wrote:  
1.57367896839 x^8-x^7-x^6+x^2-1
1.59000537390 x^7-x^5-x^4-x^3-x^2-x-1
1.60134733379 x^7-x^6-x^4-x^2-1

it seems that null polynomial coefficients are allowed, not only +1 or -1,
Is it correct?

J-F

[Bill (Smithville NJ)]

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

[Valentin Albillo]

.
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:

(11-02-2016 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.*X-1.'              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

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 J-F 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.
.

[Valentin Albillo]

 
Hi, Massimo !:

(11-02-2016 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.
 

[Valentin Albillo]

 
Hi, Jean-François !:

Thanks a lot for your kind comment and continued interest in my posts, much appreciated.

(11-02-2016 01:48 PM)J-F Garnier Wrote:  One question:
From the output examples here:

(11-02-2016 01:13 AM)Valentin Albillo Wrote:  
1.57367896839 x^8-x^7-x^6+x^2-1
1.59000537390 x^7-x^5-x^4-x^3-x^2-x-1
1.60134733379 x^7-x^6-x^4-x^2-1

it seems that null polynomial coefficients are allowed, not only +1 or -1,
Is it correct?

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" ! ... 8-D

Best regards.
V.
 

[Valentin Albillo]

 
Hi, Bill !

Thank for your kind welcome and appreciation, I'm glad you like my S&SMC's

(11-02-2016 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.
 

[Gerson W. Barbosa]

(11-02-2016 10:08 PM)Valentin Albillo Wrote:  .
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:

(11-02-2016 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.*X-1.'              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

Hello, Valentin!

Thanks for starting Season 2. I'm looking forward for the next episodes. I've always appreciated your insightful and well-thought 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 near-integers 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 built-in 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.

(11-02-2016 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 J-F 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

[J-F Garnier]

(11-02-2016 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" ! ... 8-D

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...

J-F

[Valentin Albillo]

 
Hi, J-F:

(11-06-2016 05:21 PM)J-F 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.

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.

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 ... 8-D  )

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.
 

[J-F Garnier]

(11-07-2016 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^4-x^3-1=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 HP-42S (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 :-)

J-F

[Paul Dale]

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):

• Iterate over all polynomials for the given length that satisfy the specified constraints. That the leading coefficient is 1, the constant term is ±1 and the remaining term coefficients are from {-1, 0, 1}.
• Over all lengths there are 6560 such polynomials. If the constant term can also be 0, the count is 9840.
• Find the roots of each polynomial.
• 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.
  

[J-F Garnier]

(11-08-2016 12:31 PM)Mike (Stgt) Wrote:  ... and - regarding the polynom as Sum from i=0 to <=8 of a(i)*x^i - the coefficient a(0) always -1.

I'm suspecting that a(0) is always -1 (for solutions of the given challenge), but I don't understand your argument.
J-F

[J-F Garnier]

(11-08-2016 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?
J-F

[J-F Garnier]

Here are the results of my search on HP71/Emu71:
Order 2
1.61803398875 ok x^2-x-1
Order 3
1.83928675521 ok x^3-x^2-x-1
1.46557123188 ok x^3-x^2-1
1.32471795724 ok x^3-x-1
Order 4
1.92756197548 ok x^4-x^3-x^2-x-1
1.75487766625 ok x^4-x^3-x^2-1
1.61803398875 ok x^4-x^3-x-1
1.46557123188 ok x^4-x^2-x-1
Order 5
1.61803398875 ok x^5-x^4-x^3+x^2-x-1
1.32471795724 ok x^5-x^4-1
1.32471795724 ok x^5-x^2-x-1
Order 6
1.83928675521 ok x^6-x^5-x^4-x^2-x-1
1.61803398875 ok x^6-x^5-x^4+x^2-x-1
1.46557123188 ok x^6-x^5-x^4+x^3-x^2+x-1
1.61803398875 ok x^6-x^5-x^3-x-1
1.46557123188 ok x^6-x^5-x^2-1
1.46557123188 ok x^6-x^5+x^4-x^3-x^2-x-1
1.32471795724 ok x^6-x^4-x-1
Order 7
1.83928675521 ok x^7-x^6-x^5-x^4+x^3-x^2-x-1
1.61803398875 ok x^7-x^6-x^5+x^2-x-1
...
1.32471795724 ok x^8-x^6-x^2-x-1
1.32471795724 ok x^8-x^5-x^4-x-1
884 candidates (1<root<2)
59 candidates (quasi-integer 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(D-1)
70 A(0)=1
80 A(D)=-1 ! assumed...
90 K=3^(D-1) ! numbers of coefficient combinaisons
100 FOR J=0 TO K-1
110 L=J
120 ! build the coefficients
130 FOR I=D-1 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(D-1))
230 IF IMPT(R(D-1))=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 (quasi-integer 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 D-1
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 D-I<>1 THEN DISP "^";STR$(D-I);
510 NEXT I
520 DISP STR$(A(D))
530 ! PAUSE
540 RETURN

[J-F Garnier]

(11-08-2016 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^4-x^3-x^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 :-)

J-F

[Gerson W. Barbosa]

(11-08-2016 02:03 PM)J-F Garnier Wrote:  

(11-08-2016 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?
J-F

The following is a test for 4-th 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 3-rd order polynomials).

Hopefully no typing errors since the EMU71 version I have doesn't work on Windows 10 64-bit 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
--------------

3-rd order polynomial solutions:

1.83928675521 1 -1 -1 -1
1.46557123188 1 -1 0 -1
1.32471795721 1 0 -1 -1

-------------

[Dwight Sturrock]

(11-08-2016 01:38 PM)Paul Dale Wrote:  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):

• Iterate over all polynomials for the given length that satisfy the specified constraints. That the leading coefficient is 1, the constant term is ±1 and the remaining term coefficients are from {-1, 0, 1}.
• Over all lengths there are 6560 such polynomials. If the constant term can also be 0, the count is 9840.
• Find the roots of each polynomial.
• 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.
  

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.