/**
 * Obliczanie calek rownolegle
 * Kamil Pastwa
 */

#include<pthread.h>
#include"pomiar_czasu.h"

//pthread_mutex_t mutex;

double f(double x)
{
    return(4/(1+x*x));
}

typedef struct params
{
    //thread info
    int tid;
    int nthreads;

    // data info
    double x1;  // poczatek przedziaku
    double x2;  // koniec
    int N;      // ilosc trapezow
    double S;
} params;

double S_GLOBAL;

/**
 * Calkowanie sekwencyjne
 * x1, x2 - przedzial calkowania
 * N - ilosc trapezow
 *
 * funkcja zwraca wartosc calki w zadanym zakresie
*/
double integrate1(double x1, double x2, int N)
{
    double s  = 0;                      // wartosc calki
    double dx = (x2 - x1) / N;          // odleglosc miedzy dwoma wezlami

    int i;
    for(i = 1; i < N; i++)
        s += f(x1 + i * dx);

    s = (s + (f(x1) + f(x2)) / 2) * dx;

    return s;
}


void* calc3();
/**
 * Calkowanie rownolegle. Zrownoleglenie petli
 * nthreads - liczba watkow uzytych do policzenia calki
 */
double integrate3(double x1, double x2, int N, int nthreads)
{
    S_GLOBAL=0;
    pthread_t threads_t[nthreads];      // tablica przechowujaca numery id watkow
    params *tp = (params*) malloc(nthreads * sizeof(params));   // tp - thread_params

    int i;
    for(i = 0 ; i < nthreads ; ++i)
    {
        //info o watku
        tp[i].tid = i;
        tp[i].nthreads = nthreads;
        //info o danych
        tp[i].x1 = x1;
        tp[i].x2 = x2;
        tp[i].N = N;

        pthread_create(&threads_t[i], NULL, calc3, (void*) &tp[i] );
    }

    for(i = 0 ; i < nthreads ; ++i)
        pthread_join(threads_t[i], NULL);

    double returnS = 0;
    for(i = 0 ; i < nthreads ; ++i){
   // printf("retS[%d]  = %f\n", i, tp[i].S);
        returnS += tp[i].S;

        }

    return returnS;
}

// funkcja watku dla integrate3
void* calc3(void* _p)
{
    struct params* p = (params*) _p;
    //printf("Watek nr %d. x1= %f x2=%f N=%d. Ilosc wszystkich watkow:%d\n", p->tid, p->x1, p->x2, p->N, p->nthreads);

    int tid = p->tid;
    int chunk_size = (p->N / p->nthreads);

    int start = tid * chunk_size;
    int end = start + chunk_size;
    //printf("Watek nr %d. start= %d end= %d\n", tid, start, end);

    int i;
    double dx = (p->x2 - p->x1) / p->N;

    double _s = 0;
    for(i = start ; i < end ; ++i)
        _s += f( p->x1 + i * dx);

    _s = (_s + (f(p->x1) + f(p->x2)) / 2) * dx;

    p->S = _s;
    //printf("Watek nr %d. wynik lokalny = %f\n", tid, _s);

  //  pthread_mutex_lock(&mutex);
    //S_GLOBAL += _s;
    //pthread_mutex_unlock(&mutex);
}

int main()
{
    //pthread_mutex_init(&mutex, NULL);
    const int N = 10000000;                // liczba trapezĂłw
    double x1 = 0;                      // poczatek przedzialu calkowania
    double x2 = 1;                      // koniec przedzialu calkowania
    int nthreads = 4;                   // liczba watkow
    //int n = 10000;                      // ilosc wywolan funkcji liczacej

    int i,j;
    double _s, sekwS;

    // ==== Calkowanie sekwencyjne =====
    inicjuj_czas();
    //for(i = 0 ; i < n ; ++i)
        _s = integrate1(x1,x2,N);
    printf("Wynik calkowania sekwencyjnego S= %f\n", _s);
    drukuj_czas();

    sekwS = _s;
    // ==== Calkowanie rownolegle - zrownoleglenie petli =====
    double czasCalkowityRownolegly = czas_zegara();
    double localTime;

    //for(j = 1 ; j <= nthreads ; ++j){
        //localTime=czas_zegara();
        printf("\n----\nLiczba watkow: %d\n", nthreads);
        printf("wynik= %f\n", _s);
        inicjuj_czas();
            _s = integrate3(x1,x2,N,nthreads);
        drukuj_czas();
        //localTime=czas_zegara() - localTime;

        //printf("localTime = %f\n", localTime);
    //}

    czasCalkowityRownolegly = czas_zegara() - czasCalkowityRownolegly;
    printf("\n\nWynik calkowania rownoleglego S= %f\nCzas calkowania rownoleglego: %f", _s, czasCalkowityRownolegly);

    //printf("Wynik calkowania ze zrownolegleniem petli S= %f\nWynik 2 sposobem:\t%f", S_GLOBAL, test);
    nthreads = 4;
    printf("\n\n\n===== TESTOWANIE DOKLADNOSCI dla %d watkow =======\n\n", nthreads);
    for(i = 10 ; i < 10000000 ; i*=10){
        _s = integrate3(x1,x2,i,nthreads);
        printf("Wynik dla n = %i \t=\t%f\tblad= %f\n", i, _s, _s-sekwS);
        }
    return 0;
}