HP Forums

You're currently viewing a stripped down version of our content. View the full version with proper formatting.
Pages: 1 2
Surprisingly (?) every integer N is a factor of an integer P consisting solely of 1's & 0's, for the smallest N the values of P are given here:

https://oeis.org/A004290

The challenge is to write an HP 49G programme to find minimal P/N.

Fastest unerring programme wins.
In base 10, presumably? The math is trivial for base 2!
(01-24-2020 04:06 PM)Gerald H Wrote: [ -> ]Surprisingly (?) every integer N is a factor of a integer P consisting solely of 1's & 0's, for the smallest N the values of P are given here:

https://oeis.org/A004290

The challenge is to write an HP 49G programme to find minimal P/N.

Fastest unerring programme wins.

All positive integer N that can be interpreted as valid binary gives minimal P/N = 1

Example: N=1, 10, 11, 100, 101 ...
(01-24-2020 08:48 PM)ijabbott Wrote: [ -> ]In base 10, presumably? The math is trivial for base 2!

Yes, in base 10, sorry for the omission.
(01-24-2020 11:48 PM)Albert Chan Wrote: [ -> ]
(01-24-2020 04:06 PM)Gerald H Wrote: [ -> ]Surprisingly (?) every integer N is a factor of a integer P consisting solely of 1's & 0's, for the smallest N the values of P are given here:

https://oeis.org/A004290

The challenge is to write an HP 49G programme to find minimal P/N.

Fastest unerring programme wins.

All positive integer N that can be interpreted as valid binary gives minimal P/N = 1

Example: N=1, 10, 11, 100, 101 ...

A quick solution for an infinite number of integers, leaving an infinite number still unresolved.

(01-25-2020 05:39 AM)Gerald H Wrote: [ -> ]
(01-24-2020 11:48 PM)Albert Chan Wrote: [ -> ]All positive integer N that can be interpreted as valid binary gives minimal P/N = 1

Example: N=1, 10, 11, 100, 101 ...

A quick solution for an infinite number of integers, leaving an infinite number still unresolved.

Oh, I see. "minimal P/N" meant min(P)/N, not min(P/N)

Sorry, I don't know HP-49G. Algorithm (not optimized) is like this:

Code:
f = lambda n: n if n<2 else 10*f(n>>1) + (n&1)

def p(n, mod=0, t=1, step=1):
while f(t)%n != mod: t += step
return f(t)

>>> p(85)
111010
>>> 111010 / 85
1306.0

Note: options to p(n) allowed some optimizations.
p(n that can be intepreted as valid binary) = n
Remove n trailing 0's: p(10k) = 10 * p(k) ; must be removed first
Remove n trailing 5's: p(10k + 5) = 10 * p(2k+1)
Remove even n's ﻿ ﻿ ﻿ ﻿ ﻿ ﻿: p(2k) = 10 * p(k)

What is left is n mod 10 = 1,3,7,9 → p(n) must be odd → step = 2 is ok

>>> 10 * p(17, step=2) # = p(85)
111010
Hi Gerald,

Thanks for posting this intriguing challenge.

Here is my attempt which returns a list of A004290(n) from 0 to n. It will only work on numbers up to 98 because A(99) has more than 12 digits. I opted to use approximate numbers because MOD is slow for exact integers and R->B is limited to 12 digits anyway.

The program is still not very fast, taking almost 10 minutes for an input of 98. It could be extended to larger input values with exact integers and some ListExt commands but it would be slower still.

Update: I tried an exact integer version for A(99) and it ran for over 10 minutes on EMU48 without completing. I don't think that my program is incorrect since it worked on smaller numbers.

I then ran Chai Wah Wu's Python program from the OEIS page and it returned A(99) in less than one second. I then tried A(9990) on the Python program and it took over 15 minutes to complete.

I therefor conclude that numbers above 98 are impractical on the HP49/50 without a faster algorithm.

Code:

\<< PUSH BIN I\->R \-> n
\<< 0. 1. 2. n
FOR k 1. 4095.
FOR m m R\->B \->STR 3. OVER SIZE
1. - SUB OBJ\-> DUP k MOD
IF NOT
THEN 4096. 'm' STO
ELSE DROP
END
NEXT
NEXT n 1. + \->LIST POP
\>>
\>>
(01-25-2020 05:58 PM)John Keith Wrote: [ -> ]Chai Wah Wu's Python program from the OEIS page and it returned A(99) in less than one second.
I then tried A(9990) on the Python program and it took over 15 minutes to complete.

You could try the sage code, which is also a valid Python code if "^" replaced as "**"

>>> minmulzo(999) * 10 ﻿ ﻿ ﻿ ﻿ # = p(9990)
1111111111111111111111111110

Above finished in my old Pentium 3 in under 0.1 second
Roughly, this is why above is so fast.

Since N=9990, P must also end in 0, since N divides P.
If N is odd (but not end in 5's, see my post), P must also be odd
To minimize P, we multiply by only a ten, to add a zero to P.
This optimization speedup the already fast minmulzo() by 4X

Instead of doing all the unnecessary base conversions, it build a list of mod values.
The program loops thru the combinations of sum of mods, looking for multiples of N.
Here is what each "bit" contribute to the mod value, for P(999):

>>> for i in range(21): print i, pow(10,i,999)
...
0 1
1 10
2 100
3 1
4 10
5 100
6 1
7 10
8 100
...

Above repeating pattern, P(999) required 27 1's, to have sum of mods = 9*111 = 999 = N
(01-25-2020 11:50 PM)Albert Chan Wrote: [ -> ]Here is what each "bit" contribute to the mod value, for P(999):

>>> for i in range(21): print i, pow(10,i,999)
...
0 1
1 10
2 100
3 1
4 10
5 100
6 1
7 10
8 100
...

We have the same "bit" contribution pattern for N = 333, 111
p(111) = (101*3-1)/9 = 111
p(333) = (103*3-1)/9 = 111111111
p(999) = (109*3-1)/9 = 111111111111111111111111111
We can pull this pattern out as special case, and reduce the combinatorial search.

I got an OK from Gerald for posting this fast p(n) Python code.
For my old Pentium 3, time for p(n = 0 to 9999) in 21 seconds (4 seconds for my laptop) Code:
hd = [  0,   1,  10,  11, 100, 101, 110, 111,
1000,1001,1010,1011,1100,1101,1110,1111]

def p(n, oddp=oddp, k=1):
if n<2: return n
while n%5==0: k *= 10; n //= 10-5*(n&1)
while n&1==0: k *= 10; n >>= 1
x = str(n)                  # n%10 = 1,3,7,9
if x in '139' and x.count(x) == len(x):
return 10**((ord(x)-48)*len(x)) // 9 * k
return oddp(n) * k

oddp(n) search in chunks for 3 'hex digits' (12 decimal digits) at a time.
Since returned p is odd, we only need to check 16×16×8 = 2048 cases / 3 'hex digits'

Code:
def oddp(n, t=0, td=0, mod=0):  # assumed n%10 = 1,3,7,9
h1 = [ -hd[i] % n for i in range(1,16,2)]
h2 = [10000*x % n for x in hd]
h3 = [10000*x % n for x in h2]
while 1:
for m3 in h3:           # mods of 3rd 'hex digit'
m4 = m3 + mod       # adjust for p(n) trillions
for m2 in h2:       # mods of 2nd 'hex digit'
if (m4+m2)%n not in h1: continue
p = ( hd[h1.index((m4+m2)%n)*2+1]
+ 10000*(hd[h2.index(m2)]
+ 10000*(hd[h3.index(m3)])))
if td: p += td
return p        # note: p is odd
while 1:
t += 1
td = int(bin(t)[2:]) * 10**12
mod = int(td % n)
if mod not in bad: break

>>> p(12345) ﻿ ﻿ ﻿ ﻿ ﻿ ﻿ ﻿ ﻿ ﻿ ﻿ ﻿ ﻿ ﻿ ﻿ ﻿ ﻿ # = oddp(12345 // 5) * 10
11000111010L

>>> p(123456) ﻿ ﻿ ﻿ ﻿ ﻿ ﻿ ﻿ ﻿ ﻿ ﻿ ﻿ ﻿ ﻿ ﻿ ﻿# = oddp(123456 >> 6) * 10**6
1110110100111000000L
That's fast code alright, Albert!

Now, to translate to RPL. Hmmm...
(01-27-2020 04:11 PM)Albert Chan Wrote: [ -> ]For my old Pentium 3, time for p(n = 0 to 9999) in 21 seconds (4 seconds for my laptop) To take advantage of Python's O(1) dictionary lookup, here is p(n) dict() based version.
Timings for p(n = 0 to 9999) reduced to Pentium3 = 8.1 s., Laptop = 1.5 s. (≈ 2.5X speed) Code:
hd = [0,1,10,11,100,101,110,111]
hd = hd + [1000+x for x in hd]
hd = hd + [10000+x for x in hd]
hd = hd + [100000+x for x in hd]

def bdec(x):
if x<64: return hd[x]
return hd[x&63] + bdec(x>>6) * 1000000

def oddp(n, t=0, td=0, mod=0):
m = n >> 1                  # assumed n is odd
while m % 5: m += n         # loop upto 4 times
m = m//5                    # = -1/10 (mod n)
h1 = dict(((m-x)%n, x) for x in reversed(hd))
m = 1000000%n
h2 = [m*x%n for x in hd]
while 1:                    # assumed n%10 = 1,3,7,9
for m2 in h2:
if (m2+mod)%n not in h1: continue
m = h1[(m2+mod)%n] + 1000000*hd[h2.index(m2)]
if td: m += td
return 10*m+1
while 1:
t += 1
td = bdec(t) * 10**12
mod = int(td % n)
if mod not in bad: break

def p(n, oddp=oddp, k=1):       # 2nd arg can use minmulzo
if n<2: return n
while n%5==0: k *= 10; n //= 10-5*(n&1)
while n&1==0: k *= 10; n >>= 1
x = str(n)                  # n%10 = 1,3,7,9
if x in '139' and x.count(x) == len(x):
return 10**((ord(x)-48)*len(x)) // 9 * k
return k * oddp(n)

Update: new algorithm, odd p = 10m+1, solve for m ≡ -1/10 (mod n), return 10m+1
Timings for p(n = 0 to 9999) improved to Pentium3 = 4.8 s., Laptop = 0.9 s.
Here my attempt for the 49G, incorporating insights of Albert Chan (as best I could).

Execution time remains lamentably high, improvements expected.

A non-brute force algorithm would be very useful.

Name of programme: A4290

Both

FPTR F 1A is internal SREPL

PTR 2F3A3 is internal SIZE for integers

stable from 1.19-6 to 2.10-8.

Size: 386.5

CkSum: #3E4Dh

Code:

A4290

::
CK1&Dispatch
# FF
::
" MIN K*N ONLY 1'S & 0'S CF A79339"
DispCoord1
SetDA3Temp
::
FPTR2 ^PUSHFLAGS_
DOBIN
BINT64
dostws
DUP
FPTR2 ^Z>S
BINT2
ZERO_DO
INDEX@
#>\$
NULL\$
FPTR F 1A
DROPLOOP
NULL\$?
?SEMI
FPTR2 ^ZTrialDiv2
BINT0
ROT
BEGIN
DUP
ZINT 5
FPTR2 ^ZDIVext
ZINT 0
EQUAL
WHILE
::
SWAPDROPSWAP
#1+SWAP
;
REPEAT
DROPDUP
FPTR2 ^Z>S
"9"
NULL\$
FPTR F 1A
OVER
NULL\$?
ITE
::
ROTDROP
COERCE
BINT9
#*
ZERO_DO
CHR_1
>H\$
LOOP
FPTR2 ^S>Z
;
::
2DROP
DUP
ZINT 1
EQUAL
?SEMI
ZINT 10
OVER
PTR 2F3A3
COERCE
FPTR2 ^RP#
FPTR2 ^Z>S
"# "
SWAP&\$
"b"
&\$
DOSTR>
HXS 00001 1
bit-
BINT0
BEGIN
DROP
HXS 00001 2
bit+
DUP
hxs>\$
BINT3
LAST\$
FPTR2 ^S>Z
DUP
4PICK
FPTR2 ^ZMod
ZINT 0
EQUAL
UNTIL
ROTROT2DROP
;
ZINT 10
2SWAP
#MAX
FPTR2 ^RP#
FPTR2 ^RMULText
;
FPTR2 ^POPFLAGS_
;
;
For input

123456789

a real 49G returns

1111111101

in

6.8 sec.
(01-29-2020 11:18 AM)Gerald H Wrote: [ -> ]For input

123456789

a real 49G returns

1111111101

in

6.8 sec.

Since N is odd, not divisible by 5, we know P last digit is 1

P = 123456789 M ≡ -M ≡ 1 (mod 10)
M ≡ -1 ≡ 9 (mod 10)

123456789 * 9 = 1111111101

Another example, P(17) is not as lucky, but still not too bad:

P = 17 M ≡ -3 M ≡ 1 (mod 10)
M ≡ 1/-3 ≡ -9/-3 ≡ 3 (mod 10)

17 * 3 = 51
17 * 13 = 221, 1001/17 ≈ 58.88
17 * 63 = 1071
17 * 73 = 1241, 10001/17 ≈ 588.29
17 * 593 = 10081
17 * 603 = 10251, 11001/17 ≈ 647.12
17 * 653 = 11101

We could also try P(17) that ends in 01 vs 11, and pick the smaller P:

P = 17*M ≡ 1 or 11 (mod 100)

Using my InverseMod trick:

1
2 ﻿ ﻿ ﻿ ﻿ -floor(1/2*15) = -7, 15*17=255 > 100
15
17 ﻿ ﻿ ﻿ ﻿ --floor(-7/15*100) = -47 = 17-1 (mod 100)
100

M ≡ (1 or 11)/17 ≡ (1 or 11)*(-47) ≡ 53 or 83 (mod 100)

17 * 53 = 901
17 * 83 = 1411, 10001/17 ≈ 588.29
17 * 653 = 11101

2,470,747

to see how long it takes to get computed?
p(2470747) = 110,101,011,111,111,011,101

0.027 seconds for my old Pentium 3
For

2470747

Emu 48 with 49G Rom 2.10-7 on my computer with processor

Intel(R) Core(TM) i7-4790 CPU @ 3.60GHz

returns the same answer as yours in

76 sec.
(01-29-2020 07:25 AM)Gerald H Wrote: [ -> ]A non-brute force algorithm would be very useful.

Sage code is non-brute force.
It generated a list of mod restricted possibilities. Time complexity = O(N log(P))

Example, for p(7), it produced a list, pfrem = [3, 0, 2, 1, 1, 2, 2].
This list, in modulo 7, meant this:
﻿ ﻿ ﻿
1??? ≡ 0
1 ﻿ ﻿ ﻿ ﻿≡ 1
1?? ﻿ ≡ 2
1? ﻿ ﻿ ≡ 3
1? ﻿ ﻿ ≡ 4
1?? ﻿ ≡ 5
1?? ﻿ ≡ 6

? is either 0 or 1, not yet determined. We can solve for P:

P = 1??? ≡ 0 (mod 7)
??? ≡ -1000 ≡ 1 (mod 7)

From the list, 1 ≡ 1 → ??? = 001 → P(7) = 1001

The problem is above required N item list, which may be huge, and takes time to build.
The code shines if N is small, but P is big. Example:

>>> minmulzo(142857) ﻿ ﻿ ﻿ ﻿ ﻿ ﻿ ﻿ ﻿
1010001010001110001110101110101110111111111L

Above finished on my Pentium3 in 3.6 seconds.
Brute force, oddp(142857) would take about 10 hours.

Code:
def minmulzo(n, b=10) :
if n<2 : return n
pfrem = [n] * n
pfrem = 0                # p%n=1 -> p = 1
m, i = 1, b%n
while pfrem == n:        # p%n=0 -> p < 10**n
for r in xrange(1,n) :
if pfrem[r] >= m: continue
if pfrem[r+i-n] == n: pfrem[r+i-n] = m
if pfrem[i] == n: pfrem[i] = m
i = b*i % n
m = m + 1
res = r = 0
while 1:
m = b**pfrem[r]
r = (r-m) % n
res += m
if not r: return res
(01-30-2020 12:18 AM)Albert Chan Wrote: [ -> ]The problem is above required N item list, which may be huge, and takes time to build.

That would be a problem for the HP49 due to limited memory. The largest practical list is around 5000 - 10000 items depending on size of numbers and the amount of free memory in the calculator.

I think GeraldH's time for 2470747 would be hard to beat.
A slightly improved version of A4290, smaller & faster:

Code:
Size: 358.5

CkSum: # 5EABh

::
CK1&Dispatch
# FF
::
" MIN K*N ONLY 1'S & 0'S CF A79339"
DispCoord1
SetDA3Temp
::
FPTR2 ^PUSHFLAGS_
DOBIN
BINT64
dostws
DUP
FPTR2 ^Z>S
BINT2
ZERO_DO
INDEX@
#>\$
NULL\$
FPTR F 1A
DROPLOOP
NULL\$?
casedrop
ZINT 1
FPTR2 ^ZTrialDiv2
BINT0
ROT
BEGIN
DUP
ZINT 5
FPTR2 ^ZDIVext
ZINT 0
EQUAL
WHILE
::
SWAPDROPSWAP
#1+SWAP
;
REPEAT
DROPDUP
FPTR2 ^Z>S
"9"
NULL\$
FPTR F 1A
OVER
NULL\$?
ITE
::
ROTDROP
COERCE
BINT9
#*
ZERO_DO
CHR_1
>H\$
LOOP
FPTR2 ^S>Z
;
::
2DROP
DUP
ZINT 1
EQUAL
?SEMI
%2
OVER
PTR 2F3A3
%^
%>#
HXS 00001 1
bit-
BINT0
BEGIN
DROP
HXS 00001 2
bit+
DUP
hxs>\$
BINT3
LAST\$
FPTR2 ^S>Z
DUP
4PICK
FPTR2 ^ZMod
ZINT 0
EQUAL
UNTIL
ROTROT2DROP
;
ZINT 10
2SWAP
#MAX
FPTR2 ^RP#
FPTR2 ^RMULText
;
FPTR2 ^POPFLAGS_
;
;
Pages: 1 2
Reference URL's
• HP Forums: https://www.hpmuseum.org/forum/index.php
• :