Talk:Untouchable numbers: Difference between revisions

From Rosetta Code
Content added Content deleted
Line 40: Line 40:


: Doh, what a stupid and embarrassing mistake! I've reverted the Go entry back to a limit of 1e5 until I can find time to sort out the 1e6 figure. --[[User:PureFox|PureFox]] ([[User talk:PureFox|talk]]) 17:03, 13 February 2021 (UTC)
: Doh, what a stupid and embarrassing mistake! I've reverted the Go entry back to a limit of 1e5 until I can find time to sort out the 1e6 figure. --[[User:PureFox|PureFox]] ([[User talk:PureFox|talk]]) 17:03, 13 February 2021 (UTC)

:Well, after spending some time re-examining this, I think I may have to retire hurt as I'm not going to be able to get anywhere near Julia's time for a million untouchables.

:Reverting to 100,000 untouchables where Julia was using a sieve limit of 32 million, it's time was only about 1.5 minutes on my machine. Rewriting the Go program to use a similar approach and, even after adding some optimizations, the time was about 11.5 minutes which is nearly 8 times slower! So, for a million untouchables and a sieve limit of 512 million, if Julia is taking just under an hour, the Go program will likely take north of 8 hours which is far more than I have the patience for.

:I think the difference must be due to Julia having access to a very fast factorization routine in its Primes package which is then being used to calculate the sum of a number's proper divisors and is taking most of the time here. My own fairly simple approach must be relatively pedestrian.

:FWIW here's the latest Go program if anyone has any ideas about quickening it up:
<lang go>package main

import "fmt"

func sumDivisors(n int) int {
sum := 1
k := 2
if n%2 == 0 {
k = 1
}
for i := 1 + k; i*i <= n; i += k {
if n%i == 0 {
sum += i
j := n / i
if j != i {
sum += j
}
}
}
return sum
}

func commatize(n int) string {
s := fmt.Sprintf("%d", n)
if n < 0 {
s = s[1:]
}
le := len(s)
for i := le - 3; i >= 1; i -= 3 {
s = s[0:i] + "," + s[i:]
}
if n >= 0 {
return s
}
return "-" + s
}

func main() {
const maxTarget = 100_000 // 1_000_000
const sieveLimit = 32_000_000 // 512_000_000
untouchables := make([]bool, maxTarget+1) // all false by default
primes := make([]bool, maxTarget+1) // ditto
for i := 1; i <= maxTarget; i++ {
untouchables[i] = true
}
for i := 2; i <= maxTarget; i++ {
n := sumDivisors(i)
if n <= maxTarget {
untouchables[n] = false
if n == 1 {
primes[i] = true
}
}
}
for i := maxTarget + 1; i <= sieveLimit; i++ {
n := sumDivisors(i)
if n <= maxTarget {
untouchables[n] = false
}
}
for i := 6; i <= maxTarget; i += 2 {
if untouchables[i] && (primes[i-1] || primes[i-3]) {
untouchables[i] = false
}
}

fmt.Println("List of untouchable numbers <= 2,000:")
count := 0
for i := 0; i <= 2000; i++ {
if untouchables[i] {
fmt.Printf("%6s", commatize(i))
count++
if count > 0 && count%10 == 0 {
fmt.Println()
}
}
}

fmt.Printf("\n\n%7s untouchable numbers were found <= 2,000\n", commatize(count))
p := 10
count = 0
for i, u := range untouchables {
if u {
count++
}
if i == p {
cc := commatize(count)
cp := commatize(p)
fmt.Printf("%7s untouchable numbers were found <= %9s\n", cc, cp)
p = p * 10
}
}
}</lang>
--[[User:PureFox|PureFox]] ([[User talk:PureFox|talk]]) 21:03, 13 February 2021 (UTC)


==Nice recursive solution==
==Nice recursive solution==

Revision as of 21:03, 13 February 2021

1464

Unfortunately, the factors of 14641 are {1,11,121,1331} which sum to 1464.
According to https://oeis.org/A005114/b005114.txt
There are 1212 not 1,216 untouchable numbers <= 10,000
There are 13,863 not 13,886 untouchable numbers <= 100,000 --Pete Lomax (talk) 08:50, 9 February 2021 (UTC)

Yes, I see that I need to extend the search for each range of untouchable numbers.   A fix (hopefully) will be forthcoming.   Thanks for finding those errors in the REXX's output.     -- Gerard Schildberger (talk) 09:00, 9 February 2021 (UTC)
That certainly slowed up the search when I extended the search field.     -- Gerard Schildberger (talk) 09:14, 9 February 2021 (UTC)
To match up to 1e5 I increased the limit to 18*n, a tactic I would deem rather unsound. Note added to the Phix entry. --Pete Lomax (talk) 09:17, 9 February 2021 (UTC)


reduction in the number of (counting) ranges

I've reduced the number of ranges for counting untouchable numbers.     It seems the upper limit for searching for touchable numbers was far higher than I guesstimated.     -- Gerard Schildberger (talk) 09:28, 9 February 2021 (UTC)


Number of untouchable numbers up to 1 million

According to one of the OEIS links there are 150,232 untouchable numbers up to 1e6 so it looks like the REXX result (150,233) is off by one here.

I've corrected the REXX output,   the incorrect count was computed with an overreach of   18,   not   20.     -- Gerard Schildberger (talk) 21:28, 11 February 2021 (UTC)


By running the Go program with a limit of 1e6 and a sieve factor of 20 (up from 14 for 1e5), I was in fact able to verify the figure of 150,232. Factors of 18 and 19 both gave 150,233 so there must be an 'awkward' number within this range which takes some eliminating. --PureFox (talk) 09:36, 11 February 2021 (UTC)

I've managed to isolate the 'awkward' number which is 816,422. The (first) number whose proper divisors sum to this is 19,175,641 its divisors being [1, 29, 151, 841, 4379, 22801, 126991, 661229]. This explains why we need a sieve factor between 19 and 20 to catch it since neither 816,421 or 816,419 is prime. --PureFox (talk) 11:33, 11 February 2021 (UTC)

That's mighty good detective work   ...   and a true dedicated ethic.     -- Gerard Schildberger (talk) 11:45, 11 February 2021 (UTC)
Yup, found that too, and wrote this before I saw your post: The factors of 19175641 aka 29*29*151*151 are {1,29,151,841,4379,22801,126991,661229}, aka {1,29,151,29*29,29*151,151*151,29*29*151,29*151*151}, which sum to 816422, in case that helps. Interestingly the "naughty boy" is a product of so few primes, and the original above, 14641 is 11*11*11*11. --Pete Lomax (talk) 11:50, 11 February 2021 (UTC)
Without something reliable to check against, I think that 1 million is probably as far as we can practically go. We could, of course, incrementally increase the sieve factor for higher limits until we appear to get a stable count but, judging by what it said in the Julia entry, there may still be the possibility of some outlier spoiling the party! --PureFox (talk) 14:23, 11 February 2021 (UTC)

Sorry! Both the REXX and Go entries report 150,252 not 150,232 untouchable numbers up to 1 million (ouch!).
I will quickly add that it just took the Julia entry 58mins to get the correct result, on this box. --Pete Lomax (talk) 14:45, 13 February 2021 (UTC)

Doh, what a stupid and embarrassing mistake! I've reverted the Go entry back to a limit of 1e5 until I can find time to sort out the 1e6 figure. --PureFox (talk) 17:03, 13 February 2021 (UTC)
Well, after spending some time re-examining this, I think I may have to retire hurt as I'm not going to be able to get anywhere near Julia's time for a million untouchables.
Reverting to 100,000 untouchables where Julia was using a sieve limit of 32 million, it's time was only about 1.5 minutes on my machine. Rewriting the Go program to use a similar approach and, even after adding some optimizations, the time was about 11.5 minutes which is nearly 8 times slower! So, for a million untouchables and a sieve limit of 512 million, if Julia is taking just under an hour, the Go program will likely take north of 8 hours which is far more than I have the patience for.
I think the difference must be due to Julia having access to a very fast factorization routine in its Primes package which is then being used to calculate the sum of a number's proper divisors and is taking most of the time here. My own fairly simple approach must be relatively pedestrian.
FWIW here's the latest Go program if anyone has any ideas about quickening it up:

<lang go>package main

import "fmt"

func sumDivisors(n int) int {

   sum := 1
   k := 2
   if n%2 == 0 {
       k = 1
   }
   for i := 1 + k; i*i <= n; i += k {
       if n%i == 0 {
           sum += i
           j := n / i
           if j != i {
               sum += j
           }
       }
   }
   return sum

}

func commatize(n int) string {

   s := fmt.Sprintf("%d", n)
   if n < 0 {
       s = s[1:]
   }
   le := len(s)
   for i := le - 3; i >= 1; i -= 3 {
       s = s[0:i] + "," + s[i:]
   }
   if n >= 0 {
       return s
   }
   return "-" + s

}

func main() {

   const maxTarget = 100_000                 // 1_000_000
   const sieveLimit = 32_000_000             // 512_000_000
   untouchables := make([]bool, maxTarget+1) // all false by default
   primes := make([]bool, maxTarget+1)       // ditto
   for i := 1; i <= maxTarget; i++ {
       untouchables[i] = true
   }
   for i := 2; i <= maxTarget; i++ {
       n := sumDivisors(i)
       if n <= maxTarget {
           untouchables[n] = false
           if n == 1 {
               primes[i] = true
           }
       }
   }
   for i := maxTarget + 1; i <= sieveLimit; i++ {
       n := sumDivisors(i)
       if n <= maxTarget {
           untouchables[n] = false
       }
   }
   for i := 6; i <= maxTarget; i += 2 {
       if untouchables[i] && (primes[i-1] || primes[i-3]) {
           untouchables[i] = false
       }
   }
   fmt.Println("List of untouchable numbers <= 2,000:")
   count := 0
   for i := 0; i <= 2000; i++ {
       if untouchables[i] {
           fmt.Printf("%6s", commatize(i))
           count++
           if count > 0 && count%10 == 0 {
               fmt.Println()
           }
       }
   }
   fmt.Printf("\n\n%7s untouchable numbers were found  <=     2,000\n", commatize(count))
   p := 10
   count = 0
   for i, u := range untouchables {
       if u {
           count++
       }
       if i == p {
           cc := commatize(count)
           cp := commatize(p)
           fmt.Printf("%7s untouchable numbers were found  <= %9s\n", cc, cp)
           p = p * 10
       }
   }

}</lang> --PureFox (talk) 21:03, 13 February 2021 (UTC)

Nice recursive solution

In a manner similar to Talk:Numbers_with_prime_digits_whose_sum_is_13#Nice_recursive_solution I present a scheme which requires no guessed maximums and ends when the candidate list N is empty. First I identify the primes less than the maximum value I am interested in (10 for this example) then I combine them as follows.

Important: The minimum value of Σ at level n is greater than the minimum value of Σ at level n-1. The value of Σ increase with depth along any path.

Level=root

Path=1

  Π=na
  Σ=na
Level=1 N=1..10 Untouchables=None
Path=1;2    Path=1;3    Path=1;5    Path=1;7
   Π=2         Π=3         Π=5         Π=7
   Σ=1         Σ=1         Σ=1         Σ=1
Level=2 N=2..10 Untouchables=None
Path=1;2;2  Path=1;2;3  Path=1;2;5  Path=1;2;7  Path=1;3;3  Path=1;3;5  Path=1;3;7  Path=1;5;5  Path=1;5;7  Path=1;7;7
   Π=4         Π=6         Π=10        Π=14        Π=9         Π=15        Π=21        Π=25        Π=35        Π=49
   Σ=3         Σ=6         Σ=8         Σ=10        Σ=4         Σ=9         Σ=11        Σ=6         Σ=13        Σ=8

I have not gained much yet. This elimates P+1s which I could have done by reading the task description. I have proven that 2 is untouchable. I remove 2 and the Σs found from N.

Level=3 N=5;7 Untouchables=2 (because 2 is less than minimum Σ at level 2)
Path=1;2;2;2  Path=1;2;2;3  Path=1;3;3;3
   Π=8           Π=12          Π=27
   Σ=7           Σ=14          Σ=13

Now I gain, there is no need to check i.e. Path=1;2;2;5 because this would add at least 5 to Σ[1;2;2] (3+5) which is greater than the max value in N (7).

Level=4 N=None Untouchables=5 (because 5 is less than minimum Σ at level 3)

N=None so I am done. 2 and 5 have been identified as untouchable. To find larger untouchables, probably time to turn to your poor long suffering computer. Where Talk:Numbers_with_prime_digits_whose_sum_is_13#Nice_recursive_solution may be considered 'Trees 101' This is more difficult to implement. It is described above as a breadth first search, but for for large values it may be better to use depth first. The main saving is in how well any implementation prunes the tree.--Nigel Galloway (talk) 14:55, 11 February 2021 (UTC)

I was thinking along (vaguely) similar lines, seeing it spelt out like that clarified a couple of points:
It is guaranteed to iterate through all composite numbers (as long as we carry on until Σ too big), and we have no interest in (large) primes.
Obviously it would be very silly to use trial division when we already have all the prime factors.
However, I am sure there is an easy/clever way to calculate Σ, from those primes [and/or prev Σ?], but I am ashamed to say I just can't see it.
Is that Σ=14 for Π=12 a typo, should it be 16 ie sum({1,2,3,4,6})? --Pete Lomax (talk) 08:03, 12 February 2021 (UTC)
I've only just realised there are 78,498 primes < 1 million, and if you're looking for any combination of those that sum to < 1 million, that's a pretty darn big search space... --Pete Lomax (talk) 20:57, 13 February 2021 (UTC)