(deleted post)



#29

This Message was deleted. This empty message preserves the threading when a post with followup(s) is deleted. If all followups have been removed, the original poster may delete this post again to make this placeholder disappear.


#30

Quote:
I've made lots of tests with random polynomials (of degree 10 and 20, and with coefficients -100..100) and calculated the optimal parameters for the iteration process. With the current (internal) setting of these parameters there was not one single failure for a few thousand polynomials, so I would say the program is now working quite well.

Tests using random polynomials are no good to check polynomial root finders, on the contrary you should check very specific, difficult cases which actually test the limitations and weak points of your particular implementation.

For instance, try this one and see how it fares:

x^30+30 x^29+435 x^28+4060 x^27+27405 x^26+142506 x^25+593775 x^24+2035800 x^23+5852925 x^22+14307150 x^21+30045015 x^20+54627300 x^19+86493225 x^18+119759850 x^17+145422675 x^16+155117520 x^15+145422675 x^14+119759850 x^13+86493225 x^12+54627300 x^11+30045015 x^10+14307150 x^9+5852925 x^8+2035800 x^7+593775 x^6+142506 x^5+27405 x^4+4060 x^3+435 x^2+30 x+1 = 0

then, you might consider to try this one as well:

x^20+10 x^19+55 x^18+210 x^17+615 x^16+1452 x^15+2850 x^14+4740 x^13+6765 x^12+8350 x^11+8953 x^10+8350 x^9+6765 x^8+4740 x^7+2850 x^6+1452 x^5+615 x^4+210 x^3+55 x^2+10 x+1 = 0

Best regards from V.

Edited: 5 Nov 2011, 7:36 p.m.


#31

Quote:
For instance, try this one and see how it fares:

x^30+30 x^29+435 x^28+4060 x^27+27405 x^26+142506 x^25+593775 x^24+2035800 x^23+5852925 x^22+14307150 x^21+30045015 x^20+54627300 x^19+86493225 x^18+119759850 x^17+145422675 x^16+155117520 x^15+145422675 x^14+119759850 x^13+86493225 x^12+54627300 x^11+30045015 x^10+14307150 x^9+5852925 x^8+2035800 x^7+593775 x^6+142506 x^5+27405 x^4+4060 x^3+435 x^2+30 x+1 = 0


Hi Valentin,

I must say: really very useful and 'realistic' examples! :-(

I've tried this 30-degree polynomial with 4 different 'professional' polynomial solvers on the web and all 4 programs give completely different results.

One of these programs has even several (7) different methods and each method gives also different results!

This last mentioned program can be found here:

http://www.hvks.com/Numerical/websolver.php

BTW, my program is even better than all others that I've tested: it just give NO solution (which is indeed better than wrong solutions). ;-)

So what do you really expect from a small polynomial solver made for a simple calculator???

And one other point: for EVERY program which solves any problem with approximate methods you can always find input values which make the program fail (I know from your AMx matrices that you're an expert in this) - however the question is: how useful (or realistic) are such 'strange' and 'artificial' conditions?

So what to do if EVERY approximate algorithm fails for SOME conditions? Should you avoid writing such a program at all, just because it doesn't work 100% perfectly?

Well, maybe it would be the best to remove my program again ... :-(

Franz


#32

Quote:
Well, maybe it would be the best to remove my program again ..

Please don´t ¡ Valentin's counter-examples don't signify a lack of merit from your program's part but only that the topics addressed are deeper than it appears on first (or second, or third) inspection.

As you say they are a stretch for the typical daily use of the system at hand, yet they are a wonderful challenge to the algorithms - and that's what we're also delighted to have in this forum :-)

I don't have a 34S (imagine, the 41 system quenches my thrist, how's that for (su)real-life systems :-) but otherwise I'd've saved and tried your program, thanks for sharing it.

Best,
'AM


#33

Ángel, why not try the Windows emulator?

#34

Let me add a few remarks.

  • I had to re-assemble the code to make it run on the emulator v 2.2.1789. The documentation does not mention the version required by the .dat file.

  • The documentation also says that only the numbered registers are used. This is not correct - the program uses J and K as well.

  • I do not like the way complex roots are displayed. They real part is labelled x_n while for the imaginary part x_n+1 is used. That's why I modified the program the following way:
       current        improved
    -----------------------
    x1 x1
    -1,4918 -1,4918
    [R/S]
    Re x2 Re x2,3
    -0,8058 -0,8058
    [R/S]
    Im x3 Im x2,3
    1,2229 1,2229
    [R/S]
    Re x4 Re x4,5
    0,5517 0,5517
    [R/S]
    Im x5 Im x4,5
    1,2533 1,2533
  • There are some hard-coded constants like 1E-12 and 1E-16. At Lbl 85 a value less than 1E-12 is rounded to zero. Why is this done? Step 241 tests if the current value is less than 1E-16. I wonder if this will work for any order of magnitude a solution may have. What if x is > 10^6 and 1 ULP equals 1E-9?
Dieter

#35

Quote:
Let me add a few remarks.

Fine - and let me add a few remarks to the remarks, too. ;-)
Quote:
I had to re-assemble the code to make it run on the emulator v 2.2.1789. The documentation does not mention the version required by the .dat file.

You're right - I didn't mention it because I always use the latest SVN build. And of course everyone could (or should) recompile it for their own WP34s version.
Quote:
The documentation also says that only the numbered registers are used. This is not correct - the program uses J and K as well.

No, I didn't write anywhere that ONLY the mentioned numbered registers are used, I just listed only the NUMBERED registers.

And BTW, not even your list is complete, because despite of J and K the program even uses X,Y,Z and T. ;-)
Quote:
I do not like the way complex roots are displayed. They real part is labelled x_n while for the imaginary part x_n+1 is used. That's why I modified the program the following way:

Well, I had a few different version of displaying complex solutions, and one of them was exactly the way you do it now (with n,n+1). But this needed a few more instructions, and since space is 'expensive' on the WP34s I did it the easier way.
Quote:
There are some hard-coded constants like 1E-12 and 1E-16. At Lbl 85 a value less than 1E-12 is rounded to zero. Why is this done?

Well, this 1e-12 is the precision of the calculated roots (and can be changed by the user while running the program), and thus I round
solutions to zero if their absolute value is less then this calculated accuracy. This makes results look much nicer, especially if the imaginary part is 'almost' 0 (within the calculation accuracy).
Quote:
Step 241 tests if the current value is less than 1E-16. I wonder if this will work for any order of magnitude a solution may have. What if x is > 10^6 and 1 ULP equals 1E-9?

This 1e-16 is in fact not very important and I doubt that this value is reached at all during any calculations - it's just to avoid division by zero (or overflows by VERY small denominators).
In different implementations of the Bairstow method I've even seen values for this from 1e-10 to 1e-99.

Franz


#36

Quote:
No, I didn't write anywhere that ONLY the mentioned numbered registers are used, I just listed only the NUMBERED registers.

Let's see:
Used Flags/Registers/Labels:
----------------------------
Flags: 00
Registers: R0-Rn ... Polynomial coefficients a(n) - a(0)
R(n+1)-R(2n+1) ... Polynomial backup copy
R80-R94 ... internal usage
R95 ... maxdeg (30)
R96 ... eps (1e-12)
R97 ... randfact (3)
R98 ... newinit (25)
R99 ... maxiter (500)
Labels: 11,12,13 & 80-89
This does not leave much room for interpretation. ;-) It's a simple documentation error. No big deal, but it should get fixed. Simply add J and K to the "internal usage" line. And while we're at it: the backup of the polynomial coefficents uses the registers up to R(2n+2).
Quote:
Well, I had a few different version of displaying complex solutions, and one of them was exactly the way you do it now (with n,n+1). But this needed a few more instructions, and since space is 'expensive' on the WP34s I did it the easier way.

I do not think that 10 or 15 steps will prevent the program from fitting into the available memory. If everything else fails lines 002-007 may be removed.

Finally, I am not sure if rounding everything less than 1E-12 to zero is a good idea. But I have to think harder before I can comment on this. B-)

It's a nice piece of software. Why not make it perfect?

Dieter


#37

Quote:
It's a simple documentation error.

Well, not mentioning something (or everything) isn't what I would call a 'documentation error'.
Quote:
And while we're at it: the backup of the polynomial coefficents uses the registers up to R(2n+2).

For even degrees my R(2n+1) is correct - for odd degrees (where I add one additional degree) it would even be R(2n+3). But I guess it would even confuse most users if I would describe so many internal details. ;-)
Quote:
Finally, I am not sure if rounding everything less than 1E-12 to zero is a good idea. But I have to think harder before I can comment on this. B-)

What sense would a result of x0=1e-14 make when the calculation accuracy is only 1e-12???

And I hate 'pseudo-complex' results with imaginary parts of 1e-13,14,15 ...
Quote:
It's a nice piece of software. Why not make it perfect?

Well, if you tell me what 'perfect' is, I'll think about it ... :-)

Franz.

#38

This Message was deleted. This empty message preserves the threading when a post with followup(s) is deleted. If all followups have been removed, the original poster may delete this post again to make this placeholder disappear.


#39

Quote:
an error in line 134: the SKIP 16 was one step too far,

I know I'm asking too much. ;-) This type of error can be easily avoided by using either a label (will cost a step) or the preprocessor with its symbolic labels. Both approaches would make it easier to maintain the software, especially for others.

As you might have noticed, I'm working on local registers and on a REGS command that will reduce the number of available global registers in favour of more local registers (for deeper nesting). In such a scenario it's probably better to work with lower numbered or local registers instead of the fixed high numbered ones. Data interactively entered by a user is better stored in global registers but all temporary data that is only used while the program is busy may well be put in locals.

Just a few remarks for future developments.


#40

Quote:
Just a few remarks for future developments.

Yes Marcus, you're of course right, but: ;-)

1) you know why I don't/can't use the preprocessor

2) I would really like to try all those new features like LOCL etc., but unfortunately I don't even see any of them mentioned in the new manual 3.0, so I have no clue WHAT can be done or HOW it has to be done. :-(

Franz


#41

Quote:
I don't even see any of them mentioned in the new manual 3.0, so I have no clue WHAT can be done or HOW it has to be done

Just you wait. Like many other people I've to work for a living. So I deal with the WP 34S in my spare time, where this project competes with some other activities. Be assured the documentation update is in progress and will probably be published before your Polynomial Root Solver is perfect ;-)

Edited: 6 Nov 2011, 6:26 p.m.


#42

Quote:
Just you wait. Like many other people I've to work for a living. So I deal with the WP 34S in my spare time, where this project competes with some other activities. Be assured the documentation update is in progress and will probably be published before your Polynomial Root Solver is perfect ;-)

Well, it was no complaint but just a (true) statement.

And my PolyRoots IS already perfect. ;-)


#43

When the new features become usable and stable enough I will discuss them here. If you want to play around with them, here is a short guide:

- In a subroutine, issue LOCL nn where nn denotes the highest numbered local register you want to use. Newly created registers are zeroed.

- In the same or any nested subroutine you can now address the newly created registers with a dot in the name: .00 to .nn. For direct addressing, .nn is limited to .15 (these are restrictions in the internal coding of commands.) but indirect addresses may go higher, starting at 112 for .00.

- Another call to LOCL in the same subroutine just changes the size of the local register space. LOCL in any nested subroutine will establish a new frame which hides the local registers one level up.

- Use MEM? to find out about the number of levels left. Each frame needs 2 levels + 4 for every local register. The space is taken from the return stack which extends downwards into program space. If you need many local variables, run your program off flash memory and clear the program in RAM.

- Be aware that any action that clears the return stack will remove all local variables: Manual RTN, editing your program and some more keyboard actions. This is the reason why data that is input manually by a user is better held in the normal register space.

- Each frame establishes 16 local flags .00 to .15 in addition to the registers. They can be used like any other flag.

- Use RTN from the same subroutine that has issued LOCL to get rid of the most recently allocated frame and return to the caller. Use LPOP to just delete the variables without returning.

For examples, just look at xrom.wp34s.


#44

Thanks for this explanation, Marcus!

At a first look all this seems to be a bit complicated, so I'll stay away from these new features for a while.

But now I have an other problem with the current SVN 1241: I guess it should be back-compatible to V2.2 if you don't use any of these 'local' features at all, but now with SVN 1241 almost nothing seems to work anymore. :-(

I've recompiled my yesterdays PRS (PolyRoots) program with SVN 1241, but all I get are all kinds of error messages:

1) for example a R-COPY with 0.1111 in X gives "out of range"

2) if I run my program with degree n=5, I get "domain error" after the input, and when I check the registers (00-06 should hold the coefficients) after this error stop, each of these 7 registers contain "nonnumeric".

No chance to get the program working at all.

There is something seriously broken in the latest build ...

Franz


#45

I'm in the middle of something... Thanks for your reports on errors.

Can you 'reset' the emulator from the menu and try again?

EDIT: I just did some tests.

You need to change the argument to R-COPY from .1111 to .11011. The reason is that register index values now are 3 digits in length to allow for the local registers to be addressed. The R-COPY syntax is hence: sss.ccddd for source/count/destination. The count is still limited to 2 digits. Sorry for any inconvenience.

What statement throws the "domain error"?

Edited: 7 Nov 2011, 11:19 a.m. after one or more responses were posted


#46

Quote:
Can you 'reset' the emulator from the menu and try again?

Reset? Whenever I use a new SVN I delete the file wp34s.dat before, so it automatically does a 'Reset'.

Or do you mean I should try it with the newest SVN 1842? But from what I see in the 'changes' I can't imagine that this new build would solve the problems I've mentioned above.

You must know that I can't try it here one my old internet notebook (because of Win98), so I always have to go to my other computer and copy the new version from the USB stick ...


#47

Please see my edited post above, the syntax of R-COPY has changed.


#48

Quote:
Please see my edited post above, the syntax of R-COPY has changed.

Ok, that's not nice - requires a lot of changes in almost every program! :-(

And about "Domain error": well, IIRC it was just any arithmetic operation which used one of those registers containing 'nonnumeric'.

Edited: 7 Nov 2011, 11:25 a.m.


#49

Quote:
And about "Domain error": well, IIRC it was just any arithmetic operation which used one of those registers containing 'nonnumeric'.

I have no idea where the 'nonnumeric' data might come from.

I see that your trick in line 199 no longer works. It's the art of magic digit juggling that fails here. I know that you are swapping source and destination here, but in fact you assign the count to the destination and the former source to the count, then setting the destination to zero. This makes a lot of assumptions about the values and their grouping and works only if two of the values are equal and the other is zero. This makes it almost impossible to use other index values are adopt the new alignment. A more classical approach (get the values on the stack in the proper ordering and then build the descriptor for R-COPY) would be more flexible and better to read.

<compute src>
RCL X because SRC=CNT
CLx destination
FS?C 00
x[<->] Z swap source and destination
SDR 03 move destination beyond decimal point
+ add count
SDR 02 move both count and destination to the right
+ add source
R-COPY
I know this is a bit longer but it's much more readable (and should work in the new version.)


Edited: 7 Nov 2011, 11:53 a.m.


#50

Quote:
I have no idea where the 'nonnumeric' data might come from.

I see that your trick in line 199 no longer works. It's the art of magic digit juggling that fails here.


Sorry, but what 'trick' do I use there??? It's just shifting the decimal and using LastX register.

Ok, my 2-digit shifting may not be correct anymore with the new SVN (and it's no problem for me to change this), but I don't see any reason to completely change this subroutine to your way.

And the (currently) wrong shifting should only give an error for wrong R-COPY parameters, but this does still not explain why any registers contain 'nonnumeric'.

Edit: now I see why you have problems with my R-COPY parameter calculation! I didn't write this subroutine for 'general' usage but only for this program, and here the original polynomial always starts at R00 (and both have of course the same number of elements).

Edited: 7 Nov 2011, 12:15 p.m.


#51

Quote:
Ok, my 2-digit shifting may not be correct anymore with the new SVN (and it's no problem for me to change this), but I don't see any reason to completely change this subroutine to your way.


To copy 07~13 into 00~06 you need 7.07000 in X.

To copy 00~06 into 07~13 you need 0.07007 in X.

I see no way to transform the former into the latter with a single shift operation.

Quote:
but this does still not explain why any registers contain 'nonnumeric'.

I wanted to test the situation here and was unable to do so. I need to fix the R-COPY problem first.

#52

There is a shorthand: You can leave out the count if it can be determined automatically which is the case here because source and destination are adjacent.

This should work:

RCL 94
INC X
FS?C 00
SDR 05

#53

Quote:
There is a shorthand: You can leave out the count if it can be determined automatically which is the case here because source and destination are adjacent.

LOL, so far about me using tricks! ;-)

But your method seems quite suspicious to me!? For example if X=15 (without any decimals) how should R-COPY determine how many registers to copy? And the same after SDR 05, i.e. X=0.00015.

The only thing I can imagine is that the 34s then uses 15 also for the count, but that's nowhere described in the manual (and also not really obvious).


#54

I can't speak for the manual but from what I read in the sources, R-COPY tries to determine the count from the largest non overlapping area between source and destination. This will work nicely here. :-)

I've just committed a new bugfix and ran your program. I didn't check the solutions but it returned some. :-)


#55

Quote:
I've just committed a new bugfix and ran your program. I didn't check the solutions but it returned some. :-)

Well, but here's already the next problem:

Just tried the new SVN 1843 with the changed R-COPY code, and now I see that also R-CLR doesn't work anymore.

Examples: 80.2 with R-CLR clears only R80-R89, and 80.03 clears only R80.

This change of ss.ddnn to ss.dddnn was no good idea, it has broken some other things, too.


#56

Quote:
Examples: 80.2 with R-CLR clears only R80-R89, and 80.03 clears only R80.

This change of ss.ddnn to ss.dddnn was no good idea, it has broken some other things, too.


I didn't change the parameters to R-CLR. If 80.2 does not clear 80 to 99 than that's my fault. I'll check it later. The ddd parameter is at the end of the fractional part. The length change was inevitable to support addresses 100 and above.

EDIT: It was a programming error and should be fixed now. One of the many optimizations to save precious code space.


Edited: 7 Nov 2011, 4:00 p.m.


Possibly Related Threads...
Thread Author Replies Views Last Post
  (deleted post) deleted 2 944 11-03-2013, 06:24 AM
Last Post: deleted
  trig scales on the Post Versalog slide rule Al 12 2,621 09-15-2013, 06:01 AM
Last Post: John I.
  (deleted post) deleted 19 3,057 08-31-2013, 03:46 AM
Last Post: pascal_meheut
  (deleted post) deleted 2 938 03-25-2013, 06:41 PM
Last Post: Eric Smith
  (deleted post) deleted 15 2,810 03-09-2013, 01:35 PM
Last Post: Mark Hardman
  New Blog Post: Length of a Polynomial Segment Eddie W. Shore 1 837 01-17-2013, 08:56 PM
Last Post: Namir
  (deleted post) deleted 5 1,280 01-07-2013, 09:39 AM
Last Post: Frank Boehm (Germany)
  Resuming an old post......................while crunching numberes aurelio 11 2,247 09-25-2012, 08:31 AM
Last Post: Frido Bohn
  (deleted post) deleted 7 1,563 09-18-2012, 11:11 AM
Last Post: René Franquinet
  HP41 post replacement experiment Bruce Larrabee 14 2,201 07-16-2012, 06:51 PM
Last Post: Bruce Larrabee

Forum Jump: