fmk - updating assignments for sessions 2 and 3

assignments/C-Session2/binaryFile/file3.c

@@ -0,0 +1,50 @@
+
+// program to read values from a file, each file a csv list of int and two double
+// written: fmk
+
+#include <stdio.h>
+#include <stdlib.h>
+
+int main(int argc, char **argv) {
+
+  if (argc != 3) {
+    fprintf(stdout, "ERROR correct usage appName inputFile outputBinaryFile\n");
+    return -1;
+  }
+
+  //
+  // read from ascii file
+  //
+
+  FILE *filePtr = fopen(argv[1],"r"); 
+
+  int i = 0;
+  float float1, float2;
+  int maxVectorSize = 100;
+  double *vector1 = (double *)malloc(maxVectorSize*sizeof(double));
+  double *vector2 = (double *)malloc(maxVectorSize*sizeof(double));  
+  int vectorSize = 0;
+
+  while (fscanf(filePtr,"%d, %f, %f\n", &i, &float1, &float2) != EOF) {
+    vector1[vectorSize] = float1;
+    vector2[vectorSize] = float2;
+    printf("%d, %f, %f\n",i, vector2[i], vector1[i]);
+    vectorSize++;
+
+    if (vectorSize == maxVectorSize) {
+      // some code needed here I think .. programming exercise
+    }
+  }
+
+  fclose(filePtr);
+
+  //
+  // write data to binary file
+  //
+
+  FILE *filePtrB = fopen(argv[2],"wb");
+
+  // some missing code to write vector1, followed by vector 2
+
+  fclose(filePtrB);  
+}

assignments/C-Session2/binaryFile/small.txt

@@ -0,0 +1,10 @@
+0, 0.153779, 0.560532
+1, 0.865013, 0.276724
+2, 0.895919, 0.704462
+3, 0.886472, 0.929641
+4, 0.469290, 0.350208
+5, 0.941637, 0.096535
+6, 0.457211, 0.346164
+7, 0.970019, 0.114938
+8, 0.769819, 0.341565
+9, 0.684224, 0.748597

assignments/C-Session2/matMul/CMakeLists.txt

@@ -0,0 +1,18 @@
+cmake_minimum_required (VERSION 2.6)
+set(CMAKE_CXX_STANDARD 17)
+set(CMAKE_C_STANDARD 17)
+
+project (matMUL)
+
+find_package(BLAS REQUIRED)
+
+include_directories(${PROJECT_SOURCE_DIR})
+
+add_executable(matMul matMul.c myDGEMM.c blasDGEMM.c)
+target_link_libraries(matMul m)
+target_link_libraries(matMul ${BLAS_LIBRARIES})
+
+add_executable(benchmark benchmark.cpp myDGEMM.c blasDGEMM.c)
+target_link_libraries(benchmark m)
+target_link_libraries(benchmark ${BLAS_LIBRARIES})
+

assignments/C-Session2/matMul/benchmark.cpp

@@ -0,0 +1,134 @@
+#include <algorithm>
+#include <chrono>
+#include <iomanip>
+#include <iostream>
+#include <random>
+#include <vector>
+
+#include <cmath> // For: fabs
+
+extern "C" void myDGEMM(int n, double *A, double *B, double *C);
+extern "C" void blasDGEMM(int n, double* A, double* B, double* C);
+
+void fill(double* p, int n) {
+    static std::random_device rd;
+    static std::default_random_engine gen(rd());
+    static std::uniform_real_distribution<> dis(-1.0, 1.0);
+    for (int i = 0; i < n; ++i)
+        p[i] = 2 * dis(gen) - 1;
+}
+
+/* The benchmarking program */
+int main(int argc, char** argv) {
+
+  bool test = false;
+  if (argc != 1)
+    test = true;
+
+  int ok = 0;
+
+  /*
+  int numTests = 96;
+  int testSizes[96]= {31,  32,  33,  63,  64,  65,  95,  96,  97,  127, 128, 129, 159, 160,  161,  191,
+		      192, 193, 223, 224, 225, 255, 256, 257, 287, 288, 289, 319, 320, 321,  351,  352,
+		      353, 383, 384, 385, 415, 416, 417, 447, 448, 449, 479, 480, 481, 511,  512,  513,
+		      543, 544, 545, 575, 576, 577, 607, 608, 609, 639, 640, 641, 671, 672,  673,  703,
+		      704, 705, 735, 736, 737, 767, 768, 769, 799, 800, 801, 831, 832, 833,  863,  864,
+		      865, 895, 896, 897, 927, 928, 929, 959, 960, 961, 991, 992, 993, 1023, 1024, 1025};
+  
+  int numTests = 26;
+  int testSizes[26] = {31,  32,  96,  97,  127, 128, 129, 191, 192, 229, 255, 256, 257,
+		       319, 320, 321, 417, 479, 480, 511, 512, 639, 640, 767, 768, 769};
+
+  */
+
+  int numTests = 13;
+  int testSizes[26] = {31,  32,  96,  97,  127, 128, 129, 191, 192, 229, 255, 256, 257};
+
+  double averageGFlopsBlas = 0;
+  double averageGFlopsMine = 0;
+
+  for (int i=0; i<numTests; i++) {
+
+    int n = testSizes[i];
+
+    double *A = (double *)malloc(n*n*sizeof(double));
+    double *B = (double *)malloc(n*n*sizeof(double));
+    double *C = (double *)malloc(n*n*sizeof(double));
+    double *C1 = (double *)malloc(n*n*sizeof(double));    
+
+    fill(A, n * n);
+    fill(B, n * n);
+    fill(C, n * n);
+
+    for (int i=0; i<n*n; i++)
+      C1[i]=C[i];
+
+    /* Measure performance (in Gflops/s). */
+
+    double seconds;
+    int n_iterations = 20;
+
+    /* Warm-up */
+    blasDGEMM(n, A, B, C);
+
+    /* Benchmark n_iterations runs of blasDGEMM */
+    auto start = std::chrono::steady_clock::now();
+    for (int it = 0; it < n_iterations; ++it) {
+      blasDGEMM(n, A, B, C);
+    }
+    auto end = std::chrono::steady_clock::now();
+    std::chrono::duration<double> diff = end - start;
+    seconds = diff.count();
+
+    /*  compute GFlop/s rate */
+    double GFlopsBlas = 2.e-9 * n_iterations * n * n * n / seconds;
+    averageGFlopsBlas += GFlopsBlas;
+
+
+    /* Warm-up */
+    myDGEMM(n, A, B, C1);
+
+    /* Benchmark n_iterations runs of myDGEMM */
+    start = std::chrono::steady_clock::now();
+    for (int it = 0; it < n_iterations; ++it) {
+      myDGEMM(n, A, B, C1);
+    }
+    end = std::chrono::steady_clock::now();
+    diff = end - start;
+    seconds = diff.count();
+
+    /*  compute GFlop/s rate */
+    double GFlopsMine = 2.e-9 * n_iterations * n * n * n / seconds;
+
+    averageGFlopsMine += GFlopsMine;    
+
+    // check they are the same .. take into account there will be differences
+    for (int j=0; j<n*n; j++) {
+      double diff = C1[j] - C[j];
+      double error = fabs(diff/C[j]);
+      if (error > 1e-10) {
+	ok = 1;
+	// printf("%d %d %.20g %.20g\n",i, j, C[j], C1[j]);
+	printf("%d\n", ok);
+	exit(0);
+      }
+    }
+
+    if (!test)
+      printf("%d %d %.8g %.8g \n", ok, n, GFlopsMine, GFlopsBlas);
+
+    free(A);
+    free(B);
+    free(C);
+    free(C1);    
+  }
+
+  /* Printing average percentage to screen */
+  if (!test)  
+    printf("Average GFLOP Me %.8g Blas: %.8g\n", averageGFlopsMine/numTests, averageGFlopsBlas/numTests);
+  else
+    printf("%d\n", ok);      
+
+  return 0;
+}

assignments/C-Session2/matMul/blasDGEMM.c

@@ -0,0 +1,26 @@
+void dgemm_(char* transa, char* transb, int* m, int* n, int* k,
+	    double* alpha, double* a, int* lda, double* b, int* ldb,
+	    double* beta, double* c, int* ldc);
+
+/*
+ * This routine performs a dgemm operation
+ *  C := C + A * B
+ * where A, B, and C are lda-by-lda matrices stored in column-major format.
+ * On exit, A and B maintain their input values.
+ * This function wraps a call to the BLAS-3 routine DGEMM,
+ */
+
+void blasDGEMM(int n, double* A, double* B, double* C) {
+
+  char transa = 'N';
+  char transb = 'N';
+  double alpha = 1.0;
+  double beta = 1.0;  
+  int lda = n;
+
+  dgemm_(&transa, &transb,
+	 &n, &n, &n,
+	 &alpha, A, &n,
+	 B, &n,
+	 &beta, C, &n);
+}

assignments/C-Session2/matMul/matMul.c

@@ -0,0 +1,64 @@
+#include <stdio.h>
+#include <stdlib.h>
+#include <math.h>
+
+extern void myDGEMM(int n, double *A, double *B, double *C);
+extern void blasDGEMM(int n, double* A, double* B, double* C);
+
+void fill(double* p, int n) {
+    for (int i = 0; i < n; ++i)
+        p[i] = (double)rand() / (double)RAND_MAX ;
+}
+
+/* The benchmarking program */
+int main(int argc, char** argv) {
+
+  if (argc != 2) {
+    printf("Correct usage: app matrixDimension?\n");
+    exit(0);
+  }
+
+  // get matrix size
+  int n = atoi(argv[1]);
+  n = fabs(n);
+  if (n == 0)
+    n = 10;
+
+  int result = 0;
+
+  double *A = 0;  // << SOME CODE HERE
+  double *B = 0;  // << SOME CODE HERE
+  double *C = 0;  // << SOME CODE HERE
+  double *C1 = 0; // << SOME CODE HERE
+
+  if (A == 0 || B == 0 || C == 0 || C1 == 0) {
+    printf("NO MMEORY ALLOCATED FOR ARRAYS\n");
+    exit(0);
+  }
+
+  fill(A, n * n);
+  fill(B, n * n);
+  fill(C, n * n);
+
+  for (int i=0; i<n*n; i++)
+    C1[i]=C[i];
+
+  blasDGEMM(n, A, B, C);
+
+  myDGEMM(n, A, B, C1);
+
+  // check they are the same .. take into account there will be differences due to roundoff
+  for (int j=0; j<n*n; j++) {
+    double diff = C1[j] - C[j];
+    double error = fabs(diff/C[j]);
+    if (error > 1e-10) {
+      result = 1;
+    }
+  }
+
+  // GOOD PRACTICE TO PUT 4 LINES of CODE HERE
+
+
+  printf("%d\n", result);
+  return 0;
+}

assignments/C-Session2/matMul/myDGEMM.c

@@ -0,0 +1,13 @@
+const char* dgemm_desc = "Naive, three-loop dgemm.";
+
+/*
+ * This routine performs a dgemm operation
+ *  C := C + A * B
+ * where A, B, and C are lda-by-lda matrices stored in column-major format.
+ * On exit, A and B maintain their input values.
+ *
+ * NOTE: Fortran storage: C(i,j) = C[i + j*n]
+ */
+void myDGEMM(int n, double* A, double* B, double* C) {
+  return; 
+}

assignments/C-Session3/stressTransformFile/ex2-3/CMakeLists.txt

@@ -0,0 +1,8 @@
+cmake_minimum_required (VERSION 2.6)
+set (CMAKE_C_STANDARD 99)
+
+project (Exercise2-3)
+
+include_directories(${PROJECT_SOURCE_DIR})
+add_executable(Exercise2-3 exercise2-3.c stresstransform.c)
+target_link_libraries(Exercise2-3 m)

assignments/C-Session3/stressTransformFile/ex2-3/exercise2-3.c

@@ -0,0 +1,36 @@
+#include <stdio.h>
+#include <stdlib.h>
+#include <math.h>    // need this for the constant pi = M_PI
+
+#include "stresstransform.h"
+
+
+int main(int argc, char **argv) {
+
+	// get dth from the first argument.  This is given in degrees!
+
+	// might be smart to set a default value, just in case the user
+	// forgets when calling this  program;)
+
+
+	// set the initial stress state
+	STRESS S0;
+
+	S0.sigx = 12.0;
+	S0.sigy = -5.5;
+	S0.tau  =  3.5;
+	S0.next =  0;
+
+
+	STRESS *results = &S0;
+
+	// MISSING CODE
+
+
+	STRESS *current = results;
+	while (current != NULL) {
+	  printf("sigx = %12.6f sigy' = %12.6f tau'  = %12.6f\n", current->sigx, current->sigy, current->tau);
+	  current = current->next;
+	}
+}
+

assignments/C-Session3/stressTransformFile/ex2-3/plotter.py

@@ -0,0 +1,17 @@
+import sys
+import matplotlib.pyplot as plt
+import pandas as pd
+
+try:
+    df = pd.read_csv('list.csv', header=None)
+except:
+    print("Could not find file: 'list.csv'")
+    sys.exit(-1)
+
+df.plot(x=1, y=3, label="Mohr's circle")
+
+plt.xlabel('$\sigma_x$')
+plt.ylabel('$\\tau_{xy}$')
+plt.grid(True)
+plt.axis('equal')
+plt.show()