HP 49G: Minimum Multiplier M of Integer N such that M*N Consists Only of 1's & 0's

[Gerald H]

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.

[ijabbott]

In base 10, presumably? The math is trivial for base 2!

[Albert Chan]

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

[Gerald H]

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

[Gerald H]

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

Please supply algorithm (programme) to effect your partial solution.

[Albert Chan]

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

Please supply algorithm (programme) to effect your partial solution.

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

[John Keith]

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:


@ spoiler alert!










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

[Albert Chan]

(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

[Albert Chan]

(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) = (10^1*3-1)/9 = 111
p(333) = (10^3*3-1)/9 = 111111111
p(999) = (10^9*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[0] in '139' and x.count(x[0]) == len(x):
        return 10**((ord(x[0])-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]
    bad = {}
    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
        bad[mod] = 1
        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

[John Keith]

That's fast code alright, Albert!

Now, to translate to RPL. Hmmm...

[Albert Chan]

(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]
    bad = {}
    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
        bad[mod] = 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[0] in '139' and x.count(x[0]) == len(x):
        return 10**((ord(x[0])-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.

[Gerald H]

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_
  ;
;

[Gerald H]

For input

123456789

a real 49G returns

1111111101

in

6.8 sec.

[Albert Chan]

(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

[Gerald H]

Albert, can you please try your programme on

2,470,747

to see how long it takes to get computed?

[Albert Chan]

p(2470747) = 110,101,011,111,111,011,101

0.027 seconds for my old Pentium 3

[Gerald H]

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.

[Albert Chan]

(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[1] = 0                # p%n=1 -> p = 1
    m, i = 1, b%n
    while pfrem[0] == 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

[John Keith]

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

[Gerald H]

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_
  ;
;