diff --git a/exercise7/.gitignore b/exercise7/.gitignore new file mode 100644 index 0000000..6f8f12c --- /dev/null +++ b/exercise7/.gitignore @@ -0,0 +1,4 @@ +data/ +sequential +parallel +*.png diff --git a/exercise7/Makefile b/exercise7/Makefile new file mode 100644 index 0000000..690b227 --- /dev/null +++ b/exercise7/Makefile @@ -0,0 +1,30 @@ +CC=gcc +PARALLEL_CC=nvcc +CFLAGS+= -std=c99 -O2 -Wall -Wextra +LDLIBS+= -lm +SEQUENTIAL_SRC_FILES=wave_2d_sequential.c +PARALLEL_SRC_FILES=wave_2d_parallel.cu +IMAGES=$(shell find data -type f | sed s/\\.dat/.png/g | sed s/data/images/g ) +.PHONY: all clean dirs plot movie +all: dirs ${TARGETS} +dirs: + mkdir -p data images +sequential: ${SEQUENTIAL_SRC_FILES} + $(CC) $^ $(CFLAGS) -o $@ $(LDLIBS) +parallel: ${PARALLEL_SRC_FILES} + $(PARALLEL_CC) $^ -O2 -fmad=false -o $@ $(LDLIBS) +plot: ${IMAGES} +images/%.png: data/%.dat + ./plot_image.sh $< +movie: ${IMAGES} + ffmpeg -y -an -i images/%5d.png -vcodec libx264 -pix_fmt yuv420p -profile:v baseline -level 3 -r 12 wave.mp4 +check: dirs sequential parallel + mkdir -p data_sequential + ./sequential + cp -rf ./data/* ./data_sequential + ./parallel + python3 compare.py data_sequential/00000.dat data/00000.dat + python3 compare.py data_sequential/00075.dat data/00075.dat + rm -rf data_sequential +clean: + -rm -fr sequential parallel data images wave.mp4 diff --git a/exercise7/README.md b/exercise7/README.md new file mode 100644 index 0000000..c8b89dd --- /dev/null +++ b/exercise7/README.md @@ -0,0 +1,7 @@ +* make : Creates the output folders +* make parallel : Builds the parallel version +* make sequential : builds the sequential version +* ./parallel : Runs the parallel version and fills the 'data/' directory with stored time steps +* make plot : converts saved time steps to png files under 'images/', using gnuplot. Runs faster if launched with e.g. 4 threads (make -j4 plot). +* make movie : converts collection of png files under 'images' into an mp4 movie file, using ffmpeg +* make check : builds both executeables and compares their output diff --git a/exercise7/compare.py b/exercise7/compare.py new file mode 100644 index 0000000..da09fc2 --- /dev/null +++ b/exercise7/compare.py @@ -0,0 +1,48 @@ +import sys + +import numpy as np + + +class bcolors: + OKGREEN = "\033[92m" + FAIL = "\033[91m" + ENDC = "\033[0m" + + +def main() -> None: + error_margin = 1e-4 + + if len(sys.argv) != 3: + print("Usage: python plot_difference.py file1.dat file2.dat") + sys.exit(1) + + file1, file2 = sys.argv[1], sys.argv[2] + + M = 128 + N = 128 + + data1 = load_data(file1, M, N) + data2 = load_data(file2, M, N) + + # Compute difference + diff = data2 - data1 + + if (diff > error_margin).any(): + print( + f"\n{bcolors.FAIL}Data files {file1} and {file2} differ by more than {error_margin}{bcolors.ENDC}\n" + ) + else: + print( + f"\n{bcolors.OKGREEN}Data files {file1} and {file2} are identical within the margin of {error_margin}{bcolors.ENDC}\n" + ) + + +def load_data(filename: str, M: int, N: int) -> np.ndarray: + data = np.fromfile(filename, dtype=np.float64) + if data.size != M * N: + raise ValueError(f"File {filename} does not contain M*N={M*N} entries") + return data.reshape((M, N)) + + +if __name__ == "__main__": + main() diff --git a/exercise7/plot_differences.py b/exercise7/plot_differences.py new file mode 100644 index 0000000..6653dbd --- /dev/null +++ b/exercise7/plot_differences.py @@ -0,0 +1,42 @@ +import numpy as np +import matplotlib.pyplot as plt +import sys + + +def main() -> None: + if len(sys.argv) != 3: + print("Usage: python plot_difference.py file1.dat file2.dat") + sys.exit(1) + + file1, file2 = sys.argv[1], sys.argv[2] + + M = 128 + N = 128 + + data1 = load_data(file1, M, N) + data2 = load_data(file2, M, N) + + # Compute difference + diff = data2 - data1 + + plt.figure(figsize=(6, 5)) + im = plt.imshow(diff, origin='lower', cmap='bwr', interpolation='nearest') + plt.colorbar(im, label="Difference") + plt.title(f"Difference between {file2} and {file1}") + plt.xlabel("Column") + plt.ylabel("Row") + plt.tight_layout() + + # Save figure to a file + plt.savefig("difference_plot.png", dpi=150) + print("Saved difference plot to difference_plot.png") + + +def load_data(filename:str, M:int, N:int) -> np.ndarray: + data = np.fromfile(filename, dtype=np.float64) + if data.size != M * N: + raise ValueError(f"File {filename} does not contain M*N={M*N} entries") + return data.reshape((M, N)) + +if __name__ == "__main__": + main() diff --git a/exercise7/plot_image.sh b/exercise7/plot_image.sh new file mode 100755 index 0000000..063333b --- /dev/null +++ b/exercise7/plot_image.sh @@ -0,0 +1,100 @@ +#! /usr/bin/env bash + +help() +{ + echo + echo "Plot 2D Wave Equation" + echo + echo "Syntax" + echo "--------------------------------------------------------" + echo "./plot_results.sh [-m|n|h] [data_folder] " + echo + echo "Option Description Arguments Default" + echo "--------------------------------------------------------" + echo "m y size Optional 128 " + echo "n x size Optional 128 " + echo "h Help None " + echo + echo "Example" + echo "--------------------------------------------------------" + echo "./plot_solution.sh -m 128 -n 128" + echo +} + +#----------------------------------------------------------------- +set -e + +M=128 +N=128 + +# Check if the data folder is provided +if [ $# -lt 1 ]; then + echo "Error: No data folder provided." + help + exit 1 +fi + +# Parse options and arguments +while getopts ":m:n:h" opt; do + case $opt in + m) + M=$OPTARG;; + n) + N=$OPTARG;; + h) + help + exit;; + \?) + echo "Invalid option" + help + exit;; + esac +done + +# Shift parsed options so that the remaining arguments start at $1 +shift $((OPTIND - 1)) + +# Ensure that the data folder is provided and exists +DATAFOLDER=./data +if [ ! -d "$DATAFOLDER" ]; then + echo "Error: Data folder $DATAFOLDER does not exist." + exit 1 +fi + +#----------------------------------------------------------------- +# Set up the size of the grid based on the options passed +SIZE_M=`echo $M | bc` +SIZE_N=`echo $N | bc` + +# Ensure the output directory exists +mkdir -p images + +# Loop through all .dat files in the data folder +for DATAFILE in "$DATAFOLDER"/*.dat; do + # Skip if no .dat files are found + if [ ! -f "$DATAFILE" ]; then + echo "No .dat files found in the folder." + exit 1 + fi + + # Create the corresponding output image file name + IMAGEFILE=`echo $DATAFILE | sed 's/dat$/png/' | sed 's/data/images/'` + + # Run the gnuplot command to create the plot in the background + ( + cat < +#include +#include +#include +#include +#include +#include +#include +#include + +// TASK: T1 +// Include the cooperative groups library +// BEGIN: T1 +; +// END: T1 + + +// Convert 'struct timeval' into seconds in double prec. floating point +#define WALLTIME(t) ((double)(t).tv_sec + 1e-6 * (double)(t).tv_usec) + + +// Option to change numerical precision +typedef int64_t int_t; +typedef double real_t; + +// TASK: T1b +// Variables needed for implementation +// BEGIN: T1b + +// Simulation parameters: size, step count, and how often to save the state +int_t + N = 128, + M = 128, + max_iteration = 1000000, + snapshot_freq = 1000; + +// Wave equation parameters, time step is derived from the space step +const real_t + c = 1.0, + dx = 1.0, + dy = 1.0; +real_t + dt; + +// Buffers for three time steps, indexed with 2 ghost points for the boundary +real_t + *buffers[3] = { NULL, NULL, NULL }; + +#define U_prv(i,j) buffers[0][((i)+1)*(N+2)+(j)+1] +#define U(i,j) buffers[1][((i)+1)*(N+2)+(j)+1] +#define U_nxt(i,j) buffers[2][((i)+1)*(N+2)+(j)+1] +// END: T1b + +#define cudaErrorCheck(ans) { gpuAssert((ans), __FILE__, __LINE__); } +inline void gpuAssert(cudaError_t code, const char *file, int line, bool abort=true) +{ + if (code != cudaSuccess) { + fprintf(stderr,"GPUassert: %s %s %d\n", cudaGetErrorString(code), file, line); + if (abort) exit(code); + } +} + + +// Rotate the time step buffers. +void move_buffer_window ( void ) +{ + real_t *temp = buffers[0]; + buffers[0] = buffers[1]; + buffers[1] = buffers[2]; + buffers[2] = temp; +} + + +// Save the present time step in a numbered file under 'data/' +void domain_save ( int_t step ) +{ + char filename[256]; + + // Ensure output directory exists (ignore error if it already exists) + if (mkdir("data", 0755) != 0 && errno != EEXIST) { + perror("mkdir data"); + exit(EXIT_FAILURE); + } + + snprintf(filename, sizeof(filename), "data/%05" PRId64 ".dat", step); + + FILE *out = fopen(filename, "wb"); + if (out == NULL) { + perror("fopen output file"); + fprintf(stderr, "Failed to open '%s' for writing.\n", filename); + exit(EXIT_FAILURE); + } + + for ( int_t i = 0; i < M; ++i ) { + size_t written = fwrite ( &U(i,0), sizeof(real_t), (size_t)N, out ); + if ( written != (size_t)N ) { + perror("fwrite"); + fclose(out); + exit(EXIT_FAILURE); + } + } + + if ( fclose(out) != 0 ) { + perror("fclose"); + exit(EXIT_FAILURE); + } +} + + +// TASK: T4 +// Get rid of all the memory allocations +void domain_finalize ( void ) +{ +// BEGIN: T4 + free ( buffers[0] ); + free ( buffers[1] ); + free ( buffers[2] ); +// END: T4 +} + + +// TASK: T6 +// Neumann (reflective) boundary condition +// BEGIN: T6 +void boundary_condition ( void ) +{ + for ( int_t i=0; i +#include +#include +#include +#include +#include +#include +#include +#include + + +// Convert 'struct timeval' into seconds in double prec. floating point +#define WALLTIME(t) ((double)(t).tv_sec + 1e-6 * (double)(t).tv_usec) + +// Option to change numerical precision +typedef int64_t int_t; +typedef double real_t; + +// Simulation parameters: size, step count, and how often to save the state +int_t + N = 128, + M = 128, + max_iteration = 1000000, + snapshot_freq = 1000; + +// Wave equation parameters, time step is derived from the space step +const real_t + c = 1.0, + dx = 1.0, + dy = 1.0; +real_t + dt; + +// Buffers for three time steps, indexed with 2 ghost points for the boundary +real_t + *buffers[3] = { NULL, NULL, NULL }; + +#define U_prv(i,j) buffers[0][((i)+1)*(N+2)+(j)+1] +#define U(i,j) buffers[1][((i)+1)*(N+2)+(j)+1] +#define U_nxt(i,j) buffers[2][((i)+1)*(N+2)+(j)+1] + + +// Rotate the time step buffers. +void move_buffer_window ( void ) +{ + real_t *temp = buffers[0]; + buffers[0] = buffers[1]; + buffers[1] = buffers[2]; + buffers[2] = temp; +} + + +// Save the present time step in a numbered file under 'data/' +void domain_save ( int_t step ) +{ + char filename[256]; + + // Ensure output directory exists (ignore error if it already exists) + if (mkdir("data", 0755) != 0 && errno != EEXIST) { + perror("mkdir data"); + exit(EXIT_FAILURE); + } + + snprintf(filename, sizeof(filename), "data/%05" PRId64 ".dat", step); + + FILE *out = fopen(filename, "wb"); + if (out == NULL) { + perror("fopen output file"); + fprintf(stderr, "Failed to open '%s' for writing.\n", filename); + exit(EXIT_FAILURE); + } + + for ( int_t i = 0; i < M; ++i ) { + size_t written = fwrite ( &U(i,0), sizeof(real_t), (size_t)N, out ); + if ( written != (size_t)N ) { + perror("fwrite"); + fclose(out); + exit(EXIT_FAILURE); + } + } + + if ( fclose(out) != 0 ) { + perror("fclose"); + exit(EXIT_FAILURE); + } +} + + +// Set up our three buffers, and fill two with an initial perturbation +void domain_initialize ( void ) +{ + buffers[0] = malloc ( (M+2)*(N+2)*sizeof(real_t) ); + buffers[1] = malloc ( (M+2)*(N+2)*sizeof(real_t) ); + buffers[2] = malloc ( (M+2)*(N+2)*sizeof(real_t) ); + + for ( int_t i=0; i