Programming Challenge « Next Oldest | Next Newest »

 ▼ Namir Posting Freak Posts: 2,247 Threads: 200 Joined: Jun 2005 09-25-2010, 12:18 PM Hi All folks who are NOT attending HHC 2010, Since you are not busy at the HHC, I thought you can handle the following worthwhile challenge. The Viete's formulas allow you to calculate the coefficients of a polynomial given it's roots. For example for quadratic polynomials: x1 + x2 = - b / a x1 * x2 = c / a Giving the polynomial a x^2 + b X + c. Choosing a = 1 simplifies things. In the case of a cubic polynomial you have: x1 + x2 + x3 = - b / a x1* x2 + x1 * x3 + x2 * x3 = c /a x1 * x2 * x3 = -d / a Giving the polynomial a x^3 + b x^2 + c x + d. Choosing a = 1 simplifies things. Wikipedia shows the general equations for the Viete's formulas here. All polynomials (order 2 ad higher) have the following in common: 1) The coefficient of the n-1 term (n being the polynomial order) is calculated as: sum(x(i)) = - b / a, for i = 1 to n 2) The constant term is calculated as: product(x(i) = (-1)^n * constant_term / a for i = 1 to n Your challenge is to write a program in RPN, RPL, BASIC or any other popular language (for calculator or for PC) that can calculate the coefficients of any polynomial given the array of roots. Using coefficient a as 1 simplifies things a bit. The solution requires that you handle a number of n equations to calculate the coefficients of a polynomial of order n. Namir AWOL from HHC 2010 and hiding in the south of France!! Merci Sarkozi!!! Edited: 25 Sept 2010, 12:31 p.m. ▼ Thomas Klemm Senior Member Posts: 735 Threads: 34 Joined: May 2007 09-25-2010, 05:03 PM The HP-48GX has this function built in: ```[1 2 3] PCOEF -> [1 -6 11 6] ``` You may consider it the inverse operation of PROOT. Best regards Thomas ▼ Namir Posting Freak Posts: 2,247 Threads: 200 Joined: Jun 2005 09-25-2010, 05:38 PM You are not implementing a solution using Viete's formulas. I am aware of function PCOEF and PROOT. Namir Han Senior Member Posts: 709 Threads: 104 Joined: Nov 2005 09-25-2010, 05:59 PM Are we shooting for program size, speed, or simply an implementation of the formula? =) Han Han Senior Member Posts: 709 Threads: 104 Joined: Nov 2005 09-25-2010, 06:24 PM I didn't implement Viete's formula directly, but the algorithm does (in the end) do the exact multiplication and addition used in Viete's formula. Given a list of roots: { x1 x2 x3 ... xn } it computes the list of coefficients. ```<< DUP SIZE -> r n << { 1. } 0. n 1. - FOR i 0. OVER + r n i - GET NEG * SWAP 0. + SWAP ADD NEXT >> >> ``` If anyone would like to see the math behind it, I can post that as well. ▼ Thomas Klemm Senior Member Posts: 735 Threads: 34 Joined: May 2007 09-25-2010, 07:32 PM Same idea but using STREAM: ```\<< {{1}} SWAP + \<< \-> p x \<< p 0 + 0 p + x * - \>> \>> STREAM \>> ``` Best regards Thomas ▼ Namir Posting Freak Posts: 2,247 Threads: 200 Joined: Jun 2005 09-25-2010, 08:42 PM Thomas, Your solution works great too (input has to be a list of roots). Cool!! Namir Namir Posting Freak Posts: 2,247 Threads: 200 Joined: Jun 2005 09-25-2010, 07:55 PM Sure! Please do! I am very impressed with your solution. You solution works great! Thanks! Edited: 25 Sept 2010, 8:48 p.m. ▼ Han Senior Member Posts: 709 Threads: 104 Joined: Nov 2005 09-25-2010, 09:16 PM The idea can be easily seen with the following example. Suppose our roots are 3, 5, and 7. Then the polynomial can be written in factored form is: (x-3)(x-5)(x-7) First multiply the first two factors (x-3)(x-5). (x-3)*x + (x-3)*(-5) = (x^2 - 3x) + (-5x + 15) If we use vector (or list) notation, what we are essentially doing is taking { 1 -3 } (which corresponds to the polynomial x-3) and "shifting" it by a power of x (to get x^2 - 3x). So in list form, this is { 1 -3 0 }. Then, we compute (x-3)*(-5) which gives us the list { 0 -5 15 }. We need a 3-element list though, since ADD requires same list sizes. Lastly, we add the two lists via element-wise addition. This is essentially the sum of two polynomials: one whose coefficient is { 1 -3 0 } and the other is { 0 -5 15 }. So the sum is { 1 -8 15 }, or the polynomial x^2 - 8x + 15. Now, just repeat this process with the final factor: (x-7). It's really just distribution of x and the root x_i (negated). This is best seen when using the DBUG feature of the HP48/49/50. Han Edited: 25 Sept 2010, 9:22 p.m. ▼ Han Senior Member Posts: 709 Threads: 104 Joined: Nov 2005 09-25-2010, 09:33 PM This reminds of another challenge of computing f(n) for some given number n, where f(x) is a polynomial using an RPN calculator that only has 2 levels of stack. The idea is to rewrite the polynomial in a "recursive" manner. For example, if we have 3x^2-5x+6, we can view this as: (((x * 3) - 5 ) * x) + 6 So to compute f(4), for example, we would do: 4 3 * 5 - 4 * 6 + ▼ Namir Posting Freak Posts: 2,247 Threads: 200 Joined: Jun 2005 09-26-2010, 06:14 AM I see .. so you are doing the polynomial multiplication approach. This is the easy route. I was able to write a version of PCOEF in Matlab (similar to Matlab's ROOTS function). The Viete's formulas are a bit tough to code. Hence the challenge!! :-) Namir Olivier De Smet Senior Member Posts: 275 Threads: 41 Joined: Mar 2010 09-26-2010, 11:20 AM Hi, Here is a brutal solution of Viète formula (non optimized, lot of room for optimizing I think) for HP86/7 with nice display of sums. Olivier ```10 DISP "Viete's formulas" @ OPTION BASE 0 20 DIM x(99),a(99),ik(99) 30 RESTORE 1000 @ READ n@ FOR k=1 TO n @ READ x(k)@ NEXT k @ a(n)=1 40 FOR k=1 TO n @ s=0 @ FOR i=1 TO k @ ik(i)=i @ NEXT i @ l=k @ DISP "S";n-k;" ="; 50 p=1 @ FOR i=1 TO k @ p=p*x(ik(i)) @ DISP ik(i);@ NEXT i @ s=s+p 60 ik(l)=ik(l)+1 @ IF ik(l)<= n-k+l THEN DISP "+";@ GOTO 50 70 l=l-1 @ IF l=0 THEN 200 80 IF ik(l)=n-k+l THEN 70 90 j=ik(l)+1 @ FOR i=l TO k @ ik(i)=j @ j=j+1 @ NEXT i @ l=k @ DISP "+";@ GOTO 50 200 a(n-k)=s*((-1)^k) @ DISP "" @ DISP "a(";n-k;")=";a(n-k) @ NEXT k 999 END 1000 DATA 5,1,2,3,4,5 ``` Edited: 26 Sept 2010, 11:47 a.m. ▼ Namir Posting Freak Posts: 2,247 Threads: 200 Joined: Jun 2005 09-26-2010, 01:26 PM Olivier, Quoting the famous French mathematician Bernoulli, let me say "I recognize the lion by its tale!" The quote comes from when Bernoulli challenged the mathematicians of Europe to solve a problem. At that time, Newton was sick, so Bernoulli did NOT send him a copy of the problem. After waiting several months, Bernoulli got frustrated and upon hearing that Newton's health was better, he mailed Newton a copy of the problem. Upon receipt Newton looked at the problem and mailed THE SOLUTION THE NEXT DAY!!! When Bernoulli got an envelope from England that quickly (relative to speed of mail at the time), he knew it was from Newton and said, before opening the envelope, "I recognize the lion by its tale!" I am not surprised that you solved it, since it takes a first class programmer like yourself to solve a tough problem like implementing the Viete's formulas! Thanks!! ▼ Olivier De Smet Senior Member Posts: 275 Threads: 41 Joined: Mar 2010 09-26-2010, 02:13 PM Thanks a lot for your nice comments Namir. Just a little remark ... could it be "paw" in place of "tale" in the quote of Bernouilli ? (see cite) Paw (as "patte" in french) seems better I think. Anyway thanks again for your kind comments. Olivier Edited: 27 Sept 2010, 2:07 a.m. Namir Posting Freak Posts: 2,247 Threads: 200 Joined: Jun 2005 09-26-2010, 05:50 PM I typed in your solution in the EMU71 (I have no physical 71B with me while I am on vacation). I changed the name of the array ik() to i1() (a letter and a digit) to make the code work. Other than that the code is a beauty! I am impressed on how you were able to use a two-level nested FOR loops to manage calculating the various Viete's summations (and consequently the related polynomial coefficients). I kept thinking that the solution needed more nested loops. Cool! Namir Edited: 26 Sept 2010, 5:56 p.m. Olivier De Smet Senior Member Posts: 275 Threads: 41 Joined: Mar 2010 09-27-2010, 08:03 AM Optimizing a bit (remove an unnecessary for-next) ```10 DISP "Viete's formulas" @ OPTION BASE 0@ DIM x(99),a(99),ik(99) 20 RESTORE @ READ n@ FOR k=1 TO n @ READ x(k)@ NEXT k @ a(n)=1 30 FOR k=1 TO n @ s=0 @ l=1 @ ik(1)=0 @ DISP "S";n-k;"="; 40 j=ik(l)+1 @ FOR i=l TO k @ ik(i)=j @ j=j+1 @ NEXT i @ l=k 50 p=1 @ FOR i=1 TO k @ p=p*x(ik(i)) @ DISP ik(i);@ NEXT i @ s=s+p 60 ik(l)=ik(l)+1 @ IF ik(l)<= n-k+l THEN DISP "+";@ GOTO 50 70 l=l-1 @ IF l=0 THEN 90 80 IF ik(l)=n-k+l THEN 70 ELSE DISP "+";@ GOTO 40 90 a(n-k)=s*((-1)^k) @ DISP "" @ DISP "a(";n-k;")=";a(n-k) @ NEXT k 100 END 110 DATA 5,1,2,3,4,5 ``` Olivier ▼ Namir Posting Freak Posts: 2,247 Threads: 200 Joined: Jun 2005 09-27-2010, 01:13 PM Thanks!!! I was able to translate the 71B code into Excel VBA. Had to watch for those zigzagging GOTOs!! Here is the VBA code: ```Option Explicit Option Base 0 ' Version 1.1 Sub CalcCoeff() Dim N As Integer, S As Double, P As Double Dim I As Integer, J As Integer, K As Integer, M As Integer Dim Root(99) As Double, Coeff(99) As Double, IK(99) As Integer Dim bJump As Boolean, Row As Integer, Col As Integer Dim aStr As String N = 2 Do While Cells(N, 1) <> "" Root(N - 1) = Cells(N, 1) N = N + 1 Loop N = N - 2 ' number of roots Coeff(N) = 1 ' display Coeff(N) Cells(2, 2) = Coeff(N) Range("B3:Z1000").Value = "" Row = 3 ' main loop For K = 1 To N Col = 3 S = 0 M = 1 IK(1) = 0 ' Line 40 Do J = IK(M) + 1 For I = M To K IK(I) = J J = J + 1 Next I M = K ' Line 50 Do P = 1 For I = 1 To K P = P * Root(IK(I)) aStr = Cells(Row, Col) If aStr = "" Then Cells(Row, Col) = "r" & CStr(IK(I)) Else Cells(Row, Col) = aStr & "*r" & CStr(IK(I)) End If Next I S = S + P ' L60 IK(M) = IK(M) + 1 ' move column for Viete's terms to the right If IK(M) <= N - K + M Then Col = Col + 1 Loop While IK(M) <= N - K + M bJump = False Do M = M - 1 If M = 0 Then bJump = True Exit Do End If ' move column for Viete's terms to the right If IK(M) <> (N - K + M) Then Col = Col + 1 Loop While IK(M) = N - K + M ' jump out of outer loop? If bJump Then Exit Do Loop ' Line 200 Coeff(N - K) = S * ((-1) ^ K) ' show the coefficient value Cells(Row, 2) = Coeff(N - K) Row = Row + 1 Next K End Sub ``` I replaced variables l, x() and a() with M, Root() and Coeff(), respectively. I also added the variable Row to display the polynomial coefficients in a descending power order. The program expects the roots to appear in column 1 starting with Row 2. The output polynomial coefficients appear in column 2, starting with row 2. Row 1 for columns 1 and 2 are reserved for labeling the columns. The output places the coefficients of the polynomial in a descending power order. Row 2 has 1 (the term for X^n), row 3 has the coefficient for X^(n-1), and so on. Row n+2 has the constant term. The above code also shows the different terms in Viete's formula. These symbolic terms appear in columns C, D, and up. Each cell contains a set of root products. The variable Col controls which column each root product group appears. Namir Edited: 27 Sept 2010, 11:12 p.m.

 Possibly Related Threads... Thread Author Replies Views Last Post HHC 2013 RPN Programming Challenge - Final Thought Jeff O. 7 1,290 10-17-2013, 11:02 AM Last Post: Jeff O. HHC 2012 RPN Programming Challenge Conundrum Jeff O. 15 1,870 10-08-2012, 03:34 PM Last Post: Gerson W. Barbosa Mini-challenge: HHC2012 RPL programming contest with larger input David Hayden 14 1,930 10-05-2012, 10:36 PM Last Post: David Hayden Weekend programming challenge: Euler's Totient function Allen 36 3,928 06-03-2012, 10:39 PM Last Post: Paul Dale A Sunday Programming Challenge Allen 61 5,324 03-17-2012, 06:48 PM Last Post: Allen HHC2012 programming challenge?? Pal G. 28 2,804 09-25-2011, 05:02 PM Last Post: Paul Dale HHC MMC Programming challenge inC code David Hayden 8 1,199 12-31-2010, 09:40 AM Last Post: David Hayden RE: 35s sorting routine challenge - Gene's Challenge Miguel Toro 4 850 08-01-2007, 08:36 AM Last Post: Miguel Toro HP-16C programming challenge - Hamming distance Eric Smith 9 1,199 01-27-2006, 12:53 AM Last Post: Eric Smith a simple programming challenge Don Shepherd 14 1,631 05-31-2005, 01:09 PM Last Post: Valentin Albillo

Forum Jump: