Thread: Prime Gap Theory View Single Post
2020-10-30, 15:39   #109
robert44444uk

Jun 2003
Oxford, UK

1,933 Posts

Quote:
 Originally Posted by henryzz Very interesting work on optimizing the choice of d#. Can't wait to try this new code. This has made me realize that I don't know the recommended scripts for gap searching currently. https://mersenneforum.org/showthread.php?t=20872 could do with updating for the post Thomas Nicely era and adding some links to recommended software/scripts/guides would be useful for any new participants. Is https://mersenneforum.org/showpost.p...&postcount=103 still the latest script? Is ATH's program from that thread still recommended? I have improvements to both that I can share if they aren't obsolete. I realize that some of this will probably be replaced by Seth's work.
There are two major approaches at the moment.
• Deficient primorials of the form A*p#/d form the majority of gap records, given that it is at a midpoint of a relatively smooth range 4 merits wide, and a truly efficient perl code is shown below - thanks to danaj.
• The other approach is primorial offsets, of the form A*p#+/-offset+/-midpoint where the offset is something that provides for a range if integers that are relatively smooth

Code:
#!/usr/bin/env perl
use warnings;
use strict;
#use Timer::Runtime;
use Math::GMPz;
use Math::BigFloat lib=>"GMP";
use Math::Prime::Util qw/:all/;
$|=1; my$nthreads = enter the number of threads you machine can cancel here;
my $startprim = shift || enter xth prime here; my$endprime = shift || enter yth prime here;
my $beg = shift || enter A here; my$end = shift || enter upper end A of search here;
my $mingaplen = enter the smallest gap you want to see here as a merit; my$divisor   = enter your divisor here;
my $delta = we usually put 4 here, but it will not search this further if a prime comes up in a range of merits + or - delta/2 ; my %merit; my %mindig; open(my$merits, '<', 'merits.txt') or die "Cannot open merits file\n";
while (<$merits>) { next unless /^\s*(\d+)\s+(\S+)/;$merit{$1} =$2;
$mindig{$1} = int( $1 / ($2 * 2.3025) + 1.01 );  # Slightly large
}
close($merits); my @acc = (1,5,7,11,13,17,19,23,25,29); #my @acc = (0..29); my @results :shared; foreach my$pp ($startprim ..$endprime) {
my $prim = pn_primorial($pp);
$prim = Math::BigInt->new("$prim") unless ref($prim) eq 'Math::BigInt'; my ($sn, $rem) =$prim->copy->bdiv($divisor); next unless$rem == 0;
next if ($sn *$end) < 4e18;

my $snlog = length($sn) * 2.3026;   # Approximately one merit
my $mingap = int($mingaplen * $snlog);$mingap = 1432 if $mingap < 1432; my$prev_thresh = int( $delta *$snlog );

my $nth = nth_prime($pp);
print "-$nth$pp $divisor ($beg..$end) \n"; my @threads; push @threads, threads->create('searchi',$_, $sn,$nth,$prev_thresh,$mingap) for 0 .. $nthreads-1;$_->join() for (@threads);

while (@results) {
print shift(@results);
}
}

sub merit {
my($n,$gap) = @_;
my $fgap = Math::BigFloat->new("$gap");
my $fn = Math::BigFloat->new("$n");
return sprintf("%.6lf", $fgap / log($fn));
}

sub searchi {
my($tnum,$sn,$nth,$prev_thresh,$mingap) = @_;$sn = Math::GMPz->new("$sn"); my$beg30 = int($beg/30); my$index = $tnum; while (1) { if ($tnum == 0 && @results) {
lock(@results);
while (@results) {
print shift(@results);
}
}
my $i = 30*($beg30 + int($index/10)) +$acc[$index % 10]; # my$i = 30*($beg30 + int($index/30)) + $acc[$index % 30];
$index +=$nthreads;
last if $i >$end;
my $n =$sn * $i; my($dprev,$dnext)=Math::Prime::Util::GMP::surround_primes($n,$prev_thresh); next if$dprev == 0 || $dnext == 0; # prime found inside thresh my$gap = $dprev +$dnext;
next unless $gap >$mingap;             # gap is obviously too small
my $gbeg =$n - $dprev; next if defined$mindig{$gap} && length($gbeg) > $mindig{$gap};
my $merit = merit($gbeg, $gap); next if defined$merit{$gap} &&$merit{$gap} >=$merit;
if (!is_extra_strong_lucas_pseudoprime($gbeg) || !is_extra_strong_lucas_pseudoprime($gbeg+$gap)) { print "# SPSP2 found for$i*$nth#/$divisor\n";
($dprev,$dnext) = Math::Prime::Util::GMP::surround_primes($n);$gap = $dprev +$dnext;
$gbeg =$n - $dprev;$merit = merit($gbeg,$gap);
next if defined $merit{$gap} && $merit{$gap} >= $merit; } { lock(@results); push @results, sprintf("%9.6f %8d PRP%d = %d*%d#/%d-%d\n",$merit, $gap, length($n), $i,$nth, $divisor,$dprev);
}
}
1;
}