The law is Amdahl's Law.

The formula for the execution time T(p) using p processors, given the sequential execution time T and the fraction of parallel execution f (assuming perfect speedup), is:

We use T(p) = T_s + T_p/p.

System of equations:

1. 100 = T_s + T_p/2
2. 60 = T_s + T_p/6

Subtract equation (2) from (1):

100 - 60 = (T_p/2) - (T_p/6)
40 = T_p(1/2 - 1/6) = T_p(3/6 - 1/6) = T_p(2/6) = T_p/3

T_p = 40 × 3 = 120s

Substitute T_p back into equation (1):

100 = T_s + 120/2
100 = T_s + 60
T_s = 40s

The total sequential execution time T = T_s + T_p = 40s + 120s = 160s.

The fraction of parallel execution f is T_p / T:
f = 120 / 160 = 3/4 = 0.75.

The fraction of parallel execution is 0.75 (or 75%).

According to the course summary (Section 1.1 Motivation for Parallel Computing):

1. Single-threaded performance improvements are no longer significant (due to power/heat limitations).
2. The trend in hardware architecture is towards multi-core chips, meaning performance scales with the number of cores utilized.

The two types of scaling (Section 10.1 Scaling Types) are:

1. Strong Scaling: The total problem size is kept constant. It measures how execution time decreases as the number of processors increases.
2. Weak Scaling: The workload (problem size) per processor is kept constant. It measures the ability to solve larger problems in the same amount of time as the number of processors increases.

According to the summary (Section 4.4.1 Memory Consistency), sequential consistency requires that:

1. The result of any execution is the same as if the operations of all processors were executed in some sequential order (a global interleaving exists), AND
2. The operations of each individual processor appear in this sequence in the order specified by its program (local program order is preserved).

The primary reason is scalability (Section 6.6). Maintaining cache coherence across a large number of nodes in a shared memory system (like NUMA) introduces significant complexity, communication overhead (coherence traffic), and cost, which limits the effective scaling of the system (Section 6.6.2). Distributed memory architectures avoid this overhead.

A 3x3 2D Torus:

A Data Race can occur (Section 4.1.1). This leads to undefined behavior and non-deterministic results.

A method to prevent it is using a Mutex (Mutual Exclusion Lock) (Section 2.4.1). A mutex ensures that only one thread can access the critical section (the code modifying the shared data) at a time.

#define NUM_THREADS 8
#include <thread>
#include <algorithm>

// The 'max' parameter is used here as a pointer to the storage location
// for this specific thread's local maximum.
void thread_max(int start, int end, double *data , double *max) {
    double local_max = -1.0; // Initialize based on non-negative data assumption
    for (int i = start; i < end; ++i) {
        if (data[i] > local_max) {
            local_max = data[i];
        }
    }
    *max = local_max;
}

void find_max_in_parallel(double *array, int array_size, double *array_max) {
    int local_size, local_start, local_end;
    std::thread threads[NUM_THREADS];
    double local_maxes[NUM_THREADS]; // Storage for local results

    // Calculate base chunk size and remainder for balanced distribution
    local_size = array_size / NUM_THREADS;
    int remainder = array_size % NUM_THREADS;
    int current_start = 0;

    *array_max = -1.0;

    for (int i = 0; i < NUM_THREADS; ++i) {
        local_start = current_start;
        // Distribute the remainder among the first 'remainder' threads
        local_end = local_start + local_size + (i < remainder ? 1 : 0);

        // Ensure bounds safety and handle case where array_size < NUM_THREADS
        if (local_end > array_size) local_end = array_size;
        current_start = local_end;

        // Launch the thread only if it has work
        if (local_start < local_end) {
            // Launch the thread, passing the address of its slot in local_maxes
            threads[i] = std::thread(thread_max,
                local_start, local_end, array, &local_maxes[i]
            );
        }
    }

    // Join threads (Fork/Join model, Section 2.3.2) and perform final reduction
    for (int i = 0; i < NUM_THREADS; ++i) {
        // Check if the thread was launched
        if (threads[i].joinable()) {
            threads[i].join();
            // Reduce (must happen after join to ensure data is written)
            if (local_maxes[i] > *array_max) {
                *array_max = local_maxes[i];
            }
        }
    }
}

double random_noise();

int operation(double *x, double *y, int n)
{
    int i;
    double w;
    // Initialize global minimums (assuming n > 0)
    double min_y = y[0];
    int min_y_index = 0;

    // Start parallel region. Explicitly define data sharing attributes.
    #pragma omp parallel default(none) shared(x, y, n, min_y, min_y_index) private(i, w)
    {
        // Thread-private copies for local minimum tracking
        double local_min_y = min_y;
        int local_min_y_index = min_y_index;

        // Distribute loop iterations among threads (Section 3.3.3)
        #pragma omp for
        for (i = 0; i < n - 1; ++i) {
            w = random_noise();
            x[i] = (y[i] + y[i + 1]) / 2. + w;

            // Update local minimum
            if (y[i + 1] < local_min_y) {
                local_min_y = y[i + 1];
                local_min_y_index = i;
            }
        }
        // Implicit barrier at the end of omp for

        // Update global minimum in a critical section (Section 3.5.1)
        #pragma omp critical
        {
            if (local_min_y < min_y) {
                min_y = local_min_y;
                min_y_index = local_min_y_index;
            }
        }
    } // End parallel region

    // Handle the final iteration sequentially
    if (n > 0) {
        i = n - 1;
        x[i] = y[i] + random_noise();
    }
    return min_y_index;
}

1. NUMA Effects / First Touch Policy (Section 6.5.2): Array a is initialized sequentially. The First Touch policy allocates the memory on the socket running the main thread. In the parallel loop, threads on the second socket will suffer high latency due to remote memory access.
2. False Sharing (Section 6.4): schedule(static, 1) distributes iterations cyclically. Adjacent iterations (accessing adjacent memory locations) are processed by different threads. Writes to the same cache line by different threads cause excessive cache coherence traffic (cache line bouncing).
3. Poor Spatial Locality / Cache Thrashing (Section 6.4): The cyclic distribution destroys spatial locality for individual threads, as they access scattered elements rather than contiguous blocks. This leads to inefficient cache utilization and prefetching.

(Note: Corrected initialization from 0.0f to 1.0f)

float operation(float *a, float *b, int n) {
    float prod = 1.0f; // Corrected initialization for product reduction
    #pragma omp simd reduction(*:prod)
    for (int k = 0; k < n; k++) {
        prod *= a[k] + b[k];
    }
    return prod;
}

struct Node {
    struct Node *left;
    struct Node *right;
};
void visit(struct Node *p);
void traverse(struct Node *p)
{
    if (p->left != nullptr) {
        #pragma omp task
        traverse(p->left);
    }
    if (p->right != nullptr) {
        #pragma omp task
        traverse(p->right);
    }
    // Wait for child tasks to complete before visiting the node
    #pragma omp taskwait
    visit(p);
}
void traverse_tree_in_parallel(struct Node *tree)
{
    // Start parallel region
    #pragma omp parallel
    {
        // Use a single thread to initiate the traversal (Section 4.3.3)
        #pragma omp single
        traverse(tree);
    }
}

Compilers must be conservative (Section 5.5.2). Reasons include:

1. Potential Data Dependencies / Pointer Aliasing: The compiler cannot prove that different pointers used in the loop do not point to overlapping memory regions (aliasing). To guarantee correctness, it assumes dependencies exist and prevents vectorization.
2. Complex Control Flow or Function Calls: The presence of complex branches (if/else) or calls to functions that cannot be inlined or vectorized prevents the compiler from generating efficient SIMD code (Section 5.5.4).

#include <immintrin.h>

void linear_and_avg_vec(float *a, float *b, float c, int size, float *avg_out)
{
    __m256 sum256, a256, b256, c256, prod256, new_a256;
    const int VL = 8; // Vector Length for float (32-bit) on 256-bit registers

    // Initialize sum vector to zero
    sum256 = _mm256_set1_ps(0.0f);
    // Broadcast scalar 'c' into a vector
    c256 = _mm256_set1_ps(c);

    float sum = 0.0f;
    int i = 0;

    // Main Vector Loop
    for (; i < size - (size % VL); i += VL) {
        // Load unaligned data
        a256 = _mm256_loadu_ps(&a[i]);
        b256 = _mm256_loadu_ps(&b[i]);

        // Calculate prod = c * a
        prod256 = _mm256_mul_ps(c256, a256);

        // Calculate new_a = (c * a) + b
        new_a256 = _mm256_add_ps(prod256, b256);

        // Store the result back to a (unaligned) - Required step
        _mm256_storeu_ps(&a[i], new_a256);

        // Accumulate the sum
        sum256 = _mm256_add_ps(sum256, new_a256);
    }

    // Horizontal Reduction
    float temp_sum[VL];
    _mm256_storeu_ps(temp_sum, sum256);
    for(int j=0; j<VL; ++j) {
        sum += temp_sum[j];
    }

    // Remainder Loop (Scalar) (Section 5.4.1)
    for (; i < size; ++i) {
        a[i] = c * a[i] + b[i];
        sum += a[i];
    }

    *avg_out = sum / size;
}

• Local routine: A routine whose completion depends only on the local calling process. It does not require communication or synchronization events to occur at another process to complete.
  • Example: MPI_Comm_rank, MPI_Comm_size.
• Non-local routine: A routine whose completion may depend on actions taken by other MPI processes (e.g., matching communication calls or collective participation).
  • Example: MPI_Barrier, MPI_Bcast.

• Incomplete procedure: A procedure that initiates a communication operation but may return before the operation is locally complete. Resources (like buffers) cannot be reused yet.
  • Example: MPI_Isend, MPI_Irecv.
• Completing procedure: A procedure that ensures or verifies the completion of an operation. Upon return, resources involved are available for reuse.
  • Example: MPI_Wait, MPI_Test (if completion is detected), or blocking calls like MPI_Send.

Explanation:

Passive target synchronization (Section 9.11.2, Lock/Unlock) is a mechanism in Remote Memory Access (RMA) where the target process does not actively participate in the synchronization calls. The origin process entirely controls the access epoch to the target's memory window using MPI_Win_lock (to start the epoch) and MPI_Win_unlock (to end it).

Concurrency Options (Lock Types):

1. Exclusive Lock (MPI_LOCK_EXCLUSIVE): Grants the origin process exclusive access. No other process can access it concurrently.
2. Shared Lock (MPI_LOCK_SHARED): Allows multiple processes to concurrently access the target window, provided they also use shared locks (e.g., for concurrent reads).

#include <mpi.h>
// Assume TOTAL_LENGTH, read_data(), update_top3(), print_top3() are defined.

int main(int argc, char** argv) {
    // Initialize MPI environment
    MPI_Init(&argc, &argv);
    int rank, world_size;
    MPI_Comm_rank(MPI_COMM_WORLD, &rank);
    MPI_Comm_size(MPI_COMM_WORLD, &world_size);

    const int MANAGER_RANK = world_size - 1;

    float *full_data = nullptr;
    // Initialize data only on the manager (last rank)
    if (rank == MANAGER_RANK) {
        full_data = new float[TOTAL_LENGTH];
        read_data(full_data);
    }

    int chunk_size = TOTAL_LENGTH / world_size;
    float *local_data = new float[chunk_size];

    // Distribute the data efficiently using MPI_Scatter (Section 9.7.3)
    MPI_Scatter(full_data, chunk_size, MPI_FLOAT, local_data, chunk_size, MPI_FLOAT, MANAGER_RANK, MPI_COMM_WORLD);

    float top3[3] = { -1.0f, -1.0f, -1.0f };

    if (rank == MANAGER_RANK) { // Manager
        // 1. Manager computes its own local top3
        update_top3(local_data, top3, chunk_size);

        int num_workers = world_size - 1;
        MPI_Request *requests = new MPI_Request[num_workers];
        // Buffer to hold incoming results (3 floats per worker)
        float *recv_buf = new float[3 * num_workers];

        // 2. Initiate non-blocking receives (MPI_Irecv, Section 9.5.1)
        for (int from_rank = 0; from_rank < num_workers; ++from_rank) {
            MPI_Irecv(&recv_buf[from_rank * 3], 3, MPI_FLOAT, from_rank, 0, MPI_COMM_WORLD, &requests[from_rank]);
        }

        // 3. Wait for completion and process results concurrently
        int completed = 0;
        while (completed < num_workers) {
            // Wait for any request to complete (MPI_Waitany)
            int index;
            MPI_Waitany(num_workers, requests, &index, MPI_STATUS_IGNORE);

            // Process the completed request immediately (Overlap achieved here)
            if (index != MPI_UNDEFINED) {
                update_top3(&recv_buf[index * 3], top3, 3);
                completed++;
            }
        }
        delete[] recv_buf; delete[] requests;
    } else { // Worker
        // 1. Worker computes its local top3
        update_top3 (local_data, top3, chunk_size);

        // 2. Worker sends result to the manager
        MPI_Send(top3, 3, MPI_FLOAT, MANAGER_RANK, 0, MPI_COMM_WORLD);
    }

    if (rank == MANAGER_RANK) {
        print_top3 (top3);
        delete[] full_data; // Clean up full data on manager
    }

    delete[] local_data;
    MPI_Finalize();
    return 0;
}

#include <omp.h>
#include <mpi.h>

void compute_energy(double *local_mass, double *local_velocity,
    int local_count, double *total_energy)
{
    double local_energy_sum = 0.0;

    // 1. Local computation using OpenMP parallel for with reduction (Section 3.6)
    // E = 0.5 * m * v^2
    #pragma omp parallel for reduction(+:local_energy_sum) schedule(static)
    for (int i = 0; i < local_count; ++i) {
        double m = local_mass[i];
        double v = local_velocity[i];
        local_energy_sum += 0.5 * m * v * v;
    }

    // 2. Global reduction using MPI_Allreduce (Section 9.7.4)
    // MPI_Allreduce ensures the result is available in every rank.
    MPI_Allreduce(&local_energy_sum, total_energy, 1, MPI_DOUBLE, MPI_SUM, MPI_COMM_WORLD);
}

Performance Issues:

1. Poor MPI Process Placement (Resource Imbalance): The MPI binding report shows all three ranks are bound to the first socket. Sockets 1 and 2 are idle, while Socket 0 is contended. This violates the ideal hybrid strategy of one MPI process per socket (Section 9.13.4).
2. OpenMP Thread Undersubscription: OMP_NUM_THREADS = '2'. Each socket has 6 cores. Using only 2 threads significantly underutilizes the resources.
3. Lack of OpenMP Thread Pinning (Poor Locality): OMP_PROC_BIND = 'FALSE'. Threads are not bound (Section 6.5.3) and may migrate, hurting cache locality.

Commands to Fix:

We aim to utilize the 6 physical cores per socket (ignoring SMT is often better for HPC).

1. Set OpenMP Environment Variables:

export OMP_NUM_THREADS=6       # Utilize the physical cores per socket
export OMP_PLACES=cores        # Bind to physical cores
export OMP_PROC_BIND=close     # Keep threads close together within the socket

2. Fix MPI Execution Command (assuming OpenMPI):

mpirun --map-by socket --bind-to socket -np 3 ./your_program