(34S) Prime Factors
05-30-2014, 01:53 AM (This post was last modified: 06-15-2017 01:21 PM by Gene.)
Post: #1
 Dave Britten Senior Member Posts: 1,966 Joined: Dec 2013
(34S) Prime Factors
Got my WP-34S fully built today, so naturally, that means you're getting the usual prime factoring algorithm that I've been dragging along to absolutely everything since I found it in that HP-67 program listing.

This one adds a few enhancements: factors are all obtained up front as quickly as possible, and then a simple menu is displayed allowing the user to view the list with the Up and Down arrow keys.

Usage: Enter number to factor into X, XEQ 'FAC'. Once all factors are found, use the up and down arrows to review the results.

It's pretty quick; factoring 167,699,497 down to 7 * 3,851 * 6,221 takes around 2.5 seconds.

BUG FIX EDIT: Swapped steps 88 and 89 to prevent the VWa+ from disappearing when PSE 99 times out.

Code:
001 LBL'FAC' 002 LocR 006 003 CLα 004 α'FAC' 005 α'TOR' 006 α'S:' 007 x<1? 008 GTO 99 009 FP 010 x!=0? 011 GTO 99 012 RCL L 013 STO .00 014 SQRT 015 IP 016 STO .01 017 0 018 STO .02 019 STO .03 020 2 021 XEQ 03 022 1 023 XEQ 03 024 2 025 XEQ 03 026 2 027 XEQ 03 028 RCL .00 029 x=1? 030 GTO 09 031 LBL 01 032 4 033 XEQ 03 034 2 035 XEQ 03 036 4 037 XEQ 03 038 2 039 XEQ 03 040 4 041 XEQ 03 042 6 043 XEQ 03 044 2 045 XEQ 03 046 6 047 XEQ 03 048 RCL .00 049 x=1? 050 GTO 09 051 RCL .01 052 RCL .02 053 x<=? Y 054 GTO 01 055 RCL .00 056 STO .02 057 0 058 XEQ 03 059 GTO 09 060 LBL 03 061 RCL+ .02 062 STO .02 063 RCL .00 064 x<>Y 065 MOD 066 x!=0? 067 RTN 068 RCL .00 069 RCL/ .02 070 STO .00 071 SQRT 072 IP 073 STO .01 074 INC .03 075 RCL .02 076 STO->.03 077 VWα+ .03 078 0 079 GTO 03 080 LBL 09 081 1 082 STO .04 083 LBL 10 084 CLα 085 α'⇕FC' 086 α'TR ' 087 αRC# .04 088 LBL 11 089 VWα+->.04 090 PSE 99 091 KEY? .05 092 GTO 11 093 5 094 1 095 x=? .05 096 INC .04 097 6 098 1 099 x=? .05 100 DEC .04 101 RCL .04 102 x<1? 103 1 104 x>? .03 105 RCL .03 106 STO .04 107 GTO 10 108 LBL 99 109 CLα 110 α'BAD' 111 α' IN' 112 α'PUT' 113 VWα+ X 114 RTN 115 END

LBL 'FAC', LBL 01, and LBL 03 are all just the standard prime factoring routine. LBL 09 sets up the result viewer, LBL 10 displays the currently selected result, and LBL 11 handles the keyboard for viewing other results.

CAUTION: The factors are stored in global registers starting with 01, and going up to as many as are necessary to hold all factors. This is to allow for use of the interactive factor viewer.

Suggestions and improvements are welcome. This is my first crack at a WP-34S program.
05-30-2014, 03:30 AM
Post: #2
 Thomas Klemm Senior Member Posts: 1,447 Joined: Dec 2013
RE: (WP-34S) Prime Factors
You could use NEXTP to find the next prime.

Cheers
Thomas
05-30-2014, 03:54 AM
Post: #3
 Dave Britten Senior Member Posts: 1,966 Joined: Dec 2013
RE: (WP-34S) Prime Factors
(05-30-2014 03:30 AM)Thomas Klemm Wrote:  You could use NEXTP to find the next prime.

Cheers
Thomas

Now that's interesting. I just replaced the MOD 30 loop with NEXTP to get the next candidate factor, and it took about 10 times as long to factor 167,699,497 - about 20 seconds. Knowing the amount of work NEXTP probably has to do with each call, I guess that's not a huge surprise.
05-30-2014, 08:24 AM
Post: #4
 Thomas Klemm Senior Member Posts: 1,447 Joined: Dec 2013
RE: (WP-34S) Prime Factors
(05-30-2014 03:54 AM)Dave Britten Wrote:
(05-30-2014 03:30 AM)Thomas Klemm Wrote:  You could use NEXTP to find the next prime.

Cheers
Thomas

Now that's interesting. I just replaced the MOD 30 loop with NEXTP to get the next candidate factor, and it took about 10 times as long to factor 167,699,497 - about 20 seconds. Knowing the amount of work NEXTP probably has to do with each call, I guess that's not a huge surprise.

After a 2nd thought: that was maybe a bad idea.
Well, it depends on how NEXTP is implemented.
But then have a look at the code of PF: I was not the first one to come up with this idea.
Just tested 167,699,497 on the iPhone with PF and it was reasonably fast: maybe a second to find the factor 3,851. Is the iPhone so much faster than the original WP-34s?

Cheers
Thomas
05-30-2014, 07:06 PM
Post: #5
 Thomas Klemm Senior Member Posts: 1,447 Joined: Dec 2013
RE: (WP-34S) Prime Factors
(05-30-2014 08:24 AM)Thomas Klemm Wrote:  Well, it depends on how NEXTP is implemented.

PRIME?
Code:
/* Test if a number is prime or not using a Miller-Rabin test */ #ifndef TINY_BUILD static const unsigned char primes[] = {     2, 3, 5, 7, 11, 13, 17, 19,     23, 29, 31, 37, 41, 43, 47, 53, }; #define N_PRIMES    (sizeof(primes) / sizeof(unsigned char)) #define QUICK_CHECK (59*59-1) #endif int isPrime(unsigned long long int p) { #ifndef TINY_BUILD     int i;     unsigned long long int s; #define PRIME_ITERATION 15     /* Quick check for p <= 2 and evens */     if (p < 2)  return 0;     if (p == 2) return 1;     if ((p&1) == 0) return 0;     /* We fail for numbers >= 2^63 */     if ((p & 0x8000000000000000ull) != 0) {         err(ERR_DOMAIN);         return 1;     }     /* Quick check for divisibility by small primes */     for (i=1; i<N_PRIMES; i++)         if (p == primes[i])             return 1;         else if ((p % primes[i]) == 0)             return 0;     if (p < QUICK_CHECK)         return 1;     s = p - 1;     while ((s&1) == 0)         s /= 2;     for(i=0; i<PRIME_ITERATION; i++) {         unsigned long long int temp = s;         unsigned long long int mod = expmod(primes[i], temp, p);         while (temp != p-1 && mod != 1 && mod != p-1) {             mod = mulmod(mod, mod, p);             temp += temp;         }         if(mod!=p-1 && temp%2==0)             return 0;     } #endif     return 1; }

NEXTP
Code:
/* Very minimal routine to return the next prime in sequence  */         XLBL"NEXTPRIME"             FLOOR             x[<=]1?                 JMP prime_2             PRIME?                 INC X             EVEN? prime_loop::            INC X             PRIME?                 RTN             INC X             JMP prime_loop prime_2::        CLx             _INT 2             RTN
05-30-2014, 08:07 PM
Post: #6
 Dave Britten Senior Member Posts: 1,966 Joined: Dec 2013
RE: (WP-34S) Prime Factors
My hunch is that repeatedly calling NEXTP like that causes a ton of duplicated work. I don't imagine there's any sophisticated caching going on with the limited RAM in the 30b hardware. (See: Schlemiel the Painter)

If that's the case, using it in a prime factoring algorithm is probably taking it to O(N * log N) complexity at best as the input argument grows, whereas the MOD 30 loop keeps it O(N) (even if it's not as fully optimal as it could be).
05-30-2014, 08:58 PM
Post: #7
 Thomas Klemm Senior Member Posts: 1,447 Joined: Dec 2013
RE: (WP-34S) Prime Factors
(05-30-2014 08:07 PM)Dave Britten Wrote:  My hunch is that repeatedly calling NEXTP like that causes a ton of duplicated work. I don't imagine there's any sophisticated caching going on with the limited RAM in the 30b hardware. (See: Schlemiel the Painter)

If that's the case, using it in a prime factoring algorithm is probably taking it to O(N * log N) complexity at best as the input argument grows, whereas the MOD 30 loop keeps it O(N) (even if it's not as fully optimal as it could be).

The important thing is: we don't need to know whether a factor is prime or not to figure out if it divides the number. What ever the complexity of PRIME? may be (O(k log3n) according to Wikipedia), it can't be faster than a simple division. Thus if NEXTP has to test all the next odd numbers we don't gain anything.
Using the MOD 30 loop you even gain a factor 1.875 compared to the minimal routine using only odd numbers.
How does your program compare to the built-in program PF?

Kind regards
Thomas
05-30-2014, 09:42 PM
Post: #8
 Dave Britten Senior Member Posts: 1,966 Joined: Dec 2013
RE: (WP-34S) Prime Factors
(05-30-2014 08:58 PM)Thomas Klemm Wrote:  The important thing is: we don't need to know whether a factor is prime or not to figure out if it divides the number. What ever the complexity of PRIME? may be (O(k log3n) according to Wikipedia), it can't be faster than a simple division. Thus if NEXTP has to test all the next odd numbers we don't gain anything.
Using the MOD 30 loop you even gain a factor 1.875 compared to the minimal routine using only odd numbers.
How does your program compare to the built-in program PF?

Kind regards
Thomas

Where is PF hiding? I don't see that in the catalogs or the manual.
05-30-2014, 10:57 PM
Post: #9
 Thomas Klemm Senior Member Posts: 1,447 Joined: Dec 2013
RE: (WP-34S) Prime Factors
(05-30-2014 09:42 PM)Dave Britten Wrote:  Where is PF hiding? I don't see that in the catalogs or the manual.

In the catalog:
LBL'PF'
Lib 0362

Can't you just XEQ it?
BTW: I'm using 3.3T 3606 in the actual version for the iPhone.

Cheers
Thomas
05-30-2014, 11:58 PM (This post was last modified: 05-31-2014 12:46 AM by Dave Britten.)
Post: #10
 Dave Britten Senior Member Posts: 1,966 Joined: Dec 2013
RE: (WP-34S) Prime Factors
(05-30-2014 10:57 PM)Thomas Klemm Wrote:
(05-30-2014 09:42 PM)Dave Britten Wrote:  Where is PF hiding? I don't see that in the catalogs or the manual.

In the catalog:
LBL'PF'
Lib 0362

Can't you just XEQ it?
BTW: I'm using 3.3T 3606 in the actual version for the iPhone.

Cheers
Thomas

I get "No such label" if I try to XEQ it. I'm running 3.2 3472 on a 30b. All I see in Lib are the couple of programs I put there with PSTO.

EDIT: Oh, I see what's going on. The iOS version has a larger "flash" memory, and a few preloaded programs as a result. I'll see if I can key in PF on my 'real' unit and speed test it a bit.
05-31-2014, 01:09 AM
Post: #11
 Dave Britten Senior Member Posts: 1,966 Joined: Dec 2013
RE: (WP-34S) Prime Factors
Alright, the included PF is just using the NEXTP approach with a simple square-root bailout, so it's about as slow to reach the second factor as when I attempted it with NEXTP. It's using a few BACKs/SKIPs to optimize branching, but I don't think that buys it much. It IS only 29 steps, though, so it's pretty space-efficient.

So whoever maintains the iOS version of the emulator, if you're reading this, feel free to replace PF with FAC. Not only is it faster, it's also a nice demo of KEY? to browse stored data.

I think next I'll experiment with using LocR-> to dynamically allocate local regs to avoid overwriting globals. That'll make the factor browser more complicated, since it won't be starting at reg 01, though (but that's assuming redefining LocR even works in the first place).
06-01-2014, 06:11 PM
Post: #12
 Dave Britten Senior Member Posts: 1,966 Joined: Dec 2013
RE: (WP-34S) Prime Factors
So this is funny: if you replace the NEXTP call in PF with a simple INC X, the time to find the second factor of 167,699,497 drops from 22 seconds to about 5 or 6. I didn't bother getting fancier than that, because I'd have to update a bunch of SKIPs/BACKs.
02-23-2015, 08:02 PM (This post was last modified: 02-23-2015 08:22 PM by Dirk E..)
Post: #13
 Dirk E. Junior Member Posts: 5 Joined: Jan 2014
RE: (WP-34S) Prime Factors - performance challenge
With the subroutine-driven MOD30 loop shown here I felt challenged to find some optimal between speed, resources used and rubbish left on the SBR-stack, but it's quite hard to decide...

From looking through the 'PF' in the library I decided to remain with the RPN stack registers. Since the WP 34S offers so much enhanced TEST commands, you just have to remember, when the stack lifts or drops and when not. It reminded me a little an ancient predecessor of the Rubik's cube called "Schiebefax" in german - 4x4 squares less one in a square frame which were to be put into a certain order.

Here's the very shortest (but also slowest) I found:
Code:
WP 34S (3.3 3678)   13.Feb.2015  Example: 7 x 38447 x 62119 - 1.-2.Factor 1'35.3" Find prime factors for (unknown) - Version using stack (X, Y, L) only - the shortest but slowest - trying every number from 2 usage: (unknown) XEQ 'Pr' - R/S [- R/S] ... '1' - Registers: RPN-stack -------------------------------------------------------------------------------- 001  LBL 'Pr'                             ; ( L );  X ,  Y 002  2          ;for even test            ; (   );  2 , 'X' 003  x><y                                 ; (   ); 'X', 'T' 004  SKIP 002              ;main loop: 6(7) lines per 1(1) candidates              ;i.e. 180(210) per 30(30) 005  INC Y      ;next candidate           ; ('X'); 'X/T0', 'T1' 006  x>< L                                ; (   ); 'X', 'T' 007  RCL/ Y                               ; ('X'); 'X/T', 'T' 008   x<y?      ;'X/T'<'T' i.e. 'T'^2>'X' 009   SKIP 005  ; -> end 010  FP?        ;remainder != 0 ? 011  BACK 006   ; -> next test factor 012  VIEW Y 013  STOP       ;show factor 014  BACK 007 015  x>< L                                ; (   ); 'X', 'T' 016  STOP 017  RCL/ X                               ; ('X');  1 , 'T' 018  END

Adding one line to the loop and a few for filtering out the even the required time lowers by a third:
Code:
WP 34S (3.3 3678)   13.Feb.2015  Example: 7 x 38447 x 62119 - 1.-2.Factor 53.2" Find prime factors for (unknown) - Version using stack (X, Y, L) only, testing 2 and every odd number usage: (unknown) XEQ 'Pr' - R/S [- R/S] ... '1' - Registers: RPN-stack -------------------------------------------------------------------------------- 001  LBL 'Pr'                             ; ( L );  X ,  Y 002  2          ;for even test            ; (   );  2 , 'X' 003  x><y                                 ; (   ); 'X', 'T' 004  RCL/ Y                               ; ('X'); 'X/T', 'T' 005  FP? 006  SKIP 004   ;remainder != 0 -> try '3' 007  VIEW Y 008  STOP       ;show factor 2 009  BACK 005   ;try factor again              ;main loop: 7(8) lines per 1(2) candidates              ;i.e. 105(120) per 15(30) 010  INC Y                                ; ('X'); 'X/T0', 'Tx' 011  INC Y      ;next odd candidate       ; ('X'); 'X/T0', 'T1' 012  x>< L                                ; (   ); 'X', 'T' 013  RCL/ Y                               ; ('X'); 'X/T', 'T' 014   x<y?      ;'X/T'<'T' i.e. 'T'^2>'X' 015   SKIP 005  ; -> end 016  FP?        ;remainder != 0 ? 017  BACK 007   ; -> next test factor 018  VIEW Y 019  STOP       ;show factor 2n+1 020  BACK 007 021  x>< L                                ; (   ); 'X', 'T' 022  STOP 023  RCL/ X                               ; ('X');  1 , 'T' 024  END

Still using only x, y and L, but leaving out divisions by multiples of 3, even more speed is gained, but on the cost of almost doubled line count:
Code:
WP 34S (3.3 3678)   13.Feb.2015  Example: 7 x 38447 x 62119 - 1.-2.Factor 34" Find prime factors for (unknown) - Version using stack (X, Y, L) only, trying 2, 3, excluding their multiples usage: (unknown) XEQ 'Prm' - R/S [- R/S] ... '1' - Registers: RPN-stack -------------------------------------------------------------------------------- 001  LBL 'Prm'                            ; ( L );  X ,  Y  002  2          ;for even test            ; (   );  2 , 'X' 003  x><y                                 ; (   ); 'X', 'T' 004  RCL/ Y                               ; ('X'); 'X/T', 'T' 005  FP? 006  SKIP 003   ;remainder != 0 007  VIEW Y 008  STOP       ;show factor 009  BACK 005   ;try factor again 010  CLx                                  ; ('X');  0 , 'T' 011  3                                    ; ('X');  3 , 'T' 012  x=? Y 013  SKIP 003 014  x><y                                 ; ('X'); 'T',  3 015  x>< L                                ; ('T'); 'X',  3 016  BACK 012 017  STO/ Y                               ; ('X');  3 ,  1                 ;main Loop: 13(15) lines for 2(6) candidates                 ;i.e. 65(75) for 10(30) 018  CLx                                  ; ('X');  0 , 'T' 019  4                                    ; ('X');  4 , 'T0' 020  STO+ Y                               ; ('X');  4 , 'T1' 021  x>< L                                ; ( 4 ); 'X', 'T' 022  RCL/ Y                               ; ('X'); 'X/T', 'T' 023  INT? 024  SKIP 011 025  INC Y                                ; ('X');    , 'Tx' 026  INC Y                                ; ('X');    , 'T1' 027  x>< L                                ; (   ); 'X', 'T' 028  RCL/ Y                               ; ('X'); 'X/T', 'T' 029   x<? Y     ; 'X/T'<'T' i.e. 'T'^2>'X' 030   SKIP 008  ; goto End 031  FP? 032  BACK 014 033  VIEW Y 034  STOP       ;show factor 035  BACK 007 036  VIEW Y 037  STOP       ;show factor 038  BACK 016 039  x>< L                                ; ('X/T'); 'X', 'T' 040  STOP 041  RCL/ X                               ; ('X');  1 , 'T' 042  END

Now I became curious whether the MOD 30 algorithm, i.e. leaving out the multiples of 5, too, would be worth the effort. Formerly I was always convinced, the overhead to organize the loop would eat up the speed. Now I tried - and found the way that Dave presented here is definitely the fastest.

Here it is, adapted to use the RPN stack x, y, z, t and L only. Since for a faster end condition test the root of the 'unknown' is also maintained and integer division employed, it uses a few lines and registers more:
Code:
WP 34S (3.3 3678)   17.Feb.2015  Example: 7 x 38447 x 62119 - 1.-2.Factor 28.2" Find prime factors for (unknown) - Version using stack (X, Y, Z, T, L) only,  trying mod 30 chain via a subroutine, holding sqrt(x) for end condition usage: (unknown) XEQ 'Pr' - R/S [- R/S] ... '1' - Registers: RPN-stack (ssize4 dbloff recommended) -------------------------------------------------------------------------------- 001  LBL 'Pr'                             ; ( L );  X ,  Y ,  Z ,  T 002  SKIP 021 003  LBL 10                               ; ('T0'); 'S', 'X', 'X', [sqrt('X')] 004  RCL+ L                               ; ('S'); 'T1', 'X', 'X', [sqrt('X')] 005  RMDR                                 ; ('T');  ? , 'X', [sqrt('X')] 006  x=0?      ;factor found? 007  SKIP 003 008  CLx                                  ; ('T');  0 , 'X', [sqrt('X')] 009  RCL Y     ;duplicate 'X'             ; ('T'); 'X', 'X', [sqrt('X')] 010  RTN       ;-> try next 011  x>< L                                ; ( 0 ); 'T', 'X0', [sqrt('X')] 012  STO/ Y                               ; (   ); 'T', 'X1', [sqrt('X0')] 013  STOP      ;show factor 014  XEQ 11                               ; (   ); 'T', 'X', 'X', [sqrt('X')] 015  BACK 010  ;try factor again 016  LBL 11 017  x>< Z                                ; (   ); [sqrt('X0')], 'X1', 'T' 018  RCL Y                                ; (   ); 'X1', [sqrt('X0')], 'X1', 'T' 019  STO Y     ;duplicate 'X'             ; (   ); 'X1', 'X1, 'X1', 'T' 020  SQRT()                               ; (   ); sqrt('X'), 'X', 'X', 'T' 021  IP        ;end condition value     ;(sqrt('X')); [sqrt('X')], 'X', 'X', 'T' 022  x>< T                                ; (   );, 'T', 'X', 'X', [sqrt('X')] 023  RTN 024  RCL X                                ; (   ); 'X', 'X' 025  1                                    ; (   );  1 , 'X', 'X' 026  XEQ 11                               ; (   );  1 , 'X', 'X', [sqrt('X')] 027  STO L                                ; ( 1 );  1 , 'X', 'X', [sqrt('X')] 028  XEQ 10 029  1                                    ; ('T');  1 , 'X', 'X', [sqrt('X')] 030  XEQ 10 031  2                                    ; ('T');  2 , 'X', 'X', [sqrt('X')] 032  XEQ 10 033  2                                    ; ('T');  2 , 'X', 'X', [sqrt('X')] 034  XEQ 10                               ; ('T'); 'X', 'X', [sqrt('X')] 035  ENTER^                               ; ('T'); 'X', 'X', 'X', [sqrt('X')]           ;main loop 8x 8(10) + 4 lines           ;i.e. 68 (84 w.'LBL'+skipped) for 8 ("10"/30) tests 036  CLx                                  ; ('T');  0 , 'X', 'X', [sqrt('X')] 037  4                                    ; ('T');  4 , 'X', 'X', [sqrt('X')] 038  XEQ 10 039  2                                    ; ('T');  2 , 'X', 'X', [sqrt('X')] 040  XEQ 10 041  4                                    ; ('T');  4 , 'X', 'X', [sqrt('X')] 042  XEQ 10 043  2                                    ; ('T');  2 , 'X', 'X', [sqrt('X')] 044  XEQ 10 045  4                                    ; ('T');  4 , 'X', 'X', [sqrt('X')] 046  XEQ 10 047  6                                    ; ('T');  6 , 'X', 'X', [sqrt('X')] 048  XEQ 10 049  2                                    ; ('T');  2 , 'X', 'X', [sqrt('X')] 050  XEQ 10 051  6                                    ; ('T');  6 , 'X', 'X', [sqrt('X')] 052  XEQ 10                               ; ('T'); 'X', 'X', [sqrt('X')] 053  RCL L                                ; ('T'); 'T', 'X', 'X', [sqrt('X')] 054  x<? T     ;'X'>'T'^2 -> continue 055  BACK 019 056  x>< Y                                ; ('T'); 'X', 'T', 'X', [sqrt('X')] 057  STOP 058  RCL/ X                               ; ('X');  1 , 'T', 'X', [sqrt('X')] 059  END

The time used by the equivalent routine, using global registers instead, is the same.

While testing the speed I usually terminated after the second factor and started over from the beginning - soon I found out that the WP 34S keeps it's subroutine stack, differently from HP 35s and HP 42s that purge it upon a manual XEQ, GTO or OFF.
Thus I thought it would be preferrable to avoid subroutines during which the program could be interrupted and terminated - I tried to use indirect addressing for a very short subroutine returning just the step between the division 'candidates' 1 - 2 - 2 - 4 - 2 - 4 - 2 - 4 - 6 - 2 - 6 (then continuing repeating the 8 elements starting from the first '4'). And I tried another way doing just an indirect branch, since the step value would be needed always at the same position in the main loop. Both were much slower than the following MOD 6 approach, that saves some time by employing integer arithmetics and keeping the root for the end test:
Code:
WP 34S (3.3 3678)   17.Feb.2015  Examples: 7 x 38447 x 62119 - 1.-2.Factor 29.5"    or 33013x33013 1.Factor 25.7 s Find prime factors for (unknown) - Version holding root(x) for end condition  using stack (X, Y, Z, T, L) only usage: (unknown) XEQ 'Prm' - R/S [- R/S] ... '1' - Registers: RPN-stack (ssize4 dbloff recommended) -------------------------------------------------------------------------------- 001  LBL 'Prm' 002  RCL X     ;push stack up       ; (   ), 'X', 'X' 003  2         ;test-factor 'T'     ; (   ), 'T', 'X', 'X'  004  XEQ 01                         ; (   ), 'T', 'x', 'x', sqrt('x') 005  RMDR                           ; ('T'), ? , 'X', sqrt('X') 006  x!=0? 007  SKIP 012 008  x>< L                          ; ( 0 ), 'T', 'X0', sqrt('x0") 009  STO/ Y                         ; (   ), 'T', 'x1', sqrt('x1') 010  STOP      ;show factor 2 or 3 011  BACK 007 012  LBL 01                         ; (   ), 'T', 'X1', sqrt('X0') 013  x>< Z                          ; (   ), sqrt('X0'), 'X1', 'T' 014  RCL Y                          ; (   ), 'X1', sqrt('X0'), 'X1', 'T' 015  STO Y                          ; (   ), 'X1', 'X1', 'X1', 'T' 016  SQRT()                         ; ('x1'), sqrt('x1'), 'x1', 'x1', 'T' 017  IP                             ; (sqrt'x1'), [sqrt('x1')], 'x1', 'x1', 'T' 018  x>< T                          ; (   ), 'T', 'x1', 'x1', [sqrt('x1')] 019  RTN 020  CLx                            ; ('T'), 0 , 'X', sqrt('X') 021  3                              ; ('T0'), 3 , 'X', sqrt('X') 022  x!=? L 023  BACK 019 024  STO/ L                         ; ( 1 ), 3 , 'X', sqrt('X')                ;main Loop: 14(16) lines for 2(6) candidates                ;i.e. 70(80) for 10(30) 025  CLx                            ; ('T'), 0 , 'X', sqrt('X') 026  RCL Y     ;duplicate           ; ('T'), 'X', 'X', sqrt('X') 027  4                              ; ('T0'),  4 , 'X', 'X', sqrt('X') 028  RCL+ L    ;test-factor         ; ( 4 ), 'T1', 'X', 'X', sqrt('X') 029  RMDR                           ; ('T'), ? , 'x', sqrt('x') 030  x=0? 031  SKIP 014 032  CLx                            ; ('T'), 0 , 'x', sqrt('x') 033  RCL Y     ;duplicate           ; ('T'), 'x', 'x', sqrt('x') 034  2                              ; ('T0'),  2 , 'x', 'x', sqrt('x') 035  RCL+ L    ;test-factor         ; ( 2 ), 'T1', 'X', 'X', sqrt('X') 036   x>? T     ;'T'>sqrt('X')? 037   SKIP 013  ;-> End 038  RMDR                           ; ('T'), ? , 'x', sqrt('x'), sqrt('x') 039  x!=0? 040  BACK 015 041  x>< L                          ; ( 0 ), 'T', 'X0', sqrt('X0') 042  STO/ Y                         ; ( 0 ), 'T', 'X1', sqrt('X0') 043  STOP      ;factor 6n-1 044  XEC 01                         ; (   ), 'T', 'x1', 'x1', sqrt('x1') 045  BACK 009 046  x>< L                          ; ( 0 ), 'T', 'X0', sqrt('X0') 047  STO/ Y                         ; ( 0 ), 'T', 'X1', sqrt('X0') 048  STOP      ;factor 6n+1 049  XEQ 01                         ; (   ), 'T', 'x1', 'x1', sqrt('x1') 050  BACK 021 051  x><y      ;END                 ; (   ), 'X', 'T', 'X', sqrt('X') 052  STOP 053  RCL/ X                         ; ('X'),  1 , 'T', 'X', sqrt('X') 054  END

This is almost as fast as the MOD 30 above and even 5 lines shorter. It also demonstrates that despite the main loop running one line more than the MOD 6 program with RPN stack x, y, L only cited above, the integer arithmetics RMDR and comparison with an integer gain a measurable amount of time.

For the use of indirectly addressed local registers can be reviewed here, I'd like to show here my experiment using another indirect addressing for the factor chain step.
I found it not too clear in the manual, but I understand that the global registers are indirect 00-99, the 12 stack and helper registers indirect 100-111, and local registers starting .00 become indirect 112...
Here it is, the local registers are released before showing a factor that was found, then re-allocated; anyway, if the program is interrupted, they'd need to be freed manually:
Code:
WP 34S (3.3 3678)   8.Feb.2015  Example: 7 x 38447 x 62119 - 1.-2.Factor 32" Find prime factors for (unknown) - Version using indirect local registers for MOD 30 chain  and holding root(x) for end condition usage: (unknown) XEQ 'PR' - R/S [- R/S] ... '1' - Registers: 00 - unknown 'X', 01 - tester 'T', 02 - sqrt('X'), 03 - pointer RPN-stack, Local registers .00-.10 - factor sequence steps -------------------------------------------------------------------------------- 001  LBL 'PR' 002  STO 00 003  # 111      ;loop end value 004  STO 03     ;counter/pointer 005  3 006  10^x 007  STO/ 03 008  # 122      ;loop start value 009  STO+ 03 010  2 011  STO 01     ;first test candidate 012  SKIP 004 013  RCL 01 014  PopLR      ;free local registers 015  STO/ 00    ;leave remaining nnn 016  STOP       ;show found factor 017  RCL 00 018  SQRT() 019  IP 020  STO 02     ;root(nnn) for end check 021  LocR 011   ;build factor sequence steps 022  1 023  STO .10 024  4 025  STO .07 026  STO .05 027  STO .03 028  2 029  STO .09 030  STO .08 031  STO .06 032  STO .04 033  STO .01 034  6 035  STO .02 036  STO .00 037  RCL 00     ;main loop 7x8(9)+12(14) lines 038  RCL 01     ;i.e. 68 for 8(10/30) candidates 039  RMDR 040  x=0? 041  BACK 028   ;factor found 042  RCL ->03 043  STO+ 01    ;next candidate 044  DSE 03 045  BACK 008 046  8 047  STO+ 03 048  RCL 01 049  x<=? 02    ;factor<=squrt(nnn)? 050  BACK 013   ;-> continue 051  RCL 00     ;last Factor 'X' 052  BACK 038 053  END

Maybe, this reply is a bit late, but I received my WP 34S only these weeks, and it was also my first program to do my own prime factorization - thank you for the interesting challenge with the 'PF' and the MOD 30 extention!
Hopefully someone finds something interesting in it, too!

On my somewhat slower android mobile on free42 I started exploring the WP34S's integer space and factorized 2^50+3 = 7 x 3 435 997 x 46 811 113.
This ran about 3 minutes(!) till the second factor on the mobile, while the test 7 x 38447 x 62119 ran rather 3 seconds on it. Thus I'd be a little afraid running out of batteries or patience factorizing any double precision integer on the WP 34S.
The 'PF' from the library ran more than 5 minutes factorizing 33013x33013 when I lost patience and stopped it before it found the first factor...

Cheers!
Dirk
 « Next Oldest | Next Newest »

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