/*
 * 03-parallel-primes.c  —  Parallel Prime Number Generation   [LAB 4]
 * Session 09 · 2311CSC501J Parallel Processing · Unit III — Load Balancing
 *
 * WHAT THIS SHOWS
 *   Count every prime from 2 to N. Each number is tested independently, so this
 *   looks like the easiest parallel problem in the world — "embarrassingly
 *   parallel", split the range and go.
 *
 *   It is a trap, and that is exactly why it is Lab 4.
 *
 *   The work per number is NOT uniform. Testing whether 999,983 is prime costs
 *   roughly 700x more than testing 3, because trial division runs up to sqrt(n).
 *   So if you hand thread 0 the range [2, N/4] and thread 3 the range [3N/4, N],
 *   thread 3 gets all the expensive numbers while thread 0 finishes early and
 *   sits idle. Equal COUNTS of work, wildly unequal COST of work.
 *
 *   This program runs the identical loop three ways and prints, per thread, how
 *   many numbers it tested and how long it actually spent:
 *
 *     [1] serial                  — the baseline to measure speedup against
 *     [2] schedule(static)        — equal contiguous chunks decided up front
 *     [3] schedule(dynamic, 1000) — threads grab the next chunk when they finish
 *
 *   All three MUST print the same prime count. Only the time differs. That is
 *   the whole point of Session 09: same work, same cores, same answer — and a
 *   large difference in wall-clock, decided purely by WHEN you assign the work.
 *
 * WHY TRIAL DIVISION AND NOT A SIEVE?
 *   A Sieve of Eratosthenes is a far better algorithm for counting primes. It is
 *   also the wrong teaching tool here, for two reasons: its inner steps DEPEND on
 *   each other (marking multiples), and its work per number is roughly uniform —
 *   so there is no load imbalance left to demonstrate. Trial division is slower
 *   and dumber, and it makes the imbalance impossible to miss. Use a sieve in
 *   production; use this to understand scheduling.
 *
 * THINGS TO NOTICE WHEN YOU RUN IT
 *   1. Under static, the LAST thread takes several times longer than the FIRST.
 *      Those first threads are finished and idle — that is the wasted time.
 *   2. Under dynamic, every thread finishes at almost the same moment, and the
 *      threads that got cheap numbers simply came back for more.
 *   3. Under dynamic the "numbers tested" counts drift slightly apart — threads
 *      that drew cheap numbers simply came back for more — while the TIMES all
 *      converge. That is the lesson in one line: balance the cost, not the count.
 *
 * HOW TO COMPILE & RUN
 *       gcc -O2 -fopenmp 03-parallel-primes.c -o primes -lm
 *       OMP_NUM_THREADS=4 ./primes
 *
 *   macOS (Apple clang): clang -O2 -Xpreprocessor -fopenmp -lomp 03-parallel-primes.c -o primes
 *   No compiler on your machine? Use OnlineGDB (onlinegdb.com), language C.
 *   IMPORTANT: OnlineGDB does NOT enable OpenMP by default. You must add the
 *   flag or you get "undefined reference to omp_get_thread_num" at link time:
 *       click the gear icon (next to the input panel) -> Extra Compiler Flags
 *       -> type:  -fopenmp
 *   Same on godbolt.org: add -fopenmp to the compiler options box.
 *
 *   On a real 8-core machine, raise LIMIT to 5000000 to make the gap obvious.
 */

#include <stdio.h>
#include <omp.h>

/* LIMIT controls how long the run takes. 10,000,000 gives roughly 1-2 seconds
 * serially on a modern laptop — long enough that the timings are meaningful and
 * not just noise. On a slow free online compiler, drop it to 2000000. */
#define LIMIT       10000000  /* count primes in [2, LIMIT] */
#define CHUNK       1000      /* chunk size for schedule(dynamic, CHUNK) */
#define MAX_THREADS 64

/* Per-thread statistics.
 *
 * NOTE (Session 05 callback): an array indexed by thread id is the classic
 * false-sharing bug — several threads writing neighbouring slots that share one
 * 64-byte cache line. We avoid it two ways: each entry is PADDED to a full cache
 * line, and each thread accumulates in LOCAL variables and writes here exactly
 * once, at the end. Correct code that is also fast. */
typedef struct {
    long   tested;      /* how many numbers this thread checked   */
    long   primes;      /* how many of them were prime            */
    double seconds;     /* how long this thread actually worked   */
    char   pad[64 - 2 * sizeof(long) - sizeof(double)];
} thread_stat_t;

static thread_stat_t stats[MAX_THREADS];

/* Trial division. Deliberately simple — and deliberately more expensive for
 * larger n, which is what creates the load imbalance this lab is about.
 * We test odd divisors only, up to sqrt(n), using d*d <= n to avoid sqrt(). */
static int is_prime(long n)
{
    if (n < 2)      return 0;
    if (n < 4)      return 1;      /* 2 and 3 */
    if (n % 2 == 0) return 0;

    for (long d = 3; d * d <= n; d += 2)
        if (n % d == 0)
            return 0;

    return 1;
}

/* ---------------------------------------------------------------- serial --- */

static long count_primes_serial(double *elapsed)
{
    long   count = 0;
    double t0    = omp_get_wtime();

    for (long n = 2; n <= LIMIT; n++)
        if (is_prime(n))
            count++;

    *elapsed = omp_get_wtime() - t0;
    return count;
}

/* -------------------------------------------------------------- parallel --- */

/* Runs the loop under whatever schedule was set with omp_set_schedule().
 * Fills stats[] and returns the total prime count. */
static long count_primes_parallel(double *elapsed, int nthreads)
{
    long   total = 0;
    double t0;

    for (int i = 0; i < nthreads; i++) {
        stats[i].tested  = 0;
        stats[i].primes  = 0;
        stats[i].seconds = 0.0;
    }

    t0 = omp_get_wtime();

    #pragma omp parallel num_threads(nthreads)
    {
        int    tid        = omp_get_thread_num();
        long   my_tested  = 0;
        long   my_primes  = 0;
        double my_start   = omp_get_wtime();

        /* schedule(runtime) reads its policy from omp_set_schedule() below, so
         * one function can demonstrate every schedule.
         *
         * nowait matters here: without it the implicit barrier at the end of the
         * loop would make every thread record the SAME finish time, hiding the
         * imbalance we are trying to measure. With nowait each thread stamps the
         * moment IT ran out of work. (The end of the parallel region still has a
         * barrier, so stats[] is complete before we read it.) */
        #pragma omp for schedule(runtime) nowait
        for (long n = 2; n <= LIMIT; n++) {
            my_tested++;
            if (is_prime(n))
                my_primes++;
        }

        stats[tid].seconds = omp_get_wtime() - my_start;
        stats[tid].tested  = my_tested;
        stats[tid].primes  = my_primes;
    }

    *elapsed = omp_get_wtime() - t0;

    for (int i = 0; i < nthreads; i++)
        total += stats[i].primes;

    return total;
}

/* ---------------------------------------------------------------- report --- */

static void report(const char *label, long count, long expected,
                   double elapsed, double serial_time, int nthreads)
{
    double busiest = 0.0, idlest = 1e30, worked = 0.0;
    int    balanced;

    printf("    primes found:  %ld%s\n", count,
           count == expected ? "   (correct)" : "   *** WRONG ***");
    printf("    wall-clock:    %.3f s\n", elapsed);
    printf("    speedup:       %.2fx   efficiency %.0f%%\n",
           serial_time / elapsed,
           100.0 * (serial_time / elapsed) / nthreads);

    for (int i = 0; i < nthreads; i++) {
        worked += stats[i].seconds;
        if (stats[i].seconds > busiest) busiest = stats[i].seconds;
        if (stats[i].seconds < idlest)  idlest  = stats[i].seconds;
    }

    /* "Balanced" = the fastest thread finished within 10% of the slowest, i.e.
     * nobody spent meaningful time waiting. Only call out stragglers when there
     * is genuinely a spread to call out. */
    balanced = (busiest <= 0.0) || ((busiest - idlest) / busiest < 0.10);

    printf("\n    per-thread breakdown:\n");
    for (int i = 0; i < nthreads; i++) {
        const char *note = "";
        if (balanced)
            note = "<-- all finished together";
        else if (stats[i].seconds >= busiest * 0.97)
            note = "<-- everyone waits for this one";
        else if (stats[i].seconds <= busiest * 0.75)
            note = "<-- finished early, then idle";

        printf("      thread %-2d  %9ld numbers   %6.3f s   %s\n",
               i, stats[i].tested, stats[i].seconds, note);
    }

    /* Idle time = what the cores COULD have done minus what they DID do. */
    {
        double capacity = busiest * nthreads;
        double idle     = capacity - worked;
        printf("      %s idle: %.3f core-seconds wasted (%.0f%% of the machine)\n",
               label, idle, 100.0 * idle / capacity);
    }
    printf("\n");
}

/* ------------------------------------------------------------------ main --- */

int main(void)
{
    int    nthreads = omp_get_max_threads();
    long   serial_count, static_count, dynamic_count;
    double serial_time, static_time, dynamic_time;

    if (nthreads > MAX_THREADS)
        nthreads = MAX_THREADS;

    printf("========================================================\n");
    printf(" Lab 4 — Parallel Prime Number Generation\n");
    printf(" Range: 2 .. %d      Threads available: %d\n", LIMIT, nthreads);
    printf("========================================================\n\n");

    if (nthreads < 2) {
        printf(" NOTE: only 1 thread available, so speedup will be 1.00x.\n");
        printf(" The per-thread breakdown still shows the point. Try a\n");
        printf(" multi-core machine, or OMP_NUM_THREADS=4 if you have cores.\n\n");
    }

    printf("[1] SERIAL baseline\n");
    serial_count = count_primes_serial(&serial_time);
    printf("    primes found:  %ld\n", serial_count);
    printf("    wall-clock:    %.3f s\n\n", serial_time);

    printf("[2] PARALLEL — schedule(static)\n");
    printf("    equal contiguous chunks, decided before the loop starts\n");
    omp_set_schedule(omp_sched_static, 0);
    static_count = count_primes_parallel(&static_time, nthreads);
    report("static ", static_count, serial_count, static_time, serial_time, nthreads);

    printf("[3] PARALLEL — schedule(dynamic, %d)\n", CHUNK);
    printf("    threads take the next %d numbers whenever they finish\n", CHUNK);
    omp_set_schedule(omp_sched_dynamic, CHUNK);
    dynamic_count = count_primes_parallel(&dynamic_time, nthreads);
    report("dynamic", dynamic_count, serial_count, dynamic_time, serial_time, nthreads);

    printf("========================================================\n");
    printf(" All three found %ld primes — the answer never changed.\n", serial_count);
    printf(" static  : %.3f s   (%.2fx)\n", static_time,  serial_time / static_time);
    printf(" dynamic : %.3f s   (%.2fx)   %.0f%% faster than static\n",
           dynamic_time, serial_time / dynamic_time,
           100.0 * (static_time - dynamic_time) / static_time);
    printf("\n Same code. Same cores. Same result. The only difference is\n");
    printf(" WHEN the work was assigned — up front, or on demand.\n");
    printf("========================================================\n");

    return 0;
}