main.c 2.7 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115
  1. #include <stdlib.h>
  2. #include <stdio.h>
  3. #include <omp.h>
  4. void mv_orig(int m, int n, double M[m][n], double V[n], double W[m])
  5. {
  6. int i, j;
  7. for (i = 0; i < m; i++) {
  8. W[i] = 0.0;
  9. for (j = 0; j < n; j++) {
  10. W[i] += M[i][j] * V[j];
  11. }
  12. }
  13. }
  14. double reduce(int n, double M_i[n], double V[n], int actual);
  15. void mv(int m, int n, double M[m][n], double V[n], double W[m])
  16. {
  17. int i, j;
  18. for (i = 0; i < m; i++) {
  19. W[i] = reduce(n, M[i], V, 0);
  20. }
  21. }
  22. double reduce(int n, double* M_i, double* V, int actual) {
  23. // recursion tail, if indices coincide computation is trivial
  24. if (n == 1) {
  25. return M_i[0] * V[0];
  26. }
  27. // halve workload using a recursive call, with parallel for O(log n) time steps
  28. double sum[2];
  29. for (int i = 0; i < 2; i++) {
  30. sum[i] = reduce(
  31. n / 2 + (i * (n % 2)),
  32. &M_i[(n / 2) * i],
  33. &V[(n / 2) * i],
  34. actual + (n / 2) * i
  35. );
  36. }
  37. return sum[0] + sum[1];
  38. }
  39. void funky_sort(int n, int a[n]) {
  40. int b[n];
  41. #pragma omp parallel for shared (a, b)
  42. for (int i = 0; i < n; i++) {
  43. int j;
  44. int count = 0;
  45. for (j = 0; j < i; j++) {
  46. if (a[j] <= a[i]) count++;
  47. }
  48. j++;
  49. for (; j < n; j++) {
  50. if (a[j] < a[i]) count++;
  51. }
  52. b[count] = a[i];
  53. }
  54. #pragma omp parallel for shared (a, b)
  55. for (int i=0; i<n; i++)
  56. a[i] = b[i];
  57. }
  58. //int main(int argc, char *argv[]) {
  59. int main(void) {
  60. int max_threads = omp_get_max_threads();
  61. printf("max number of threads is %d\n", max_threads);
  62. #pragma omp parallel for
  63. for (int i = 0; i < max_threads * 2; i++) {
  64. printf("hello from thread %d\n", omp_get_thread_num());
  65. }
  66. #pragma omp parallel
  67. {
  68. printf("hello %d\n", omp_get_thread_num());
  69. printf("again %d\n", omp_get_thread_num());
  70. }
  71. double M[9][9] = {
  72. { 1, 2, 3, 4, 5, 6, 7, 8, 9 },
  73. { 1, 2, 3, 4, 5, 6, 7, 8, 9 },
  74. { 1, 2, 3, 4, 5, 6, 7, 8, 9 },
  75. { 1, 2, 3, 4, 5, 6, 7, 8, 9 },
  76. { 1, 2, 3, 4, 5, 6, 7, 8, 9 },
  77. { 1, 2, 3, 4, 5, 6, 7, 8, 9 },
  78. { 1, 2, 3, 4, 5, 6, 7, 8, 9 },
  79. { 1, 2, 3, 4, 5, 6, 7, 8, 9 },
  80. { 1, 2, 3, 4, 5, 6, 7, 8, 9 },
  81. };
  82. double V[9] = { 1, 2, 3, 4, 5, 6, 7, 8, 9 };
  83. double W[9];
  84. double W2[9];
  85. mv(9, 9, M, V, W);
  86. mv_orig(9, 9, M, V, W2);
  87. for (int i = 0; i < 9; i++) {
  88. printf("%d %f %f\n", i, W[i], W2[i]);
  89. }
  90. int funky_sort_test[9] = { 123, 1328, 132, 19325, 538, 13589, 248, 24187, 248 };
  91. funky_sort(9, funky_sort_test);
  92. for (int i = 0; i < 9; i++)
  93. printf("%d %d\n", i, funky_sort_test[i]);
  94. return 0;
  95. }