how to compute sum of n/m Gregory-Leibniz terms in C language

Question

get the two values named m & n from the command line arguments and convert them into integers. now after that create m threads and each thread computes the sum of n/m terms in Gregory-Leibniz Series.

pi = 4 * (1 - 1/3 + 1/5 - 1/7 + 1/9 - ...)

Now when thread finishes its computation, print its partial sum and atomically add it to a shared global variable.

& how to check that all of the m computational threads have done the atomic additions?

I share my source code, what I tried

#include<stdio.h>
#include<pthread.h>
#include <stdlib.h>
#include<math.h>

pthread_barrier_t barrier;
int count;
long int term;
// int* int_arr;
double total;

void *thread_function(void *vargp)
{
    int thread_rank = *(int *)vargp;
    // printf("waiting for barrier... \n");
    
    pthread_barrier_wait(&barrier);
    // printf("we passed the barrier... \n");
    
    double sum = 0.0;
    int n = count * term;
    int start = n - term;
    
    // printf("start %d & end %d \n\n", start, n);
    for(int i = start; i < n; i++) 
    {
        sum += pow(-1, i) / (2*i+1);
        // v +=  1 / i - 1 / (i + 2);
    }
    total += sum; 
    // int_arr[count] = sum;
    count++;
    
    printf("thr %d : %lf \n", thread_rank, sum);
    
    return NULL;
}

int main(int argc,char *argv[])
{
    if (argc <= 2) {
        printf("missing arguments. please pass two num. in arguments\n");
        exit(-1);
    }
    
    int m = atoi(argv[1]); // get value of first argument
    int n = atoi(argv[2]); // get value of second argument
    // int_arr = (int*) calloc(m, sizeof(int));
    
    count = 1;
    term = n / m;

    pthread_t thread_id[m];
    int i, ret;
    double pi;
    
    /* Initialize the barrier. */
    pthread_barrier_init(&barrier, NULL, m);

    for(i = 0; i < m; i++)
    {
        ret = pthread_create(&thread_id[i], NULL , &thread_function, (void *)&i);
        if (ret) {
            printf("unable to create thread! \n");
            exit(-1);
        } 
    }
    
    for(i = 0; i < m; i++)
    {
        if(pthread_join(thread_id[i], NULL) != 0) {
            perror("Failed to join thread");
        }
    }
    
    pi = 4 * total;
    printf("%lf ", pi);
    
    pthread_barrier_destroy(&barrier);


    return 0;
}

what I need :- create M thread & each thread computes the sum of n/m terms in the Gregory-Leibniz Series.

first thread computes the sum of term 1 to n/m , the second thread computes the sum of the terms from (n/m + 1) to 2n/m etc.

when all the thread finishes its computation than print its partial sum and Value of Pi.

I tried a lot, but I can't achieve exact what I want. I got wrong output value of PI

for example : m = 16 and n = 1024

then it sometimes return 3.125969, sometimes 12.503874 , 15.629843, sometimes 6.251937 as a output of Pi value

please help me

---

Edited Source Code :

#include <inttypes.h>
#include <math.h>
#include <pthread.h>
#include <string.h>
#include <stdlib.h>
#include <stdio.h>

struct args {
    uint64_t thread_id;
    struct {
        uint64_t start;
        uint64_t end;
    } range;
};

pthread_mutex_t mutex = PTHREAD_MUTEX_INITIALIZER;
pthread_barrier_t barrier;
long double total = 0;
uint64_t total_iterations = 0;

void *partial_sum(void *arg)
{
    struct args *args = arg;
    long double sum = 0;
    
    printf("waiting for barrier in thread -> %" PRId64 "\n", args->thread_id);
    pthread_barrier_wait(&barrier);
    // printf("we passed the barrier... \n");
    
    for (uint64_t n = args->range.start; n < args->range.end; n++) 
        sum += pow(-1.0, n) / (1 + n * 2);
    
    if (pthread_mutex_lock(&mutex)) {
        perror("pthread_mutex_lock");
        exit(EXIT_FAILURE);
    }
    
    total += sum;
    total_iterations += args->range.end - args->range.start;
    
    if (pthread_mutex_unlock(&mutex)) {
        perror("pthread_mutex_unlock");
        exit(EXIT_FAILURE);
    }
    
    printf("thr %" PRId64 " : %.20Lf\n", args->thread_id, sum);
    
    return NULL;
}

int main(int argc,char *argv[])
{
    if (argc <= 2) {
        fprintf(stderr, "usage: %s THREADS TERMS.\tPlease pass two num. in arguments\n", *argv);
        return EXIT_FAILURE;
    }
    
    int m = atoi(argv[1]); // get value of first argument  & converted into int
    int n = atoi(argv[2]); // get value of second argument & converted into int

    if (!m || !n) {
        fprintf(stderr, "Argument is zero.\n");
        return EXIT_FAILURE;
    }
    
    uint64_t threads = m;
    uint64_t terms = n;
    
    uint64_t range = terms / threads;
    uint64_t excess = terms - range * threads;
    
    pthread_t thread_id[threads];
    struct args arguments[threads];
    
    int ret;
    
    /* Initialize the barrier. */
    ret = pthread_barrier_init(&barrier, NULL, m);
    if (ret) {
        perror("pthread_barrier_init");
        return EXIT_FAILURE;
    }
    
    for (uint64_t i = 0; i < threads; i++) {
        arguments[i].thread_id = i;
        arguments[i].range.start = i * range;
        arguments[i].range.end = arguments[i].range.start + range;
        
        if (threads - 1 == i)
            arguments[i].range.end += excess;
            
        printf("In main: creating thread %ld\n", i);
        ret = pthread_create(thread_id + i, NULL, partial_sum, arguments + i);
        if (ret) {
            perror("pthread_create");
            return EXIT_FAILURE;
        }
    }
    
    for (uint64_t i = 0; i < threads; i++)
        if (pthread_join(thread_id[i], NULL))
            perror("pthread_join");

    pthread_barrier_destroy(&barrier);

    printf("Pi value is : %.10Lf\n", 4 * total);
    printf("COMPLETE? (%s)\n", total_iterations == terms ? "YES" : "NO");
 
    return 0;   
}

Answer 1

In each thread, the count variable is expected to be of a steadily increasing value in this expression

int n = count * term;

being one larger than it was in the "previous" thread, but count is only increased later on in each thread.

Even if you were to "immediately" increase count, there is nothing that guards against two or more threads attempting to read from and write to the variable at the same time.

The same issue exists for total.

The unpredictability of these reads and writes will lead to indeterminate results.

When sharing resources between threads, you must take care to avoid these race conditions. The POSIX threads library does not contain any atomics for fundamental integral operations.

You should protect your critical data against a read/write race condition by using a lock to restrict access to a single thread at a time.

The POSIX threads library includes a pthread_mutex_t type for this purpose. See:

• pthread_mutex_init / pthread_mutex_destroy
• pthread_mutex_lock / pthread_mutex_unlock

---

Additionally, as pointed out by @Craig Estey, using (void *) &i as the argument to the thread functions introduces a race condition where the value of i may change before any given thread executes *(int *) vargp;.

The suggestion is to pass the value of i directly, storing it intermediately as a pointer, but you should use the appropriate type of intptr_t or uintptr_t, which are well defined for this purpose.

pthread_create(&thread_id[i], NULL , thread_function, (intptr_t) i)

int thread_rank = (intptr_t) vargp;

---

How to check that all of the m computational threads have done the atomic additions?

Sum up the number of terms processed by each thread, and ensure it is equal to the expected number of terms. This can also naturally be assumed to be the case if all possible errors are accounted for (ensuring all threads run to completion and assuming the algorithm used is correct).

---

A moderately complete example program:

#define _POSIX_C_SOURCE 200809L
#include <inttypes.h>
#include <math.h>
#include <pthread.h>
#include <stdint.h>
#include <stdio.h>
#include <stdlib.h>

struct args {
    uint64_t thread_id;
    struct {
        uint64_t start;
        uint64_t end;
    } range;
};

pthread_mutex_t mutex;
long double total = 0;
uint64_t total_iterations = 0;

void *partial_sum(void *arg)
{
    struct args *args = arg;
    long double sum = 0;

    for (uint64_t n = args->range.start; n < args->range.end; n++)
        sum += pow(-1.0, n) / (1 + n * 2);

    if (pthread_mutex_lock(&mutex)) {
        perror("pthread_mutex_lock");
        exit(EXIT_FAILURE);
    }

    total += sum;
    total_iterations += args->range.end - args->range.start;

    if (pthread_mutex_unlock(&mutex)) {
        perror("pthread_mutex_unlock");
        exit(EXIT_FAILURE);
    }

    printf("thread(%" PRId64 ") Partial sum: %.20Lf\n", args->thread_id, sum);

    return NULL;
}

int main(int argc,char **argv)
{
    if (argc < 3) {
        fprintf(stderr, "usage: %s THREADS TERMS\n", *argv);
        return EXIT_FAILURE;
    }

    uint64_t threads = strtoull(argv[1], NULL, 10);
    uint64_t terms = strtoull(argv[2], NULL, 10);

    if (!threads || !terms) {
        fprintf(stderr, "Argument is zero.\n");
        return EXIT_FAILURE;
    }

    uint64_t range = terms / threads;
    uint64_t excess = terms - range * threads;

    pthread_t thread_id[threads];
    struct args arguments[threads];

    if (pthread_mutex_init(&mutex, NULL)) {
        perror("pthread_mutex_init");
        return EXIT_FAILURE;
    }

    for (uint64_t i = 0; i < threads; i++) {
        arguments[i].thread_id = i;
        arguments[i].range.start = i * range;
        arguments[i].range.end = arguments[i].range.start + range;

        if (threads - 1 == i)
            arguments[i].range.end += excess;

        int ret = pthread_create(thread_id + i, NULL , partial_sum, arguments + i);

        if (ret) {
            perror("pthread_create");
            return EXIT_FAILURE;
        }
    }

    for (uint64_t i = 0; i < threads; i++)
        if (pthread_join(thread_id[i], NULL))
            perror("pthread_join");

    pthread_mutex_destroy(&mutex);

    printf("%.10Lf\n", 4 * total);
    printf("COMPLETE? (%s)\n", total_iterations == terms ? "YES" : "NO");
}

Using 16 threads to process 1 billion terms:

$ ./a.out 16 10000000000
thread(14) Partial sum: 0.00000000000190476190
thread(10) Partial sum: 0.00000000000363636364
thread(2) Partial sum: 0.00000000006666666667
thread(1) Partial sum: 0.00000000020000000000
thread(8) Partial sum: 0.00000000000555555556
thread(15) Partial sum: 0.00000000000166666667
thread(0) Partial sum: 0.78539816299744868408
thread(3) Partial sum: 0.00000000003333333333
thread(13) Partial sum: 0.00000000000219780220
thread(11) Partial sum: 0.00000000000303030303
thread(4) Partial sum: 0.00000000002000000000
thread(5) Partial sum: 0.00000000001333333333
thread(7) Partial sum: 0.00000000000714285714
thread(6) Partial sum: 0.00000000000952380952
thread(12) Partial sum: 0.00000000000256410256
thread(9) Partial sum: 0.00000000000444444444
3.1415926535
COMPLETE? (YES)