meditationatae

Just another WordPress.com site

Source code for superabundant numbers computation program in C

/**********************************************************

Source code for:  superabun1848a.c

$ date
Thu Dec 4 13:14:45 EST 2014

********************************************************/

#include <stdio.h>
#include <math.h>
#include <stdlib.h>
long double splogs[100000];
long double theta[100000];
long double prod100k[100000];
long ptr_exponents[100000];
long pts_primes[100000];
long ptr_isstep[100000];
int main(void)
{
long theIndex;
long double euler = 0.57721566490153286060651L;
long j;
long numprimes;
long num_new_primes;
long file_offset;
long jj;
long kij;
long numBigPrimesUsed;
long double sum_log_cofactor;
long double sum_theta_Cheb;
long double puiss3sigmalogs[50];
long numsteps;
long expon;
long sigma;
long * ptr_primes;
long list_steps[1000];
long power;
long old_prime;
long powerb;
int nownow;
int go_skip;
int answer;
long k;
long l;
long numloops;
long vintageloops;
long tmp_exponent;
long m;
long p;
long sum;
int not_done;
int junk;
int leave_loop;
long ljunk;
int skip_base_case;
long current_step;
long now_step;
long double RecordRobinQuotient;
long double RecordQuotient;
long double tmp_now_n;
long double log2nprime;
long double robin_bound;
long double base_tmp_now_n;
long RecordStep;
long double tmp_now_robin;
long double constant;
long double prod_sigma_r;
long double prod_sigma_rone;
long double base_tmp_now_sigma_r;
long double prod_sigma_r_cofactor;
long double tmp_now_sigma_r;
long double tmp_now_quotient;
long double delta;
long double delta_log;
long double delta_estimate;
long double variation_delta;
long double record_low_variation_delta;
long double record_low_log2n;
long double tmp_el_now_robin;
long double now_n;
long current_step_number;
long double sum_logs_primes;
long double sum_logs_all;
long double ratio;
long double now_robin;
FILE *in;
record_low_variation_delta = (long double) 100;
record_low_log2n = (long double) 0;
ptr_primes = (long *) calloc(80000000, sizeof(long));
in = fopen(“/home/david2/eratosthenes/primes100meg01a_tmp.txt”, “r”);
for(j=0; j< 100000 ; j++)
{
fscanf(in, “%ld”, &theIndex);
fscanf(in, “%ld”, &pts_primes[j]);
}
fclose(in);

constant = ((long double) 323336)/((long double) 1000000);

in = fopen(“/home/david2/eratosthenes/prout100k”, “r”);
for(j=0; j<100000; j++)
{
fscanf(in, “%Lf”, &prod100k[j]);
}
fclose(in);
in = fopen(“/home/david2/eratosthenes/LargePrimes_to_1e10a.txt”, “r”);
for(j=0; j<79920000 ; j++)
{
fscanf(in, “%ld”, &ljunk);
}
for(j=0; j< 80000000 ; j++)
{
fscanf(in, “%ld”, &ptr_primes[j]);
}

fclose(in);
in = fopen(“/home/david2/eratosthenes/puiss3sigmalogs”, “r”);

for(j=0; j< 50 ; j++)
{
fscanf(in, “%d”, &junk);
fscanf(in, “%Lf”, &puiss3sigmalogs[j]);
}
fclose(in);
for(j=0; j< 100000 ; j++)
{
ptr_exponents[j] = (long) 0;
}
prod_sigma_r = (long double) 1;

prod_sigma_rone = (long double) 1;

for(j=0; j< 100000 ; j++)
{
ptr_isstep[j] = (long) 0;
}
in = fopen(“/home/david2/eratosthenes/exponents_out802a.txt”, “r”);
for(j=0; j< 100000 ; j++)
{
fscanf(in, “%ld”, &theIndex);
fscanf(in, “%ld”, &ptr_exponents[j]);
}
fclose(in);
numprimes = (long) 80000000;

numBigPrimesUsed = (long) 79000000;

leave_loop = 0;
in = fopen(“/home/david2/eratosthenes/primes_out802a.txt”, “r”);

sum_logs_primes = (long double) 0;
sum_logs_all = (long double) 0;
for(j=0; j<80000000; j++)
{
fscanf(in, “%ld”, &ljunk);
fscanf(in, “%ld”, &old_prime);
sum_logs_primes = sum_logs_primes + logl((long double)old_prime);

if(j < 100000)
{
sum_logs_all = sum_logs_all + ((long double)ptr_exponents[j])*logl((long double)old_prime);
}
else
{
sum_logs_all = sum_logs_all + logl((long double)old_prime);
}
}

fclose(in);

for(j=0; j<100000; j++)
{
splogs[j] = logl((long double) pts_primes[j]);
}

sum_theta_Cheb = (long double) 0;
for(j=0; j<100000; j++)
{
sum_theta_Cheb = sum_theta_Cheb + logl((long double) pts_primes[j]) ;
theta[j] = sum_theta_Cheb;
}
skip_base_case = (int) 0;
for(j=0; j< 100000 ; j++)
{

p = pts_primes[j];

expon = ptr_exponents[j];

ratio = (long double) 1;

for(jj = 0; jj < 1; jj++)
{
ratio = ((long double)1) + (ratio/((long double) p) );
}

if(expon == 1)
{
prod_sigma_rone = prod_sigma_rone*ratio;
}
prod_sigma_r = prod_sigma_r*ratio;
}

in = fopen(“/home/david2/eratosthenes/primes_out802a.txt”, “r”);

for(j=0 ; j<numprimes; j++)
{
fscanf(in, “%ld”, &ljunk);
fscanf(in, “%ld”, &old_prime);

if( j > 99999 )
{
prod_sigma_r = prod_sigma_r + (prod_sigma_r/((long double)old_prime)) ;
prod_sigma_rone = prod_sigma_rone + (prod_sigma_rone/((long double)old_prime)) ;
}
}

fclose(in);

// printf(“sum_logs_all = %.16Lf\n\n”, sum_logs_all);
// printf(“log of sum_logs_all or log log n is: %.16Lf\n\n”, logl(sum_logs_all) );
for(numloops = 0; numloops < 80002284 ; numloops++)
{

ptr_isstep[0] = (long) 1;
list_steps[0] = (long) 0;
numsteps = (long) 1;
not_done = (int) 1;
j = (long) 0;

while(not_done == 1)
{
j++;

if(ptr_exponents[j]< ptr_exponents[j-1])
{
ptr_isstep[j] = (long) 1;
list_steps[numsteps] = j;
numsteps++;
}
else
{
ptr_isstep[j] = (long) 0;
}

if(ptr_exponents[j] == ((long) 1))
{
not_done = (int) 0;
}

}

list_steps[numsteps] = numprimes;
numsteps++;
if( (2283 ==(numloops%64000)) || (((numloops – 80002283) < 6)&&((80002283-numloops) < 6)) )
{
// for(j=0; j< numsteps; j++)
// {
// printf(“step at %ld\n”, list_steps[j]);
// }
// printf(“\nNumber of big primes used in making current n is: %ld\n”, numBigPrimesUsed);
// printf(“numloops = %ld\n”, numloops);
// printf(“\n”);

if(numloops == 80002283)
{
// printf(“numloops = %ld\n”, numloops);
// printf(“break from current loop, ok?\n”);
// printf(“enter 1 for yes, 0 for no:\n”);

if(1 == 1)
{
leave_loop = 1;
}
else
{
leave_loop = 1;
}
}
}

if(leave_loop == 0)
{
if((2283 ==(numloops%64000)) || (((numloops – 80002283) < 6)&&((80002283-numloops) < 6)) )
{
// printf(“number of steps: %ld\n”, numsteps);
}
RecordRobinQuotient = (long double) 0;
if(skip_base_case == 0)
{

base_tmp_now_n = (long double) 0;

base_tmp_now_sigma_r = (long double) 1;

now_step = (long) 0;

for(j=0; j<= list_steps[numsteps-2]; j++)
{
if(j == list_steps[now_step] )
{
now_step++;
}
else
{
base_tmp_now_n = base_tmp_now_n + ((long double)ptr_exponents[j])*splogs[j];
p = pts_primes[j];
power = (long) 1;
for(l=0; l<ptr_exponents[j]+((long)1); l++)
{
power = power*p;
}
powerb = power/p;
sigma = (power – ((long)1) )/(p – ((long)1));
ratio = ((long double)sigma)/((long double)powerb);
base_tmp_now_sigma_r = base_tmp_now_sigma_r*ratio;
}
}
}

sum_log_cofactor = sum_logs_primes;

prod_sigma_r_cofactor = prod_sigma_r;
sum_log_cofactor = sum_log_cofactor – theta[list_steps[numsteps-2]];

prod_sigma_r_cofactor = prod_sigma_r_cofactor/prod100k[list_steps[numsteps-2]];
for(m=0; m< numsteps ; m++)
{
current_step_number = m;

current_step = list_steps[m];

tmp_now_n = sum_log_cofactor + base_tmp_now_n;
tmp_now_sigma_r = base_tmp_now_sigma_r*prod_sigma_r_cofactor;
if(m < numsteps-((long)1) )
{
for(j=0; j<= numsteps- ((long)2) ; j++)
{
if(j == m)
{
tmp_exponent = ptr_exponents[list_steps[j]] + ((long)1);
}
else
{
tmp_exponent = ptr_exponents[list_steps[j]];
}
p = pts_primes[list_steps[j]];
tmp_now_n = tmp_now_n + ((long double)tmp_exponent)*splogs[list_steps[j]];
power = (long) 1;
for(k=0; k<(tmp_exponent+((long)1) ); k++)
{
power = power*p;
}
powerb = power/p;
sigma = (power – ((long)1) )/(p – ((long)1));

ratio = ((long double)sigma)/((long double)powerb);
tmp_now_sigma_r = tmp_now_sigma_r*ratio;
}

tmp_now_quotient = tmp_now_sigma_r;

tmp_now_robin = tmp_now_quotient/expl(euler);

tmp_el_now_robin = tmp_now_robin/logl(tmp_now_n);

tmp_now_robin = tmp_el_now_robin;

if(tmp_now_robin > RecordRobinQuotient)
{
RecordRobinQuotient = tmp_now_robin;
RecordStep = current_step;
}
if((2283 ==(numloops%64000)) || (((numloops – 80002283) < 6)&&((80002283-numloops) < 6)) )
{

// printf(“current step number = %ld “, current_step_number);

// printf(“quotient de Robin: %.16Lf “, tmp_now_robin);

// printf(“record step is at: %ld\n”, RecordStep);

}
}
else
{
for(j=0; j<= numsteps- ((long)2) ; j++)
{
tmp_exponent = ptr_exponents[list_steps[j]];

p = pts_primes[list_steps[j]];
if(list_steps[j] != 1)
{
tmp_now_n = tmp_now_n + ((long double)tmp_exponent)*splogs[list_steps[j]];
power = (long) 1;
for(k=0; k<(tmp_exponent+1); k++)
{
power = power*p;
}

powerb = power/p;
sigma = (power – ((long)1) )/(p – ((long)1));

ratio = ((long double)sigma)/((long double)powerb);

tmp_now_sigma_r = tmp_now_sigma_r*ratio;
}
else
{
tmp_now_n = tmp_now_n + ((long double)tmp_exponent)*splogs[1];
ratio = expl( puiss3sigmalogs[tmp_exponent]-(((long double)tmp_exponent)*splogs[1]) );

tmp_now_sigma_r = tmp_now_sigma_r*ratio;
}
}

tmp_now_n = tmp_now_n + logl((long double) ptr_primes[current_step-80000000]);
ratio = ((long double)1)/( (long double)ptr_primes[current_step-80000000] );

tmp_now_sigma_r = tmp_now_sigma_r + (ratio*tmp_now_sigma_r) ;

tmp_now_quotient = tmp_now_sigma_r;
tmp_now_robin = tmp_now_quotient/expl(euler);

tmp_el_now_robin = tmp_now_robin/logl(tmp_now_n);

tmp_now_robin = tmp_el_now_robin;
if(tmp_now_robin > RecordRobinQuotient)
{
RecordRobinQuotient = tmp_now_robin;
RecordStep = current_step;
}

if((2283 ==(numloops%64000)) || (((numloops – 80002283) < 6)&&((80002283-numloops) < 6)) )
{
// printf(“current step number = %ld “, current_step_number);

// printf(“quotient de Robin: %.16Lf “, tmp_now_robin);

// printf(“record step is at: %ld\n”, RecordStep);
}
}

}

if(RecordStep < numprimes)
{
ptr_exponents[RecordStep]++;
}
if(RecordStep == numprimes)
{
sum_logs_primes = sum_logs_primes + logl( (long double) ptr_primes[numprimes-80000000]);

ratio = ((long double)1)/((long double)ptr_primes[numprimes-80000000]);

if( RecordStep > 999999 )
{
numBigPrimesUsed++;
}

prod_sigma_r = prod_sigma_r + (prod_sigma_r*ratio) ;

numprimes++;
skip_base_case = (int) 1;
}
else
{
skip_base_case = (int) 0;
}

if( 2283 ==(numloops%64000) )
{
for(j=0; j< 200 ; j++)
{
// printf(“%ld “, ptr_exponents[j]);
}

// printf(“\n\n”);
}
}
/*** if leave loop == 0 ****/
if( leave_loop == 1 )
{
break;
}
}

/**************** next 80 million batch … *****************/
/*************************************************************
prime 160,000,000 is:

3340200037

or

3,340,200,037 .

:159920000p
3340200037
3340200037
************************************************************/

in = fopen(“/home/david2/eratosthenes/LargePrimes_to_1e11ofNov22a.txt”, “r”);
for(j=0; j< 159920000 ; j++)
{
fscanf(in, “%ld”, &old_prime);
}
for(j=0; j< 80000000 ; j++)
{
fscanf(in, “%ld”, &ptr_primes[j]);
}
/**************************************

==============================================
next: loops ….

******************************************************/
file_offset = (long) 160000000;
num_new_primes = (long) 0;

numloops = (long) 0;
while( num_new_primes < ((long)80000000) )
{
numloops++;

ptr_isstep[0] = (long) 1;
list_steps[0] = (long) 0;
numsteps = (long) 1;
not_done = (int) 1;
j = (long) 0;

while(not_done == 1)
{
j++;

if(ptr_exponents[j]< ptr_exponents[j-1])
{
ptr_isstep[j] = (long) 1;
list_steps[numsteps] = j;
numsteps++;
}
else
{
ptr_isstep[j] = (long) 0;
}

if(ptr_exponents[j] == ((long) 1))
{
not_done = (int) 0;
}

}

list_steps[numsteps] = numprimes;
numsteps++;

if( 9999 == (num_new_primes%10000) )
{
for(j=0; j< numsteps; j++)
{
// printf(“step at %ld\n”, list_steps[j]);
}
// printf(“\nNumber of big primes used in making current n is: %ld\n”, numBigPrimesUsed);
// printf(“numloops = %ld\n”, numloops);
// printf(“\n”);
}
if( 9999 == (num_new_primes%10000) )
{
// printf(“number of steps: %ld\n”, numsteps);
}
RecordRobinQuotient = (long double) 0;
if(skip_base_case == 0)
{

base_tmp_now_n = (long double) 0;

base_tmp_now_sigma_r = (long double) 1;

now_step = (long) 0;

for(j=0; j<= list_steps[numsteps-2]; j++)
{
if(j == list_steps[now_step] )
{
now_step++;
}
else
{
base_tmp_now_n = base_tmp_now_n + ((long double)ptr_exponents[j])*splogs[j];
p = pts_primes[j];
power = (long) 1;
for(l=0; l<ptr_exponents[j]+((long)1); l++)
{
power = power*p;
}
powerb = power/p;
sigma = (power – ((long)1) )/(p – ((long)1));
ratio = ((long double)sigma)/((long double)powerb);
base_tmp_now_sigma_r = base_tmp_now_sigma_r*ratio;
}
}
}

sum_log_cofactor = sum_logs_primes;

prod_sigma_r_cofactor = prod_sigma_r;
sum_log_cofactor = sum_log_cofactor – theta[list_steps[numsteps-2]];

prod_sigma_r_cofactor = prod_sigma_r_cofactor/prod100k[list_steps[numsteps-2]];

for(m=0; m< numsteps ; m++)
{
current_step_number = m;

current_step = list_steps[m];

tmp_now_n = sum_log_cofactor + base_tmp_now_n;
tmp_now_sigma_r = base_tmp_now_sigma_r*prod_sigma_r_cofactor;
if(m < numsteps-((long)1) )
{
for(j=0; j<= numsteps- ((long)2) ; j++)
{
if(j == m)
{
tmp_exponent = ptr_exponents[list_steps[j]] + ((long)1);
}
else
{
tmp_exponent = ptr_exponents[list_steps[j]];
}
p = pts_primes[list_steps[j]];
tmp_now_n = tmp_now_n + ((long double)tmp_exponent)*splogs[list_steps[j]];
power = (long) 1;
for(k=0; k<(tmp_exponent+((long)1) ); k++)
{
power = power*p;
}
powerb = power/p;
sigma = (power – ((long)1) )/(p – ((long)1));

ratio = ((long double)sigma)/((long double)powerb);
tmp_now_sigma_r = tmp_now_sigma_r*ratio;
}

tmp_now_quotient = tmp_now_sigma_r;

tmp_now_robin = tmp_now_quotient/expl(euler);

tmp_el_now_robin = tmp_now_robin/logl(tmp_now_n);

tmp_now_robin = tmp_el_now_robin;

if(tmp_now_robin > RecordRobinQuotient)
{
RecordRobinQuotient = tmp_now_robin;
RecordStep = current_step;
}
if( 9999 == (num_new_primes%10000) )
{

// printf(“current step number = %ld “, current_step_number);

// printf(“quotient de Robin: %.16Lf “, tmp_now_robin);

// printf(“record step is at: %ld\n”, RecordStep);

}
}
else
{
for(j=0; j<= numsteps- ((long)2) ; j++)
{
tmp_exponent = ptr_exponents[list_steps[j]];

p = pts_primes[list_steps[j]];
if(list_steps[j] != 1)
{
tmp_now_n = tmp_now_n + ((long double)tmp_exponent)*splogs[list_steps[j]];
power = (long) 1;
for(k=0; k<(tmp_exponent+1); k++)
{
power = power*p;
}

powerb = power/p;
sigma = (power – ((long)1) )/(p – ((long)1));

ratio = ((long double)sigma)/((long double)powerb);

tmp_now_sigma_r = tmp_now_sigma_r*ratio;
}
else
{
tmp_now_n = tmp_now_n + ((long double)tmp_exponent)*splogs[1];
ratio = expl( puiss3sigmalogs[tmp_exponent]-(((long double)tmp_exponent)*splogs[1]) );

tmp_now_sigma_r = tmp_now_sigma_r*ratio;
}
}

tmp_now_n = tmp_now_n + logl((long double) ptr_primes[current_step- file_offset]);
ratio = ((long double)1)/( (long double)ptr_primes[current_step- file_offset] );

tmp_now_sigma_r = tmp_now_sigma_r + (ratio*tmp_now_sigma_r) ;

tmp_now_quotient = tmp_now_sigma_r;
tmp_now_robin = tmp_now_quotient/expl(euler);

tmp_el_now_robin = tmp_now_robin/logl(tmp_now_n);

tmp_now_robin = tmp_el_now_robin;
if(tmp_now_robin > RecordRobinQuotient)
{
RecordRobinQuotient = tmp_now_robin;
RecordStep = current_step;
}

if( 9999 == (num_new_primes%10000) )
{
// printf(“current step number = %ld “, current_step_number);

// printf(“quotient de Robin: %.16Lf “, tmp_now_robin);

// printf(“record step is at: %ld\n”, RecordStep);
}
}

}

if(RecordStep < numprimes)
{
ptr_exponents[RecordStep]++;
}
if(RecordStep == numprimes)
{
sum_logs_primes = sum_logs_primes + logl( (long double) ptr_primes[numprimes- file_offset]);

ratio = ((long double)1)/((long double)ptr_primes[numprimes-file_offset]);

if( RecordStep > 999999 )
{
numBigPrimesUsed++;
}

num_new_primes++;

prod_sigma_r = prod_sigma_r + (prod_sigma_r*ratio) ;

numprimes++;
skip_base_case = (int) 1;
}
else
{
skip_base_case = (int) 0;
}

if( 9999 == (num_new_primes%10000) )
{
for(j=0; j< 200 ; j++)
{
// printf(“%ld “, ptr_exponents[j]);
}

// printf(“\n\n”);
}
}
/********************** Next batch 80 M and loop *******************/

for(kij = 1; kij < 49; kij++)
{

for(j=0; j< 80000000 ; j++)
{
fscanf(in, “%ld”, &ptr_primes[j]);
}

file_offset = file_offset + ((long)80000000) ;

/*** 240 000 000 for kij=1 ***/
num_new_primes = (long) 0;

numloops = (long) 0;

while( num_new_primes < ((long)80000000) )
{
numloops++;

ptr_isstep[0] = (long) 1;
list_steps[0] = (long) 0;
numsteps = (long) 1;
not_done = (int) 1;
j = (long) 0;

while(not_done == 1)
{
j++;

if(ptr_exponents[j]< ptr_exponents[j-1])
{
ptr_isstep[j] = (long) 1;
list_steps[numsteps] = j;
numsteps++;
}
else
{
ptr_isstep[j] = (long) 0;
}

if(ptr_exponents[j] == ((long) 1))
{
not_done = (int) 0;
}

}

list_steps[numsteps] = numprimes;
numsteps++;

if( 9999 == (numloops%10000) )
{
for(j=0; j< numsteps; j++)
{
// printf(“step at %ld\n”, list_steps[j]);
}
// printf(“\nNumber of big primes used in making current n is: %ld\n”, numBigPrimesUsed);
// printf(“numloops = %ld\n”, numloops);
// printf(“\n”);
}
if( 9999 == (numloops%10000) )
{
// printf(“number of steps: %ld\n”, numsteps);
}
if( 9999 == (numloops%10000) )
{
for(j=0; j< 200 ; j++)
{
// printf(“%ld “, ptr_exponents[j]);
}

// printf(“\n”);
}

RecordRobinQuotient = (long double) 0;
if(skip_base_case == 0)
{

base_tmp_now_n = (long double) 0;

base_tmp_now_sigma_r = (long double) 1;

now_step = (long) 0;

for(j=0; j<= list_steps[numsteps-2]; j++)
{
if(j == list_steps[now_step] )
{
now_step++;
}
else
{
base_tmp_now_n = base_tmp_now_n + ((long double)ptr_exponents[j])*splogs[j];
p = pts_primes[j];
power = (long) 1;
for(l=0; l<ptr_exponents[j]+((long)1); l++)
{
power = power*p;
}
powerb = power/p;
sigma = (power – ((long)1) )/(p – ((long)1));
ratio = ((long double)sigma)/((long double)powerb);
base_tmp_now_sigma_r = base_tmp_now_sigma_r*ratio;
}
}
}

sum_log_cofactor = sum_logs_primes;

prod_sigma_r_cofactor = prod_sigma_r;
sum_log_cofactor = sum_log_cofactor – theta[list_steps[numsteps-2]];

prod_sigma_r_cofactor = prod_sigma_r_cofactor/prod100k[list_steps[numsteps-2]];

for(m=0; m< numsteps ; m++)
{
current_step_number = m;

current_step = list_steps[m];

tmp_now_n = sum_log_cofactor + base_tmp_now_n;
tmp_now_sigma_r = base_tmp_now_sigma_r*prod_sigma_r_cofactor;
if(m < numsteps-((long)1) )
{
for(j=0; j<= numsteps- ((long)2) ; j++)
{
if(j == m)
{
tmp_exponent = ptr_exponents[list_steps[j]] + ((long)1);
}
else
{
tmp_exponent = ptr_exponents[list_steps[j]];
}
p = pts_primes[list_steps[j]];
tmp_now_n = tmp_now_n + ((long double)tmp_exponent)*splogs[list_steps[j]];
power = (long) 1;
for(k=0; k<(tmp_exponent+((long)1) ); k++)
{
power = power*p;
}
powerb = power/p;
sigma = (power – ((long)1) )/(p – ((long)1));

ratio = ((long double)sigma)/((long double)powerb);
tmp_now_sigma_r = tmp_now_sigma_r*ratio;
}

tmp_now_quotient = tmp_now_sigma_r;

tmp_now_robin = tmp_now_quotient/expl(euler);

log2nprime = logl(tmp_now_n);

robin_bound = expl(euler)*log2nprime;

tmp_el_now_robin = tmp_now_robin/logl(tmp_now_n);

tmp_now_robin = tmp_el_now_robin;

if(tmp_now_robin > RecordRobinQuotient)
{
RecordRobinQuotient = tmp_now_robin;
RecordStep = current_step;
}
if( 9999 == (numloops%10000) )
{

// printf(“current step number = %ld “, current_step_number);

// printf(“quotient de Robin: %.16Lf “, tmp_now_robin);

// printf(“record step is at: %ld\n”, RecordStep);

// printf(“log(log(n’)) = %.16Lf “, log2nprime);

// printf(“Robin upper bound = %.16Lf “, robin_bound);

// printf(“sigma(n’)/n’ = %.16Lf \n”, tmp_now_quotient);
}
}
else
{
for(j=0; j<= numsteps- ((long)2) ; j++)
{
tmp_exponent = ptr_exponents[list_steps[j]];

p = pts_primes[list_steps[j]];
if(list_steps[j] != 1)
{
tmp_now_n = tmp_now_n + ((long double)tmp_exponent)*splogs[list_steps[j]];
power = (long) 1;
for(k=0; k<(tmp_exponent+1); k++)
{
power = power*p;
}

powerb = power/p;
sigma = (power – ((long)1) )/(p – ((long)1));

ratio = ((long double)sigma)/((long double)powerb);

tmp_now_sigma_r = tmp_now_sigma_r*ratio;
}
else
{
tmp_now_n = tmp_now_n + ((long double)tmp_exponent)*splogs[1];
ratio = expl( puiss3sigmalogs[tmp_exponent]-(((long double)tmp_exponent)*splogs[1]) );

tmp_now_sigma_r = tmp_now_sigma_r*ratio;
}
}

tmp_now_n = tmp_now_n + logl((long double) ptr_primes[current_step- file_offset]);
ratio = ((long double)1)/( (long double)ptr_primes[current_step- file_offset] );

tmp_now_sigma_r = tmp_now_sigma_r + (ratio*tmp_now_sigma_r) ;

tmp_now_quotient = tmp_now_sigma_r;
tmp_now_robin = tmp_now_quotient/expl(euler);
log2nprime = logl(tmp_now_n);
robin_bound = expl(euler)*log2nprime;
tmp_el_now_robin = tmp_now_robin/logl(tmp_now_n);

tmp_now_robin = tmp_el_now_robin;
if(tmp_now_robin > RecordRobinQuotient)
{
RecordRobinQuotient = tmp_now_robin;
RecordStep = current_step;
}

if( 9999 == (numloops%10000) )
{
// printf(“current step number = %ld “, current_step_number);

// printf(“quotient de Robin: %.16Lf “, tmp_now_robin);

// printf(“record step is at: %ld\n”, RecordStep);

// printf(“log(log(n’)) = %.16Lf “, log2nprime);

// printf(“Robin upper bound = %.16Lf “, robin_bound);

// printf(“sigma(n’)/n’ = %.16Lf \n”, tmp_now_quotient);

}
}

}

if(RecordStep < numprimes)
{
ptr_exponents[RecordStep]++;
}
if(RecordStep == numprimes)
{

if(num_new_primes < 79999500)
{

for(j=0; j< 100 ; j++)
{

sum_logs_primes = sum_logs_primes + logl( (long double) ptr_primes[numprimes- file_offset]);

ratio = ((long double)1)/((long double)ptr_primes[numprimes-file_offset]);

if( RecordStep > 999999 )
{
numBigPrimesUsed++;
}

num_new_primes++;

prod_sigma_r = prod_sigma_r + (prod_sigma_r*ratio) ;

numprimes++;
}
}
else
{

sum_logs_primes = sum_logs_primes + logl( (long double) ptr_primes[numprimes- file_offset]);

ratio = ((long double)1)/((long double)ptr_primes[numprimes-file_offset]);

if( RecordStep > 999999 )
{
numBigPrimesUsed++;
}

num_new_primes++;

prod_sigma_r = prod_sigma_r + (prod_sigma_r*ratio) ;

numprimes++;
}
skip_base_case = (int) 1;
}
else
{
skip_base_case = (int) 0;
}

if( 9999 == (numloops%10000) )
{
// printf(“=================================\n”);
}
}
}

fclose(in);
/********************** Next batch 80 M and loop *******************/

/*****************************************************************
? primepi(99036319010)

%37 = 4080000000

4,080,000,000 primes included in `n’
in = fopen(“/home/david2/eratosthenes/LargePrimes_to_1e11ofNov22a.txt”, “r”);
**********************************************************************/
in = fopen(“/home/david2/eratosthenes2/LargePrimes_to_2e11ofNov26a.txt”, “r”);
for(kij = 1; kij < 50; kij++)
{

for(j=0; j< 80000000 ; j++)
{
fscanf(in, “%ld”, &ptr_primes[j]);
}

file_offset = file_offset + ((long)80000000) ;
num_new_primes = (long) 0;

numloops = (long) 0;

while( num_new_primes < ((long)80000000) )
{
numloops++;

ptr_isstep[0] = (long) 1;
list_steps[0] = (long) 0;
numsteps = (long) 1;
not_done = (int) 1;
j = (long) 0;

while(not_done == 1)
{
j++;

if(ptr_exponents[j]< ptr_exponents[j-1])
{
ptr_isstep[j] = (long) 1;
list_steps[numsteps] = j;
numsteps++;
}
else
{
ptr_isstep[j] = (long) 0;
}

if(ptr_exponents[j] == ((long) 1))
{
not_done = (int) 0;
}

}

list_steps[numsteps] = numprimes;
numsteps++;

if( 9999 == (numloops%10000) )
{
for(j=0; j< numsteps; j++)
{
// printf(“step at %ld\n”, list_steps[j]+1);
}
// printf(“\nNumber of big primes used in making current n is: %ld\n”, numBigPrimesUsed);
// printf(“numloops = %ld\n”, numloops);
// printf(“\n”);
}
if( 9999 == (numloops%10000) )
{
// printf(“number of steps: %ld\n”, numsteps);
}
if( 9999 == (numloops%10000) )
{
for(j=0; j< 200 ; j++)
{
// printf(“%ld “, ptr_exponents[j]);
}

// printf(“\n”);
}

if( 799999 == (numloops%800000) )
{
printf(“record_low_variation_delta = %.16Lf “, record_low_variation_delta);
printf(“record_low_log2n = %.16Lf\n\n”, record_low_log2n);
fflush(stdout);

record_low_variation_delta = (long double) 100;
record_low_log2n = (long double) 0;
}

RecordRobinQuotient = (long double) 0;
if(skip_base_case == 0)
{

base_tmp_now_n = (long double) 0;

base_tmp_now_sigma_r = (long double) 1;

now_step = (long) 0;

for(j=0; j<= list_steps[numsteps-2]; j++)
{
if(j == list_steps[now_step] )
{
now_step++;
}
else
{
base_tmp_now_n = base_tmp_now_n + ((long double)ptr_exponents[j])*splogs[j];
p = pts_primes[j];
power = (long) 1;
for(l=0; l<ptr_exponents[j]+((long)1); l++)
{
power = power*p;
}
powerb = power/p;
sigma = (power – ((long)1) )/(p – ((long)1));
ratio = ((long double)sigma)/((long double)powerb);
base_tmp_now_sigma_r = base_tmp_now_sigma_r*ratio;
}
}
}

sum_log_cofactor = sum_logs_primes;

prod_sigma_r_cofactor = prod_sigma_r;
sum_log_cofactor = sum_log_cofactor – theta[list_steps[numsteps-2]];

prod_sigma_r_cofactor = prod_sigma_r_cofactor/prod100k[list_steps[numsteps-2]];

for(m=0; m< numsteps ; m++)
{
current_step_number = m;

current_step = list_steps[m];

tmp_now_n = sum_log_cofactor + base_tmp_now_n;
tmp_now_sigma_r = base_tmp_now_sigma_r*prod_sigma_r_cofactor;
if(m < numsteps-((long)1) )
{
for(j=0; j<= numsteps- ((long)2) ; j++)
{
if(j == m)
{
tmp_exponent = ptr_exponents[list_steps[j]] + ((long)1);
}
else
{
tmp_exponent = ptr_exponents[list_steps[j]];
}
p = pts_primes[list_steps[j]];
tmp_now_n = tmp_now_n + ((long double)tmp_exponent)*splogs[list_steps[j]];
power = (long) 1;
for(k=0; k<(tmp_exponent+((long)1) ); k++)
{
power = power*p;
}
powerb = power/p;
sigma = (power – ((long)1) )/(p – ((long)1));

ratio = ((long double)sigma)/((long double)powerb);
tmp_now_sigma_r = tmp_now_sigma_r*ratio;
}

tmp_now_quotient = tmp_now_sigma_r;

tmp_now_robin = tmp_now_quotient/expl(euler);

log2nprime = logl(tmp_now_n);

robin_bound = expl(euler)*log2nprime;

delta = robin_bound – tmp_now_quotient;

delta_log = logl(delta);

delta_estimate = constant – log2nprime/((long double) 2);

variation_delta = delta_log – delta_estimate;

if( variation_delta < record_low_variation_delta)
{
record_low_variation_delta = variation_delta;
record_low_log2n = log2nprime;
}

tmp_el_now_robin = tmp_now_robin/logl(tmp_now_n);

tmp_now_robin = tmp_el_now_robin;

if(tmp_now_robin > RecordRobinQuotient)
{
RecordRobinQuotient = tmp_now_robin;
RecordStep = current_step;
}
if( 9999 == (numloops%10000) )
{

// printf(“current step number = %ld “, current_step_number+1);

// printf(“quotient de Robin: %.16Lf “, tmp_now_robin);

// printf(“record step is at: %ld\n”, RecordStep+1);

// printf(“log(log(n’)) = %.16Lf “, log2nprime);

// printf(“Robin upper bound = %.16Lf “, robin_bound);

// printf(“sigma(n’)/n’ = %.16Lf \n”, tmp_now_quotient);
}
}
else
{
for(j=0; j<= numsteps- ((long)2) ; j++)
{
tmp_exponent = ptr_exponents[list_steps[j]];

p = pts_primes[list_steps[j]];
if(list_steps[j] != 1)
{
tmp_now_n = tmp_now_n + ((long double)tmp_exponent)*splogs[list_steps[j]];
power = (long) 1;
for(k=0; k<(tmp_exponent+1); k++)
{
power = power*p;
}

powerb = power/p;
sigma = (power – ((long)1) )/(p – ((long)1));

ratio = ((long double)sigma)/((long double)powerb);

tmp_now_sigma_r = tmp_now_sigma_r*ratio;
}
else
{
tmp_now_n = tmp_now_n + ((long double)tmp_exponent)*splogs[1];
ratio = expl( puiss3sigmalogs[tmp_exponent]-(((long double)tmp_exponent)*splogs[1]) );

tmp_now_sigma_r = tmp_now_sigma_r*ratio;
}
}

tmp_now_n = tmp_now_n + logl((long double) ptr_primes[current_step- file_offset]);
ratio = ((long double)1)/( (long double)ptr_primes[current_step- file_offset] );

tmp_now_sigma_r = tmp_now_sigma_r + (ratio*tmp_now_sigma_r) ;

tmp_now_quotient = tmp_now_sigma_r;
tmp_now_robin = tmp_now_quotient/expl(euler);
log2nprime = logl(tmp_now_n);
robin_bound = expl(euler)*log2nprime;

delta = robin_bound – tmp_now_quotient;

delta_log = logl(delta);

delta_estimate = constant – log2nprime/((long double) 2);

variation_delta = delta_log – delta_estimate;

if( variation_delta < record_low_variation_delta)
{
record_low_variation_delta = variation_delta;
record_low_log2n = log2nprime;
}

tmp_el_now_robin = tmp_now_robin/logl(tmp_now_n);

tmp_now_robin = tmp_el_now_robin;
if(tmp_now_robin > RecordRobinQuotient)
{
RecordRobinQuotient = tmp_now_robin;
RecordStep = current_step;
}

if( 9999 == (numloops%10000) )
{
// printf(“current step number = %ld “, current_step_number+1);

// printf(“quotient de Robin: %.16Lf “, tmp_now_robin);

// printf(“record step is at: %ld\n”, RecordStep+1);

// printf(“log(log(n’)) = %.16Lf “, log2nprime);

// printf(“Robin upper bound = %.16Lf “, robin_bound);

// printf(“sigma(n’)/n’ = %.16Lf \n”, tmp_now_quotient);

}
}

}

if(RecordStep < numprimes)
{
ptr_exponents[RecordStep]++;
}
if(RecordStep == numprimes)
{

if(num_new_primes < 79999500)
{

for(j=0; j< 100 ; j++)
{

sum_logs_primes = sum_logs_primes + logl( (long double) ptr_primes[numprimes- file_offset]);

ratio = ((long double)1)/((long double)ptr_primes[numprimes-file_offset]);

if( RecordStep > 999999 )
{
numBigPrimesUsed++;
}

num_new_primes++;

prod_sigma_r = prod_sigma_r + (prod_sigma_r*ratio) ;

numprimes++;
}
}
else
{

sum_logs_primes = sum_logs_primes + logl( (long double) ptr_primes[numprimes- file_offset]);

ratio = ((long double)1)/((long double)ptr_primes[numprimes-file_offset]);

if( RecordStep > 999999 )
{
numBigPrimesUsed++;
}

num_new_primes++;

prod_sigma_r = prod_sigma_r + (prod_sigma_r*ratio) ;

numprimes++;
}
skip_base_case = (int) 1;
}
else
{
skip_base_case = (int) 0;
}

if( 9999 == (numloops%10000) )
{
// printf(“=================================\n”);
}
}
}

fclose(in);

in = fopen(“/home/david2/eratosthenes2/LargePrimes_to_5e11ofNov29a.txt”, “r”);
for(kij = 1; kij < 142; kij++)
{

for(j=0; j< 80000000 ; j++)
{
fscanf(in, “%ld”, &ptr_primes[j]);
}

file_offset = file_offset + ((long)80000000) ;
num_new_primes = (long) 0;

numloops = (long) 0;

while( num_new_primes < ((long)80000000) )
{
numloops++;

ptr_isstep[0] = (long) 1;
list_steps[0] = (long) 0;
numsteps = (long) 1;
not_done = (int) 1;
j = (long) 0;

while(not_done == 1)
{
j++;

if(ptr_exponents[j]< ptr_exponents[j-1])
{
ptr_isstep[j] = (long) 1;
list_steps[numsteps] = j;
numsteps++;
}
else
{
ptr_isstep[j] = (long) 0;
}

if(ptr_exponents[j] == ((long) 1))
{
not_done = (int) 0;
}

}

list_steps[numsteps] = numprimes;
numsteps++;

if( 9999 == (numloops%10000) )
{
for(j=0; j< numsteps; j++)
{
// printf(“step at %ld\n”, list_steps[j]+1);
}
// printf(“\nNumber of big primes used in making current n is: %ld\n”, numBigPrimesUsed);
// printf(“numloops = %ld\n”, numloops);
// printf(“\n”);
}
if( 9999 == (numloops%10000) )
{
// printf(“number of steps: %ld\n”, numsteps);
}
if( 9999 == (numloops%10000) )
{
for(j=0; j< 200 ; j++)
{
// printf(“%ld “, ptr_exponents[j]);
}

// printf(“\n”);
}

if( 799999 == (numloops%800000) )
{
printf(“record_low_variation_delta = %.16Lf “, record_low_variation_delta);
printf(“record_low_log2n = %.16Lf\n\n”, record_low_log2n);
fflush(stdout);

record_low_variation_delta = (long double) 100;
record_low_log2n = (long double) 0;
}

RecordRobinQuotient = (long double) 0;
if(skip_base_case == 0)
{

base_tmp_now_n = (long double) 0;

base_tmp_now_sigma_r = (long double) 1;

now_step = (long) 0;

for(j=0; j<= list_steps[numsteps-2]; j++)
{
if(j == list_steps[now_step] )
{
now_step++;
}
else
{
base_tmp_now_n = base_tmp_now_n + ((long double)ptr_exponents[j])*splogs[j];
p = pts_primes[j];
power = (long) 1;
for(l=0; l<ptr_exponents[j]+((long)1); l++)
{
power = power*p;
}
powerb = power/p;
sigma = (power – ((long)1) )/(p – ((long)1));
ratio = ((long double)sigma)/((long double)powerb);
base_tmp_now_sigma_r = base_tmp_now_sigma_r*ratio;
}
}
}

sum_log_cofactor = sum_logs_primes;

prod_sigma_r_cofactor = prod_sigma_r;
sum_log_cofactor = sum_log_cofactor – theta[list_steps[numsteps-2]];

prod_sigma_r_cofactor = prod_sigma_r_cofactor/prod100k[list_steps[numsteps-2]];

for(m=0; m< numsteps ; m++)
{
current_step_number = m;

current_step = list_steps[m];

tmp_now_n = sum_log_cofactor + base_tmp_now_n;
tmp_now_sigma_r = base_tmp_now_sigma_r*prod_sigma_r_cofactor;
if(m < numsteps-((long)1) )
{
for(j=0; j<= numsteps- ((long)2) ; j++)
{
if(j == m)
{
tmp_exponent = ptr_exponents[list_steps[j]] + ((long)1);
}
else
{
tmp_exponent = ptr_exponents[list_steps[j]];
}
p = pts_primes[list_steps[j]];
tmp_now_n = tmp_now_n + ((long double)tmp_exponent)*splogs[list_steps[j]];
power = (long) 1;
for(k=0; k<(tmp_exponent+((long)1) ); k++)
{
power = power*p;
}
powerb = power/p;
sigma = (power – ((long)1) )/(p – ((long)1));

ratio = ((long double)sigma)/((long double)powerb);
tmp_now_sigma_r = tmp_now_sigma_r*ratio;
}

tmp_now_quotient = tmp_now_sigma_r;

tmp_now_robin = tmp_now_quotient/expl(euler);

log2nprime = logl(tmp_now_n);

robin_bound = expl(euler)*log2nprime;

delta = robin_bound – tmp_now_quotient;

delta_log = logl(delta);

delta_estimate = constant – log2nprime/((long double) 2);

variation_delta = delta_log – delta_estimate;

if( variation_delta < record_low_variation_delta)
{
record_low_variation_delta = variation_delta;
record_low_log2n = log2nprime;
}

tmp_el_now_robin = tmp_now_robin/logl(tmp_now_n);

tmp_now_robin = tmp_el_now_robin;

if(tmp_now_robin > RecordRobinQuotient)
{
RecordRobinQuotient = tmp_now_robin;
RecordStep = current_step;
}
if( 9999 == (numloops%10000) )
{

// printf(“current step number = %ld “, current_step_number+1);

// printf(“quotient de Robin: %.16Lf “, tmp_now_robin);

// printf(“record step is at: %ld\n”, RecordStep+1);

// printf(“log(log(n’)) = %.16Lf “, log2nprime);

// printf(“Robin upper bound = %.16Lf “, robin_bound);

// printf(“sigma(n’)/n’ = %.16Lf \n”, tmp_now_quotient);
}
}
else
{
for(j=0; j<= numsteps- ((long)2) ; j++)
{
tmp_exponent = ptr_exponents[list_steps[j]];

p = pts_primes[list_steps[j]];
if(list_steps[j] != 1)
{
tmp_now_n = tmp_now_n + ((long double)tmp_exponent)*splogs[list_steps[j]];
power = (long) 1;
for(k=0; k<(tmp_exponent+1); k++)
{
power = power*p;
}

powerb = power/p;
sigma = (power – ((long)1) )/(p – ((long)1));

ratio = ((long double)sigma)/((long double)powerb);

tmp_now_sigma_r = tmp_now_sigma_r*ratio;
}
else
{
tmp_now_n = tmp_now_n + ((long double)tmp_exponent)*splogs[1];
ratio = expl( puiss3sigmalogs[tmp_exponent]-(((long double)tmp_exponent)*splogs[1]) );

tmp_now_sigma_r = tmp_now_sigma_r*ratio;
}
}

tmp_now_n = tmp_now_n + logl((long double) ptr_primes[current_step- file_offset]);
ratio = ((long double)1)/( (long double)ptr_primes[current_step- file_offset] );

tmp_now_sigma_r = tmp_now_sigma_r + (ratio*tmp_now_sigma_r) ;

tmp_now_quotient = tmp_now_sigma_r;
tmp_now_robin = tmp_now_quotient/expl(euler);
log2nprime = logl(tmp_now_n);
robin_bound = expl(euler)*log2nprime;

delta = robin_bound – tmp_now_quotient;

delta_log = logl(delta);

delta_estimate = constant – log2nprime/((long double) 2);

variation_delta = delta_log – delta_estimate;

if( variation_delta < record_low_variation_delta)
{
record_low_variation_delta = variation_delta;
record_low_log2n = log2nprime;
}

tmp_el_now_robin = tmp_now_robin/logl(tmp_now_n);

tmp_now_robin = tmp_el_now_robin;
if(tmp_now_robin > RecordRobinQuotient)
{
RecordRobinQuotient = tmp_now_robin;
RecordStep = current_step;
}

if( 9999 == (numloops%10000) )
{
// printf(“current step number = %ld “, current_step_number+1);

// printf(“quotient de Robin: %.16Lf “, tmp_now_robin);

// printf(“record step is at: %ld\n”, RecordStep+1);

// printf(“log(log(n’)) = %.16Lf “, log2nprime);

// printf(“Robin upper bound = %.16Lf “, robin_bound);

// printf(“sigma(n’)/n’ = %.16Lf \n”, tmp_now_quotient);

}
}

}

if(RecordStep < numprimes)
{
ptr_exponents[RecordStep]++;
}
if(RecordStep == numprimes)
{

if(num_new_primes < 79999500)
{

for(j=0; j< 100 ; j++)
{

sum_logs_primes = sum_logs_primes + logl( (long double) ptr_primes[numprimes- file_offset]);

ratio = ((long double)1)/((long double)ptr_primes[numprimes-file_offset]);

if( RecordStep > 999999 )
{
numBigPrimesUsed++;
}

num_new_primes++;

prod_sigma_r = prod_sigma_r + (prod_sigma_r*ratio) ;

numprimes++;
}
}
else
{

sum_logs_primes = sum_logs_primes + logl( (long double) ptr_primes[numprimes- file_offset]);

ratio = ((long double)1)/((long double)ptr_primes[numprimes-file_offset]);

if( RecordStep > 999999 )
{
numBigPrimesUsed++;
}

num_new_primes++;

prod_sigma_r = prod_sigma_r + (prod_sigma_r*ratio) ;

numprimes++;
}
skip_base_case = (int) 1;
}
else
{
skip_base_case = (int) 0;
}

if( 9999 == (numloops%10000) )
{
// printf(“=================================\n”);
}
}
}

fclose(in);
free(ptr_primes);
return 0;

}

Advertisements

Written by meditationatae

December 4, 2014 at 6:26 pm

%d bloggers like this: