From fc83484c1928ab71d3f28502318415800b5f5f3d Mon Sep 17 00:00:00 2001 From: "lisa.pizzo" Date: Fri, 9 Jan 2026 16:21:15 +0100 Subject: [PATCH] . --- Sheet7/Ex_9to13/accu_template/uv.txt | 711 ++++++++++++++++++ Sheet7/Ex_9to13/accu_template/vdop.cpp | 135 ++++ Sheet7/Ex_9to13/accu_template/vdop.h | 166 ++++ .../accu_template/visualize_results.m | 20 + 4 files changed, 1032 insertions(+) create mode 100644 Sheet7/Ex_9to13/accu_template/uv.txt create mode 100644 Sheet7/Ex_9to13/accu_template/vdop.cpp create mode 100644 Sheet7/Ex_9to13/accu_template/vdop.h create mode 100644 Sheet7/Ex_9to13/accu_template/visualize_results.m diff --git a/Sheet7/Ex_9to13/accu_template/uv.txt b/Sheet7/Ex_9to13/accu_template/uv.txt new file mode 100644 index 0000000..b975e97 --- /dev/null +++ b/Sheet7/Ex_9to13/accu_template/uv.txt @@ -0,0 +1,711 @@ +189 2 332 3 +0 0 +1 0 +1 1 +0 1 +0.1 0 +0.2 0 +0.3 0 +0.4 0 +0.5 0 +0.6 0 +0.7 0 +0.8 0 +0.9 0 +1 0.1 +1 0.2 +1 0.3 +1 0.4 +1 0.5 +1 0.6 +1 0.7 +1 0.8 +1 0.9 +0.9 1 +0.8 1 +0.7 1 +0.6 1 +0.5 1 +0.4 1 +0.3 1 +0.2 1 +0.1 1 +0 0.9 +0 0.8 +0 0.7 +0 0.6 +0 0.5 +0 0.4 +0 0.3 +0 0.2 +0 0.1 +0.479684 0.513096 +0.0481001 0.0438139 +0.965456 0.0484346 +0.0451033 0.961897 +0.951568 0.955138 +0.754828 0.354966 +0.360392 0.264508 +0.267427 0.637367 +0.656472 0.728531 +0.0895934 0.160441 +0.160541 0.910786 +0.911386 0.844141 +0.841625 0.0900995 +0.451559 0.770058 +0.557953 0.236427 +0.229318 0.453048 +0.687177 0.534082 +0.814618 0.688238 +0.258637 0.157666 +0.852304 0.258371 +0.158345 0.741149 +0.745635 0.841799 +0.734225 0.0924118 +0.266818 0.909939 +0.0900392 0.266867 +0.845905 0.489216 +0.340275 0.504784 +0.482979 0.335446 +0.532712 0.650471 +0.609668 0.857218 +0.390369 0.135169 +0.137448 0.608978 +0.711137 0.231036 +0.311061 0.779084 +0.220383 0.310202 +0.542551 0.0753476 +0.442762 0.912034 +0.0870955 0.443563 +0.718892 0.663674 +0.911061 0.643693 +0.910252 0.745257 +0.175659 0.104224 +0.825832 0.898216 +0.104444 0.824712 +0.899601 0.174158 +0.605353 0.422553 +0.396664 0.624666 +0.371373 0.387899 +0.824315 0.584766 +0.825103 0.373429 +0.524605 0.742235 +0.459558 0.258909 +0.261176 0.536511 +0.606947 0.570757 +0.662867 0.13484 +0.141276 0.353324 +0.353007 0.858336 +0.824377 0.795396 +0.180506 0.217545 +0.21789 0.819436 +0.7886 0.183202 +0.920144 0.449301 +0.445028 0.0729222 +0.552079 0.922008 +0.0761811 0.552606 +0.645993 0.929394 +0.0688885 0.645589 +0.353835 0.0678044 +0.33465 0.18669 +0.267936 0.0754627 +0.930563 0.255545 +0.188302 0.664576 +0.0758913 0.731951 +0.73311 0.923739 +0.67084 0.805934 +0.927247 0.0922677 +0.0756899 0.909417 +0.90786 0.925626 +0.0924729 0.0754185 +0.568222 0.332532 +0.700102 0.44222 +0.916963 0.544333 +0.738472 0.754101 +0.388486 0.711163 +0.287351 0.384911 +0.629828 0.654561 +0.900814 0.346453 +0.41132 0.456786 +0.465434 0.593118 +0.514559 0.831524 +0.165281 0.517019 +0.484085 0.160379 +0.765138 0.511037 +0.772882 0.273336 +0.241305 0.724548 +0.274741 0.24004 +0.469165 0.706886 +0.334436 0.584173 +0.409135 0.327415 +0.506546 0.420168 +0.291267 0.463861 +0.560193 0.497192 +0.405521 0.200759 +0.20339 0.593316 +0.593794 0.779935 +0.708431 0.167068 +0.286585 0.842293 +0.645259 0.197405 +0.157539 0.286478 +0.0634282 0.35253 +0.638015 0.0744131 +0.352258 0.936361 +0.429133 0.833293 +0.1657 0.43037 +0.7748 0.630926 +0.320743 0.124262 +0.125329 0.678845 +0.681069 0.871955 +0.856774 0.417367 +0.141959 0.0447961 +0.956726 0.141515 +0.858417 0.95601 +0.0449026 0.858609 +0.51729 0.290331 +0.791421 0.430793 +0.642597 0.273685 +0.663496 0.353015 +0.67685 0.601987 +0.308736 0.701675 +0.297005 0.306339 +0.632153 0.492996 +0.404572 0.546354 +0.436385 0.385368 +0.741249 0.587512 +0.208229 0.386035 +0.386067 0.791133 +0.821222 0.320652 +0.338641 0.652462 +0.450215 0.657713 +0.340125 0.439662 +0.344719 0.333947 +0.588452 0.145922 +0.58736 0.711465 +0.529103 0.564926 +0.709038 0.296367 +1 0.95 +0.05 1 +1 0.05 +0.95 1 +13 2 43 +85 15 111 +82 6 110 +7 8 108 +43 2 188 +42 5 119 +119 5 160 +6 7 110 +92 47 143 +9 10 76 +76 10 151 +111 16 127 +186 3 189 +17 18 102 +102 18 122 +12 13 53 +53 13 116 +43 14 116 +15 16 111 +40 1 42 +134 46 185 +21 22 52 +83 24 114 +25 26 106 +52 22 118 +45 23 118 +24 25 114 +126 49 183 +27 28 77 +77 28 152 +30 31 51 +51 31 117 +44 32 117 +84 33 113 +33 34 113 +34 35 107 +1 5 42 +50 40 119 +118 23 162 +4 32 44 +93 48 144 +36 37 78 +78 37 150 +39 40 50 +116 14 161 +117 32 163 +44 31 187 +45 22 186 +19 20 80 +128 41 172 +166 73 185 +16 17 127 +129 41 184 +8 9 103 +138 67 172 +35 36 105 +40 42 119 +170 47 181 +31 44 117 +169 48 178 +22 45 118 +26 27 104 +13 43 116 +115 49 123 +135 48 169 +29 30 64 +136 47 170 +11 12 63 +79 49 126 +38 39 65 +92 55 164 +94 57 168 +80 20 81 +168 57 174 +109 47 136 +39 50 65 +90 46 177 +12 53 63 +112 48 135 +30 51 64 +20 21 81 +79 58 123 +63 53 101 +139 68 173 +64 51 100 +91 54 137 +65 50 99 +93 56 141 +21 52 81 +121 46 165 +124 54 176 +125 56 175 +148 55 182 +128 88 173 +121 57 171 +93 67 138 +28 29 152 +115 62 158 +10 11 151 +109 59 156 +37 38 150 +112 61 157 +85 60 101 +120 55 166 +84 61 100 +97 74 176 +82 59 99 +96 75 175 +120 68 164 +11 63 151 +129 69 179 +29 64 152 +128 67 180 +38 65 150 +80 58 89 +91 69 183 +81 52 98 +89 58 155 +98 62 123 +58 80 81 +96 65 149 +108 71 156 +83 62 98 +106 70 158 +97 64 147 +107 72 157 +95 63 146 +134 60 177 +159 66 165 +120 86 140 +94 69 184 +124 87 179 +92 68 139 +125 88 180 +58 79 155 +18 19 122 +101 60 134 +89 66 122 +115 70 145 +130 91 145 +99 59 136 +109 71 143 +100 61 135 +112 72 144 +155 79 174 +140 86 142 +101 73 146 +146 73 148 +99 75 149 +131 78 154 +100 74 147 +130 77 153 +58 81 98 +52 83 98 +132 92 143 +50 82 99 +131 93 144 +51 84 100 +148 73 166 +53 85 101 +19 80 122 +102 66 159 +103 76 132 +9 76 103 +104 77 130 +27 77 104 +105 78 131 +36 78 105 +104 70 106 +26 104 106 +105 72 107 +35 105 107 +103 71 108 +8 103 108 +125 75 170 +110 108 156 +59 82 110 +7 108 110 +60 85 111 +46 90 165 +124 74 169 +113 107 157 +61 84 113 +34 107 113 +62 83 114 +25 106 114 +49 79 123 +114 106 158 +85 53 116 +14 15 161 +84 51 117 +32 33 163 +83 52 118 +23 24 162 +82 50 119 +5 6 160 +132 76 182 +121 86 167 +142 86 171 +133 121 165 +80 89 122 +66 102 122 +58 98 123 +62 115 123 +77 97 153 +138 87 178 +78 96 154 +139 88 181 +145 91 183 +69 94 126 +17 102 127 +60 111 127 +137 124 179 +141 125 180 +142 94 184 +129 87 172 +54 91 130 +70 104 130 +56 93 131 +72 105 131 +55 92 132 +71 103 132 +66 89 133 +57 121 133 +127 90 177 +73 101 134 +74 100 135 +61 112 135 +75 99 136 +59 109 136 +69 91 137 +54 124 137 +67 128 172 +48 93 138 +140 128 173 +47 92 139 +68 120 140 +41 128 140 +67 93 141 +56 125 141 +57 94 171 +41 140 142 +47 109 143 +71 132 143 +48 112 144 +72 131 144 +49 115 145 +70 130 145 +151 95 182 +63 101 146 +74 97 147 +64 100 147 +86 120 167 +95 146 148 +75 96 149 +65 99 149 +96 78 150 +65 96 150 +55 132 182 +63 95 151 +97 77 152 +64 97 152 +153 97 176 +54 130 153 +154 96 175 +56 131 154 +133 89 174 +57 133 174 +71 109 156 +59 110 156 +72 112 157 +61 113 157 +62 114 158 +70 115 158 +127 102 159 +90 127 159 +6 82 160 +82 119 160 +15 85 161 +85 116 161 +24 83 162 +83 118 162 +33 84 163 +84 117 163 +68 92 164 +55 120 164 +66 133 165 +90 159 165 +73 134 185 +55 148 166 +46 121 167 +120 166 167 +126 94 168 +79 126 168 +87 124 178 +74 135 169 +88 125 181 +75 136 170 +86 121 171 +94 142 171 +41 129 172 +87 138 172 +88 139 173 +68 140 173 +89 155 174 +79 168 174 +75 125 175 +56 154 175 +74 124 176 +54 153 176 +60 127 177 +46 134 177 +48 138 178 +124 169 178 +87 129 179 +69 137 179 +88 128 180 +67 141 180 +47 139 181 +125 170 181 +95 148 182 +76 151 182 +69 126 183 +49 145 183 +69 129 184 +41 142 184 +167 166 185 +46 167 185 +4 44 187 +14 43 188 +23 45 189 +45 186 189 +0 +0 +4 +0 +0 +0 +0 +0 +0 +0 +0 +0 +0 +1.41421 +2 +1.41421 +2.44929e-16 +-1.41421 +-2 +-1.41421 +-4.89859e-16 +1.41421 +1.62 +1.28 +0.98 +0.72 +0.5 +0.32 +0.18 +0.08 +0.02 +0 +-0 +-0 +-0 +-0 +0 +0 +0 +0 +-0.178548 +0.000780529 +0.346087 +0.0019439 +0.849854 +0.197346 +0.113567 +-0.0684595 +-0.229399 +0.00764267 +0.0197013 +0.282226 +0.460446 +-0.0475098 +0.298657 +-0.0212814 +-0.410328 +-0.510522 +0.0632297 +0.651413 +-0.0111809 +0.17926 +0.35781 +0.0541124 +0.00701453 +-0.461354 +-0.0848919 +0.113266 +-0.261777 +0.16147 +0.133057 +-0.0188449 +0.490766 +-0.0158238 +0.0314844 +0.164206 +0.151083 +-0.00254504 +-0.453517 +-0.781636 +-0.345363 +0.0225303 +0.475441 +0.00210396 +0.792672 +-0.0645703 +-0.154399 +0.0130884 +-0.674637 +0.141046 +-0.120619 +0.188989 +-0.0599069 +-0.358711 +0.383099 +0.00715405 +0.0551179 +-0.0245705 +0.0322734 +0.00721919 +0.616485 +-0.319708 +0.107329 +0.249376 +-0.00540613 +0.354771 +-0.00444466 +0.0635661 +0.111379 +0.0401008 +0.784844 +-0.0309942 +-0.00293374 +0.443894 +0.0209654 +0.569942 +0.00433926 +0.68755 +0.00477416 +0.163194 +-0.159566 +-0.761732 +-0.192361 +-0.0969628 +0.00976245 +-0.360816 +0.331295 +-0.0729785 +-0.216312 +0.0648866 +-0.0217185 +0.223083 +-0.448251 +0.500969 +-0.0325213 +0.071781 +-0.147004 +-0.110984 +0.0903415 +-0.0404745 +-0.040789 +-0.216954 +0.164444 +-0.0413103 +-0.055336 +0.485181 +0.0267824 +0.416273 +0.0193108 +0.00146543 +0.224589 +0.108906 +0.0476062 +-0.00648705 +-0.582694 +0.0852059 +-0.0127907 +0.248406 +-0.0998178 +0.00694483 +0.820448 +0.693334 +0.000895672 +0.203023 +-0.150011 +0.34569 +0.158791 +-0.458069 +-0.0665075 +0.0591938 +-0.266608 +-0.149364 +0.0218359 +-0.54681 +0.00474609 +-0.0103716 +0.393607 +-0.10508 +-0.182224 +-0.0354569 +0.0589187 +0.315509 +-0.221021 +-0.269395 +0.365485 +1.84776 +0.005 +0.765367 +1.805 diff --git a/Sheet7/Ex_9to13/accu_template/vdop.cpp b/Sheet7/Ex_9to13/accu_template/vdop.cpp new file mode 100644 index 0000000..902ef94 --- /dev/null +++ b/Sheet7/Ex_9to13/accu_template/vdop.cpp @@ -0,0 +1,135 @@ +#include "vdop.h" +#include // assert() +#include +#include +#include +using namespace std; + + +void vddiv(vector & x, vector const& y, + vector const& z) +{ + assert( x.size()==y.size() && y.size()==z.size() ); + size_t n = x.size(); + +#pragma omp parallel for + for (size_t k = 0; k < n; ++k) + { + x[k] = y[k] / z[k]; + } + return; +} + +//****************************************************************************** + +void vdaxpy(std::vector & x, std::vector const& y, + double alpha, std::vector const& z ) +{ + assert( x.size()==y.size() && y.size()==z.size() ); + size_t n = x.size(); + +#pragma omp parallel for + for (size_t k = 0; k < n; ++k) + { + x[k] = y[k] + alpha * z[k]; + } + return; +} +//****************************************************************************** + +double dscapr(std::vector const& x, std::vector const& y) +{ + assert( x.size()==y.size()); + size_t n = x.size(); + + double s = 0.0; +//#pragma omp parallel for reduction(+:s) + for (size_t k = 0; k < n; ++k) + { + s += x[k] * y[k]; + } + + return s; +} + +//****************************************************************************** +//void DebugVector(vector const &v) +//{ + //cout << "\nVector (nnode = " << v.size() << ")\n"; + //for (size_t j = 0; j < v.size(); ++j) + //{ + //cout.setf(ios::right, ios::adjustfield); + //cout << v[j] << " "; + //} + //cout << endl; + + //return; +//} +//****************************************************************************** +bool CompareVectors(std::vector const& x, int const n, double const y[], double const eps) +{ + bool bn = (static_cast(x.size())==n); + if (!bn) + { + cout << "######### Error: " << "number of elements" << endl; + } + //bool bv = equal(x.cbegin(),x.cend(),y); + bool bv = equal(x.cbegin(),x.cend(),y, + [eps](double a, double b) -> bool + { return std::abs(a-b)(x.size())==n); + cout << "######### Error: " << "values" << endl; + } + return bn && bv; +} + +//****************************************************************************** +double par_scalar(vector const &x, vector const &y, MPI_Comm const& icomm) +{ + const double s = dscapr(x,y); + double sg; + MPI_Allreduce(&s,&sg,1,MPI_DOUBLE,MPI_SUM,icomm); + + return(sg); +} + +//****************************************************************************** +void ExchangeAll(vector const &xin, vector &yout, MPI_Comm const &icomm) +{ + int myrank, numprocs,ierr(-1); + MPI_Comm_rank(icomm, &myrank); // my MPI-rank + MPI_Comm_size(icomm, &numprocs); + int const N=xin.size(); + int const sendcount = N/numprocs; // equal sized junks + assert(sendcount*numprocs==N); // really all junk sized? + assert(xin.size()==yout.size()); + + auto sendbuf = xin.data(); + auto recvbuf = yout.data(); + ierr = MPI_Alltoall(sendbuf, sendcount, MPI_DOUBLE, + recvbuf, sendcount, MPI_DOUBLE, icomm); + assert(0==ierr); + + return; +} + +//****************************************************************************** +void ExchangeAllInPlace(vector &xin, MPI_Comm const &icomm) +{ + int myrank, numprocs,ierr(-1); + MPI_Comm_rank(icomm, &myrank); // my MPI-rank + MPI_Comm_size(icomm, &numprocs); + int const N=xin.size(); + int const sendcount = N/numprocs; // equal sized junks + assert(sendcount*numprocs==N); // really all junk sized? + + auto sendbuf = xin.data(); + ierr = MPI_Alltoall(MPI_IN_PLACE, sendcount, MPI_DOUBLE, + sendbuf, sendcount, MPI_DOUBLE, icomm); + assert(0==ierr); + + return; +} diff --git a/Sheet7/Ex_9to13/accu_template/vdop.h b/Sheet7/Ex_9to13/accu_template/vdop.h new file mode 100644 index 0000000..a609aee --- /dev/null +++ b/Sheet7/Ex_9to13/accu_template/vdop.h @@ -0,0 +1,166 @@ +#ifndef VDOP_FILE +#define VDOP_FILE +#include +#include // MPI +#include +#include + +/** @brief Element-wise vector divison x_k = y_k/z_k. + * + * @param[out] x target vector + * @param[in] y source vector + * @param[in] z source vector + * + */ +void vddiv(std::vector &x, std::vector const &y, + std::vector const &z); + +/** @brief Element-wise daxpy operation x(k) = y(k) + alpha*z(k). + * + * @param[out] x target vector + * @param[in] y source vector + * @param[in] alpha scalar + * @param[in] z source vector + * + */ +void vdaxpy(std::vector &x, std::vector const &y, + double alpha, std::vector const &z ); + + +/** @brief Calculates the Euclidean inner product of two vectors. + * + * @param[in] x vector + * @param[in] y vector + * @return Euclidean inner product @f$\langle x,y \rangle@f$ + * + */ +double dscapr(std::vector const &x, std::vector const &y); + + +inline +double L2_scapr(std::vector const &x, std::vector const &y) +{ + return dscapr(x, y) / x.size(); +} + + +/** Parallel inner product + @param[in] x vector + @param[in] y vector + @param[in] icomm MPI communicator + @return resulting Euclidian inner product +*/ +double par_scalar(std::vector const &x, std::vector const &y, + MPI_Comm const& icomm=MPI_COMM_WORLD); + + + +/* ReadId : Input and broadcast of an integer */ +inline +int ReadIn(std::string const &ss = std::string(), MPI_Comm const &icomm = MPI_COMM_WORLD) +{ + MPI_Barrier(icomm); + int myrank; /* my rank number */ + MPI_Comm_rank(icomm, &myrank); + int id; + + if (myrank == 0) { + std::cout << "\n\n " << ss << " : Which process do you want to debug ? \n"; + std::cin >> id; + } + MPI_Bcast(&id, 1, MPI_INT, 0, icomm); + + return id; +} + +/** + * Print entries of a vector to standard output. + * + * @param[in] v vector values + * @param[in] ss string containing the vector name + * @param[in] icomm communicator group for MPI + * +*/ +//void DebugVector(std::vector const &v); +template +void DebugVector(std::vector const &v, std::string const &ss = std::string(), MPI_Comm const &icomm = MPI_COMM_WORLD) +{ + MPI_Barrier(icomm); + int numprocs; /* # processes */ + MPI_Comm_size(icomm, &numprocs); + int myrank; /* my rank number */ + MPI_Comm_rank(icomm, &myrank); + + int readid = ReadIn(ss); /* Read readid */ + + while ( (0 <= readid) && (readid < numprocs) ) { + if (myrank == readid) { + std::cout << "\n\n process " << readid; + std::cout << "\n .... " << ss << " (nnode = " << v.size() << ")\n"; + for (size_t j = 0; j < v.size(); ++j) { + std::cout.setf(std::ios::right, std::ios::adjustfield); + std::cout << v[j] << " "; + } + std::cout << std::endl; + fflush(stdout); + } + + readid = ReadIn(ss, icomm); /* Read readid */ + } + MPI_Barrier(icomm); + return; +} + +/** @brief Compares an STL vector with POD vector. + * + * The accuracy criteria @f$ |x_k-y_k| < \varepsilon \left({1+0.5(|x_k|+|y_k|)}\right) @f$ + * follows the book by + * Stoyan/Baran, p.8. + * + * @param[in] x STL vector + * @param[in] n length of POD vector + * @param[in] y POD vector + * @param[in] eps relative accuracy criteria (default := 0.0). + * @return true iff pairwise vector elements are relatively close to each other. + * + */ +bool CompareVectors(std::vector const &x, int n, double const y[], double const eps = 0.0); + + +/** Output operator for vector + * @param[in,out] s output stream, e.g. @p cout + * @param[in] v vector + * + * @return output stream +*/ +template +std::ostream &operator<<(std::ostream &s, std::vector const &v) +{ + for (auto vp : v) { + s << vp << " "; + } + return s; +} + +/** Exchanges equal size partions of vector @p xin with all MPI processes. + * The received data are return in vector @p yout . + * + * @param[in] xin input vector + * @param[out] yout output vector + * @param[in] icomm MPI communicator + * +*/ +void ExchangeAll(std::vector const &xin, std::vector &yout, MPI_Comm const &icomm = MPI_COMM_WORLD); + +/** Exchanges equal size partions of vector @p xin with all MPI processes. + * The received data are return in vector @p xin . + * + * @param[in,out] xin input/output vector + * @param[in] icomm MPI communicator + * +*/ +void ExchangeAllInPlace(std::vector &xin, MPI_Comm const &icomm = MPI_COMM_WORLD); + + + +#endif diff --git a/Sheet7/Ex_9to13/accu_template/visualize_results.m b/Sheet7/Ex_9to13/accu_template/visualize_results.m new file mode 100644 index 0000000..e40da61 --- /dev/null +++ b/Sheet7/Ex_9to13/accu_template/visualize_results.m @@ -0,0 +1,20 @@ +%% Visualize results +% +% flatpak run org.octave.Octave +% or +% octave --no-window-system --no-gui -qf +% +% or +% matlab -nosplash < + +clear all +clc + +%% +fname = 'uv.txt'; + +[xc,ia,v] = ascii_read_meshvector(fname); + +h = trisurf(ia, xc(:,1), xc(:,2), v); + +waitfor(h) % wait for closing the figure \ No newline at end of file