Re: Challenge: Chasing Knuth's Conjecture

[hv]

The code below finds 1..10 in 2 seconds on this machine, having calculated nothing greater than 201!; 1..37 takes 8s, and calculates 366!. Beyond that we start to find some more difficult, but at 37 mins and 24MB it has so far found all but 8 numbers (59, 64, 72, 86, 87, 92, 97, 99) out of 1..99.

Update: still missing (92, 97, 99) after 274 mins (process size 67MB).

The basic ideas are: a) keep a sorted list (@try) of the numbers we've reached but not yet calculated a factorial for; b) at each iteration, take the smallest untried number and find it's factorial, and repeatedly take the square root of the result until we get to a number we've seen before; c) use a binary chop to insert new pending numbers.

#!/usr/bin/perl
use strict;
use bigint; # optionally with eg C< lib => 'Pari' >

my $max = shift || 10;
my @try;

# deal with '1' explicitly
my %seen = (1 => 's');
my $waiting = $max - 1;
print "1 => s\n";

insert(3, '');

while (@try) {
  my $n = shift @try;
  my $s = 'f' . $seen{$n};
  $n->bfac; # $n = ($n)!
  $s = "f$s";
  while (!defined $seen{$n}) { 
    insert($n, $s);
    $n = sqrt($n);
    $s = "s$s";
  } 
}

sub insert {
  my($n, $s) = @_;
  if ($n <= $max) { 
    print "$n => $s\n";
    exit 0 unless --$waiting;
  } 
  $seen{$n} = $s;
  my($min, $max) = (0, scalar @try);
  while ($min + 1 < $max) { 
    my $new = ($min + $max) >> 1;
    if ($try[$new] > $n) {
      $max = $new; 
    } else {
      $min = $new; 
    }
  } 
  ++$min if $min < @try && $try[$min] < $n;
  splice @try, $min, 0, $n;
}

In the output, eg "5 => ssff" means that 5 = int(sqrt(int(sqrt(fact(fact(3)))))).

Hugo

[tall_man]

Here is a variation of that one that deals with log factorials, so it needs only regular arithmetic. It uses the gammln recently posted here by ewijaya.

It is extremely fast, solving 1..99 in about 2 seconds. (There may be roundoff error issues. I have not checked all the answers.)

Update: I added an independent confirmation with Math::Pari for every number I insert in the queue. They all come out correct. With checking added, the solution to 1..99 takes a minute and 23 seconds. The hardest to get is 92, which passes through 12985 factorial.

Read more... (3 kB)

[hv]

Interesting, I would have expected rounding errors large enough to skew the results, but your results agree with mine as far as the final value I've discerned so far in 1..99:

87 => 
ssssssssfsssssssssssfsssssfssssssssssfss
ssssssssssfsssssssssssssfssssfsssssssssf
sssssssssssfsssssssfssssssssfssssssssssf
ssssfsssssssfsssssssssssfssssssssfssssss
ssssfsssssssssssfsssssssfsssssssssfsssss
sfsssssssssfsssssssssfsssssssfssssssfsss
fsssssssssfssssssfssssssfssssfssssssssfs
fssfssssssfsssssfssfsfssff

Hugo

Update 2023-11-27: broke long line

[Roy Johnson]

I took your ideas for avoiding re-testing things and revamped my sieve. This generates all the answers that don't require factorial on anything higher than whatever you set (350 gets me all solutions for 1..30). It runs in about 5 seconds on my 1 GHz PC with that setting.
Read more... (1484 Bytes)
Update: pregenerated factorials to save time; also put in a limit on how high to print, so the huge numbers don't show up. Even taking factorials as high as 999 (takes 30-40 seconds), there's no solution for 38.

[tall_man]

I computed an answer for 38 with my modified version of hv's program, using logs. The result is:

"sssssssssfsssssssssssfssssssfssssssfssssssfssssssssssfssssssssfssssss
+ssfssssssfssssssfssssfssssssssfsfssfssssssfsssssfssfsfssff"

I confirmed this answer with Math::Pari, and it is correct. The highest factorial required is 1861.

[Roy Johnson]

My code looks like it agrees with yours; I set the factorial limit to 5000, and it ran for 20 or 30 minutes, I think. 64 is the first skip at that setting.

I did a simple translation of your Inline C sub into Perl, and it is still very fast. Running it for 64 took 33 seconds, and it popped out:

64 => sssssssssfssssssfsssssssssssssfsssssssssssfssssssssssssfssssssss
+ssfssssssf
ssssssssfsssssssssfsssssssssssfsssssssfssssssfsssssssssssfssssssssfsss
+ssssssssfs
sssssfssssssfssssssfssssssssssfssssssssfssssssssfssssssfssssssfssssfss
+ssssssfsfs
sfssssssfsssssfssfsfssff
biggest tried was 5071

Code follows.

Read more... (3 kB)