The Last Prime Digit

April 16, 2019

We don’t want to store lots of prime numbers, only to throw them away, so we generate the primes and calculate the last digits on the fly, trading time for space:

(define (grimes n)
  (let ((ps (primegen)) (counts (make-matrix 4 4 0)))
    (ps) (ps) (ps) ; discard 2,3,5
    (let loop ((prev (ps)) (n (- n 4)))
      (if (zero? n) counts
        (let ((curr (modulo (ps) 10)))
          (matrix-set! counts (quotient prev 3) (quotient curr 3)
            (+ 1 (matrix-ref counts (quotient prev 3) (quotient curr 3))))
          (loop curr (- n 1)))))))

This calculation took a few minutes:

> (grimes 10000000)
#(#(446808 756071 769923 526953)
  #(593195 422302 714795 769915)
  #(639384 681759 422289 756851)
  #(820368 640076 593275 446032))

You can run the program at https://ideone.com/IcT9ib.

Pages: 1 2

4 Responses to “The Last Prime Digit”

  1. Zack said

    Here is my take in Julia (1.0.2):

    lastdigit(n::Int64) = n % 10

    for i = 2:N
    ld1 = lastdigit(P[i-1])
    ld2 = lastdigit(P[i])
    a = b = 0

    if ld1 == 1
        a = 1
    elseif ld1 == 3
        a = 2
    elseif ld1 == 7
        a = 3
    elseif ld1 == 9
        a = 4
    end
    
    if ld2 == 1
        b = 1
    elseif ld2 == 3
        b = 2
    elseif ld2 == 7
        b = 3
    elseif ld2 == 9
        b = 4
    end
    
    if a*b != 0
        Z[a,b] += 1
    end
    

    end

    Z /= N

    4×4 Array{Float64,2}:
    0.0462304 0.0742944 0.0750461 0.0544234
    0.0601098 0.0444256 0.0704369 0.075029
    0.0637398 0.0675519 0.0443935 0.0743187
    0.0799143 0.0637294 0.0601274 0.0462292

    On another note, I wonder if Dr. Grimes changed his surname so that it rhymes with “primes” because he’s so intrigued by them…

  2. Paul said

    Sourcecode in Python is on ideone.
    Here are some results:

    primes < 1_000_000
     0.041   0.080   0.083   0.045
     0.056   0.036   0.074   0.085
     0.065   0.069   0.037   0.079
     0.089   0.065   0.056   0.040
    
    10_000_000
     0.043   0.078   0.080   0.049
     0.058   0.039   0.073   0.080
     0.064   0.069   0.039   0.078
     0.085   0.065   0.058   0.042
    
    100_000_000
     0.044   0.076   0.078   0.052
     0.059   0.042   0.072   0.078
     0.064   0.068   0.042   0.076
     0.083   0.064   0.059   0.044
    
    1_000_000_000
     0.046   0.075   0.076   0.054
     0.060   0.044   0.071   0.076
     0.064   0.068   0.044   0.075
     0.080   0.064   0.060   0.046
    
    10_000_000_000
     0.047   0.074   0.074   0.055
     0.060   0.046   0.070   0.074
     0.064   0.067   0.046   0.074
     0.079   0.064   0.060   0.047
    
     350_000_000_000 15045.708864189168
     0.049   0.072   0.072   0.057
     0.061   0.048   0.069   0.072
     0.063   0.067   0.048   0.072
     0.077   0.063   0.061   0.049
    
    1_140_000_000_000 32903.831111384796
     0.049   0.072   0.072   0.057
     0.061   0.048   0.069   0.072
     0.063   0.066   0.048   0.072
     0.076   0.063   0.061   0.049
    
     3_130_000_000_000 88038.33787216355
     0.050   0.072   0.071   0.057
     0.061   0.049   0.069   0.071
     0.063   0.066   0.049   0.072
     0.076   0.063   0.061   0.050
    
     4_240_000_000_000 118131.35975553434
     0.050   0.072   0.071   0.057
     0.061   0.049   0.068   0.071
     0.063   0.066   0.049   0.072
     0.076   0.063   0.061   0.050
    
  3. Paul said

    I let it run longer and with primes upto 10_900_000_000_000 the result only changes a little.

     10_900_000_000_000 298241 seconds
     0.050   0.071   0.071   0.058
     0.061   0.049   0.068   0.071
     0.063   0.066   0.049   0.071
     0.075   0.063   0.061   0.050
    
  4. Daniel said

    Here’s a solution in C++. It uses the non-segmented Sieve of Eratosthenes, limiting how large the primes can get before the system runs out of memory.

    #include <cstdlib>
    #include <cstring>
    #include <iomanip>
    #include <iostream>
    #include <string>
    #include <vector>
    
    using namespace std;
    
    typedef unsigned long long ull;
    
    int main(int argc, char* argv[]) {
      if (argc != 2) return EXIT_FAILURE;
      ull n = stoull(argv[1], NULL, 10);
      vector<bool> sieve(n + 1, true);
      sieve[0] = false;
      sieve[1] = false;
      for (ull i = 2; i * i <= n; ++i) {
        for (ull j = i * 2; j <= n; j += i) {
          sieve[j] = false;
        }
      }
      int idxs[10] = {-1, 0, -1, 1, -1, -1, -1, 2, -1, 3};
      ull total = 0;
      ull counts[4][4];
      memset(counts, 0, sizeof(counts));
      // ignore 2 and 5 by initializing 'last' to 3 and
      // starting iteration at 7.
      int last = 3;
      for (ull i = 7; i <= n; ++i) {
        if (!sieve[i]) continue;
        int current = i % 10;
        ++counts[idxs[last]][idxs[current]];
        ++total;
        last = current;
      }
      cout << fixed;
      cout << setprecision(3);
      for (int i = 0; i < 4; ++i) {
        for (int j = 0; j < 4; ++j) {
          if (j > 0) cout << " ";
          cout << (float)counts[i][j] / total;
        }
        cout << endl;
      }
      return EXIT_SUCCESS;
    }
    

    Example Usage:

    $ ./praxis 1000000
    0.041 0.080 0.083 0.045
    0.056 0.036 0.074 0.085
    0.065 0.069 0.037 0.079
    0.089 0.065 0.056 0.040
    

Leave a Reply

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out /  Change )

Google photo

You are commenting using your Google account. Log Out /  Change )

Twitter picture

You are commenting using your Twitter account. Log Out /  Change )

Facebook photo

You are commenting using your Facebook account. Log Out /  Change )

Connecting to %s

%d bloggers like this: