added pseudo random generator examples

CMakeLists.txt

@@ -23,6 +23,7 @@ set(DWARF_PREFIX 7_montecarlo) # The prefix of the name of the binaries produced
 
 # Add the examples
 add_subdirectory(pi)
+add_subdirectory(prng)
 add_subdirectory(integral_basic)
 
 # ==================================================================================================

prng/CMakeLists.txt

+# Packages are optional: if they are not present, certain code samples are not compiled
+cmake_minimum_required(VERSION 2.8.10 FATAL_ERROR)
+
+find_package(OpenMP)   # Built-in in CMake
+find_package(MPI)      # Built-in in CMake
+
+include(${CMAKE_CURRENT_SOURCE_DIR}/../../cmake/common.cmake)
+
+# ==================================================================================================
+
+if ("${DWARF_PREFIX}" STREQUAL "")
+  set(DWARF_PREFIX 7_montecarlo)
+endif()
+set(NAME_OMP ${DWARF_PREFIX}_prng_omp)
+set(NAME_SERIAL ${DWARF_PREFIX}_prng_serial)
+set(NAME_MPI ${DWARF_PREFIX}_prng_mpi)
+
+
+# C compiler settings
+
+find_package(Common)
+
+if (OPENMP_FOUND)
+    set (CMAKE_C_FLAGS "${CMAKE_C_FLAGS} ${OpenMP_C_FLAGS}")
+    add_executable(${NAME_OMP} prng_omp.c prng_shared.c)
+    install(TARGETS ${NAME_OMP} DESTINATION bin)
+    message("** Enabling '${NAME_OMP}': with OpenMP")
+else()
+     message("## Skipping 'prng_omp': no OpenMP support found")
+#     dummy_install(${NAME} "OpenMP")
+endif()
+
+if (MPI_FOUND)
+    cmake_policy(SET CMP0003 OLD)
+    set (CMAKE_C_FLAGS "${CMAKE_C_FLAGS}")
+    include_directories(${MPI_INCLUDE_PATH})
+    add_executable(${NAME_MPI} prng_mpi.c prng_shared.c)
+    target_link_libraries(${NAME_MPI} ${MPI_LIBRARIES} stdc++)
+    install(TARGETS ${NAME_MPI} DESTINATION bin)
+    message("** Enabling '${NAME_MPI}': with MPI")
+else()
+     message("## Skipping '${NAME_MPI}': no OpenMP support found")
+#     dummy_install(${NAME} "OpenMP")
+endif()
+
+add_executable(${NAME_SERIAL} prng_serial.c prng_shared.c)
+message("** Enabling '${NAME_SERIAL}': serial version")
+
+install(TARGETS ${NAME_SERIAL} DESTINATION bin)
+
+set(CMAKE_C_FLAGS "${CMAKE_C_FLAGS} ${C_FLAGS}")
+
+unset(NAME)
+# ==================================================================================================

prng/prng_mpi.c

+#include <stdio.h> // printf
+#include <stdlib.h> // EXIT_SUCCESS
+#include <mpi.h>
+
+#include "prng_shared.h"
+
+int main(int argc, char* argv[]) {
+    int seed;
+    int N;
+    int i;
+    int number_of_workers;
+    int N_per_worker;
+    int my_process_id;
+    double *seq;
+    double *my_seq;
+
+    // initialize MPI sybsystem
+    MPI_Init(&argc,&argv);
+    
+    // get number of available processes, and a unique id for this process
+    MPI_Comm_rank(MPI_COMM_WORLD, &my_process_id);
+    MPI_Comm_size(MPI_COMM_WORLD, &number_of_workers);
+
+    // get parameters from command line arguments
+    read_arguments(&seed, &N, argc, argv);
+
+    // get number of workers and this workers id
+    MPI_Comm_rank(MPI_COMM_WORLD, &my_process_id);
+    MPI_Comm_size(MPI_COMM_WORLD, &number_of_workers);
+
+    // exit if work cannot be divided equally among workers
+    if ( (N % number_of_workers) != 0 ) {
+        printf("Cannot equally divide %d numbers by %d workers \n", N, number_of_workers);
+        exit(EXIT_FAILURE);
+    }
+
+    // calculate numbers to be generated by each process
+    N_per_worker = N / number_of_workers;
+
+    // allocate and array to hold N_per_worker random numbers
+    my_seq = malloc(N_per_worker * sizeof(double));
+    // exit if allocation fails
+    if (my_seq == NULL) {
+        printf("Cannot allocate space for %d numbers\n", N_per_worker);
+        exit(EXIT_FAILURE);
+    }
+
+    // adjust seed so that each process creates a different number sequence
+    seed = seed + my_process_id;
+
+    // generate N_per_worker random numbers and save them in the array
+    // please remember that each process will have its own memory context
+    for (i=0; i<N_per_worker; i++) {
+        my_seq[i] = (double)rand_r(&seed) / (double)RAND_MAX;
+    }
+
+    // in the following step, sequences from all workers will be gathered on
+    // the main process. in order to do that, memory space to accommodate all
+    // the numbers have to be allocated
+    if (my_process_id == 0) {
+        seq = malloc(N * sizeof(double));
+        // exit if allocation fails
+        if (seq == NULL) {
+            printf("Cannot allocate space for %d numbers\n", N);
+            exit(EXIT_FAILURE);
+        }
+    }
+
+    // everyone (including the main process) will send N_per_worker numbers from
+    // their my_seq buffer to the seq buffer in the main process, using different offsets
+    MPI_Gather(my_seq, N_per_worker, MPI_DOUBLE, seq, N_per_worker, MPI_DOUBLE, 0, MPI_COMM_WORLD);
+
+    // print the generated random numbers, showing 8 digits after the decimal point
+    // this step will only be executed in the main process
+    if(my_process_id == 0) {
+        for (i=0; i<N; i++) {
+            printf("%10.8lf\n", seq[i]);
+        }
+    }
+
+    // cleanly terminate MPI subsystem
+    MPI_Finalize();
+
+    return EXIT_SUCCESS;
+
+}
+

prng/prng_omp.c

+#include <stdio.h> // printf
+#include <stdlib.h> // EXIT_SUCCESS
+#include <omp.h>
+
+#include "prng_shared.h"
+
+int main(int argc, char* argv[]) {
+    int seed;
+    int N;
+    int i;
+    int number_of_workers;
+    int my_seed;
+    int my_i;
+    int my_thread_id;
+    double *seq;
+
+
+    // get parameters from command line arguments
+    read_arguments(&seed, &N, argc, argv);
+
+    // allocate and array to hold N random numbers
+    seq = malloc(N * sizeof(double));
+    // exit if allocation fails
+    if (seq == NULL) {
+        printf("Cannot allocate space for %d numbers\n", N);
+        exit(EXIT_FAILURE);
+    }
+
+    // create OMP_NUM_THREADS threads
+    // every thread will have its own my_i, my_thread_id and  my_seed variables
+    // the computer will run number_of_workers different for loops concurrently,
+    // one for each thread, and threads will alter the same seq array.
+    // however, due to the striding access patterns, areas affected by different
+    // threads will not overlap, and race conditions are avoided
+#pragma omp parallel private (my_i, my_thread_id, my_seed)
+    {
+        // get number of threads
+        number_of_workers = omp_get_num_threads();
+    
+        // exit if work cannot be divided equally among workers
+        if ( (N % number_of_workers) != 0 ) {
+            printf("Cannot equally divide %d numbers by %d workers \n", N, number_of_workers);
+            exit(EXIT_FAILURE);
+        }
+    
+        my_thread_id = omp_get_thread_num();
+        my_seed = seed + omp_get_thread_num();
+    
+        for (my_i = my_thread_id; my_i < N; my_i += number_of_workers) {
+            seq[my_i] = (double)rand_r(&my_seed) / (double)RAND_MAX;
+        }
+
+    }
+
+    // print the generated random numbers, showing 8 digits after the decimal point
+    for (i=0; i<N; i++) {
+        printf("%10.8lf\n", seq[i]);
+    }
+
+    return EXIT_SUCCESS;
+
+}
+

prng/prng_serial.c

+#include <stdio.h> // printf
+#include <stdlib.h> // EXIT_SUCCESS
+
+#include "prng_shared.h"
+
+int main(int argc, char* argv[]) {
+    int seed;
+    int N;
+    int i;
+    double *seq;
+
+
+    // get parameters from command line arguments
+    read_arguments(&seed, &N, argc, argv);
+
+    // allocate and array to hold N random numbers
+    seq = malloc(N * sizeof(double));
+    // exit if allocation fails
+    if (seq == NULL) {
+        printf("Cannot allocate space for %d numbers\n", N);
+        exit(EXIT_FAILURE);
+    }
+
+    // generate N random numbers and save them in the array
+    for (i=0; i<N; i++) {
+        seq[i] = (double)rand_r(&seed) / (double)RAND_MAX;
+    }
+
+    // print the generated random numbers, showing 8 digits after the decimal point
+    for (i=0; i<N; i++) {
+        printf("%10.8lf\n", seq[i]);
+    }
+
+    return EXIT_SUCCESS;
+
+}
+

prng/prng_shared.c

+#include <stdlib.h> // rand_r, RAND_MAX, EXIT_FAILURE, atoi
+#include <stdio.h> // printf
+#include <time.h> // time
+
+void read_arguments(int *seed, int *N, int argc, char *argv[]) {
+    /*
+     * parse command line arguments to obtain
+     * - seed: seed for the pseudo random number generator
+     * - N: number of random numbers to generate
+     */
+
+    // check if enough arguments are given
+    if (argc != 3) { 
+        // print usage information
+        printf("USAGE %s <seed> <N>\n", argv[0]);
+        printf("  seed (integer): seed for the pseudo random number generator\n");
+        printf("                  a negative value indicates current timestamp\n");
+        printf("  N (integer)   : number of random numbers to generate\n");
+        // do not continue if number of arguments is invalid
+        exit(EXIT_FAILURE);
+    }
+
+    // read seed for the pseudo random number generator
+    // WARNING: according to the man page, atoi does not detect errors
+    *seed = atoi(argv[1]);
+    if (*seed < 0) {
+        // use current timestamp as seed
+        *seed = time(NULL);
+    }
+
+    // read number of random samples
+    // WARNING: according to the man page, atoi does not detect errors
+    *N = atoi(argv[2]);
+
+}
+

prng/prng_shared.h

+/* prng_shared.c */
+void read_arguments(int *seed, int *N, int argc, char *argv[]);