meditationatae

Just another WordPress.com site

PARI/gp code to find abundant products of primorials

Abundant here refers to the the sum of divisors function sigma(n) relative to the size of n, a large positive integer. J. Lagarias found a simple criterion equivalent to the Riemann Hypothesis involving “extremely large” values of sigma(n) relative to n.  I’d say the criterion is in the same vein as Gronwall’s Theorem as well as Robin’s inequality and is known as Robin’s Theorem (after Guy Robin, 1984).

In computer emperiments using PARI/gp, it seems that products of primorials are good candidates for extreme values of sigma(n).  The challenge is in finding suitable n.  Today, I’ll post just the PARI/gp command-line, which runs over about six lines, plus the output so far.  I’m approaching using the primes to the first 50,000 after which I’ll have to use larger arrays than 50000.  For the primorials data, there are currently only 22 non-trivial primorials in the product of primorials, so that takes up very little space.  

 

 

    ? for(YYY=1,200,for(ZZZZ=1,5,thebest=0;www3=vector(100);for(Y=1,100,www3[Y]=www[Y]);stdwww3=vector(50000);for(Y=1,50000,stdwww3[Y]=sum(Z=1,22,www3[Z]>(Y-1)));cst=ratio(stdwww3);print(”          “,cst);for(WW=1,1,for(X=1,22,www2=vector(100);for(Y=1,100,www2[Y]=www[Y]);www2[X]=www2[X]+1;stdwww2=vector(50000);for(Y=1,50000,stdwww2[Y]=sum(Z=1,22,www2[Z]>(Y-1)));newr2=ratio(stdwww2);if(newr2>cst,print(X,”         “,newr2));if(newr2>cst,cst=newr2;thebest=X)));print(“Best:  “,thebest,”  Score:  “,cst);if(thebest>0,www[thebest]=www[thebest]+1));print(www))
          9999.1894989261621907312811884519039384
1         9999.1895061636623047308531394503822126
Best:  1  Score:  9999.1895061636623047308531394503822126
          9999.1895061636623047308531394503822126
1         9999.1895137712762868807522356383356695
Best:  1  Score:  9999.1895137712762868807522356383356695
          9999.1895137712762868807522356383356695
1         9999.1895208275164828986454332589360438
Best:  1  Score:  9999.1895208275164828986454332589360438
          9999.1895208275164828986454332589360438
1         9999.1895275298966623386738974278834245
Best:  1  Score:  9999.1895275298966623386738974278834245
          9999.1895275298966623386738974278834245
1         9999.1895341416878657804050573572990624
Best:  1  Score:  9999.1895341416878657804050573572990624
[46949, 175, 30, 12, 7, 5, 4, 3, 3, 2, 2, 2, 2, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]
 
 
 
 
          9999.1957308870381038582045326605722729
1         9999.1957422461064549342801633081826736
Best:  1  Score:  9999.1957422461064549342801633081826736
          9999.1957422461064549342801633081826736
1         9999.1957538370670765035292311754897536
Best:  1  Score:  9999.1957538370670765035292311754897536
          9999.1957538370670765035292311754897536
1         9999.1957656599046314379409343968226769
Best:  1  Score:  9999.1957656599046314379409343968226769
          9999.1957656599046314379409343968226769
1         9999.1957775867461696658468560724706557
Best:  1  Score:  9999.1957775867461696658468560724706557
[47588, 176, 30, 12, 7, 5, 4, 3, 3, 2, 2, 2, 2, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]
          9999.1957775867461696658468560724706557
1         9999.1957893618826497594344654203795770
Best:  1  Score:  9999.1957893618826497594344654203795770
          9999.1957893618826497594344654203795770
1         9999.1958006018159804827750178847894967
2         9999.1958009928724465797876465647760308
Best:  2  Score:  9999.1958009928724465797876465647760308
          9999.1958009928724465797876465647760308
1         9999.1958124548874291084631469937406003
Best:  1  Score:  9999.1958124548874291084631469937406003
          9999.1958124548874291084631469937406003
1         9999.1958242765714072965695512855941375
Best:  1  Score:  9999.1958242765714072965695512855941375
          9999.1958242765714072965695512855941375
1         9999.1958350518377948890924384068961494
Best:  1  Score:  9999.1958350518377948890924384068961494
[47592, 177, 30, 12, 7, 5, 4, 3, 3, 2, 2, 2, 2, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]
          9999.1958350518377948890924384068961494
1         9999.1958461867485862439632965225049802
Best:  1  Score:  9999.1958461867485862439632965225049802
          9999.1958461867485862439632965225049802
1         9999.1958573617577219934809842008260444
Best:  1  Score:  9999.1958573617577219934809842008260444
          9999.1958573617577219934809842008260444
1         9999.1958676822814811026506173924886487
Best:  1  Score:  9999.1958676822814811026506173924886487
          9999.1958676822814811026506173924886487
1         9999.1958774679169080712113587517544242
Best:  1  Score:  9999.1958774679169080712113587517544242
          9999.1958774679169080712113587517544242
1         9999.1958874853672691525439166438792715
Best:  1  Score:  9999.1958874853672691525439166438792715
[47597, 177, 30, 12, 7, 5, 4, 3, 3, 2, 2, 2, 2, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]
          9999.1958874853672691525439166438792715
1         9999.1958974790750359025587175474328713
Best:  1  Score:  9999.1958974790750359025587175474328713
          9999.1958974790750359025587175474328713
1         9999.1959077684565907996537266820226472

Advertisements

Written by meditationatae

July 10, 2013 at 5:33 am

Posted in History

3 Responses

Subscribe to comments with RSS.

  1. Function definitions (harmonic and ratio):

    ? \u
    harmonic =
    (Z)->Euler+psi(Z+1)

    ratio =
    (ZZZ)->base*(prod(Y=1,50000,sum(W=0,ZZZ[Y],prime(Y)^W))/(harmonic(prod(Y=1,50000,prime(Y)^ZZZ[Y]))+log(harmonic(prod(Y=1,50000,prime(Y)^ZZZ[Y])))*exp(harmonic(prod(Y=1,50000,prime(Y)^ZZZ[Y])))))

    meditationatae

    July 10, 2013 at 5:40 am

  2. I’m now verifying the data (computation) surrounding:
    [47597, 177, 30, 12, 7, 5, 4, 3, 3, 2, 2, 2, 2, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, ETC.]
    9999.1958874853672691525439166438792715

    As explained in the post, to avoid having to factor numbers each time the sum of divisors function sigma() is called, I now work with representations of integers by a vector of the exponents of the primes 2, 3, 5, 7, …. p_50000 that appear in the prime factorization of an integer. In the previous line, p_50000 is the 50,000th prime. The vectors of exponents used are therefore of length 50,0000. This means that I compute sigma(.) without calling/invoking PARI/gp ‘s built-in sigma function.

    N.B.: 9999.1958 etc. is the result of multiplying a quotient that was 0.99991958 etc. by 10,000.

    For that reason, I wish to check my home-grown sigma-related function using PARI-gp’s own sigma(.) built-in.

    There’s a product of 22 not always distinct primorials that is represented by [47597, 177, 30, 12, 7, 5, 4, 3, 3, 2, 2, 2, 2, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, ETC.] :

    a40 :=
    primorial(47597)*
    primorial(177)*
    primorial(30)*
    primorial(12)*
    primorial(7)*
    primorial(5)*
    primorial(4)*
    30*30*
    (6^4)*
    (2^9);

    I wrote ’30’ here instead of primorial(2) and ‘2’ instead of primorial(1).

    First, the function definitions (user-defined, by me):

    ? \u

    Qr =
    (W)->sigma(W)/(harmonic(W)+log(harmonic(W))*exp(harmonic(W)))

    harmonic =
    (Z)->Euler+psi(Z+1)

    primorial =
    (W)->prod(X=1,W,prime(X))
    ———–

    So,
    Qr(n):= sigma(n)/[ H_n + log(H_n)*exp(H_n) ],

    where H_n is the n’th harmonic number, or
    H_n := sum_{k = 1 … n} 1/k
    (partial sums of the harmonic series).

    For complenetess, psi(Z+1) in my definition of harmonic(.) is the digamma function, at Z+1, and can be used together with the Euler-Mascheroni constant, ‘Euler’ in PARI/gp, to define
    H_n using digamma at n+1 and the Euler-Mascheroni constant.

    We see that the user-define Qr uses PARI/gp’s built-in sigma sum of divisors.

    I’m waiting on the evaluation of Qr(a40) to finish,

    ? Qr(a40)

    ==========

    The hoped-for result is about

    0.9999195887485367269 etc

    [not done yet … ]

    meditationatae

    July 10, 2013 at 7:25 am

    • Now done:

      ? Qr(a40)
      %67 = 0.99991958874853672691525439166438792716
      ?

      vs.

      9999.1958874853672691525439166438792715

      Ok. Pass.

      Note:

      ratio(ZZZ):=

      base*(prod(Y=1,50000,sum(W=0,ZZZ[Y],prime(Y)^W))
      /(harmonic(prod(Y=1,50000,prime(Y)^ZZZ[Y]))
      +log(harmonic(prod(Y=1,50000,
      prime(Y)^ZZZ[Y])))*exp(harmonic(prod(Y=1,50000,
      prime(Y)^ZZZ[Y])))))

      with base = 10000.

      meditationatae

      July 10, 2013 at 7:51 am


Comments are closed.

%d bloggers like this: