#include <pthread.h>
#include <stdio.h>
#include <stdlib.h>
#include <time.h>
#include <primesieve.h>
/*
SIZE%(NTHREADS*2) must be == 0
*/
#define FILENAME "/dev/stdout"
#define SIZE 100000000 // number of fractions to compute
#define NTHREADS 4 // number of threads to utilize
/*digests per thread. please, choose reasonable number,
not smth like 123. checks are not performed to guarantee accuracy*/
#define RESOLUTION 10
long double sum;
pthread_mutex_t mutex = PTHREAD_MUTEX_INITIALIZER;
pthread_barrier_t barrier;
unsigned long int primes[NTHREADS] = { 0 };
void *job (void *arg)
{
unsigned long int index, n, p, size, end;
long double localsum;
primesieve_iterator it;
time_t t;
struct tm *time_p;
FILE *f;
char buf[BUFSIZ];
pthread_mutex_lock(&mutex);
index = *((unsigned long int*)arg);
size = SIZE/NTHREADS;
n = size*index;
if (index == NTHREADS - 1)
{
end = SIZE - 1;
}
else
{
end = n + size - 1;
}
f = fopen(FILENAME, "w");
setvbuf(f, buf, _IOLBF, BUFSIZ);
fprintf(f, "%ld: from %ld to %ld\n", index, n + 1, end + 1);
time(&t);
time_p = localtime(&t);
fprintf(f, "[%04d-%02d-%02d %02d:%02d:%02d] %ld: skipping...\n",
time_p->tm_year + 1900, time_p->tm_mon + 1, time_p->tm_mday,
time_p->tm_hour, time_p->tm_min, time_p->tm_sec, index);
if (primes[index])
{
primesieve_init(&it);
p = primes[index];
primesieve_skipto(&it, p - 1, p + size);
}
else
{
primesieve_init(&it);
p = primesieve_nth_prime(n + 1, 0);
printf("%ld\n", p);
primesieve_skipto(&it, p - 1, p + size);
}
pthread_mutex_unlock(&mutex);
pthread_barrier_wait(&barrier);
time(&t);
time_p = localtime(&t);
fprintf(f, "[%04d-%02d-%02d %02d:%02d:%02d] %ld: calculating...\n",
time_p->tm_year + 1900, time_p->tm_mon + 1, time_p->tm_mday,
time_p->tm_hour, time_p->tm_min, time_p->tm_sec, index);
while (n <= end)
{
localsum -= (++n)/(long double)primesieve_next_prime(&it);
time(&t);
time_p = localtime(&t);
fprintf(f, "[%04d-%02d-%02d %02d:%02d:%02d] %ld: %ld\n",
time_p->tm_year + 1900, time_p->tm_mon + 1, time_p->tm_mday,
time_p->tm_hour, time_p->tm_min, time_p->tm_sec, index, n);
localsum += (++n)/(long double)primesieve_next_prime(&it);
time(&t);
time_p = localtime(&t);
fprintf(f, "[%04d-%02d-%02d %02d:%02d:%02d] %ld: %ld\n",
time_p->tm_year + 1900, time_p->tm_mon + 1, time_p->tm_mday,
time_p->tm_hour, time_p->tm_min, time_p->tm_sec, index, n);
for (unsigned long int i = 0; i < size/RESOLUTION - 1; i++)
{
localsum -= (++n)/(long double)primesieve_next_prime(&it);
localsum += (++n)/(long double)primesieve_next_prime(&it);
}
}
primesieve_free_iterator(&it);
time(&t);
time_p = localtime(&t);
fprintf(f, "[%04d-%02d-%02d %02d:%02d:%02d] %ld: calculated! partial result: %.20Lf\n",
time_p->tm_year + 1900, time_p->tm_mon + 1, time_p->tm_mday,
time_p->tm_hour, time_p->tm_min, time_p->tm_sec, index, localsum);
fclose(f);
sum += localsum;
return NULL;
}
int main (int argc, char *argv[])
{
pthread_t threads[NTHREADS];
unsigned long int args[NTHREADS];
int i;
FILE *f;
char buf[BUFSIZ];
time_t t;
struct tm *time_p;
sum = 0.0L;
pthread_barrier_init(&barrier, NULL, NTHREADS);
f = fopen(FILENAME, "w");
setvbuf(f, buf, _IOLBF, BUFSIZ);
time(&t);
time_p = localtime(&t);
fprintf(f, "[%04d-%02d-%02d %02d:%02d:%02d] M: started\n\n",
time_p->tm_year + 1900, time_p->tm_mon + 1, time_p->tm_mday,
time_p->tm_hour, time_p->tm_min, time_p->tm_sec);
for (i = 0; i < argc - 1; i++)
{
primes[i] = strtoul(argv[i + 1], NULL, 10);
printf("%ld\n", primes[i]);
}
for (i = 0; i < NTHREADS; i++)
{
args[i] = i;
pthread_create(&threads[i], NULL, job, &args[i]);
}
for (i = 0; i < NTHREADS; i++)
{
pthread_join(threads[i], NULL);
}
fprintf(f, "\nM: total result: %.20Lf\n", sum);
fclose(f);
return EXIT_SUCCESS;
}
Comments