mirror of
https://github.com/fverdugo/XM_40017.git
synced 2026-11-15 01:47:17 +01:00
5 files changed
+11059
-1653
No files matched your search
+2
-2
@@ -124,8 +124,8 @@ makedocs(;
|
||||
"MPI (point-to-point)" => "julia_mpi.md",
|
||||
"MPI (collectives)" => "mpi_collectives.md",
|
||||
"Jacobi method" => "jacobi_method.md",
|
||||
#"All pairs of shortest paths" => "asp.md",
|
||||
#"Gaussian elimination" => "LEQ.md",
|
||||
"All pairs of shortest paths" => "asp.md",
|
||||
"Gaussian elimination" => "LEQ.md",
|
||||
#"Traveling salesperson problem" => "tsp.md",
|
||||
#"Partial differential equations" => "pdes.md",
|
||||
],
|
||||
|
||||
+295
-160
@@ -22,7 +22,7 @@
|
||||
"In this notebook, we will learn\n",
|
||||
"\n",
|
||||
"- How to parallelize Gaussian elimination\n",
|
||||
"- How to fix static load imbalance"
|
||||
"- What is load imbalance, and how to fix it in this algorithm"
|
||||
]
|
||||
},
|
||||
{
|
||||
@@ -45,24 +45,69 @@
|
||||
"using Printf\n",
|
||||
"function answer_checker(answer,solution)\n",
|
||||
" if answer == solution\n",
|
||||
" \"🥳 Well done! \"\n",
|
||||
" \"🥳 Well done!\"\n",
|
||||
" else\n",
|
||||
" \"It's not correct. Keep trying! 💪\"\n",
|
||||
" end |> println\n",
|
||||
"end\n",
|
||||
"function ge_par_why()\n",
|
||||
" msg = \"The outer loop of the algorithm is not parallelizable, since the iterations depend on the results of the previous iterations. However, we can extract parallelism from the inner loops.\"\n",
|
||||
" println(msg)\n",
|
||||
"end\n",
|
||||
"ge_par_check(answer) = answer_checker(answer, \"a\")\n",
|
||||
"ge_dep_check(answer) = answer_checker(answer, \"b\")"
|
||||
"ge_dep_check(answer) = answer_checker(answer, \"b\")\n",
|
||||
"function ge_lb_answer()\n",
|
||||
" msg = \"It is a form of static load balancing. We know in advance the load distribution and the partition strategy does not depend on the actual values of the input matrix\"\n",
|
||||
" println(msg)\n",
|
||||
"end\n",
|
||||
"println(\"🥳 Well done!\")"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"id": "8dcee319",
|
||||
"id": "8ad95aa5",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"## Gaussian elimination\n",
|
||||
"\n",
|
||||
"\n",
|
||||
"[Gaussian elimination](https://en.wikipedia.org/wiki/Gaussian_elimination) is a method to solve systems of linear equations, e.g.\n",
|
||||
"[Gaussian elimination](https://en.wikipedia.org/wiki/Gaussian_elimination) is a method to solve systems of linear equations, like this one:\n",
|
||||
"\n",
|
||||
"$$\n",
|
||||
"\\left[\n",
|
||||
"\\begin{matrix}\n",
|
||||
"1 & 3 & 1 \\\\\n",
|
||||
"1 & 2 & -1 \\\\\n",
|
||||
"3 & 11 & 5 \\\\\n",
|
||||
"\\end{matrix}\n",
|
||||
"\\right]\n",
|
||||
"\\left[\n",
|
||||
"\\begin{matrix}\n",
|
||||
"x \\\\\n",
|
||||
"y \\\\\n",
|
||||
"z \\\\\n",
|
||||
"\\end{matrix}\n",
|
||||
"\\right]\n",
|
||||
"=\n",
|
||||
"\\left[\n",
|
||||
"\\begin{matrix}\n",
|
||||
"9 \\\\\n",
|
||||
"1 \\\\\n",
|
||||
"35 \\\\\n",
|
||||
"\\end{matrix}\n",
|
||||
"\\right].\n",
|
||||
"$$\n",
|
||||
"\n",
|
||||
"This is just a small example with three unknowns, but practical applications need to solve linear equations with large number of unknowns. Parallel processing is needed in these cases.\n",
|
||||
"\n",
|
||||
"### Problem statement\n",
|
||||
"\n",
|
||||
"Let us consider a system of linear equations written in matrix form $Ax=b$, where $A$ is a nonsingular square matrix, and $x$ and $b$ are vectors. $A$ and $b$ are given, and $x$ is unknown. The goal of Gaussian elimination is to transform the system $Ax=b$, into a new system $Ux=c$ such that\n",
|
||||
"- both system have the same solution vector $x$,\n",
|
||||
"- the matrix $U$ of the new system is *upper triangular* with unit diagonal, namely $U_{ii} = 1$ and $U_{ij} = 0$ for $i>j$.\n",
|
||||
"\n",
|
||||
"\n",
|
||||
"For the particular system shown above, the transformed one is the following:\n",
|
||||
"\n",
|
||||
"$$\n",
|
||||
"\\left[\n",
|
||||
@@ -87,51 +132,97 @@
|
||||
"35 \\\\\n",
|
||||
"\\end{matrix}\n",
|
||||
"\\right]\n",
|
||||
"\\longrightarrow\n",
|
||||
"\\left[\n",
|
||||
"\\begin{matrix}\n",
|
||||
"1 & 3 & 1 \\\\\n",
|
||||
"0 & 1 & 2 \\\\\n",
|
||||
"0 & 0 & 1 \\\\\n",
|
||||
"\\end{matrix}\n",
|
||||
"\\right]\n",
|
||||
"\\left[\n",
|
||||
"\\begin{matrix}\n",
|
||||
"x \\\\\n",
|
||||
"y \\\\\n",
|
||||
"z \\\\\n",
|
||||
"\\end{matrix}\n",
|
||||
"\\right]\n",
|
||||
"=\n",
|
||||
"\\left[\n",
|
||||
"\\begin{matrix}\n",
|
||||
"9 \\\\\n",
|
||||
"8 \\\\\n",
|
||||
"4 \\\\\n",
|
||||
"\\end{matrix}\n",
|
||||
"\\right]\n",
|
||||
"$$\n",
|
||||
"\n",
|
||||
"The steps of the Gaussian elimination will transform the system into an upper triangular matrix. The system of linear equations can now easily be solved by backward substitution. \n",
|
||||
"The most challenging part of solving a system of linear equations is to transform it to upper triangular form. Afterwards, the solution vector can be obtained easily with a backward substitution. For this reason, we will study here the triangulation step only."
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"id": "980cb36f",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"\n",
|
||||
"\n",
|
||||
"### Augmented system matrix\n",
|
||||
"\n",
|
||||
"In practice, vector $b$ is added as an additional column to A forming the so-called *augmented* matrix $A^* = [A b]$. The augmented matrix in the example above is\n",
|
||||
"\n",
|
||||
"$$\n",
|
||||
"\\left[\n",
|
||||
"\\begin{matrix}\n",
|
||||
"1 & 3 & 1 & 9 \\\\\n",
|
||||
"1 & 3 & 1 \\\\\n",
|
||||
"1 & 2 & -1 \\\\\n",
|
||||
"3 & 11 & 5 \\\\\n",
|
||||
"\\end{matrix}\n",
|
||||
"\\right]\n",
|
||||
"\\left[\n",
|
||||
"\\begin{matrix}\n",
|
||||
"x \\\\\n",
|
||||
"y \\\\\n",
|
||||
"z \\\\\n",
|
||||
"\\end{matrix}\n",
|
||||
"\\right]\n",
|
||||
"=\n",
|
||||
"\\left[\n",
|
||||
"\\begin{matrix}\n",
|
||||
"9 \\\\\n",
|
||||
"1 \\\\\n",
|
||||
"35 \\\\\n",
|
||||
"\\end{matrix}\n",
|
||||
"\\right]\\longrightarrow\n",
|
||||
"A^*=\n",
|
||||
"\\left[\n",
|
||||
"\\begin{matrix}\n",
|
||||
"1 & 3 & 1 & 9 \\\\\n",
|
||||
"1 & 2 & -1 & 1 \\\\\n",
|
||||
"3 & 11 & 5 & 35 \\\\\n",
|
||||
"3 & 11 & 5 & 35\\\\\n",
|
||||
"\\end{matrix}\n",
|
||||
"\\right]\n",
|
||||
"\\rightarrow\n",
|
||||
"\\right].\n",
|
||||
"$$\n",
|
||||
"\n",
|
||||
"With this new notation, the goal of Gaussian elimination is to find the augmented matrix containing $U$ and $c$, namely $U^*= [U c]$, given the augmented matrix $A^* = [A b]$. These are $A^*$ and $U^*$ in our example:\n",
|
||||
"\n",
|
||||
"$$\n",
|
||||
"A^*=\n",
|
||||
"\\left[\n",
|
||||
"\\begin{matrix}\n",
|
||||
"1 & 3 & 1 & 9 \\\\\n",
|
||||
"1 & 2 & -1 & 1 \\\\\n",
|
||||
"3 & 11 & 5 & 35\\\\\n",
|
||||
"\\end{matrix}\n",
|
||||
"\\right]\\longrightarrow\n",
|
||||
"U^*=\n",
|
||||
"\\left[\n",
|
||||
"\\begin{matrix}\n",
|
||||
"1 & 3 & 1 & 9 \\\\\n",
|
||||
"0 & -1 & -2 & -8 \\\\\n",
|
||||
"0 & 2 & 2 & 8 \\\\\n",
|
||||
"0 & 1 & 2 & 8\\\\\n",
|
||||
"0 & 0 & 1 & 4\\\\\n",
|
||||
"\\end{matrix}\n",
|
||||
"\\right]\n",
|
||||
"\\rightarrow\n",
|
||||
"\\left[\n",
|
||||
"\\begin{matrix}\n",
|
||||
"1 & 3 & 1 & 9 \\\\\n",
|
||||
"0 & 1 & 2 & 8 \\\\\n",
|
||||
"0 & 2 & 2 & 8 \\\\\n",
|
||||
"\\end{matrix}\n",
|
||||
"\\right]\n",
|
||||
"\\rightarrow\n",
|
||||
"\\left[\n",
|
||||
"\\begin{matrix}\n",
|
||||
"1 & 3 & 1 & 9 \\\\\n",
|
||||
"0 & 1 & 2 & 8 \\\\\n",
|
||||
"0 & 0 & -2 & -8 \\\\\n",
|
||||
"\\end{matrix}\n",
|
||||
"\\right]\n",
|
||||
"\\rightarrow\n",
|
||||
"\\left[\n",
|
||||
"\\begin{matrix}\n",
|
||||
"1 & 3 & 1 & 9 \\\\\n",
|
||||
"0 & 1 & 2 & 8 \\\\\n",
|
||||
"0 & 0 & 1 & 4 \\\\\n",
|
||||
"\\end{matrix}\n",
|
||||
"\\right]\n",
|
||||
"\\right].\n",
|
||||
"$$\n",
|
||||
"\n"
|
||||
]
|
||||
@@ -142,9 +233,13 @@
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"### Serial implementation\n",
|
||||
"The following algorithm computes the Gaussian elimination on a matrix which represents a system of linear equations.\n",
|
||||
"- The first inner loop in line 4 divides the current row by the value of the diagonal entry, thus transforming the diagonal to contain only ones. \n",
|
||||
"- The second inner loop beginning in line 8 substracts the rows from one another such that all entries below the diagonal become zero. "
|
||||
"\n",
|
||||
"\n",
|
||||
"The following algorithm computes the Gaussian elimination on a given augmented matrix `B`, representing a system of linear equations. The result is given by overwriting `B`, avoiding the allocation of an additional matrix.\n",
|
||||
"\n",
|
||||
"- The outer loop is a loop over rows.\n",
|
||||
"- The first inner loop in line 4 divides the current row by the value of the diagonal entry, thus transforming the diagonal to contain only ones. Note that we skip the first entries in the row, as we know that these values are zero at this point. The cells updated in this loop at iteration $k$ are depicted in red in the figure below.\n",
|
||||
"- The second inner loop beginning in line 8 subtracts the rows from one another such that all entries below the diagonal become zero. The entries updated in this loop are depicted in blue in the figure below."
|
||||
]
|
||||
},
|
||||
{
|
||||
@@ -172,13 +267,28 @@
|
||||
"end"
|
||||
]
|
||||
},
|
||||
{
|
||||
"attachments": {
|
||||
"g27851.png": {
|
||||
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAzAAAAJsCAYAAAAiM0YjAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAAewgAAHsIBbtB1PgAAABl0RVh0U29mdHdhcmUAd3d3Lmlua3NjYXBlLm9yZ5vuPBoAACAASURBVHic7N15nBx1nf/xV/X0zOS+LxIIwYRDDkFuFBFUUFbE+wQPFBRXXQVcL1xss3ggrgfeeKKuq7grIh4/BeUSOQxBbgEFQgIkAZJJyDVX1++PmoGeTnfPVKe7q2rm9Xw8+pGu6u936t09k5n+9Lfq+w2QJI12eeCsku2/A5cO3O8EXgG8GdgD2BkoAvcDvwMuAB6u8nWfBZwKPB94xsDXWgvcBPwQ+GUjn4QkSZKksWE8EJbcfjqwf1/g7rLHym8bgdeUfb1xwNeB/mH6/hqY0KTnJEmSJGmUqlTAHAxspnYBMnjrA44b+FrjgD+OsF8I/BYImvz8JEmSJI0i5QXMH4AHBu5vA75NNMpyBHAM8FFgdVmfFQNf57sl+24hOjXtRcBhwOuITk0rL2JObv5TlCRJkjRalBcwg7d7gd2r9Jk78Hhp+x8P/FsE/h3IVen74bJ+tzTiSUiSJEkaGyoVMOuAXYbpd3yFfiGwdJh+AXBzWZ8968wuSZIkaYypVMB8YAT9cmx/KtlKoH0EfT9Q1u9NsVNLklRBteF/SdLo1Q/8YATtisDysn0/AnpH0PevZduOwEiSGsICRpLGnjuArhG2XVW2fV2d/aaNsJ8kSTVZwEjS2FNtYcpKNpVtPzLCfpvLtifHOKYkSVVZwEjS2FNelNRSrLNveT/XgpEkNYQFjCSNPWHSASRJqpcFjCRJkqTMsICRJEmSlBkWMJIkSZIywwJGkiRJUmZYwEiSJEnKDAsYSZIkSZlhASNJkiQpMyxgJEmSJGWGBYwkSZKkzLCAkSRJkpQZ+aQDSJKaLgTWl2xvjtF3a1nf/jqPuSXGMSVJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkqTWGA+cDFwGtCecRZIkxTcTuBJ4M9HfdUkadQLgKOASYDWwDfjfRBNJkqQdcS3R3/M1wKXA84FcookkqQGWAF8AVgKPAeHAbQ2wc4K5JEnSjjmQoX/b1wGPAl8B9kgwlyTFNh14L3An0WhLH0//chu83ZFYOkmS1CgPsf3f+CLwBHAv8AFgRmLpNOoFSQdQprUDxwPnEH3qMo7a17f8GfhRC3JJkqTm+QDwzBqP9wFbgduBzwG/BXpbkEtjhAWM6nU48AeigmVcwlkkSVI69RBdM3MocE/CWTRKWMCoXjngGOBM4LnAJKBtmD73E30aI0mSsuswYN4wbXqB9cBvgAuAvzU7lCTFMQ54PfAXYC3RL63yc2ND4LqkAkqSpIa5j8p/57cSXdB/PfA6XDZBUkbsDCwFHmToLCUh0cX9ExNLJkmSdtR8or/ng3/b+4lmGb0DeDcwNblokrTjDia6aH81sIHoHNgPJ5pIkiTtiG8SFS2PEy2XcD7wjEQTSVITtAMvB64mmmLRIWVJkrJnJtF1LT8nuv7Va6kljQmdSQeQJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSpOFMB9aV3L6WbBxJktIhn3QASVJFAVERM2hCUkEkSUqTXNIBJEmSJGmkLGAkSZIkZYYFjCRJkqTMsICRJEmSlBkWMJIkSZIywwJGkiRJUmY4jbIkjT47A3PL9t0BdCeQRZIkSdIYMAMIS27fH2G/Y4ENJf2KQKEJ+SRJkiTpKfUUMG8Dekr6dANvalI+SZIkSXpKnAImIBplKW2/Hji6qQklSZIkacBIC5g8cGFZ24eB/VuQUZIkSZKAkRUwk4DflrW7HdilRRklSZIkCRi+gNkJWF7W5gpgagszSpIkSRJQu4DZB1hR9vhFQHuLM0qSJEkSUL2AOYboAv3Sx75MdCG/JEmSJCWiUgHzWmBryb4+4N1JBZQkSZKkQeUFzKNEi1IObm8CXppYOkmSJEkqUV7AlN8uBzoSSydJkiRJJcoLmOuIRmFK9/0GGJdUQEmSJEkaVOkamN2A+8v2XwVMTiaiJEmSJEWqzUK2ELiv7LFrgSkJZJQkSZIkoPY6MPOA28seXwbMbHFGSZIkSQJqFzAAc4C/lbW5E9iphRklSZIkCRi+gAGYDtxQ1u5uYEGLMkqSJEkSMLICBmAq8OeytvcTXfAvSZIkSS0x0gIGYALRujCl7VcAS5qcUZIkSZKAeAUMQCfwy7I+jwL7NjGjJEmSJAHxCxiADuB/y/qtAfZvUkZJkiRJAuorYADagB+U9V0HHNb4iJIktV6QdABJUkU5YFHJ9pPAYyPs2wbsWrZvC7B6x2NJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJktKrHZiRdAhJkiRJ6ZRPOkCZdwGvAF6UdBCNVLgbsEfSKSRJUqr9GYLNSYfQ6BAkHaDMXcBcYAmwPuEsGpHww8Bnk04hSZJSbR8I7ko6hEaHXNIBSiwBZgETgZMTzqJhheMhnAK/PCbpJJIkKe0+enz0viFsTzqJsi9NIzBfAM4YuP8A8IwEs2hY4aXAiUmnkCRJmXIaBN9JOoSyLS0jMAHw2pLtqXhdhSRJkqQyabmI/3nA+JLt6UQX9J+VTBzVYTnwlqRDbO/A38DCi+GXFyWdpLIXfRQ27gs3nZR0ksq+sRC+/gu4/eCkk1T3rBvh9W+Fs/+edJLKDvsBTHgQriwkHKSK170O7nsH3PLipJNUt9+tcPv+Saeo7oAr4Jlfg/+5JOkklT3/09A9G244LekklZ29P1z2Lbjt8KSTVLffzfAfx8Pr1iadpLJDfgrTl8EfPp90kjKTgeuTDiE1yyVAWHZbRbpOcdMQ4aUQhiW3a5NOVMW9wJlJh6jhAtL9y30JkPZZYzYCByQdooargAuTDlHD6USn7aZZb9IBhvEQcErSIWq4CLg86RA1HA50JR1iGFuBBUmHqOEm4PykQ2wvnFb2XiGE8NSkUyn70nAK2TjgiAr7JwBHtjiLJEmSpBRLQwHzGqKZx8pNBz7Q4iySJEmSUiwNBcz7gUlVHnsu0UiMJEmSJCVewMwFdqnx+CTglS3KIkmSJCnlkp6F7FSiU8WqmUg0QvPfrYmjhM0F9gbaBrbvAVYmF0eSJElpk3QBcxrRKNAaopmOdgY6gPuJipdJwO7APGB1QhnVPM8FjgEOBg4i+v6Xeh/w1VaHkiRJUnoleQrZvkSFyReB/YHFwGMDj70W2BV4L9APpHSNDO2gTwD/Cbyc7YsXSZIkaTtJjsDcD+wGPFrl8W7gB8DFwNQWZVIyeoDbgJuJCteXJBtHkiRJaZXkCMwWqhcv9bRT9nye6NSxScAhRAvq3ZBoIkmSJKVa0tfAaGz7Q9IBJEmSlC1JT6MsSZIkSSNmASNJkiQpMyxgJEmSJGWGBYwkSZKkzLCAkSRJkpQZFjCSJEmSMsMCRpIkSVJmWMBIkiRJygwLGEmSJEmZYQEjSZIkKTMsYCRJkiRlhgWMJEmSpMywgJEkSZKUGRYwkiRJkjLDAkZJejkQlt0KJY9/pcLjR7Y2oiRJktLEAkaSJElSZuSTDqAx7S7gIzH7PNiEHJIkScoICxgl6T7gvKRDSJIkKTs8hUySJElSZljASJIkScoMCxhJkiRJmWEBI0mSJCkzLGAkSZIkZYYFjCRJkqTMsICRJEmSlBmuA6MGWTUPOD3pFBXMAJ4PbEk6SBXPBhaQztcOYC7QTnrzAXQCrwcOTzpIFQuBcaT3NXwBMJX05oPow7Y055sCHEv0s5hG+5Du7/ESotcurfkger/0VmBd0kGq2Ak4hNS9hruPj5Z8kxorSDpAmVVEb+YOApYnnEU1hZcCJz69fdM2OGx1YnGq2wXYBKxPOkgVs4j+cD+cdJAq2oH5wIqkg9SwK7Aa6E46SBXzgT5gbdJBqpg6cHso6SA1LAIeTDhDLQuBLmBj0kGqmAu0AY8kHaSKcUQZ0/x7ZhHRe5S+hHNUs4Dod+DjSQcZakYOnlhYtvM0CL6TSBypSVYBIXBg0kE0nPBSCMOS27VJJ6riXuDMpEPUcAFwfdIhalgCbE46xDA2AgckHaKGq4ALkw5Rw+nAA0mHGEZv0gGG8RBwStIhargIuDzpEDUcTlQAptlWoiIhrW4Czk86xPbCaWXvFUIIT006lbLPa2AkSZIkZYYFjCRJkqTMsICRJEmSlBkWMJIkSZIyw2mUpeHtDOzB0wX/DUQzm0mSJKnFLGCk7b2UaD79g4mm9J5X9vgBwK2tDiVJkiQLGKmS/yVal0CSJEkp4zUwUnVbiU4X+xouJSxJkpQKjsBI23sX8DfgLp5edflyYPfEEkmSJAmwgJEq+WHSASRJklSZp5BJkiRJygwLGEmSJEmZYQEjSZIkKTMsYCRJkiRlhgWMJEmSpMywgJEkSZKUGRYwkiRJkjLDAkaSJElSZljASJIkScoMCxhJkiRJmWEBI0mSJCkzLGAkSZIkZYYFjCRJkqTMsICRJEmSlBn5pANIKbQUeG/Zvkkl968B+ku2Q2Bms0NJkiTJAkaqZDwwvcbjU8q2wyZmkSRJUgkLGGl7VzJ0hGU4FjCSJEktYgEjbe+3AzdJkiSljBfxS5IkScoMCxhJkiRJmWEBI0mSJCkzLGAkSZIkZYYFLine truncated
|
||||
}
|
||||
},
|
||||
"cell_type": "markdown",
|
||||
"id": "cc6ce65e",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"<div>\n",
|
||||
"<img src=\"attachment:g27851.png\" align=\"left\" width=\"300\"/>\n",
|
||||
"</div>"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"id": "3763b000",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"<div class=\"alert alert-block alert-info\">\n",
|
||||
"<b>Note:</b> This algorithm is not correct for all matrices: if any diagonal element <code>B[k,k]</code> is zero, the computation in the first inner loop fails. To get around this problem, another step can be added to the algorithm that swaps the rows until the diagonal entry of the current row is not zero. This process of finding a nonzero value is called <b>pivoting</b>. \n",
|
||||
"<b>Note:</b> This algorithm is not correct for all nonsingular matrices: if any diagonal element <code>B[k,k]</code> is zero, the computation in the first inner loop fails. To get around this problem, another step can be added to the algorithm that swaps the rows until the diagonal entry of the current row is not zero. This process of finding a nonzero value is called <b>pivoting</b>. We are not going to consider pivoting here for simplicity. \n",
|
||||
"</div>"
|
||||
]
|
||||
},
|
||||
@@ -203,6 +313,16 @@
|
||||
"gaussian_elimination!(B)"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"id": "a7d6f90a",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"### Complexity of the algorithm\n",
|
||||
"\n",
|
||||
"The number of operations of the algorithm is $O(N^3)$, where $N$ is the number of unknowns in the system. Intuitively, this makes sense as there are three nested loops. However, the length of the loops is not equal to $N$. In any case, it can be proven that the total number of operations is still $O(N^3)$. The actual proof is a bit challenging and it is not discussed here."
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"id": "39f2e8ef",
|
||||
@@ -220,17 +340,15 @@
|
||||
"\n",
|
||||
"```julia\n",
|
||||
"n,m = size(B)\n",
|
||||
"for k in 1:n\n",
|
||||
" for t in (k+1):m\n",
|
||||
" B[k,t] = B[k,t]/B[k,k]\n",
|
||||
" end\n",
|
||||
" B[k,k] = 1\n",
|
||||
" for i in (k+1):n \n",
|
||||
" for j in (k+1):m\n",
|
||||
" B[i,j] = B[i,j] - B[i,k]*B[k,j]\n",
|
||||
" end\n",
|
||||
" B[i,k] = 0\n",
|
||||
"for t in (k+1):m\n",
|
||||
" B[k,t] = B[k,t]/B[k,k]\n",
|
||||
"end\n",
|
||||
"B[k,k] = 1\n",
|
||||
"for i in (k+1):n \n",
|
||||
" for j in (k+1):m\n",
|
||||
" B[i,j] = B[i,j] - B[i,k]*B[k,j]\n",
|
||||
" end\n",
|
||||
" B[i,k] = 0\n",
|
||||
"end\n",
|
||||
"```"
|
||||
]
|
||||
@@ -262,72 +380,50 @@
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"id": "14d57c52",
|
||||
"cell_type": "code",
|
||||
"execution_count": null,
|
||||
"id": "1169c86e",
|
||||
"metadata": {},
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"### Two possible data partitions"
|
||||
"ge_par_why()"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"id": "6b17aee4",
|
||||
"id": "74491f64",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"The outer loop of the algorithm is not parallelizable, since the iterations depend on the results of the previous iterations. However, we can extract parallelism from the inner loops. Let's have a look at two different parallelization schemes. \n",
|
||||
"### Data partition\n",
|
||||
"\n",
|
||||
"1. **Block-wise partitioning**: Each processor gets a block of subsequent rows. \n",
|
||||
"2. **Cyclic partitioning**: The rows are alternately assigned to different processors. "
|
||||
"Let start considering a row-wise block partition, as we did in previous algorithms.\n",
|
||||
"\n",
|
||||
"In the figure below, we use different colors to illustrate which entries are assigned to a CPU. All entries with the same color are assigned to the same CPU."
|
||||
]
|
||||
},
|
||||
{
|
||||
"attachments": {
|
||||
"fig-asp-1d-partition.png": {
|
||||
"image/png": "iVBORw0KGgoAAAANSUhEUgAABL4AAAI0CAYAAAATGeCIAAAACXBIWXMAAB7CAAAewgFu0HU+AAAAGXRFWHRTb2Z0d2FyZQB3d3cuaW5rc2NhcGUub3Jnm+48GgAAIABJREFUeJzs3XmcXFWZ//HvU1XdnU5nIyEhCYRVGBYBWRwQERFk2JRVUBxEXAYXXEcUB3QE11FREcF12BQFQWRRGJBNVFBQdhAwkIXsezq9d1Xd5/fHrfzsVN1Oeq17cvN5v155Kae2p+69ffv09557jrm7AGy+zOxsSVdXNX/N3S+oYw3Nkjqrmue4+y71qmFzYmbtklr6NLW5+4S06gEAAEA2mdmnJX2jqvk8d/9WGvUAaSikXQC2XGZ2uqRJfZpud/elI/Teh2tgx3e7pF5JayV1S+p097UjUQMAAACyzcxOkbR1n6Y73X3hCL33oZLGDOCpHZJ6FPdnexT3Z9eMRA3YcpnZREn/lfDQbe7+53rXAwwHwRdSYWY5SVdpw1Evz0kakeBL0i3aMFQbMDNbJ+kFSX+X9Lyk+yQ97gyPBAAAwIZ+KGlqn/8+RtKIBF+SfiFp1lBeaGZtkl7UP/uz90v6K/1ZDMJ4SecntC+URPCFzQrBF9JyqDYMvUIyQdK/Vv6tt8DMbpf0fXf/ezplAQAAIBRmtr82DL1CMl7SgZV/6y02s99I+oG7P5VOWQBQf7m0C8CWx8yaJH0t7ToGaZakcyU9Y2a/NLM90y4IAAAA6TCzBtXOmxS6mZI+IOkJM7vFzPZNuyAAqAeCL9SVmc2U9FtJh6RdyxDlJJ0u6XEz+5SZWdoFAQAAoH7MbKqkX0s6Mu1ahsgknSTpb2Z2oZnl0y4IAEYTtzpi1JjZBMWTfU6TtI+kNyn+JTuQSTpHw/skJQ3rbpHUKGmypB0k7SLp1ZIOUv8/I02SLpH0ejM7w917Rr5cAAAApMnMxkuaImkbxf3DwyWdrPSm7PiokudXGqu4f7qV/tmf3UvSwZIa+nmvgqQvSzrUzE519+oVurFl65D044T25+pdCDBcBF8YFWY2W9Kr0q6jyovu/thAn1xZyeTNkj6k/q/onSzpDjN7q7t3jUCNAAAACICZPSHpNWnXUWX2IPuz4xX3Y89RPPF+0t0Kx0j6nZkd6+5tI1MmNneVlUE/kHYdwEjgVkeMls3+2HL3Vne/2d3fLOm1ild3THKkpCu57REAACBTstCfbXP3W939OEn7Sbqzn6e+XtJ1lZXXASBTOLEBA+Duf5N0lKSPSColPOUMSZ+oa1EAAADAALn7U+5+vOLpP5Km6ThB0oX1rQoARh+3OmK03KB4fq/+nKR47q/Nhru7pCvMbK7iCU2bqp7yZTO71d3n1r+6gTOzaZK2U7yyT0HSMknz3H1JqoX1ozLh6s6Kj5epklzSCklLFdcd1bGOrSv/xkhaK2k5twQAAJBZN0v6y0YeP17StnWqZcS4+1WV/uxvFc8N1tfnzOxX7v58CqVtkpmNUdwXmyqpLGm1pBVZnG/XzJr1z75nJr+rmTXqn98xr39+R+abw4gi+MKocPeNXi0ys320mQVf67n7nWb2QUlXVz00VtK3JJ1S/6o2rhJ2fVjxyLTd+nnO3yT9QtIP3L27juUl1VJQHI6erHjeicn9PHWZmd2hOIi8sxJOjmQdhyju1B6t+PaAmlGyZjZPcaf4Pkm3uPuqkaxhE/Vto3ikYfVttt2SvkGnAQCAoXP3L27scTO7R5th8CVJ7v6AmZ0t6ZfasB/RKOm7kv4tjbqqVYKuExTXc7Tii7fVSmb2tOJJ/++QdK+7FxPe61TF05f09aK7V/fph1Ln8ZLeUNX8F3e/dRDv0awNv2vSsVXs811/K+n+pO86UszsI6rd5l9z99Yhvl9e8Xc7VvH3TPq7JDKz5yU9Iun/FPfx6dNiWAi+gCFw92vM7C2STq166EQz+xd3fzGNuqqZ2ThJn5J0nqRxm3j6gZV/55jZe909acWgUWdmR0u6XANbHGEbSe+t/HvSzD7g7o+OQA2vl/QVSW8cwNN3rPx7h6QfmNndkr7r7vcMt46NMbOtJd2reIWpvroknUQHAQAAbIy731QJbN5d9dBRZra/uz+eRl2SZGYNivt3n1Ny2NVXQdL+lX/nSlptZtcp7o/N6fO8LknnV722s3LHxpph1JqTdJniOxT66m9xrOrXNyq+/fRCbTpIbZB0QOXfRyStMrOfSbpslO46ebfivw/6ulzSoIKvylzIp0r6oqQ9NvH0nOIVSfdSfAx0mNmvJF3q7k8O5nOB9ZjjCxi6/5RUfYUlJ+ljKdRSw8wOlvS0pIu06dCrr90l3W9mR4xGXRtjZp9RfGVnKCuCvkbSn8zsXcP4/LyZXSrpTxpY6FWtoHiE2O/M7KKh1rEpZraVpHuUHHqd6O6/G63PBgAAmXK+4v5DtY/Xu5D1zGxbxaPpf6hNh15JJivuj//DzPbt036XpJernjtW0tlD+Iy+jlFt6PWcpAc29UIzm6V4ZNP3NbTRg1MUj/6fbWbV/cIgmNkExaPTbtKmQ68kLYoDuCfM7OSRrA1bDoIvYIjc/RXFt9hVOymAFR5PlvRHSTsN8fVjJN1mZnVbwtvMviHp60peZluK5zZYqnhOsv5uaWyQdK2ZDXqhgcpy37dr5Dp6LSP0Phsws4mS7lbt8uqdkk4Y7ZFmAAAgO9x9meKpLqq9tTL1RF2Z2f6SHlU8emu48pKa1/9HZV7YHyY874PD7Lufm9B2+aam4DCzAxWHXiPR384r7r8Hxcx2kPSQpONG6C0HczEf+P8IvoDh+XlC20yNzC/r4dhdtbcyr5V0naTPSnqXpLMq//9+xaFStXGSflKPZa3N7CxJn054qFvSpZJeJ6nR3We4+3TFAdfhkn6i2lU2TdK3zOyoQXx+QdKvlPxLuaz4CtW/S9pT8VXEqZL2lnS64rnelg/0s4ajcuvqnaqdn2J96HVvPeoAAACZktSf3UrS6+tZhJntJul3ivvS1ZZL+rbieaF2kjRB8WiwgxQHT3creaXKalepdoTbbpLePMSad1Y84quvtZJ+tonX7a74u85IeHiZ4nmDj1I8ncb673qwpI9WXhf8BPdmNkXJdyhIUrukHyue02w3xd9xhuJbON+jeGGJ9vpUii0Bc3wBw/N7xbc7NlS1HyzpsbpXk+yvkr6qeGLI3oTHv14ZBn6dan8xHSjp/Yp/MY0KM9tJ0vcSHnpa0inuXj0kXe5elvSgpAfN7PuSbpW0Q5+n5CRdY2b7DHCy+c8peRLXRyS9392fTXhspaRnJd1UmYfi3xUHif8ygM8bNDMbq3iY+CFVD3VKequ73z8anwsAADLvIUkdqh2tfrDi/taoq8xz9UvFt+5Vu0zSfydMqN4maZHiEWLfryz685+SPiRpfNLnuPtqM7te8dxRfX1YcUgzWB9S7WCSq9y9o78XmFmTpBsVh4vVLpX0BXdfV9W+/rs+IulyM5uueB7fDyrcUVBXSto1of1Xkj7u7our2tsU393xuOJ+/HjF++U/tZkuioZwMOILGAZ3b5OUNJH9nvWuJcFixZOc/6u739pP6CVJcvenJB0hKWlSzNGe4+Griq/y9PWCpCOTQq9qlUkuD5e0pOqhmYonCd0oM3uVpAsSHvqdpCP6Cb2qayi6+zWK9/snlTxXxpBVVjS6TbXzjnVIeguhFwAAGKpKH/GZhIfq2Z/9mJJv+fuku398IKsIuvsydz9f8XxbN27kqVcktL2lMt/WgFVWYXxPVXOkeL6ujfmk4jsH+nJJH3P3TyaEXjXcfam7f1rSLopHRwXFzE6UdGLCQ5dLOj0h9Krh7m3u/nXFI/y+q3jbAkNC8AUM3+yEtqSleevpAUmvdvfbBvoCd1+h+KpVtT0r8y2MuMrkpdUrY0aKR1mtHOj7uPs8JS8q8N7K1aKN+ZxqR+wtlPT2wa6M6O6Ru18qaT/FV+SGrXIF9GbVDsHvkHS8u29y4lQAAIBNSOrPJo3WGXGVAOmzCQ9dU+lXDYq7r3T3t0t6u6Sakf+V1SqrVy8vSDpnkB/1DtWOULtzYxduzaxF0mcSHrrS3ZPugNgod1/u7m+T9E5JQ16ZchRcnND2R8UjvTY691k1d+9090+o/4v0wCYRfAHDtyyhbVLdq9jQX4ayLLO7363kWzTfNvySEp2j2tDpN+7+0GDfyN1/pfi2zr4mKp7PLJGZTVLcKar2CXdfO9ga+tTyYqWeYancQvlL1c491i7pOHevy+0HAAAg85L6s0m34o2GU1UbIK1VvFrhkLn7je6eFOhJyaOy3l+54DhQH05o21R4dZpqt+tqxbctDpm7Xz+QOyXqobKy/L5VzesvbA951Ja7P+jufxpWcdhiEXwBw5c08eKmRhmFLGmlygNG6bOS5tW6bhjvl/TajU1yf7hqV8BZovi2wlSZWV7StZJOqnqoQ/FE9n+of1UAACCjkuakqld/9tiEtp8P5PbGYbhRtWHfdEmnDOTFlXDnwKrmf2jT84QlfdefDeT2xs1I9WT/knSPu/+j7pUAFQRfwPAlrYi4OS8ckbQy4Egss7yByrxV+1U1u+JVeYbqjoS2QzayRHX1nFmS9Et3r14psq4qK2leI+mMqofaJB3D7Y0AAGCEJfV96tWfPSyh7YbR/MDKvGZXJjyUNO1HkqTnXTaA2/jq/l1TkPQdr697FUAfBF/A8FWvgCNJm/NVmxcS2qZVliQeSa+R1FTVNreyYMBQzVEcDvU1TfGkmEkOTWh7dBifP2yVkO4nks6seqhN0rEM8QYAAKMglf5sZXXv7aqaS4pX9httP1LtBezDzKx6lfMNmNlUSadXNa+T9NNNvG5XxaPK+uqV9MSmS908VKbpODjhoerpSIC6IvgChi9p/oPRHJo9qipDrVckPDTS8zzMTGhLCt0GrHKVLWmVzW37eUl1R0uSnhpODcNRCb2uUO0S2+sUj/Qa9NxnAAAAA5BWfzapLine truncated
|
||||
"g31767.png": {
|
||||
"image/png": "iVBORw0KGgoAAAANSUhEUgAAArUAAAHoCAYAAABXWBf1AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAAewgAAHsIBbtB1PgAAABl0RVh0U29mdHdhcmUAd3d3Lmlua3NjYXBlLm9yZ5vuPBoAACAASURBVHic7N13nB11vf/x1/ecs9l0UiCEhC5FShAISBGkiHAFERW7KAq2a7lYsF69ev1ZLhYUwYIoYi+ABRVFQAWU3gnSpJsQCKRnk83uOd/fH7PBzdnZzfaZb3g9H4/zgP3Ome98dmZ289453/lO4GnxauAAJEmS8p0D4W1FFyHlqXT9dyIsm1hoJZIkqeQemgFUi65CyhMg/h2WTYLa7jAhFF2QJEkqq0URpt4AtXshvKHoaqTuAsRYdBGSJCkp8yDMKboIqbtac8PWu8EWzyqilL7d9If/iXNf9OnSXkm+7bIv1HfY94TqhE1mFV1Krntv/AGTp23HzO0PLrqUXIsX3M7CB/5e3/Wg/yztx1plPwdv//OX6tvt+YrqpGnbFl1Krvtv+Rmt4zdtbLnzEZUNv3v0LXvyPh658+LGnENOKWV9kJ2De//Hp0II5Szxjr9+pb71bi+ubrLZjkWXkuuhO35FqLQ0ttntxaXcgWvXLGfelV9r7H3kx0tX3wO3wKJHi65C6luPULvFDjDn0AIq2YBbLzuLOYd9uugyenXX389p7LjvkdXps8sZah9/+CJmbncQux9SzlD74G33snzxTxpzDitvqL310jNLfQ7ec+13G8/a+8Dq5tttW3QpuRY/9gcmT98hzjnsiKJLybXgvkdY9MgPGnMOK2+ovfXSs5hz6Ccpaablvhu/19huzznV2TuVM9SuWHwZ1ZbxjTmHlTPUti1byT3XfjPOOezjRZfSw5KFhlqVXyl/sCVJkqSBMNRKkiQpeYZaSZIkJc9QK0mSpOQZaiVJkpQ8Q60kSZKSZ6iVJElS8gy1kiRJSp6hVpIkSckz1EqSJCl5hlpJkiQlz1ArSZKk5BlqJUmSlDxDrSRJkpJnqJUkSVLyDLWSJElKnqFWkiRJyTPUSpIkKXmGWkmSJCXPUCtJkqTkGWolSZKUPEOtJEmSkmeolSRJUvIMtZIkSUqeoVaSJEnJM9RKkiQpeYZaSZIkJc9QK0mSpOQZaiVJkpQ8Q60kSZKSZ6iVJElS8gy1kiRJSp6hVpIkSckz1EqSJCl5hlpJkiQlz1ArSZKk5AWIsXtDreV/YrXl60XV06vOtatCbcyEuOF3FqPR6AwhVGII5fw7oXNtW6hUW6hUW0q5Dxv1DhqNjlBrGV/K+gA6O1aHWsu40tYX6x2BSiWGUC26lFydHatDpVKlUh1Tyn3YqHfQqK8t9e+Zemd7qNZaS1tfo94RQqgQKtVS1ljvWE0I1VCplfMcJEZirIdQqZWuvs61Z4V652u7N82DMKeoeqQ8PULtfsfBnEMLqqYPP/zvafENn10ciq6jN7/4zI4dL3jz+S3TZ+9ZdCm5LjvvFczc7iB2P+S9RZeS68HbLmDelWd0HPueq1qKrqU3P/zY1PiGzy0p7Tl4wed36Tj4Nd9t2Xy7A4suJdcVPzmRydN3qO911CdKmboX3Hc51130oc6XfeCmWtG19OaHH5sWT/jMohAqpdyF/PKLe3Tsd9zpLbN3OqLoUnJdfeG7qLaM79zvJV8s5TFuW7aAX5++b/11/zu/dAf4yp/Cvdev12SoVemU87KiJEmSNACGWkmSJCWvlB/BSJIkDU4cDxwFHAzsBWwHzADGAZ3Akq7XP4Frul5XQVg7gG0cBWzTzzevBJYBS4E7ICzv/3ae3t5cYG5T468hPDHwvtbr9zBgx6bGH0BYM7R+ByvOBvYBNm9acBeEqza0tqFWkiRtBOIuwAeA1wATenlTDdis67UTcHRX+yKIPwLOgXBXPzb2TuAlgyiyAfEe4LfAtyHc38/1jgU+2dR2OzDEUMtJwAlNbb8CRiHUxi3Igvo+3f47s5c3fxvYYKh1+IEkSUpYnAjxa8AdwMn0Hmj7shnwPmAexG9DnDGcFXZTAXYBPgTcC/GbEAdTb8LisRDnAwvIwv0ngRfTe6DtN0OtJElKVNwZuB54D9DXrBGdwELgIaCtj/dVgLeSBc5xw1RkX9t6B3ADxM1GeFtlsgUwayQ6dviBJElKUJwDXEY2XrbZCuAC4BLgCggLm9adCOwKvIjsKuFcoPuUjZvQd0jO81Egb47hiWRjRPcH8qZB2wX4PcTnFzeWtVTagdbBrGiolSRJiYmzgT/QM9B2At8C/l/fN1GFlWRXeK8H/hfifsBngRcMoajTIGzgwRlxLvB1YL+mBfuSjdM9fQjbT9Ea4Dbgpq7XjcCjwOLBdGaolSRJCYkV4OfA7KYFbcArIVw88D7DdcAREP8D+D75V3+HQbgpuyLL74AXNi38CMSzBjYLQ5JuBd5CFmLvhNCx/uI4abAdO6ZWkiSl5F3A85ra2oEjBxdouwt/JLsL//oNvXMI21hLNuPAqqYFm5ENUdjIheshfBfCrT0D7dAYaiVJUiLiJOATOQs+CeHvw7ON8CjwfGD18PSXu40nyK42Nzts5La58TPUSpKkVJxIdkWzu5uBLw3vZkI7hPrw9tnDFTltW47wNjdqhlpJkpSKk3LazhiFADoSFuS0bTrqVWxEDLWSJCkBcSuyx952txg4v4BiVEKGWkmSlIK86bYuhzCCY19HVPPsDQBPjnoVGxFDrSRJSsHeOW03jHoVw+f5OW3/GvUqNiKGWkmSlIJn5bQlGmrjDODVOQv+OsqFbFQMtZIkKQXTctrmj3oVQxZbgB8CE5oWPAVcM/r1bDwMtZIkKQXTc9qWjXoVQxL3Aa4CjsxZ+H/PgKeJjSgfkytJklKQd6W2TKH2wxBjTvt4YCZwADCnl3VvBr4+UoU9UxhqJUlSqvJCZFE+P8j17gOOSXgWh9Jw+IEkSUrB4py2TUa9iuETge8CcyEsLLqYjYFXaiVJUgryQu0UYNFoFzIEEbgX+B3wbQj3FlzPRsVQK0mSUpAXarcg+/i+DN5B/nCINrKxv0uAOyCUaRxwKLqA4WSolSRJKXggp21f4MrRLqQX34YwUmN86zlt44ah37w+Ooeh30I4plaSJKXg5py25456FcVYmtM2cRj6ndT0dQNYPgz9FsJQK0mSUnB5TtvhEMeOeiWjLy/U5s3bO1CbNn29HELeVeEkGGolSVICwsPA7U2NmwKvKKCY0fZETttuQ+syVoCd+7GdZBhqJUlSKs7NaXtPV0DbmN1Iz5vQ9hpin88mezBEdzcMsc9CbewngSRJ2nicBzzV1PZc4JTh3UwcU66gHJ4C7mlqfD7E5uEDA/HynLarh9Bf4Up0wCRJkvoSlgGfzVnwWYj7DM824izgz/S8ilm0PzZ93QKcNLiuYivw5qbGBnDJ4PorB0OtJElKydeA65raxgGXQ3zB0LqOLwBuAp43tH5GxFn0nNrrExC3HURfHwW2b2r7LYT7B1NYWRhqJUlSQkIdeDXweNOCycDFEL8IcdrA+ox7Q7wYuAyYORxVDr9wP3BBU+NE4K8Qd+lfHzFA/BDwP00LGsAXhlhg4Xz4giRJSkx4GOIxwJ+A7gF2DHAq8BaIPyf7OP2vEJasv34cS3bn/zFdrwNI4+la7ya7irxlt7ZtgFshngv8GLgewtr1V4ubA0cA7wf2zun3NAijOJ42ngLs2svClpy2gyGe3UeHX4Nwp6FWkiQlKNwE8RDgV8AOTQunAG/vegFxNfAk2SNrN+9a3pfF5D/Fq2DhSYivIhtfO7nbgjFkj+l9B7AG4kKyG+payeaz3aKPTi+h55XbkXY0cOQA3r9L16s3vwLudPiBJElKVJgHzAXOJvsIvTfjgK3Irs72FWg7ga8DO0FYPVxVDq9wDbAfPWdDWGcssC3ZftmdvgPtOcBxEJJ9NG53hlpJkpSwsBzCO4A9gR8CawbRyRPAl4HdILy7awqtEgt3A/sAHwQeG0QHfwIOgfA2CO3DWlqBHH4gSZI2AuEO4I0Q3wX8B3AQ2QMKtie7WrnuQt5SsuEF9wLXAteQjbtd29xjH74B/L5p+80PRxhhYSXwJYhnko2XPZjse96KbJzxRLKr14u7XncAVwKXQ7hzdGvt4QzgwmHs706AAHG9gzBu0un1sRO/V7pxJKuWPDpmwtStBnLCjaq2ZQtaWsdPa1RbxpZu3wG0LV/YUmsZy5hxUzqKriVPR/uKSkf7yur4yVuUsj5I4BxcvqBl7PhNG5XamFKeg6uW/Ks2ZtwmoWXspFIe4441Kyod7Stq4zeZVdpjvGrJo2MmTNmqg9DjyUKl0LZsQcuYcVMatTHjS3kOrlmxqBYqldA6YXopz8HYaFRWL3+sOn7K7NLVt2rpZ2prVx/X/dPdeRDmFFbQoMQa0IDQ1zCFjUhsgVC6c2kk9bhSu+2cl1S323PPahHF9OWyc1/OAS87Y0zRdfTmyp+e1Njj8A9VJ2/6rNLtO4CbL/k007bYnW33eHkp9+HC+6/iodsvbOz/sq+Wsj6Ay753fLnPwZ+d3Njt4PdUp8zctbTn4OTp28Ud9nlDKffhokeu597rvtco8zG+7HvHx/1f9pWWEMo5cuzv5/9n/dkHvq06ffZepTwH77zya1RrYxvPPvBtpTzG7W2LuebC98QynoO3/3kX/nV30VUM1cYxbrT/nlmBFnJC7eTNdmDWjs03ERYvVCpx1o6Hl3a6jUqlpb7Z1vtWps/es+hScv3j799g8qY7MGvHw4suJVd722Lm33dZfdaOh5fzX2uyjzXKfA7WamPr07fcu7L5dgcWXUqu+274PpOmbd+YtePhpQw8EHng1vMb5T4HA7N2OIxQKecurLa0NqbNek61rL9nHrr9Qqot40t7jNuWLYBQKeXPyD9vLLoCacNK+YMtSZIkDYShVpIkSckz1EqSJCl5hlpJkiQlz1ArSZKk5BlqJUmSlDxDLine truncated
|
||||
}
|
||||
},
|
||||
"cell_type": "markdown",
|
||||
"id": "c518f863",
|
||||
"id": "7460b0af",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"<div>\n",
|
||||
"<img src=\"attachment:fig-asp-1d-partition.png\" align=\"left\" width=\"400\"/>\n",
|
||||
"<img src=\"attachment:g31767.png\" align=\"left\" width=\"250\"/>\n",
|
||||
"</div>"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"id": "102d6fa2",
|
||||
"id": "3ebddda1",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"## What is the work per process at iteration k?\n",
|
||||
"To evaluate the efficiency of both partitioning schemes, consider how much work the processors do in the following example. \n",
|
||||
"In any iteration k, which part of the matrix is updated in the inner loops? \n",
|
||||
"### Load imbalance\n",
|
||||
"\n",
|
||||
"### Block-wise partition"
|
||||
]
|
||||
},
|
||||
{
|
||||
"attachments": {
|
||||
"fig-asp-data-updated.png": {
|
||||
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAtcAAAJNCAYAAAD6c1l4AAAACXBIWXMAAB7CAAAewgFu0HU+AAAAGXRFWHRTb2Z0d2FyZQB3d3cuaW5rc2NhcGUub3Jnm+48GgAAIABJREFUeJzs3XeYJFXVx/Hv6ZnZvLC7LGFZwrLkDIpECRJEEMVXl/waEQyor6IoBhAFEwZEBUFUVEREEQwgIkEyIkGQKJkl7i5h8+ykPu8ft4btrqmeTtVTNTO/z/P0A32n6tbZ6urqU7du3Wvujohkx8w6gXHR2053n5BlPKOFmc0GHi8pusvdd8ggjo2B24GpJcW3APu6+4qhjkdkuDCz04AvlhQd5+5nZxVP3pjZGGBtYDxQBJYBL7p7bw3rfhU4KVb8Pnf/ZeqBjkDtWQdQjZl9B5g8yCLdhANmafT/c4GHgYf1wyQieWZmqwB/pDyxvhk4sNL5y8zWAy6NFX/C3W9tTZTSKmb2e2CDkqJL3f3rWcUjw5uZjQPeARwI7AasD7TFFltuZr90948OVpe7n2xmHcCJJcXnmNmD7n5HmnGPRLlProF3A2s0sF7RzO4B/kQ4Yd2fbljpMrMdgXeWFN3g7ldmFY+Mbmb2JmD/kqK/uvuNWcUzgp0DbFHy/r/A2919ySDrjANeHytbtdLCZjaZ8ta9+e7+vXoDlXJmtivw9pKia9396jqr2YLyz//OpgOTUcfMJgL/BxwPrFZl8QnAZrXU6+6fN7OZhDwMwrnnEjPbxt0XNRrvaDAckutGFYDXRa+vmNkVwCnunteT13bA52JlSq4lK7tSfjy+BCi5TpGZHQIcUVK0kJBYv5rypiZR/lk+CCi5bt7rKd+vXUC9ybVIU8xsF+DXwOwWbeJYYBNgp+j9eoTzx9Et2t6IMByT6yVAaX+hVRh42yPJW4EDzOx04GR372lFcCIi1ZjZNOCsWPHH3P2RLOIRkeHHzOYQEuuxsT8VgbuAfwILCF1n1wVmAXvWsw13X2FmRwL3sLKL7gfM7HfuflXj0Y9swzG5PsDdbyktMDMDpgBrARsRWhTeDOxIeeJdIPQf2tbM5rj78qEJWUSkzMnA6iXvL3H3C2tc91lgv1jZv1OJSobaB4GJJe+fyyoQGV7MbF/gN0BHSXEn4aL92+4+v8J6HcDG9WzL3Z8ws+OB80qKv2dm29bycORoNByT6wE8DHnyavR6CPgLcIqZbQR8BvgA5QfgAcAfzexAHRgiMpTMbBOg9GGi5cAna10/ahS4Ju24ZOi5+21ZxyDDj5mtxcDE+ingIHd/YLB1o7v2Dzaw2Z8RuoLsHL3fgtBlRKOzJChkHUAruftj7v5h4A2EpLvUfsC3hj4qERnlTqT8R/EMd1eLpYjU6kzK73w9A+xaLbFuRtSI+ZlY8YlRS7jEjIiW62rc/V4z2w34O1A6ju0nzez37v7PWusys3ZC36X1CE/ljgH6CK1P8whDAC5NLfgmmNkkwjBPMwl9pdoIsb4CPA086e592UWYPTPbENgSWIewb54HbnP3l1KqfzPCFf46QA/h7sp/CMdJMYX6OwjH4rrAdMJ3uo8wNOWL0XY6m91OmsysDdieMEzUmoSL/JeBJwhjTTd9NynqKrY94UGcGYTv5wLg3+7+ZLP1NxHX2sCRJUULgdMzCicVZrYusA3hx34K4RhfQPgepf1wZt2ic3bpd6SD8B1Zxspz9rLsImyemW0AbEX4DCYT9v+zwL9aNSRttM2tCfvVCeeb29z9hVZsTwIz2wk4pKSoB5gzFPvd3W8xs8uBg6KidYHDgQtave3hZlQk1wDu/qqZHUzom9g/tF8B+C5hPMiKoocG9iK0gG9HSKgHW/4e4DLgR+7+yiDLTQRuiN5Oj/35vVGfqkqucvfS4bUws9WA90dxvoHy8VOTdJrZ1YQvxh98CGYUMrNbWbn/utx90H0fW/cA4NSSot+6+3cGWf47hM+t37HufneUkB4NHEMYTSau18yuAk5w9/gdj1ribAc+RLj1v0WFxZ43s4uBHzRQ/3sIx+wbCD+og7UcFM3sLuD3wI8Hu/Azs3UIYy5DSEhL/Z+ZHT7Idgb9LKL6NyIMCXcQA4/3fovM7BLgNHd/arD6KmxjInAC4XuwXoVlHiE8BPS3eutPwbGUP3z0C3dfXG8lZnYl5S1Xu7l7V2yZnxCO7/jxMdvMBhs1abm771Fl++OA44D3EBLrJH1m9k/gW+7+l8HqK6n357H65pQeB9Ht8P8Bdif8sI8jjGRzm7t/tWS5w4E9CN+RbRn8O+Jm9m/C2OFnufvCQeKbysoRQeJDxB5jZm8dZDuXu/spsfpOBOaUFH2x1ofEoglCjiN0e9yqwmLLomPlVHf/Ty31RnX/kdAg0G9vd18cndveTzi/xYeDhLAvrwc+l/exkM3sy8DbSoqcEPd1GYVUq88AVvL+Z+7+ryHc/pmsTK4BPoGS64HcPdcvQsuCl7x2a7K+98fqc2CXKus8kLBOLa+XCT8OlepdpcF6HfhNQn07N1HfLcDaQ/B5dpZsc3md6/5vLOYzqyx/WWz5PQkPuj5W4z7pBP63zhg3A+6rY7/3EJ7srnmfAPMb/IyfB/YbpN4Nmzh+Kn4WhDsmpxMmeaq1vi7gA3Xu+zcRWutq3UY8njuH4Ph/pGR7fcBGDdbzXCz28QnLXN/gZ7mkyrb3I9z5qqfOPwOTavh33RJbb7OofE3Cg1qVjqG7Y/U83uC/fT5w8CDxrd5gvQ6cn1DfWbFljqjx89+B2s9jThhh6/tAW431x+ueBuwTO36rfX8PG4Lv02mx7X60xvW+HFuvSJjdsaXxpvDvnRb7DvQC6wxxDEbos126/zbOet/k7TWi+1xXcAHhB7jUUS3a1jTgt2Z2aIvqT9OuwI1R6/dI9XlCa+WGNS4/Dji/SmvUa8zs9cBNVG5FStJOeStEK80A/lLljkiqopbkywityfHWw2eBewnPQ8QnJBgD/MzMPlHjdt5BGBd+Zh3hDWlfQTPbgfKn9G9298eGMoZmmdn7gSsYeFdgEeFzvIeB51cILYRXm1l8yLBatnkoYXKdj9L6z2x1wiQZB7d4Ow2LJni6juTzmBO6QMW1ESYZ+X0jnwHwVUK3ylpHmRgD/NLMtm5gWy1lZp8DTikpcsIwmPGhMfPorZR/B25y96TvW8t4yLDjU6AfNpQxDAejpltIP3fvNbPfUt4x/8AaV19OaFm5mdA6+QRh3O0xhOmL1yXManckITGDcFI7y8yudfeXY/V1sfKhym2Bt5T87XZCy1Ml91aJ9RngH8C/CC3v86LtTSIk/VsRLip2LllnwyieD1ape7gqnXGwmzDiwnWEZKBA+PcfSuhH2K8d+KmZbeqD3L6Pxi3+IwO7O9xIeKr7fkL/2tUIx8mmhBNlUreUWi2J6r8lqv9pQl/rsYTjcQPCsX0IK0/IY4HzzGxLHzgU5UJWHo+7EG6r9+s/liqpNPX2eZTfel1I6Ir1S3d/pr8w6oe9I6HbSOnFzLfN7FYfZPInM9sCuJDy7hZO+DwuJbTCLSFcXMwgdO16O2HYzqH05tj7mrpKNOEiwji3kwjdB/q9RHjyv5LupEIz2x/4KeUPwl9K+Dxv95LnN6K+5R8kTLIyISremXB81TwyCmGyigPqWD5uGeF83X/Oforyc/b6rDxn93dXaydM83y9D5yFbjkrvyOvB0ovVPt/Gyq5q+F/RcTM1ifs88klxYuB7xC6fj3q7n1mtgrhTt3/EVqc+/0P8E3gU3VuuvT46SV0jbmG8DvTRvguHUZ5w8JY4NuU/65lysw+Sfj393NCi/WPq6zXRjhWWqnTq/eb3jH2vuwcEnUV2pxwkTiVcNH7KvCgp/tcwZ8p349vJtxFkH5ZN53XcAsi1W4hUZ1vYeBtrOmDLP8V4GBgbI31b0N4oKS0/o9XWefY2PLfbODftTHwNWCbOtb5Umy73cCUFn6eWXYLccKP7enA6hXWKRASvPh6p1XZ1sWx5Tup4bYoIdErvc1XS7eQ/h+sjhr3266EpLs0vsOrrBPfB59p4LM+OlbHf4ENq6xjhGStdL2rBlm+jdBaWrr8i9XOE9F25sTWa2m3EEKLb+n2NmmirqrdQkqWnRFb9oEGtrcW5d2ReoH317DeTtF3rn+9HmD9QZaPdwspfT1F+AHfi3AhvD4hmTsIOCpWz9cIF3W1nrNfR3jQu3R7x1ZZ5+Ox5b/SwH6tq1sIA7v6PEyVbgGEu0al6xSBfaqsk9TlpAv4ETCjwjptwI9j6/QN9nmn8J2quVsI8BHKu+AVgQ/XuJ3VBjku03pdU0McN8fW2ZfQkPceQkNRZ4W6+wh3VT9AjV2Daojl0ZL6O2v9ro2WV+YB1PABtiK5Xjfh4Ns95bg/GKv/r1WWbzq5biLWW2PbPqSF28oyub6FGvuGJfxIPFvppES4mIqftN82FPukjm18Pvbv+VWV5ZtKrgmtVs+XrL+QKol1ybrtlCfMRWDzCsseFotzGbBVjduZHVu3Zck1IZl/tWRbzzVZ31An199p9HgAPh1b91uDLJuUXM8DPgy0t+rzibZ9XGy7f6yy/JAm14SW6NJllwKzatzOj2Lr3lBl+XhyfRewZQ3biX/vHXh3Cz+zmpJrQlIZP0d/qI7tZJ5cExp+lsTW+TL1PWfihHNrTefIKvGcF6t30GfXRttrNPa5hnAwxodBWzPlbfwx9r7ayB1Z+lPs/awsghgCX3D3R2tc9mSgdAirmYQRCpLEn96+0GscHWEIxY/HWS3e3pGUjzjyNXd/vJYVPQzF98OSImPgjIT9Toi9P83d7685yqGzBmGYun7VunXlhpmtShhZp9/9wBl1VHEu5f2A491jBnMV4U7cOd76Cb/yfh78WOz9d732EXW+QOg+0m8PM9u2jm3v7zWMoexhxJrLY8WVRpMZEmZ2FCER7D9HO6HF+tzsomrIaoQuXqVOob7nTCB0Qb3JzAYdEagG8XPYpk3WN6KMuj7XLine truncated
|
||||
}
|
||||
},
|
||||
"cell_type": "markdown",
|
||||
"id": "a67e0aad",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"<div>\n",
|
||||
"<img src=\"attachment:fig-asp-data-updated.png\" align=\"left\" width=\"250\"/>\n",
|
||||
"</div>"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"id": "d9d29899",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"It is clear from the code that at any given iteration `k`, the matrix is updated from row `k` to `n` and from column `k` to `m`. If we look at how that reflects the distribution of work over the processes, we can see that CPU 1 does not have any work, whereas CPU 2 does a little work and CPU 3 and 4 do a lot of work. "
|
||||
"Each CPU will be on charge of updating the entries stored locally. However, this algorithm is a bit different than the previous ones. At iteration $k$ only a subset of the entries in the matrix `B` are updated, i.e. only the entries in the bottom right block `B[k:end,k:end]` are updated. This block is illustrated in a darker color in the next figure for $k=6$. "
|
||||
]
|
||||
},
|
||||
{
|
||||
@@ -347,34 +443,16 @@
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"id": "6409890d",
|
||||
"id": "e341fda6",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"### Load imbalance\n",
|
||||
"In this particular example CPU 1 will not have any work to do at iteration $k=6$ and it will be waiting, while other CPUs are computing. In general, CPUs containing rows previous to row $k$ will have less work to do than CPUs containing rows after row $k$. This will worsen as the value of $k$ increases. At some point, only the last CPU will have work to do and the others will idle. \n",
|
||||
"\n",
|
||||
"The block-wise partitioning scheme leads to load imbalance across the processes: CPUs with rows $<k$ are idle during any iteration $k$. The bad load balance leads to bad speedups, as some CPUs are waiting instead of doing useful work. \n",
|
||||
" \n",
|
||||
"### Data dependencies\n",
|
||||
" \n",
|
||||
"<div class=\"alert alert-block alert-success\">\n",
|
||||
"<b>Question:</b> What are the data dependencies in the block-wise partitioning?\n",
|
||||
"<div class=\"alert alert-block alert-info\">\n",
|
||||
"<b>Definition:</b> *Load imbalance*: is the problem when work is not equally distributed over all processes and consequently some processes do more work than others.\n",
|
||||
"</div>\n",
|
||||
"\n",
|
||||
" a) CPUs with rows >k need all rows <=k\n",
|
||||
" b) CPUs with rows >k need part of row k\n",
|
||||
" c) All CPUs need row k \n",
|
||||
" d) CPUs with row k needs all rows >k \n"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": null,
|
||||
"id": "e0565e92",
|
||||
"metadata": {},
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"answer = \"x\" # replace x with a, b, c, or d \n",
|
||||
"ge_dep_check(answer)"
|
||||
"Having processors waiting for others is a waist of computational resources and affects negatively parallel speedups. The optimal speedup (speedup equal to the number of processors) assumes that the work is perfectly parallel and that it is evenly distributed. If there is load imbalance, the last assumption is not true anymore and the speedup will be suboptimal.\n"
|
||||
]
|
||||
},
|
||||
{
|
||||
@@ -382,9 +460,11 @@
|
||||
"id": "51498a44",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"### Cyclic partition\n",
|
||||
"### Fixing load imbalance\n",
|
||||
"\n",
|
||||
"In contrast, if we look at how the work is balanced for the same example and cyclic partitioning, we find that the processes have similar work load. "
|
||||
"In this application, is relatively easy to fix the load imbalance problem. We know in advance which data is going to be processes at each CPU and we can design a more clever data partition.\n",
|
||||
"\n",
|
||||
"We can consider row-wise cyclic partition to fix the problem. See figure below. In this case, the CPUs will have less work as the value of $k$ increases, but amount of work will be better distributed than with the previous row block partition."
|
||||
]
|
||||
},
|
||||
{
|
||||
@@ -404,58 +484,121 @@
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"id": "866824c6",
|
||||
"id": "1c9e4c8a",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"## Conclusion\n",
|
||||
"Cyclic partitioning tends to work well in problems with predictable load imbalance. It is a form of **static load balancing** which means using a pre-defined load schedule based on prior information about the algorithm (as opposed to **dynamic load balancing** which can schedule loads flexibly during runtime). The data dependencies are the same as for the 1d block partitioning.\n",
|
||||
"### Static vs dynamic load balancing\n",
|
||||
"\n",
|
||||
"At the same time, cyclic partitioning is not suitable for all communication patterns. For example, it can lead to a large communication overhead in the parallel Jacobi method, since the computation of each value depends on its neighbouring elements."
|
||||
"Load balancing is the process of distributing work to processors uniformly with the aim to efficiently exploit a parallel system. There are two key forms of load balancing: Static and dynamic load balancing.\n",
|
||||
"\n",
|
||||
"- **Static load balancing** the work distribution strategy is based on prior information of the algorithm and it does not depend on runtime values.\n",
|
||||
"- **Dynamic load balancing** the work distribution strategy is based on runtime values.\n",
|
||||
"\n",
|
||||
"Static load balancing is often used in algorithms for which the load distribution is known in advance and it does not depend on runtime values. On the other hand, dynamic load balancing is often needed in problems in which the work distribution cannot be predicted in advance and depends on runtime values.\n",
|
||||
"\n",
|
||||
"\n",
|
||||
"\n",
|
||||
"<div class=\"alert alert-block alert-success\">\n",
|
||||
"<b>Question:</b> Using a cyclic partition for Gaussian elimination, is a form of static or dynamic load balancing?\n",
|
||||
"</div>\n",
|
||||
"\n"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": null,
|
||||
"id": "a6741a25",
|
||||
"metadata": {},
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"ge_lb_answer()"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"id": "20982b04",
|
||||
"id": "e151fb8d",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"## Exercise\n",
|
||||
"The actual implementation of the parallel algorithm is left as an exercise. Implement both 1d block and 1d cyclic partitioning and compare their performance. The implementation is closely related to that of Floyd's algorithm. To test your algorithms, generate input matrices with the function below (a random matrix is not enough, we need a non singular matrix that does not require pivoting). "
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": null,
|
||||
"id": "a65cf8e6",
|
||||
"metadata": {},
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"function tridiagonal_matrix(n)\n",
|
||||
" C = zeros(n,n)\n",
|
||||
" stencil = [(-1,2,-1),(-1,0,1)]\n",
|
||||
" for i in 1:n\n",
|
||||
" for (coeff,o) in zip((-1,2,-1),(-1,0,1))\n",
|
||||
" j = i+o\n",
|
||||
" if j in 1:n\n",
|
||||
" C[i,j] = coeff\n",
|
||||
" end\n",
|
||||
" end\n",
|
||||
"### Data dependencies\n",
|
||||
"\n",
|
||||
"Using a cyclic partition, we managed to distribute the work uniformly. But we still need to study the data dependencies in order to implement this algorithm in parallel efficiently.\n",
|
||||
"\n",
|
||||
"Look again to the algorithm\n",
|
||||
"\n",
|
||||
"```julia\n",
|
||||
"n,m = size(B)\n",
|
||||
"for t in (k+1):m\n",
|
||||
" B[k,t] = B[k,t]/B[k,k]\n",
|
||||
"end\n",
|
||||
"B[k,k] = 1\n",
|
||||
"for i in (k+1):n \n",
|
||||
" for j in (k+1):m\n",
|
||||
" B[i,j] = B[i,j] - B[i,k]*B[k,j]\n",
|
||||
" end\n",
|
||||
" C\n",
|
||||
"end"
|
||||
" B[i,k] = 0\n",
|
||||
"end\n",
|
||||
"```\n",
|
||||
"\n",
|
||||
"Note that all updates on the loop over i and j we do the following: \n",
|
||||
"\n",
|
||||
"```julia\n",
|
||||
"B[i,j] = B[i,j] - B[i,k]*B[k,j]\n",
|
||||
"```\n",
|
||||
"\n",
|
||||
"As we are using row-wise partitions, the CPU that updates `B[i,j]` will also have entry `B[i,k]` in memory (both are in the same row). However, `B[k,j]` is in another row and it might be located on another processor. We might need to communicate `B[k,j]` for `j=(k+1):m`. This corresponds to the cells marked in red in the figure below. These red entries are the data dependencies of this algorithm. The owner of these entries has to send them to the other processors. This is very similar to the communication pattern seen in previous notebook for Floyd's algorithm. There is a key difference however. In the current case, we do not need to send the full row, only the entries beyond column $k$ (the red cells in the figure).\n"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": null,
|
||||
"id": "31d8586a",
|
||||
"attachments": {
|
||||
"g27851.png": {
|
||||
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAzAAAAJsCAYAAAAiM0YjAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAAewgAAHsIBbtB1PgAAABl0RVh0U29mdHdhcmUAd3d3Lmlua3NjYXBlLm9yZ5vuPBoAACAASURBVHic7N15nBx1nf/xV/X0zOS+LxIIwYRDDkFuFBFUUFbE+wQPFBRXXQVcL1xss3ggrgfeeKKuq7grIh4/BeUSOQxBbgEFQgIkAZJJyDVX1++PmoGeTnfPVKe7q2rm9Xw8+pGu6u936t09k5n+9Lfq+w2QJI12eeCsku2/A5cO3O8EXgG8GdgD2BkoAvcDvwMuAB6u8nWfBZwKPB94xsDXWgvcBPwQ+GUjn4QkSZKksWE8EJbcfjqwf1/g7rLHym8bgdeUfb1xwNeB/mH6/hqY0KTnJEmSJGmUqlTAHAxspnYBMnjrA44b+FrjgD+OsF8I/BYImvz8JEmSJI0i5QXMH4AHBu5vA75NNMpyBHAM8FFgdVmfFQNf57sl+24hOjXtRcBhwOuITk0rL2JObv5TlCRJkjRalBcwg7d7gd2r9Jk78Hhp+x8P/FsE/h3IVen74bJ+tzTiSUiSJEkaGyoVMOuAXYbpd3yFfiGwdJh+AXBzWZ8968wuSZIkaYypVMB8YAT9cmx/KtlKoH0EfT9Q1u9NsVNLklRBteF/SdLo1Q/8YATtisDysn0/AnpH0PevZduOwEiSGsICRpLGnjuArhG2XVW2fV2d/aaNsJ8kSTVZwEjS2FNtYcpKNpVtPzLCfpvLtifHOKYkSVVZwEjS2FNelNRSrLNveT/XgpEkNYQFjCSNPWHSASRJqpcFjCRJkqTMsICRJEmSlBkWMJIkSZIywwJGkiRJUmZYwEiSJEnKDAsYSZIkSZlhASNJkiQpMyxgJEmSJGWGBYwkSZKkzLCAkSRJkpQZ+aQDSJKaLgTWl2xvjtF3a1nf/jqPuSXGMSVJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkqTWGA+cDFwGtCecRZIkxTcTuBJ4M9HfdUkadQLgKOASYDWwDfjfRBNJkqQdcS3R3/M1wKXA84FcookkqQGWAF8AVgKPAeHAbQ2wc4K5JEnSjjmQoX/b1wGPAl8B9kgwlyTFNh14L3An0WhLH0//chu83ZFYOkmS1CgPsf3f+CLwBHAv8AFgRmLpNOoFSQdQprUDxwPnEH3qMo7a17f8GfhRC3JJkqTm+QDwzBqP9wFbgduBzwG/BXpbkEtjhAWM6nU48AeigmVcwlkkSVI69RBdM3MocE/CWTRKWMCoXjngGOBM4LnAJKBtmD73E30aI0mSsuswYN4wbXqB9cBvgAuAvzU7lCTFMQ54PfAXYC3RL63yc2ND4LqkAkqSpIa5j8p/57cSXdB/PfA6XDZBUkbsDCwFHmToLCUh0cX9ExNLJkmSdtR8or/ng3/b+4lmGb0DeDcwNblokrTjDia6aH81sIHoHNgPJ5pIkiTtiG8SFS2PEy2XcD7wjEQTSVITtAMvB64mmmLRIWVJkrJnJtF1LT8nuv7Va6kljQmdSQeQJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSJEmSpOFMB9aV3L6WbBxJktIhn3QASVJFAVERM2hCUkEkSUqTXNIBJEmSJGmkLGAkSZIkZYYFjCRJkqTMsICRJEmSlBkWMJIkSZIywwJGkiRJUmY4jbIkjT47A3PL9t0BdCeQRZIkSdIYMAMIS27fH2G/Y4ENJf2KQKEJ+SRJkiTpKfUUMG8Dekr6dANvalI+SZIkSXpKnAImIBplKW2/Hji6qQklSZIkacBIC5g8cGFZ24eB/VuQUZIkSZKAkRUwk4DflrW7HdilRRklSZIkCRi+gNkJWF7W5gpgagszSpIkSRJQu4DZB1hR9vhFQHuLM0qSJEkSUL2AOYboAv3Sx75MdCG/JEmSJCWiUgHzWmBryb4+4N1JBZQkSZKkQeUFzKNEi1IObm8CXppYOkmSJEkqUV7AlN8uBzoSSydJkiRJJcoLmOuIRmFK9/0GGJdUQEmSJEkaVOkamN2A+8v2XwVMTiaiJEmSJEWqzUK2ELiv7LFrgSkJZJQkSZIkoPY6MPOA28seXwbMbHFGSZIkSQJqFzAAc4C/lbW5E9iphRklSZIkCRi+gAGYDtxQ1u5uYEGLMkqSJEkSMLICBmAq8OeytvcTXfAvSZIkSS0x0gIGYALRujCl7VcAS5qcUZIkSZKAeAUMQCfwy7I+jwL7NjGjJEmSJAHxCxiADuB/y/qtAfZvUkZJkiRJAuorYADagB+U9V0HHNb4iJIktV6QdABJUkU5YFHJ9pPAYyPs2wbsWrZvC7B6x2NJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJkiRJktKrHZiRdAhJkiRJ6ZRPOkCZdwGvAF6UdBCNVLgbsEfSKSRJUqr9GYLNSYfQ6BAkHaDMXcBcYAmwPuEsGpHww8Bnk04hSZJSbR8I7ko6hEaHXNIBSiwBZgETgZMTzqJhheMhnAK/PCbpJJIkKe0+enz0viFsTzqJsi9NIzBfAM4YuP8A8IwEs2hY4aXAiUmnkCRJmXIaBN9JOoSyLS0jMAHw2pLtqXhdhSRJkqQyabmI/3nA+JLt6UQX9J+VTBzVYTnwlqRDbO/A38DCi+GXFyWdpLIXfRQ27gs3nZR0ksq+sRC+/gu4/eCkk1T3rBvh9W+Fs/+edJLKDvsBTHgQriwkHKSK170O7nsH3PLipJNUt9+tcPv+Saeo7oAr4Jlfg/+5JOkklT3/09A9G244LekklZ29P1z2Lbjt8KSTVLffzfAfx8Pr1iadpLJDfgrTl8EfPp90kjKTgeuTDiE1yyVAWHZbRbpOcdMQ4aUQhiW3a5NOVMW9wJlJh6jhAtL9y30JkPZZYzYCByQdooargAuTDlHD6USn7aZZb9IBhvEQcErSIWq4CLg86RA1HA50JR1iGFuBBUmHqOEm4PykQ2wvnFb2XiGE8NSkUyn70nAK2TjgiAr7JwBHtjiLJEmSpBRLQwHzGqKZx8pNBz7Q4iySJEmSUiwNBcz7gUlVHnsu0UiMJEmSJCVewMwFdqnx+CTglS3KIkmSJCnlkp6F7FSiU8WqmUg0QvPfrYmjhM0F9gbaBrbvAVYmF0eSJElpk3QBcxrRKNAaopmOdgY6gPuJipdJwO7APGB1QhnVPM8FjgEOBg4i+v6Xeh/w1VaHkiRJUnoleQrZvkSFyReB/YHFwGMDj70W2BV4L9APpHSNDO2gTwD/Cbyc7YsXSZIkaTtJjsDcD+wGPFrl8W7gB8DFwNQWZVIyeoDbgJuJCteXJBtHkiRJaZXkCMwWqhcv9bRT9nye6NSxScAhRAvq3ZBoIkmSJKVa0tfAaGz7Q9IBJEmSlC1JT6MsSZIkSSNmASNJkiQpMyxgJEmSJGWGBYwkSZKkzLCAkSRJkpQZFjCSJEmSMsMCRpIkSVJmWMBIkiRJygwLGEmSJEmZYQEjSZIkKTMsYCRJkiRlhgWMJEmSpMywgJEkSZKUGRYwkiRJkjLDAkZJejkQlt0KJY9/pcLjR7Y2oiRJktLEAkaSJElSZuSTDqAx7S7gIzH7PNiEHJIkScoICxgl6T7gvKRDSJIkKTs8hUySJElSZljASJIkScoMCxhJkiRJmWEBI0mSJCkzLGAkSZIkZYYFjCRJkqTMsICRJEmSlBmuA6MGWTUPOD3pFBXMAJ4PbEk6SBXPBhaQztcOYC7QTnrzAXQCrwcOTzpIFQuBcaT3NXwBMJX05oPow7Y055sCHEv0s5hG+5Du7/ESotcurfkger/0VmBd0kGq2Ak4hNS9hruPj5Z8kxorSDpAmVVEb+YOApYnnEU1hZcCJz69fdM2OGx1YnGq2wXYBKxPOkgVs4j+cD+cdJAq2oH5wIqkg9SwK7Aa6E46SBXzgT5gbdJBqpg6cHso6SA1LAIeTDhDLQuBLmBj0kGqmAu0AY8kHaSKcUQZ0/x7ZhHRe5S+hHNUs4Dod+DjSQcZakYOnlhYtvM0CL6TSBypSVYBIXBg0kE0nPBSCMOS27VJJ6riXuDMpEPUcAFwfdIhalgCbE46xDA2AgckHaKGq4ALkw5Rw+nAA0mHGEZv0gGG8RBwStIhargIuDzpEDUcTlQAptlWoiIhrW4Czk86xPbCaWXvFUIIT006lbLPa2AkSZIkZYYFjCRJkqTMsICRJEmSlBkWMJIkSZIyw2mUpeHtDOzB0wX/DUQzm0mSJKnFLGCk7b2UaD79g4mm9J5X9vgBwK2tDiVJkiQLGKmS/yVal0CSJEkp4zUwUnVbiU4X+xouJSxJkpQKjsBI23sX8DfgLp5edflyYPfEEkmSJAmwgJEq+WHSASRJklSZp5BJkiRJygwLGEmSJEmZYQEjSZIkKTMsYCRJkiRlhgWMJEmSpMywgJEkSZKUGRYwkiRJkjLDAkaSJElSZljASJIkScoMCxhJkiRJmWEBI0mSJCkzLGAkSZIkZYYFjCRJkqTMsICRJEmSlBn5pANIKbQUeG/Zvkkl968B+ku2Q2Bms0NJkiTJAkaqZDwwvcbjU8q2wyZmkSRJUgkLGGl7VzJ0hGU4FjCSJEktYgEjbe+3AzdJkiSljBfxS5IkScoMCxhJkiRJmWEBI0mSJCkzLGAkSZIkZYYFLine truncated
|
||||
}
|
||||
},
|
||||
"cell_type": "markdown",
|
||||
"id": "0050b6c4",
|
||||
"metadata": {},
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"n = 12\n",
|
||||
"C = tridiagonal_matrix(n)\n",
|
||||
"b = ones(n)\n",
|
||||
"B = [C b]\n",
|
||||
"gaussian_elimination!(B)"
|
||||
"<div>\n",
|
||||
"<img src=\"attachment:g27851.png\" align=\"left\" width=\"300\"/>\n",
|
||||
"</div>"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"id": "cf189776",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"### Parallel implementation\n",
|
||||
"\n",
|
||||
"\n",
|
||||
"At iteration $k$,\n",
|
||||
"\n",
|
||||
"1. The CPU owning row $k$ does the loop over $t$ to update row $k$.\n",
|
||||
"2. The CPU owning row $k$ sets $B_{kk} = 1$.\n",
|
||||
"2. This CPU sends the red cells in figure above to the other processors.\n",
|
||||
"3. All processors receive the updated values in row $k$ and do the loop over i and j locally (blue cells).\n",
|
||||
"\n",
|
||||
"\n",
|
||||
"As you probably see, the parallel implementation of this method is closely related to Floyd's algorithm. But there are some differences. \n",
|
||||
"\n",
|
||||
"1. The process that owns row $k$ updates its values before sending them.\n",
|
||||
"2. We do not send the full row $k$, only the entries beyond column $k$.\n",
|
||||
"3. We need a cyclic partition to balance the load properly.\n",
|
||||
"\n",
|
||||
"A key similarity between the two algorithms, however, is that they both suffer from synchronization problems. We need to make sure that the rows arrive in the right order. The strategies discussed for Floyd's algorithm also apply in this current case.\n",
|
||||
"\n",
|
||||
"The actual implementation of the parallel algorithm is left as an open exercise."
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"id": "866824c6",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"## Conclusion\n",
|
||||
"\n",
|
||||
"\n",
|
||||
"We studied the parallelization of an algorithm that leads to load imbalance. We fixed the problem using a cyclic partition. This is a form of static load balancing since we were able to distribute the load in advance without using runtime values. This is opposed to dynamic load balancing which can schedule loads flexibly during runtime. Note however that cyclic partitioning is not suitable for all communication patterns. For example, it can lead to a large communication overhead in the parallel Jacobi method, since the computation of each value depends on its neighboring elements.\n"
|
||||
]
|
||||
},
|
||||
{
|
||||
@@ -469,14 +612,6 @@
|
||||
"\n",
|
||||
"This notebook is part of the course [Programming Large Scale Parallel Systems](https://www.francescverdugo.com/XM_40017) at Vrije Universiteit Amsterdam and may be used under a [CC BY 4.0](https://creativecommons.org/licenses/by/4.0/) license."
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": null,
|
||||
"id": "8ab22f67",
|
||||
"metadata": {},
|
||||
"outputs": [],
|
||||
"source": []
|
||||
}
|
||||
],
|
||||
"metadata": {
|
||||
|
||||
+355
-275
@@ -38,18 +38,10 @@
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": 17,
|
||||
"execution_count": null,
|
||||
"id": "1dc78750",
|
||||
"metadata": {},
|
||||
"outputs": [
|
||||
{
|
||||
"name": "stdout",
|
||||
"output_type": "stream",
|
||||
"text": [
|
||||
"🥳 Well done!\n"
|
||||
]
|
||||
}
|
||||
],
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"using Printf\n",
|
||||
"\n",
|
||||
@@ -65,11 +57,19 @@
|
||||
"function q1_answer(bool)\n",
|
||||
" bool || return\n",
|
||||
" msg = \"\"\"\n",
|
||||
" It is not needed to start the loop over `j` with `j=k` as none of the entries of the matrix will be updated\n",
|
||||
" in this particular case. Rememebr: `C[i,j] = min(C[i,j],C[i,k]+C[k,j])` if we substitute `j=k`, we get:\n",
|
||||
" `C[i,k] = min(C[i,k],C[i,k]+C[k,k])`. Rememeber that `C[k,k]=0, thus, `C[i,k] = min(C[i,k],C[i,k])`, or\n",
|
||||
" `C[i,k] = C[i,k]`. In other words, the new value of `C[i,k]` will correspond to the old value.\n",
|
||||
" The same is true for `i=k`.\n",
|
||||
" The we can change the loop order over i and j without changing the result. Rememebr:\n",
|
||||
" \n",
|
||||
" C[i,j] = min(C[i,j],C[i,k]+C[k,j])\n",
|
||||
" \n",
|
||||
" if we substitute j=k, we get\n",
|
||||
" \n",
|
||||
" C[i,k] = min(C[i,k],C[i,k]+C[k,k]).\n",
|
||||
" \n",
|
||||
" Since C[k,k]=0, thus, C[i,k] = min(C[i,k],C[i,k]), and C[i,k] = C[i,k].\n",
|
||||
" \n",
|
||||
" In other words, the value of C[i,k] will not be updated at iteration k.\n",
|
||||
" \n",
|
||||
" The same is true for i=k.\n",
|
||||
" \"\"\"\n",
|
||||
" println(msg)\n",
|
||||
"end\n",
|
||||
@@ -97,8 +97,8 @@
|
||||
"id": "1faddbfa",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"We represent the distance table as a matrix $C$, where $C_{ij}$ is the distance from node $i$ to node $j$ via a direct connection (a single hop in the graph). If there is no direct connection from $i$ to $j$, this is represented using a very large value in $C_{ij}$ (e.g. the largest possible floating point number, `inf`). \n",
|
||||
"The next figure shows a simple directed graph with 4 nodes an its corresponding distance matrix (matrix labeled as \"input\")."
|
||||
"We represent the distance table as a matrix $C$, where $C_{ij}$ is the distance from node $i$ to node $j$ via a direct connection (a single hop in the graph). If there is no direct connection from $i$ to $j$, this is represented using a large value in $C_{ij}$ representing infinity. \n",
|
||||
"The next figure shows a simple directed graph with 4 nodes an its corresponding distance matrix (labeled as \"input\")."
|
||||
]
|
||||
},
|
||||
{
|
||||
@@ -121,7 +121,7 @@
|
||||
"id": "ade31d26",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"The ASP problem consists in computing the minimum distance between any two pair of nodes $i$ and $j$ in the graph. All these values can be also represented as a matrix (matrix labeled as \"output\" in the figure above). For instance, the minimum distance from node 2 to node 3 is 8 as highlighted in the figure."
|
||||
"The ASP problem consists in computing the minimum distance between any two pair of nodes $i$ and $j$ in the graph. All these values can be also represented as a matrix (labeled as \"output\" in the figure above). For instance, the minimum distance from node 2 to node 3 is 8 as highlighted in the figure. You can understand both input and output matrices as distance tables. The key difference is that the input contains the distance using direct connections only, whereas the output contains the (minimum) distance allowing indirect connections."
|
||||
]
|
||||
},
|
||||
{
|
||||
@@ -136,21 +136,10 @@
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": 4,
|
||||
"execution_count": null,
|
||||
"id": "4fe447c5",
|
||||
"metadata": {},
|
||||
"outputs": [
|
||||
{
|
||||
"data": {
|
||||
"text/plain": [
|
||||
"floyd! (generic function with 1 method)"
|
||||
]
|
||||
},
|
||||
"execution_count": 4,
|
||||
"metadata": {},
|
||||
"output_type": "execute_result"
|
||||
}
|
||||
],
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"function floyd!(C)\n",
|
||||
" n = size(C,1)\n",
|
||||
@@ -176,25 +165,10 @@
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": 6,
|
||||
"execution_count": null,
|
||||
"id": "860e537c",
|
||||
"metadata": {},
|
||||
"outputs": [
|
||||
{
|
||||
"data": {
|
||||
"text/plain": [
|
||||
"4×4 Matrix{Int64}:\n",
|
||||
" 0 9 6 1\n",
|
||||
" 2 0 8 3\n",
|
||||
" 5 3 0 6\n",
|
||||
" 10 8 5 0"
|
||||
]
|
||||
},
|
||||
"execution_count": 6,
|
||||
"metadata": {},
|
||||
"output_type": "execute_result"
|
||||
}
|
||||
],
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"inf = 1000\n",
|
||||
"C = [\n",
|
||||
@@ -228,7 +202,7 @@
|
||||
"source": [
|
||||
"### The algorithm explained\n",
|
||||
"\n",
|
||||
"The main idea of the algorithm is to perform as many iterations as nodes in the graph. At iteration $k$, we update the distance matrix $C$ by finding the shortest paths between each pair of nodes, allowing indirect paths via nodes from 1 to $k$. At the last iteration, it is allowed to visit all nodes and, thus, the distance table will contain the minimum possible distances, i.e. the solution of the ASP problem.\n",
|
||||
"The main idea of the Floyd–Warshall algorithm is to perform as many iterations as nodes in the graph. At iteration $k$, we update the distance matrix $C$ by finding the shortest paths between each pair of nodes, allowing indirect paths via nodes from 1 to $k$. At the last iteration, it is allowed to visit all nodes and, thus, the distance table will contain the minimum possible distances, i.e. the solution of the ASP problem.\n",
|
||||
"\n",
|
||||
"This process is cleverly done with three nested loops:\n",
|
||||
"\n",
|
||||
@@ -243,7 +217,7 @@
|
||||
"end\n",
|
||||
"```\n",
|
||||
"\n",
|
||||
"At each outer iteration $k$, we do a loop over the distance matrix $C$. For each pair of nodes $i$ and $j$ we compare the current distance $C_{ij}$ against the distance via node $k$, namely $C_{ik}+C_{kj}$, and update $C_{ij}$ with the minimum."
|
||||
"At each outer iteration $k$, we do a loop over the distance matrix $C$. For each pair of nodes $i$ and $j$ we compare the current distance $C_{ij}$ against the distance via node $k$, namely $C_{ik}+C_{kj}$, and update $C_{ij}$ with the minimum. I.e., at iteration $k$ one checks if it is beneficial to visit node $k$ to reduce the distance between nodes $i$ and $j$."
|
||||
]
|
||||
},
|
||||
{
|
||||
@@ -266,7 +240,7 @@
|
||||
"id": "722e330c",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"The update of the distance matrix at each iteration is illustrated in the next figure for a small ASP presented above. We highlight in green the distances that are updated in each iteration."
|
||||
"The update of the distance matrix at each iteration is illustrated in the next figure for the small ASP presented above. We highlight in green the distances that are updated in each iteration. Note that some distances that were initially infinity (i.e., no connection) will be updated with a finite value in this process. You can understand this as adding new edges in the graph as illustrated in the figures below."
|
||||
]
|
||||
},
|
||||
{
|
||||
@@ -352,28 +326,17 @@
|
||||
"source": [
|
||||
"### Serial performance\n",
|
||||
"\n",
|
||||
"This algorithm is memory bound, meaning that the main cost is in getting and setting data from the input matrix `C`. In this situation, the order in which we traverse the entries of matrix `C` has a significant performance impact.\n",
|
||||
"Before starting to parallelize this code, we want to make sure that we have an efficient sequential implementation. In this algorithm, the order in which we traverse the entries of matrix `C` has a significant performance impact.\n",
|
||||
"\n",
|
||||
"The following function computes the same result as for the previous function `floyd!`, but the nesting of loops over i and j is changed.\n"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": 8,
|
||||
"execution_count": null,
|
||||
"id": "75cac17e",
|
||||
"metadata": {},
|
||||
"outputs": [
|
||||
{
|
||||
"data": {
|
||||
"text/plain": [
|
||||
"floyd2! (generic function with 1 method)"
|
||||
]
|
||||
},
|
||||
"execution_count": 8,
|
||||
"metadata": {},
|
||||
"output_type": "execute_result"
|
||||
}
|
||||
],
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"function floyd2!(C)\n",
|
||||
" n = size(C,1)\n",
|
||||
@@ -399,19 +362,10 @@
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": 10,
|
||||
"execution_count": null,
|
||||
"id": "907bc8c9",
|
||||
"metadata": {},
|
||||
"outputs": [
|
||||
{
|
||||
"name": "stdout",
|
||||
"output_type": "stream",
|
||||
"text": [
|
||||
" 1.777208 seconds\n",
|
||||
" 2.967238 seconds\n"
|
||||
]
|
||||
}
|
||||
],
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"n = 1000\n",
|
||||
"C = rand(n,n)\n",
|
||||
@@ -457,32 +411,19 @@
|
||||
"\n",
|
||||
"\n",
|
||||
"<div class=\"alert alert-block alert-success\">\n",
|
||||
"<b>Question (hard):</b> Can we really parallelize the loops over `i` and `j` ? To compute `C[i,j]` at iteration `k`, we first need to compute `C[i,k]` and `C[k,j]`. In order words, it seems that the order of the loops over `i` and `j` cannot be arbitrary. The first value of `i` and `j` in the loops should be `k`. However, this is not really necessary, why?\n",
|
||||
"<b>Question (hard):</b> Can we really parallelize the loops over `i` and `j` ? To compute `C[i,j]` at iteration `k`, we first need to compute `C[i,k]` and `C[k,j]`. In order words, it seems that the order of the loops over `i` and `j` cannot be arbitrary. However, this is not really true, why?\n",
|
||||
"</div>"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": 18,
|
||||
"execution_count": null,
|
||||
"id": "8d05b686",
|
||||
"metadata": {},
|
||||
"outputs": [
|
||||
{
|
||||
"name": "stdout",
|
||||
"output_type": "stream",
|
||||
"text": [
|
||||
"It is not needed to start the loop over `j` with `j=k` as none of the entries of the matrix will be updated\n",
|
||||
"in this particular case. Rememebr: `C[i,j] = min(C[i,j],C[i,k]+C[k,j])` if we substitute `j=k`, we get:\n",
|
||||
"`C[i,k] = min(C[i,k],C[i,k]+C[k,k])`. Rememeber that `C[k,k]=0, thus, `C[i,k] = min(C[i,k],C[i,k])`, or\n",
|
||||
"`C[i,k] = C[i,k]`. In other words, the new value of `C[i,k]` will correspond to the old value.\n",
|
||||
"The same is true for `i=k`.\n",
|
||||
"\n"
|
||||
]
|
||||
}
|
||||
],
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"uncover = true\n",
|
||||
"q1_answer(uncover)"
|
||||
"uncover = false\n",
|
||||
"q1_answer(true)"
|
||||
]
|
||||
},
|
||||
{
|
||||
@@ -616,7 +557,7 @@
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"<div class=\"alert alert-block alert-success\">\n",
|
||||
"<b>Question:</b> How much data is communicated in each iteration in this parallel algorithm?\n",
|
||||
"<b>Question:</b> How much data is send from the owner of row k in each iteration in this parallel algorithm?\n",
|
||||
"</div>\n",
|
||||
"\n",
|
||||
" a) O(N²/P)\n",
|
||||
@@ -627,18 +568,10 @@
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": 25,
|
||||
"execution_count": null,
|
||||
"id": "1bf4de56",
|
||||
"metadata": {},
|
||||
"outputs": [
|
||||
{
|
||||
"name": "stdout",
|
||||
"output_type": "stream",
|
||||
"text": [
|
||||
"It's not correct. Keep trying! 💪\n"
|
||||
]
|
||||
}
|
||||
],
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"answer = \"x\" # replace x with a, b, c or d\n",
|
||||
"floyd_check(answer)"
|
||||
@@ -659,7 +592,7 @@
|
||||
"source": [
|
||||
"### Computation complexity\n",
|
||||
"\n",
|
||||
"Each process updates $N^2/P$ entries per iteration. The computation complexity is $O(N^2/P)$."
|
||||
"Each process updates $N^2/P$ entries per iteration. The computation complexity per iteration is $O(N^2/P)$."
|
||||
]
|
||||
},
|
||||
{
|
||||
@@ -699,7 +632,7 @@
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"<div>\n",
|
||||
"<img src=\"attachment:fig-asp-efficiency-comm-2.png\" align=\"left\" width=\"300\"/>\n",
|
||||
"<img src=\"attachment:fig-asp-efficiency-comm-2.png\" align=\"left\" width=\"250\"/>\n",
|
||||
"</div>"
|
||||
]
|
||||
},
|
||||
@@ -714,7 +647,7 @@
|
||||
"- On the receive side $O(N)/O(N^2/P) = O(P/N)$\n",
|
||||
"\n",
|
||||
"\n",
|
||||
"In summary, the send/computation ratio is $O(P^2/N)$ and the receive/computation ratio is $O(P/N)$. The algorithm is potentially scalable if $P^2<<N$."
|
||||
"In summary, the send/computation ratio is $O(P^2/N)$ and the receive/computation ratio is $O(P/N)$. The algorithm is potentially scalable if $P^2<<N$. Note that this is worse than for matrix-matrix multiplication, which is scalable for $P<<N$. I.e., you need a larger problem size in the current algorithm than in matrix-matrix multiplication."
|
||||
]
|
||||
},
|
||||
{
|
||||
@@ -749,159 +682,339 @@
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"id": "680e56cf",
|
||||
"id": "4c2e0010",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"### Code\n",
|
||||
"We are going to implement the method in a function with the following signature\n",
|
||||
"\n",
|
||||
"We split the code into two functions. The first function is called on the main process (the process running this notebook). It splits the input matrix into blocks of rows. Then, we use a remotecall to compute Floyd's algorithm in each worker with its corresponding block of rows.\n"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": null,
|
||||
"id": "22234b2f",
|
||||
"metadata": {},
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"function floyd_dist!(C)\n",
|
||||
" m = size(C,1)\n",
|
||||
" @assert mod(m,nworkers()) == 0\n",
|
||||
" nrows_w = div(m,nworkers())\n",
|
||||
" @sync for (i,w) in enumerate(workers())\n",
|
||||
" rows_w = (1:nrows_w) .+ (i-1)*nrows_w\n",
|
||||
" Cw = C[rows_w,:]\n",
|
||||
" @async C[rows_w,:] = remotecall_fetch(floyd_worker!,w,Cw,rows_w)\n",
|
||||
" end\n",
|
||||
" C\n",
|
||||
"end"
|
||||
"```julia\n",
|
||||
" floyd_mpi!(C,comm)\n",
|
||||
"```\n",
|
||||
"\n",
|
||||
"The signature is similar to the sequential function `floyd!(C)`, but there are some important differences\n",
|
||||
"\n",
|
||||
"- This function is going to be called on several MPI ranks, but only one rank (rank 0 in the implementation below) will receive the input matrix, the other ranks will receive an empty matrix. It would not be memory efficient if all ranks receive a copy of the input matrix as this one could be very large.\n",
|
||||
"\n",
|
||||
"- On output, only the value of C will be correct on rank 0. It will not be efficient to try to write the result on all ranks.\n",
|
||||
"\n",
|
||||
"- The parallel function takes an MPI communicator object. This allows the user to decide which communicator to use. For instance, `MPI.COMM_WORLD` directly or a duplicate of this one (which is the recommended approach)."
|
||||
]
|
||||
},
|
||||
{
|
||||
"attachments": {
|
||||
"g3642.png": {
|
||||
"image/png": "iVBORw0KGgoAAAANSUhEUgAABRQAAAMuCAYAAABo148ZAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAAewgAAHsIBbtB1PgAAABl0RVh0U29mdHdhcmUAd3d3Lmlua3NjYXBlLm9yZ5vuPBoAACAASURBVHic7N15fFTV/f/xV0IChMiusijIorK4oagoQusC1rbgtyp8tdaqSLUFQRCqBCWTZMKqYLW41VYK9VurFtpfv9J+bWMVC9a6gLsiqICyibIq+zK/Pz6MTiaT5NzJ3Nnyfj4e82i5nHPvBySfuffccz4HRDLLQ0Ao4jM3pdHUTWOgC9ARaJjiWEQkOeZSOYc9mNJoRETi979UzmdTUhtOXJoDnYEeQDugILXhiIjUC9Oo/P3x/1IbjsQrL9UBiNQjucAPgB8D/YAjo35/PbAImAP8M6mRiYi4yQdOBs4Eukb93h6gNNkBiYg4ygPOAb5z+H97UfVeDGAF8DLwZ2zQ9GCyAhQRkbRyPPATD+2/BLYD7wOvHf7/6aAB0B27f++GjUtEKgN2JzsokVTI1BmKvYBlVI69ps+zwHEpiVRE/DSXzJmh2AA4CbgWuA9YAuyi+ry1NTVhikiKZMIMxRzsJe4jwBbc78PCnzXAVUmPWkQku2XKDMUBeP/eCH8OAf8GrgMaJTHmXGzW/TXAvcBi4KtaYm0R78U0Q1HEf98F5gNNPPS5CHgVuBh4w4+gRERq8TLQO9VBiIjUwY+Ax+rQvyPwB+BS4HpgXwJiEhGR7JcDnHv483NslWIynuv/CFyehOsAGlAU8duZwAKq1uTZBfwV+Aj7OewGXIItJww7Cvg79kC/1vdIRUQq0z2CiGS6BjX83kFgFbAJuy9rDvQECmO0/SFW7/pKtARaRES8ORmbrXgxtuLHT0m9f9fDgoh/GmJvtaMHEx8DxgFfRB0/Fvg1NrAYdjS2RHKAPyGKiHh2AHgXWx6h0gwikkk2AY9jL2xfxOpdRcoHvofVg+0V9XtXACOB2f6GKCIiaewO4B8xjjcFWmOTgYZi9RcjFQB/wb5bPvUzwGocAj7AlkR3S8H1RdJCJtVQHEXV+gQPYdOfq5OPJZroft/xNVIRSZa5ZE4NxTeA/djg4e+AMVgtsnD5ht+iGooi9Vkm1FC8DovtZeAyKq8EqUkjLO/FqhXbPPFhiojUK5lcQ9Glrm4ucAuwN0b/3/gS6Tf+gs2kfx+byDQW6I8NeIK9MEtYDUWRTJMpA4qFwAYqx7oStzqKbbDZi5F9X6PmgUgRyQxzyZwBxZ5UnWEdSQOKIvVbJgwo9gcGxdk3D1uaFv3gdXViQhMRqbeyfUAx7IYY/fcDRyY2zEpO4JvBw1hKY8QU94Bi9HbRIpIYg4G2UccmYzV6avMZ8MuoY72BMxIQl4iIq/eA3akOQkSkDhYDC+PsewBb2hbtkhjHREREos3BaidGysM2YPXLSqqW8/CNaiimTg5VR4J3ULnQcz5wHtAH6AS0PNzmXeC+Ws5/JDYI1RPoDByBLd/YevjzFnaTtTHO+POoOvIdPTulAPth6QMcg9UU/AJYAzyDTcNNheZUHUzfQ2IfnC+N+vU2bMclV3OAYir/jP4XsLSOcYkkSi5Vl31tx+pzhDXEZoecxTc5bBvwJrXPzDsKy2EnHe57xOHzbcFyzZtYDtsUZ/zROSx0OLZITbA3k32A9of7fIEV8X8GWBHntesqVg7bjeUxEfGuAdAs6tg2LC+ENQa+jW221hHLZ5uxpbRzazl/O765J+uI5bM8vslnS7F6fpvjjD//8DnDDmH5OFJTrBj8mYfjycXy2YdYPvs4zmvXVQuqrsDYhS0TSwdLsP9GLSOOdUxRLCKSnlpG/fpL7IVEWAPgHKAv9lzeEstzK4HptZy7BXYf3QPoiuX6Ar55pn8PeAFYF2fsLt9/DYELsN2KO2Dfh19gdQD/gY0rpEJTqo5n7cVtAk8y/Q/23z7SOcCTKYhFskgrqk41DRfHzANuBdbHaBPCBuRi6Q3cBbyN3UzG6hv9+RdWfNqrPjHOFX7AzQcmUnXZbvTnxcPn8aKuS56LqPp3s5SqswnrIh9L8JHX8DKYGPafqHOkKlmLxHIMVX+mOxz+vYbABGy2bayf/epeJpwNzMJemrjmsH8CA+OI//yo80S+UGgEBKj6cxz9WYT3mcNzo87hZclzDrZMIfrv5j/YAGyyacmzZItTqPrzHR6gK8B+7jbHaBOi6swDsPuhftjmHR9U0y/6cwj46+F+Xl0eda7IB8tC7IH1y1qu/3/YgKcXdVnynAv8IkYcz5J+NQrfoHKMb6c2HBFJI7lUzWN9In5vBPbsHivvRr/4CTsFy6fLsMlGLt8h4TqxXlegdotxrlaHf68BVj97Yy3XXoq9cPOirkueR2CDtpHneI/EbxZY1yXPEHvc5PcJjNGr0hjxqIZiBqpuQPEo7Oa0ph/aWLsC/bSWPrV9fk/NtbKiVTeg2AZLaK7X3Ye3H8p4BxQbYA/usW6gj6ihXzxOjXGdW+M4z8yocxyi5noIIslU3YBiO6zmZ00/9x/EON/YWvrU9nkUGwh0dX5U//CA4jHA6x6uuwf4gYfrzo3q7zqgmI/NXI6+/l9wq83qh99GxaIBRclU1Q0odsYeUGrKAS/FOF9ZLX1q+/wCb6uIqhtQ7Aos93DdL4ELPVw33gHFxsD8GNd/DHshlW4+onKc/0ptOCKSRqobUGyOvSCpLedGu7KWPrV9/oK358XqBhRbAs97uO5B4EYP1413QDEnRt9wXm5VQ794JWJA8YQY5/hbAmP0qjRGPKqhmCUKsQGuc2tpF2tzjpp2rdvPN8tqqnM19ma8QS3XrknB4XOc7aFPPjYN+Jw6XLc2TYA/Y28yIs3Bah1+leDrdY9x7NU4zhPdJ6eac4uki2bY0ofetbTzmsP2YTkseklypBuwn/O6fK8dgeXgXh76NAKe8tjHq6bA08CwqOMPYgMJ6ba0QyQbHI3NgO5RSzuv+Wwvls+qm5kC9oJlXi3XrU1roIJvVr+4OALLNV3reO2atMYesq+IOj4NuBbL9+nkCKoucf4kFYGISMZoiL1wqa1OntfvjwN880wfqqbNpdh3l5eX7LFi+DP28t1VLvArrLSGXxpi4wZFUcf/ePi6W3y8dl3EmryU6PGHlNGAYnqZyjcP4gewH5ghwMnYzJ/TgR8Dz9VwjkPYUuLbsSTQEvvha42N2jfGBu/upurD+QVY3b54/Toi/o3AnYd/3QwbqGwP/Ah4J6pfA2z7dD/+PR6NvV0ZHHW8DPgJletbJEqsh4/VcZwnVp/aHmxEUmkWlq/AXmTMxZZfnITlsDOA67Dah9U5iL1lHA98C3tj1gjLYS2xFxd9gXuo+mb3u8Btccaegz3An3L41+uwG5ZeWA7Lw2YvXovN+ImUj+U/P7THauN8J+JYCCsrcTOV6+6KSOI8iM1QBJuJ/Ah2L9ETy2dnYfcRNb0wPIANnt2C1cRuit2HtcZyWyGW5x6kah3nq/E22yNSHvaiIxz/KmAclp+bYjmrw+Hzr4rq2wT/dp/vjN2jnhdx7CD2wvcOqn9ATqX/pups0YpUBCIiGeNOLLeDPZv/EZvVdgqWe087/OuaZqmFgFew3DgA+97I55tn+kZYTdzJWJmxSGcBM+oQ/y/4ZgnzZuy5uQ/2vZWLlQobiq1IipSDDSr6MdO8BVbv9+oYsV5FetcRjzUh6POkRyFZJ9aS5/DnQ755qHU1DLgXb4Wij6Hq8updWKKqTawlz+HPfGpeRtwQewMe3c9ltyMvS55PwP4uI9vvx2Yy+emxqGvuI76Zn+2o+ncUTFCMInUVa8lz+PM+3mfT/gxb5n+Mhz4dqbq8egdVi0vHcj7Vx/97al5GXIA9UEb3iy64HMvcqD41Pbj3xF4sRLbfi72YSQe/pXJsWvIsmSrWkufwZynfDMy5uh0ox15qujqRqsurN+D2YBa95Dk6x9R0jmbYsu3ofi453MuS5zOpWodrJ1Vf+KaTPKx+dfR98pGpDEpE0kqsJc/hz1pscM+LIcDDwPEe+hxF1eXV+/mmtnlNYi15Dn+eoealsA2Ax2P0i56BHouXJc8dsNq1ke0PYrP5/ZaIJc9PxDjHzxIYo1elMeJRDcUMVN2A4md4e6Cuq+bY2+nIGG536FfdgOLfcJtp2Iqqm7Y86tDPdUDxXGzkP7Ltl8AlDteoq+gb7Hh33cql6uYLte3uLZIs1Q0oforVUk2WVtjPWGQMox36nU/s+P9M7CUo0drwzS544c/9Dv3mRvWpbkDx29jSjci22/BW38xvGlCUbFHdgOJyknuTfQxVf+5dXiBUN6A4x/G6XbAZkpF9XV5gug4ofh9b3hXZdhPeN+ZLtglU/Tuty6wfEck+1Q0obsPboGBdNcFWAUbGMNWhX3UDikuoefl1WCF27x/Zd4FDP9cBxdOwgdnItruxgddkqOuA4llU3TwmRHL/bUQrjRGPaihmkbHEPwAVj+3YNOZI8ez6DPbWdhg2CFabLdhMvkheai/W5DKsdkTkG+SN2AP6Mwm6Rk0Ko34dvYzJ1SGqTt/WpiyS7kZhL0aSZQtVH2K/G+e5dgDDsS/W2nyGvXGMlKgcdhXwd2yJd9haoD81l7wQkcS6kZprtybaOqycQ6R489kG3F6uAHyMrRyJlKh8dhP2oBh5b/QhNqP75QRdww+9gZKoY5uwh2ARkdpMxHJdsuwCJkUdi/eZ/gBW4me/Q9udWOmySIn6/hiIlUGKnGy15fDx+Qm6hp96Y99/0SsV/0Vy/234SgOK6WUtVvMm2aJvIs/E2+6CYY/hbSDhH1G/7k7dCsiC3Tz/kco7Vi/HZiwuq+O5XUUv965LTYfowchE70gtkkgfYjNWki06h/XBbZZhtDl4K+gcncNOoe7fq7dhy0cic+HbWA57u47nFhF3r1JzvVe/ROezeDetewh70HMVnc9Oi/O6YTlYba9fUfmeLine truncated
|
||||
}
|
||||
},
|
||||
"cell_type": "markdown",
|
||||
"id": "9fc3be11",
|
||||
"id": "99554811",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"The second function is the one that runs on the workers. Note that we use MPI for communication in this case."
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": null,
|
||||
"id": "d5a8d3a8",
|
||||
"metadata": {},
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"@everywhere function floyd_worker!(Cw,rows_w)\n",
|
||||
" comm = MPI.Comm_dup(MPI.COMM_WORLD)\n",
|
||||
" rank = MPI.Comm_rank(comm)\n",
|
||||
" nranks = MPI.Comm_size(comm)\n",
|
||||
" m,n = size(Cw)\n",
|
||||
" C_k = similar(Cw,n)\n",
|
||||
" for k in 1:n\n",
|
||||
" if k in rows_w\n",
|
||||
" # Send row k to other workers if I have it\n",
|
||||
" myk = (k-first(rows_w))+1\n",
|
||||
" C_k .= view(Cw,myk,:)\n",
|
||||
" for proc in 0:(nranks-1)\n",
|
||||
" if rank == proc\n",
|
||||
" continue\n",
|
||||
" end\n",
|
||||
" MPI.Send(C_k,comm;dest=proc,tag=0)\n",
|
||||
" end\n",
|
||||
" else\n",
|
||||
" # Wait until row k is received\n",
|
||||
" MPI.Recv!(C_k,comm,source=MPI.ANY_SOURCE,tag=0)\n",
|
||||
" end\n",
|
||||
" for j in 1:n\n",
|
||||
" for i in 1:m\n",
|
||||
" @inbounds Cw[i,j] = min(Cw[i,j],Cw[i,k]+C_k[j])\n",
|
||||
" end\n",
|
||||
" end\n",
|
||||
" end\n",
|
||||
" Cw\n",
|
||||
"end"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"id": "bc81a83d",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"<div class=\"alert alert-block alert-success\">\n",
|
||||
"<b>Question:</b> Which of the following statements is true?\n",
|
||||
"<div>\n",
|
||||
"<img src=\"attachment:g3642.png\" align=\"left\" width=\"300\"/>\n",
|
||||
"</div>"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"id": "4ec6718c",
|
||||
"id": "075ee09d",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
" a) The processes are synchronized in each iteration due to the blocking send and receive of row k.\n",
|
||||
" b) Receiving processes may overwrite the data in row k, which can lead to incorrect behavior.\n",
|
||||
" c) The sending process can only continue the computation after the data are received in every other process.\n",
|
||||
" d) The receiving process does not know the source of the received data."
|
||||
"The implementation of the parallel function is given in the following code snipped:"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": null,
|
||||
"id": "4f4a57de",
|
||||
"id": "37863d46",
|
||||
"metadata": {},
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"answer = \"x\" # replace x with a, b, c or d\n",
|
||||
"floyd_impl_check(answer)"
|
||||
"code1 = quote\n",
|
||||
" function floyd_mpi!(C,comm)\n",
|
||||
" myC = distribute_input(C,comm)\n",
|
||||
" floyd_iterations!(myC,comm)\n",
|
||||
" collect_result!(C,myC,comm)\n",
|
||||
" end\n",
|
||||
"end;"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"id": "c624722a",
|
||||
"id": "c6cdae12",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"### Testing the implementation"
|
||||
"In the following cells, we discuss the implementation of the helper functions `distribute_input`, `floyd_iterations!`, and `collect_result!`."
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"id": "627378c3",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"### Distributing the input matrix\n",
|
||||
"\n",
|
||||
"Since only rank 0 receives the input matrix C, this rank needs to split it row-wise and\n",
|
||||
"send the pieces to all other ranks. This is done in the function below. We start by communicating the problem size N to all ranks (at the start only rank 0 knows the problem size). This is trivially done with an `MPI.Bcast!`. Once all ranks know the problem size, they can allocate space for their local part of `C` , called `myC` in the code. After this, rank 0 sends the pieces to all other ranks. We do it here with `MPI.Send` and `MPI.Recv!`. This can also be done with `MPI.Scatter!`, but it is more challenging since we are using a row partition and Julia stores the matrices in column major order. Note that this algorithm can also be implemented using a column partition. In this case, using `MPI.Scatter!` would probably be the best option.\n"
|
||||
]
|
||||
},
|
||||
{
|
||||
"attachments": {
|
||||
"g5958.png": {
|
||||
"image/png": "iVBORw0KGgoAAAANSUhEUgAABSkAAATJCAYAAAA7PJMlAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAAewgAAHsIBbtB1PgAAABl0RVh0U29mdHdhcmUAd3d3Lmlua3NjYXBlLm9yZ5vuPBoAACAASURBVHic7N15eJTlof7x70wSdgFxAXFDLSouVQ9WRQhbCbYq2ipUjwdBxZ/WDaypqEcTIJFqtWkLlVo90IpaTy3QcyxqjwYECWrdagWVKlWxyC4oKGvIvL8/RiOEQN6EzLxJ5vu5rvkj7zzPzA0ql94+C0iSJEmSJEmK2hFAUOXVKdJEaZQddQBJkiRJkiRJTcKBwKnAKcA+Vd77C/D87iZaUkqSJEmSJEmqrXbAiUD3HV7dgNhuxn+KJaUkSZIkSZKkevIdkisj6028Pj9MkiRJkiRJUpNX7wsfLSklSZIkSZIk1Zc1wLzaTnK7tyRJkiRJkqS6WA+8Bby+w+sd4CBgWW0+yJJSkiRJkiRJUm28BHQBPqqvD7SklCRJkiRJklQba7981RtLSkmSJEmSJNW31kCzHX7eCmyqMuYIoD9wHLA/ybtTVgIPAov38NnZwDeBE4Gjv5zbDvgC+BRYCrwA/B2oqGP+fdi5N9sMbKkyphswADgG2A/4HFhNcpXhHHb99aZDc6BVNc8/TXcQSZIkSZIkKWqPAMEOr/t2eO9UYFaV93d8fb+az+sAXAk8DWzYw9wdX58AY7+cW1vzqnzW6B3eO62a96u+1gO3kSwNwzqims/pVIv5RwHvVZlfTvL3Ld06s+uv5ZYIckiSJEmSJCmD7a6kvAnYzp4Lvguq+bxNNczZ02s50KOW+XdXUt4YIv+Or9lAm5DfuTcl5enAqipzPwfODjm/vtW6pIynIZQkSZIkSZJ0LVACZNUwrrq+KmcP49cD60huKa/OQSS3X/esKWANrgB+Qc35d9Qf+O1efm9NzgOeAw7c4dlKoC/JlaeSJEmSJElSRqq6knIOyTMdv/p5IfBjklunDwO6kiz0fgHkVfN55Xy9KvLXwGDgG0CsyrjDgeHAX9l1Jd9HwL4h81ddSTljh/wJ4I/A90gWoNkkz+A8DZi0Q9YdX+eE+M66rKS8ll1Xdi768rOi5HZvSZIkSZIkRa5qSfnVK0Fy63RtViMCzAUGEX5XcIzk1uyKKt9/R8j5uztzci3Qp4a5ZwPbqsz7vxDfWZuSMgbcXc34+SQv8YmaJaUkSZIkSVItHUL1NyKr7nZXUuanOUdBle//mJ1v7d6d6krKrSQv/Qnjp1XmVrDzduzqhC0pmwOPVTN2OtAiZL5Us6SUJEmSJEkK6UKSW2O/KqDup+YiSeFUV1K+QPrvR2kGrKiSo3uIedWVlONq8b2HsOsqzpousQlTUrYnuXW+6rhf0rDunrGklCRJkiRJqkE3kttvq1vp9wUwluRqNdVddSXl4IiyPFQlx7Uh5lQtKbcAB9TyexdW+YzbahhfU0l5GPBWlfcrSN6Y3tB4u7ckSZIkSdJudAAmkCyPztrNmNbAmC/HDElTrkywFZgZ0Xd/WOXn4+vwGa8Ca2o55+0qP9e25NzRycBL7Jx9C3Ax8PO9+NwGI8wefEmSJEmSpMYsB7gcGA/sH3JOV5I3OM8meQHLW6mJljHeIllU1pc2wDeBE0lurW795as6J1X5OewN3zt6ow5z1lX5uV0dPgNgIMnzJvep8tnfA8rq+JkNjiWlJEmSJElqygYAvwBOqOP8b5MsqH4L3A58Uk+5Mk3V1Yx1dTbw/4DvUvct+e3rMGdtHeZsqvJzyzp8xuUkz8LM2eHZEpK/D4vq8HkNltu9JUmSJElSU3QM8BRQSt0Lyq9kA1eRLIWuAbL28vMy0fq9nH8Ayb+eT5FcQbg3Z4bWZe7mOswJ6jCnqp+wc0H5d+BMmlhBCa6klCRJkiRJTUt7oAC4gZ3Lnd069dRTeffdd/n8889rGro/8Gvgh8CPgOf2Imem2bIXcw8keaP1cbt5fwXwMcki9AugvMr7x5DcGt4U7A+0JflrliRJkiRJUgMTB4YBK6n+1u5dXp07dw4mTJgQLF++PHj77beDESNGBFlZWaHmfvmaCRyZvl9io1L1du/79uKz/sCuv/fzgB8Q7jKa0VXmzg4xp+rt3qNrnRrurfIZj9UwvrrbvccDiSrPVrD3q4NTrda3e0uSJEmSJDV2fUieGxmqXGzVqlWQn58ffPjhh8Hy5ct3ej377LPBGWecUZuicivJG8Pbpu+X2yjUV0l5Mrv+nv8EiNXiM+6sMr8xlZSdgCuBiirP1wD/VodM6VLrktIzKSVJkiRJUmN1CPAwya3AJ9c0OBaLMWjQIObNm0d+fj7Nm+96NOEJJ5zAn/70J6ZOncqhhx4aJkMzYCTwD5LnVtq11K8Lqvy8ALiD2p33eGD9xYnEZJKrhLfv8Gx/kmXrGZEkSgH/wZEkSZIkSY1NK2As8B5wKSFW1Z188sn8+c9/5oEHHqBz5841fkFeXh7PP/88t99+O23atAmT6SDgAeAVoFeYCQrlW1V+nkFy+3NtnFpPWaL0e+AiYNsOz9oDz5JcSdzoWVJKkiRJkqTGIgYMAd4BxgAta5rQsWNHJkyYwFNPPUX37t1r9WUtWrTguuuuY968eQwePJhYLNQO4+4ktwr/ETi8Vl+o6lRdBfmvWs4/iKZzac6fSK4s3fESon2Ap4GBkSSqR5aUkiRJkiSpMegOlBGy/PuqYCwrK2PIkCFhC8ZqderUiYkTJ9am6NyxTB1LiDJVu5Vd5ecWtZz/QyCrnrI0BE8B5wIbd3jWCvgzMCiSRPXEklKSJEmSJDVkO26j7hlmQl5eHnPnzq3NVu1QvtoyPmHCBA48MNQxh61Irvh8j+SZgnVvSjPXmio/n16Lud8A8usxS0MxG/gOsGGHZ81JboUfEkmiemBJKUmSJEmSGqIcYBS1uJBmx0tvDjvssJSEisViDBkyhBdffHG3l+9U4xBgKiEv+NFOXq3y88Uky8eaHEiytGtd74kahvnAAGDdDs9ygP8meU5ro1N1yawkSZIkSVLUBgG/AI4KM3jfffflRz/6EZdffjlZWenZ2duqVSvy8/O58MILueuuu5g5c2aYaX2A10legnIzsCqVGZuI/wFu3eHnFsBfSJ7NuHA3cwYAv+Hrv3++AOpvSW3D8SrQn+TlOV8t7c0CHiJ5xMCDKf7+09j90Qv7VvPsRPaw0tOSUpIkSZIkNRTHkiwnvxNmcE5ODsOGDePmm2+mbdu2qU22G126dOGBBx7g0ksvZcyYMSxatKimKXGSK92+D5QAdwFbUxyzMXsFeJLkOYxf+QbwBsmyci6wguSKyS7Ad4FTdhj7AcnVhbenPmok3gT6ArOAr66tj5MsaVsCE1L43ddTu1Wb//Hlq1pu95YkSZIkSVHrQLJMWUjIgjI3N5dnn32W4uLiyArKHfXq1YtnnnmGe+65h/322y/MlDYkz6tcSCM+RzBNRgBLqjzLIllc/ozkytQHgf9k54JyJXA2sD71ESO1COgNfLTDsxjwS3ZehdqgWVJKkiRJkqSoZJM8b/JdYCQhdnweeeSRPPzwwzz++OMcc8wxqc5XK9nZ2QwdOpT58+czYsQIsrNDbWDtSvLG8lLghJQGbLxWA7nAC7WYM4fkduR3U5Ko4XmfZFH5zyrP7wLGpT9O7XmrlCRJkiRJisK3SW7tPjHM4LZt23LDDTdw1VVXkZOTk9pk9eT9999n7NixzJ49O+yUcuB3JLcmf5KyYOlxETuvanwJeGIvPzMOnE9yZWVvYJ8d3guA5cDzwMMkz2kMvnyvN8kVlV95H/ivGr7rhyS3j3/lL19+dm2c/eV3f+VNklvPd6cDMLrKs/HA57X4zoNIFv5VO7+pJFdc1qeqf433iiWlJEmSJElKp28APyHkFud4PM4FF1xAYWEh+++/f2qTpUhZWRl33HEHixcvDjtlHVAE3AdUpCxY45YFHEDygpYNJC8h2h5pIu0VS0pJkiRJkpQOrUneaH0r0DzMhDPPPJOioiKOO+64lAZLh/LycqZOncq9997L55+HXhi3CPgR8EzqkkkNgyWlJEmSJElKpTgwFLgH6BhmQufOnbnlllsYMqTp3Sezbt06fvGLX/DQQw9RURF6keSTwCiSN1VLTZIlpSRJkiRJSpXTSN7afUaYwS1btuTaa6/l+uuvp3nzUIstG62FCxdSWFjIyy+/HHbKNuA3wB3U7oxCqVGwpJQkSZIkSfXtYJK3Cg8lRPcQi8U499xzKSws5OCDD055uIaktLSUO+64g6VLl4adspzkbc2TgUTKgklpZkkpSZIkSZLqSyvgBpKr/dqEmXDSSSdRXFzMqaeemtJgDdnmzZv57W9/yy9/+Us2btwYdtprJLeAv5i6ZFL6WFJKkiRJkqT6MAiYCHQJM7hjx47k5+dzySWXEI/HUxqssVi5ciU/+clPmDFjBkEQhJkSANOBHwP/Smk4KcUsKSVJkiRJ0t74N+CXQG6YwTk5OQwbNoxbbrmFNm1CLbbMOG+88QaFhYW8/vrrYadsAu4F7ga2pCyYlEJZUQeQJEmSJEmN0v4kz538L0KunszLy+Phhx/m+9//Ps2aNUtltkbtoIMO4uKLL+bwww/nb3/7W5gt4DlAX2A4sBZ4M8URpXrnSkpJkiRJklQbOcC1QBHQNsyErl27Mm7cOPr27ZvKXE3Spk2buP/++/nVr37Ftm3bwk6bC9yIZaUaEUtKSZIkSZIU1gCS5052CzO4ffv23HTTTVx++eVkZbmZc28sWbKEu+66i5kzZ4adkgB+T/K8ytUpCybVE0tKSZIkSZJUk2OAnwNnhxmck5PDRRddxK233kqHDh1SmyzDzJ8/nzFjxrBo0aKwUz4jeVblL4DQSzGldLOklCRJkiRJu7MvcAvwIyDUIZK5ubmMGzeOY489NqXBMtn27dv5wx/+wE9/+lPWrl0bdtp7QD7wLine truncated
|
||||
}
|
||||
},
|
||||
"cell_type": "markdown",
|
||||
"id": "d1ddc881",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"<div>\n",
|
||||
"<img src=\"attachment:g5958.png\" align=\"left\" width=\"300\"/>\n",
|
||||
"</div>"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": null,
|
||||
"id": "09937668",
|
||||
"id": "7efe441c",
|
||||
"metadata": {},
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"function rand_distance_table(n)\n",
|
||||
" threshold = 0.4\n",
|
||||
" mincost = 3\n",
|
||||
" maxcost = 10\n",
|
||||
" infinity = 10000*maxcost\n",
|
||||
" C = fill(infinity,n,n)\n",
|
||||
" for j in 1:n\n",
|
||||
" for i in 1:n\n",
|
||||
" if rand() > threshold\n",
|
||||
" C[i,j] = rand(mincost:maxcost)\n",
|
||||
" end\n",
|
||||
"code2 = quote\n",
|
||||
" function distribute_input(C,comm)\n",
|
||||
" rank = MPI.Comm_rank(comm)\n",
|
||||
" P = MPI.Comm_size(comm)\n",
|
||||
" # Communicate problem size\n",
|
||||
" if rank == 0\n",
|
||||
" N = size(C,1)\n",
|
||||
" if mod(N,P) !=0\n",
|
||||
" println(\"N not multplie of P\")\n",
|
||||
" MPI.Abort(comm,-1)\n",
|
||||
" end\n",
|
||||
" Nref = Ref(N)\n",
|
||||
" else\n",
|
||||
" Nref = Ref(0)\n",
|
||||
" end\n",
|
||||
" MPI.Bcast!(Nref,comm;root=0)\n",
|
||||
" N = Nref[]\n",
|
||||
" # Distribute C row-wise\n",
|
||||
" L = div(N,P)\n",
|
||||
" myC = similar(C,L,N)\n",
|
||||
" if rank == 0\n",
|
||||
" lb = L*rank+1\n",
|
||||
" ub = L*(rank+1)\n",
|
||||
" myC[:,:] = view(C,lb:ub,:)\n",
|
||||
" for dest in 1:(P-1)\n",
|
||||
" lb = L*dest+1\n",
|
||||
" ub = L*(dest+1)\n",
|
||||
" MPI.Send(view(C,lb:ub,:),comm;dest)\n",
|
||||
" end\n",
|
||||
" else\n",
|
||||
" source = 0\n",
|
||||
" MPI.Recv!(myC,comm;source)\n",
|
||||
" end\n",
|
||||
" return myC\n",
|
||||
" end\n",
|
||||
" C[j,j] = 0\n",
|
||||
" end\n",
|
||||
" C\n",
|
||||
"end"
|
||||
"end;"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"id": "72d82138",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"### Running Floyd's updates in parallel\n",
|
||||
"\n",
|
||||
"As discussed above, we need to communicate row $k$ of matrix $C$ in order to perform iteration $k$ of the algorithm. The function below is similar to the sequential function `floyd!`, but there is a key difference. At the start of iteration $k$, the owner of row $k$ sends it to the other processors. Once row $k$ is available, all ranks can do the Floyd update locally on their portion of matrix $C$."
|
||||
]
|
||||
},
|
||||
{
|
||||
"attachments": {
|
||||
"fig-asp-efficiency-comm-2.png": {
|
||||
"image/png": "iVBORw0KGgoAAAANSUhEUgAAA1AAAAMbCAYAAABDhHSeAAAACXBIWXMAAB7CAAAewgFu0HU+AAAAGXRFWHRTb2Z0d2FyZQB3d3cuaW5rc2NhcGUub3Jnm+48GgAAIABJREFUeJzs3Xt81dWV9/HPSggJIJcqNxEUW8TitYpYFFCsYlWIM9MiUYeZWqVFGwUqGrk8WquDIIpVaFrQUeplShEfnWnUto/ai0FEW7x3ULxFCYpEBCn3QNbzxz4hh3BOrifnknzfrxevnvM7+7d/C6EkK3vvtczdERERERGRzGRmxwP/BIwC+gPdAQO+AN4ESoEH3X1dA+aaD+TWM2wHsC0y/9vAq+7+ST3zTgaOibxd6u5/rC+WqHv7AD+JunSlJziJMbPpwG2Rt3nuvivuWCVQIiIiIiKZx8y+CswBxhISprrsAe4HbnD3L+uYcyvQqQnh/BUocvc/x5n3D8C5kbeT3X1+Qyc2s2OBt6IuZbt7VRNijDd/NvC/wMDIpToTqKxEPVhERERERJLDzM4EVgEXUZM8OWFF6C/As8A7QHWi0Q6YCLxqZke1QEhDgGfMbFwLzN3S5lKTPNWrXQsGIiIiIiIiCWZmI4HfU7PVbiswG/hV7a10ZtYD+AEwg7CydCTwFzP7hrtvqOdRPwLejfNZe+AwYARwMZBDyC0eMLMX3X1tY39fyWZmQ4AbgfzG3KcESkREREQkQ5hZd+ARapKnD4Bz3P3DWOPdvQK4zcz+m7AqdWjkVzFh9aouL7v7qnrG3GdmDwJPE5KqTsD1wKQG/HaSzsxmA6cBXwd6NWUObeETEREREckctxFWfgC+BL4VL3mK5u7/C3w/6tJ3zOxriQjI3Z8DHoi69N1EzNtCzoz8ik6enm/MBEqgREREREQygJn1BP4t6tIN7v5RQ+939z8AqyNvs2jk1rV6/N+o133M7JAEzp1IfwUeBxYCUwiVAb/dmAm0hU9EREREJDNcCuRFXm8EHmzCHAuB4ZHXOxIRVETtEumdCDGmFXefXPuameXFGhuPEigRERERkcwwIur1/7j7zsZOECkf3uAS4o3Qo9b7tEueEkVb+EREREREMsPpUa9fSlkUsV0Y9foDd9+WskhamBIoEREREZE0F2n2Gl344J1UxVJbpKx6YdSl/05RKEmhLXwiIiIiIunvK9Q0zIVQga+ldTazr8S43gnoAgwC/olwNis78tkm4PYkxJYySqBERERERNJfx1rvE1kAIp4/NXL8VuDCBjTozWjawiciIiIikv4213rfJSVRxLYXeAI4wd2XpzqYlqYVKBERERGR9PcPYA81378no8/Sy8CWGNerCFsI1wGvAc+6e+0y5rXtiXqdHXdUbNE5i0d+pYwSKBERERGRGMysH3ADcCqwEri9AYlCi3B3N7N3CeeOAE4Eft/Cj/2Ru69K0Fxbo153auS9B0W93unuSqBERERERNKFmR0CXA9MAjpELg8BJprZr4Cb3P2zFIS2nJoEalhTJjCzY4A+kbefuvvfExFYA0QnUI1dPYsen4ziGXXSGSgREREREcDMupjZT4EPCStPHWoNaQ/8EFhjZjeZWeckh/hc1OvzzOzQJsxxN/BM5NelCYmqYT6Ken1MI++NHl+egFiaRQmUiIiIiLRpZtbBzKYC7wM3AfUlRl2AnwLvm9kUM8tt6RgjHiecOwLIAW5szM1m1hUYHnXpbwmKqyHeino91MxqJ6d1OSvOPCmhBEpERERE2iQzyzGzHwJrgDuB7rXHdOjQgTFjxpCXlxdrih7AzwgrUleYWYsej3H3ykic1a40swsbMcVN1KyqbQJ+l6jYGuBlaoo/dAEua8hNZnYscHbUpb8kNqzGUwIlIiIiIm2KBRcRVjMWAX1rj8nJyWH8+PGsWLGCe++9l5dffpnCwkLat28fa8rDgf8E3jKzfzezlvweez41/ZkMWGZm3zczi3eDmWWZWRHw46jLt7v7zhaMcz/uXg78IerSHWZ2QV33mNlRwDJqqvZtBv5vy0TYcJbiIhYiIiIiIkljZucAtwMnx/o8KyuL0aNHM336dPr373/A5+Xl5cyfP58lS5awd+/eeI95C7jF3ZclKOz9mFkf4I/A0VGXXwR+BTxP2OZnwGGELXsTgcFRY/8IfNvdo0uLV8+9lZoqeacksAofZvZ1wkpU9RZJJ5zrehx4k1BoIhcYAIwCLo68r/ZDd78vUfFExZXH/o2J89x9V9zxSqBEREREpLUzs9OA24CR8caMGDGCm266iWOPPbbe+d555x3uuusunnzySer4fnoFMMPdE77tLFIp8AlgRCNv/S1wibtvjzNviyVQkfnPIcR9UH1ja5nl7v8nkbFUa2wCpS18IiIiItJqmdlxZvYoIZkZGWvMkCFDePzxx1m6dGmDkieAo48+mkWLFvHb3/6W4cOHxxt2OvBnM3vGzGKueDWVu28k/H7GA+814JY3CInTP8VLnpLB3Z8FTgGeJjTkrc8bwOiWSp6aQitQIiIiItLqmFl/YDowgTiLBt/4xjf48Y9/zKhRo5r9vNLSUmbNmsUbb7wRb4gDjwEz3f3dZj+wFjM7DjgHOIJQ3MKAL4C3gRfc/bUGznMWNWeOXnb3LYmONepZfQkFIo4jFPDoBPwD+DTy63l3f7Olnh8VRxbwrahLf3T3uMmdEigRERERaTUi35TfCFwOxKyKN2DAAK6//nrGjBlDHbUXGs3defbZZ5kzZw6rV6+ON6wSWEw4I7Uu3iBJX0qgRERERCTjRc4EXQ9M4sAGuAD06dOHKVOmcPHFF9OuXctVHK+qquKpp57itttu46OPPoo3bDeh6MNN7v5ZiwUjCacESkREREQylpkdBBQStut1jTXm4IMP5qqrrmLChAnk5iar5y1UVlaydOlS7rzzTjZs2BBv2FagGJjt7l8mLThpMiVQIiIiIpJxzKw9oRnrLUCvWGMOOuggvve97zFp0iQ6d+4ca0hSbN++ncWLF7NgwQK2bIl7pGgjcAdwTzL7M0njKYESERERkYxhZu2AS4GbgSNjjenQoQOXXnopkydPpnv37skMr06bN2+muLiY+++/n5074+ZIa4H/AB6I1adJUk8JlIiIiIikPQvVHsYSkouBscbk5ORQUFDA1KlT6dUr5qJUWli/fj3FxcU8/PDD7N69O96wdwh9qx6pqyKcJJ8SKBERERFJa5Hmq7cDMXspZWVlMXr0aKZPn07//v2TGltzlJeXM3/+fJYsWcLevXvjDXuLULFvWRJDkzoogRIRERGRtGRmpxFWYUbGGzNixAhuuummBjfATUfvvPMOd911F08++SR1fG++Apjh7n9JYmgSgxIoEREREUkrkaawNwEXxRszZMgQpk+fztChQ5MXWAtbtWoVc+bM4YUXXqhr2LNAkbu/mqSwpBYlUCIiIiKSFsysP6Ec+QQgK9aYb3zjG/z4xz9m1KhRSYwsuUpLS5k1axZvvPFGvCEOPAbMdPd3kxeZgBIoEREREUkxM+sL3AhcDsTscDtgwACuv/56xowZQ6gn0bq5O88++yxz5sxh9erV8YZVAosJZ6TWJS+6tk0JlIiIiIikhJkdAlwPTAI6xBrTp08fpkyZwsUXX0y7djFzq1atqqqKp556ittuu42PPvoo3rDdwK+AG909bsdeSQwlUCIiIiKSVGZ2EFBI2K7XNdaYgw8+mKuuuooJEyaQm5ub1PjSUWVlJUuXLuXOO+9kw4a4OdJWoBi4zd3jduyV5lECJSIiIiJJYWbtgcuAW4CYjZo6derEZZddxqRJk+jcuXMyw8sI27dvZ/HixSxYsIAtW+LmSJ8DdwL3uHvcjr3SNEqgRERERKRFmVkW8F1CL6cjY43p0KEDl156KZMnT6Z79+5JjS8Tbd68meLiYu6//3527oybI60lNB5+wN33JC+61k0JlIiIiIi0CAvVHsYSvokfGGtMTk4OBQUFTJ06lV69Yi5KSR3Wr19PcXExDz/8MLt374437B1CP61H3L0qedG1TkqgRERERCThzOwcworTybE+z8rKYvTo0UyfPp3+/fsnNbbWaO3atSxYsIAlS5awd+/eeMPeBG5192VJDK3VUQIlIiIiIgljZqcRVjtGxhszYsQIbrrpJo499tikxdVWvPPOO9x1112UlJTUNWwFMMPd/5KksFoVJVAiIiIi0mxmdhxwE3BRvDFDhgxh+vTpDB06NHmBtVGrVq1izpw5vPDCC3UNexYocvdXkxRWq6AESkRERESazMz6E8qRTwCyYo0ZNGgQU6ZMIT8/P4mRCUBpaSmzZs3ijTfeiDfEgceAme7+bvIiy1xKoERERESk0cysL3AjcDkQs8PtgAEDuP766xkzZgyhnoSkgrvz5JNPMnfuXN5///14wyqBxcAt7r4uedFlHiVQIiIiItJgZnYIcD0wCegQa0yfPn2YMmUKF198Me3axcytJAWqqqp46qmnmDVrFh9//HG8YbuBXwE3unvcjr1tmRIoEREREamXmR0EFBK263WNNebggw/mqquuYsKECeTm5iY1Pmm4yspKli5dyp133smGDXFzpK1AMXCbu8ft2NsWKYESERERkbjMrD1wGXALELNRU6dOnbjsssuYNGkSnTt3TmZ40gzbt29n8eLFLFiwgC1b4uZInwN3Ave4e9yOvW2JEigREREROYCZZQHfJfRyOjLWmA4dOnDppZcyefJkunfvntT4JHE2b95McXEx999/Pzt3xs2R1hIaIj/g7nuSF136UQIlIiIiIvtYqPYwlvDN8sBYY3JycigoKGDq1Kn06hVzUUoy0Pr16ykuLuahhx6isrIy3rC3CeXqH/M2mkgogRIRERERAMzsHMKK08mxPs/KymL06NFMnz6d/v37JzU2SZ61a9eyYMEClixZwt69e+MNexO41d2XJTG0tKAESkRERKSNM7PTgFnAWfHGjBgxghtvvJHjjjsueYFJSr3zzjvcddddlJSULine truncated
|
||||
}
|
||||
},
|
||||
"cell_type": "markdown",
|
||||
"id": "6abd8c34",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"<div>\n",
|
||||
"<img src=\"attachment:fig-asp-efficiency-comm-2.png\" align=\"left\" width=\"250\"/>\n",
|
||||
"</div>"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": null,
|
||||
"id": "dd77ee3d",
|
||||
"id": "52961b52",
|
||||
"metadata": {},
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"using Test\n",
|
||||
"load = 10\n",
|
||||
"n = nworkers()*load\n",
|
||||
"C = rand_distance_table(n)\n",
|
||||
"C_seq = floyd!(copy(C))\n",
|
||||
"C_par = floyd_dist!(copy(C))\n",
|
||||
"@test C_seq == C_par"
|
||||
"code3 = quote\n",
|
||||
" function floyd_iterations!(myC,comm)\n",
|
||||
" L = size(myC,1)\n",
|
||||
" N = size(myC,2)\n",
|
||||
" rank = MPI.Comm_rank(comm)\n",
|
||||
" P = MPI.Comm_size(comm)\n",
|
||||
" lb = L*rank+1\n",
|
||||
" ub = L*(rank+1)\n",
|
||||
" C_k = similar(C,N)\n",
|
||||
" for k in 1:N\n",
|
||||
" if (lb<=k) && (k<=ub)\n",
|
||||
" # Send row k to other workers if I have it\n",
|
||||
" myk = (k-lb)+1\n",
|
||||
" C_k[:] = view(myC,myk,:)\n",
|
||||
" for dest in 0:(P-1)\n",
|
||||
" if rank == dest\n",
|
||||
" continue\n",
|
||||
" end\n",
|
||||
" MPI.Send(C_k,comm;dest)\n",
|
||||
" end\n",
|
||||
" else\n",
|
||||
" # Wait until row k is received\n",
|
||||
" MPI.Recv!(C_k,comm,source=MPI.ANY_SOURCE)\n",
|
||||
" end\n",
|
||||
" # Now, we have the data dependencies and\n",
|
||||
" # we can do the updates locally\n",
|
||||
" for j in 1:N\n",
|
||||
" for i in 1:L\n",
|
||||
" myC[i,j] = min(myC[i,j],myC[i,k]+C_k[j])\n",
|
||||
" end\n",
|
||||
" end\n",
|
||||
" end\n",
|
||||
" myC\n",
|
||||
" end\n",
|
||||
"end;"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"id": "54081b85",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"### Collecting back the results\n",
|
||||
"\n",
|
||||
"At this point, we have solved the ASP problem, but the solution is cut in different pieces, each one stored on a different MPI rank. It is often useful to gather the solution into a single matrix, e.g., to compare it against the sequential algorithm.The following function collects all pieces and stores them in $C$ on rank 0. Again, we implement this with `MPI.Send` and `MPI.Recv!` as it is easier as we are working with a row partition. However, we could do it also with `MPI.Gather!`."
|
||||
]
|
||||
},
|
||||
{
|
||||
"attachments": {
|
||||
"g7043.png": {
|
||||
"image/png": "iVBORw0KGgoAAAANSUhEUgAABSkAAATJCAYAAAA7PJMlAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAAewgAAHsIBbtB1PgAAABl0RVh0U29mdHdhcmUAd3d3Lmlua3NjYXBlLm9yZ5vuPBoAACAASURBVHic7N1peFXlofbx/95JmAXEAcQ6F5GqVYsTQphKsFXRVkF6PAoqHm1FwTYV8WhCTKRaa1pBqdUDrai1tUDPsah9a0BQUOtUFVSqVMUis6CgjCF7vR+2RgiBrITsvZLs/++68iFrP8/aN7b10rvPAJIkSZIkSZKidgQQVPnpFGmiNMqOOoAkSZIkSZKkJuFA4GTgJGCfKp/9FXhmdxMtKSVJkiRJkiTVVjvgeKD7Dj/dgNhuxn+CJaUkSZIkSZKkevIdkisj6028Pl8mSZIkSZIkqcmr94WPlpSSJEmSJEmS6ssa4NnaTnK7tyRJkiRJkqS6WA+8Cby6w8/bwEHAstq8yJJSkiRJkiRJUm28ABwOfFhfL7SklCRJkiRJklQba7/4qTeWlJIkSZIkSapvrYFmO/y+FdhUZcwRQH/gG8D+JO9OWQncDyzew7uzgW8CxwNHfzG3HfA58AmwFHgOeB2oqGP+fdi5N9sMbKkyphswAOgK7Ad8BqwmucpwDrv+edOhOdCqmuefpDuIJEmSJEmSFLWHgGCHn3t2+OxkYFaVz3f8+X417+sAXAE8AWzYw9wdfz4Gir6YW1vPVnnXmB0+O7Waz6v+rAduJFkahnVENe/pVIv5RwHvVplfTvKvW7p1Ztc/yw0R5JAkSZIkSVIG211J+RNgO3su+M6v5n0ba5izp5/lQI9a5t9dSXldiPw7/swG2oT8zr0pKU8DVlWZ+xlwVsj59a3WJWU8DaEkSZIkSZKkq4FSIKuGcdX1VTl7GL8eWEdyS3l1DiK5/bpnTQFrcDnwK2rOv6P+wG/38ntrci7wNHDgDs9WAn2BJ1P83ZIkSZIkSVKDVXUl5RySZzp++ftC4Kckt04fCnQhWej9Csir5n3b+GpV5K+BwcDXgViVcYcBw4G/s+tKvg+BfUPmr7qScsYO+RPAn4DvkSxAs0mewXkqMInkFuuq3312iO+sy0rKq9l1ZeeiL94VJbd7S5IkSZIkKXJVS8ovfxIkt07XZjUiJEvOQYTfFRwjuTW7osr33xxy/u7OnFwL9Klh7ll8Vap++fP/QnxnbUrKGHB7NePnk7zEJ2qWlJIkSZIkSbV0CMmVcKo/uyspf5rmHAVVvv8jdr61e3eqKym3krz0J4yfV5lbwc7bsasTtqRsDjxSzdjpQIuQ+VLNklKSJEmSJCmkwcA/+aqAug/oGGmipqO6kvJ50n8/SjNgRZUc3UPMq66kvKUW3/s1dl3FWdMlNmFKyvYkV5VWHXcXDevuGUtKSZIkSZKkGnQjuf22upV+nwNFJFerqe6qKymHRJTlgSo5rg4xp2pJuQXYv5bfu7DKO26sYXxNJeWhwJtVPq8geWN6Q+Pt3pIkSZIkSbvRAZhAsjw6czdjWgPjSJZBUZVqTdFW4C8RffcHVX4/tg7veBn4uJZz3qry+wF1+N4vnQi8wM7ZtwA/AH65F+9tMMLswZckSZIkSWrMckiunisiuV02jK+TvMF5NvBjksWm6u5NkkVlfWkDfBM4nuQKw9ZAy92MPaHK72Fv+N7Ra3WYs67K7+3q8A6AgSTPm9ynyru/B8yr4zsbHEtKSZIkSZLUlA0geV5fXVbPAXwb+AfwW5I3Q6+pp1yZpupqxrr6LnAl8B3qfklM2KJ6R2vrMGdTld93V6LuyWUkz8LM2eHZEpLnWy6qw/saLLd7S5IkSZKkpqgr8CRQRt0Lyi9lkyzGFpFckZm1l+/LROv3cv4BwBMk/zP9Hnt3i3VdzhvdXIc5QR3mVPUzdi4oXwfOoIkVlOBKSkmSJEmS1LS0BwqBa9i53Nmtk08+mXfffZcNGzbUNHQ/YBLwQ+A64Om9yJlptuzF3ANJ3mj9jd18vgL4iGQR+jlQXuXzriS3hjcF+wNtSf6ZJUmSJEmS1MDEgWHASqq/tXuXn8MPPzy47777guXLlwdvvfVWMGLEiCArKyvU3C9+ZgJHpe+P2KhUvd37nr141x/Z9a/9s8CFhLuMZkyVubNDzKl6u/eYWqeGX1R5xyM1jK/udu/xQKLKsxXAcXXIk061vt1bkiRJkiSpsetDchtsqHKxVatWQX5+fvDBBx8Ey5cv3+nnqaeeCk4//fTaFJXbSN4Y3jZ9f9xGob5KyhPZ9a/5z4BYLd5xa5X5jamk7ARcAVRUeb4G+FYdMqVLrUtKz6SUJEmSJEmN1deAB0luBa56g/MuYrEYgwcP5vnnnyc/P5/mzXc9mvC4447jz3/+M1OnTuXQQw8NkyEHGAX8k+S5lZ5XWb/Or/L7ApIXGNXmvMcD6y9OJCaTXCW8fYdn+5MsW0+PJFEKWFJKkiRJkqTGphVQBLwLXEKIVXUnnngif/nLX5g4cSIHHlhzZ5WXl8fcuXO56aabaNOmTZhMBwH3AS8CvcJMUCinVPl9Bsntz7Vxcj1lidLvgaEkV+5+qT3wFMmVxI2eJaUkSZIkSWosYsAQ4G1gHNCypgmdOnViwoQJPPHEE3Tv3r1WX9aiRQtGjhzJs88+y+DBg4nFQu0w7k5yq/CfgMNq9YWqTtVG+d+1nH8QTefSnD+TXFm64yVE+5C88XxgJInqkSWlJEmSJElqDLoD8whZ/u1YMA4ZMiRswVitTp06MXHixNoUnTuWqUWEKFO1W9lVfm9Ry/k/pGltwX8COAfYuMOzVsBfgEGRJKonlpSSJEmSJKkh+3Ib9UtAzzAT6rBVO5Qvt4xPmDAh1JZxkuXROJLb0odRu8telLSmyu+n1WLu14H8eszSUMwGvgNs2OFZc5Jb4YdEkqgeWFJKkiRJkqSGqBkwmq8upKmxw6jDpTe1FovFGDJkyB4v36nG14CpwAvUrmQTvFzl9x+QLB9rciDJ0q51vSdqGOYDA4B1OzzLAf5A8pzWRqfqkllJkiRJkqSoDQLuAo4MM3jfffflxz/+MZdddhlZWenZ2duqVSvy8/O54IILuO2225g5c2aYaacBz5O8BOV6YFUqMzYR/wuM3eH3FsBfSZ7NuHA3cwYAvwGO+uL3z4H6W1LbcLwM9Cd5ec6XS3uzgAdIHjFwf4q//1R2f/TCvtU8O549rPS0pJQkSZIkSQ3FMcCvSG5lrVFOTg7Dhg3j+uuvp23btqlNthuHH3449913H8OGDaOwsJBFixbVNCVOcqXb94FS4DZga4pjNmYvAY+TPIfxS18HXiNZVs4FVpBcMXk48F3gpB3Gvk9ydeFNqY8aiTeAvsAsoPMXz+IkS9qWwIQUfvdIkscYhPWfX/xUy+3ekiRJkiQpah1IlikLCVlQ5ubmUlZWRklJSWQF5Y569uxJWVkZEyZMYL/99gszpQ3J8yoX0ojPEUyTEcCSKs+ySBaXd5JcmXo/8N/sXFCuBM4C1qc+YqQWAb2BD3d4FiO5GnlstTMaIEtKSZIkSZIUlWyS502+A4wixI7PI488kgcffJBHH32Uo48+OtX5aiUejzNkyBDmz5/PyJEjycnJCTOtC8kby8uA41IasPFaDeQCz9VizhyS25HfSUmihuc9kkXlv6o8vw24Jf1xas9bpSRJkiRJUhS+TXKlV6hirm3btlx77bVceeWVYcu/yL333nsUFRUxe/bssFO2A78luTX545QFS4+h7Lyq8QXgsb18Zxw4j+TKyt7APjt8FgDLgWeAB0me0xh88Vlvkisqv/Qe8D81fNcPSW4f/9Jfv3h3bZz1xXd/6Q2SW893pwMwpsqz8cBntfjOg0gW/lU7v6kkV1zWpwuBb9XXyywpJUmSJElSOnUhWbyE2uIcj8c5//zzKSwsZP/9909tshSZN28eBQUFvPvuu2GnrAOKgXuAipQFa9yygANIXtCygeQlRNsjTaS9YkkpSZIkSZLSoQ3wU5Jn5DUPM6Fnz54UFxfTrVu3lAZLh/LycqZOncqdd97Jhg0bwk5bBPwY+FvqkkkNgyWlJEmSJElKpThwMXAH0DHMhM6dO3PDDTcwZEjTu0/mk08+4Ze//CUPPPAAFRWhF0k+DowmeVO11CRZUkqSJEmSpFQ5jeS5k6eHGdyqVSt+9KMfcc0119C8eajFlo3WwoULKSws5MUXXww7ZRvwG6CA5PZmqUmxpJQkSZIkSfXta8DPSK6grLF7iMVinHPOORQWFnLwwQenPFxDUlZWxs0338zSpUvDTllO8rbmyUAiZcGkNLOklCRJkiRJ9aUVcC1wM8kzKGt0wgknUFJSwsknn5zSYA3Zli1bmDJlCnfddRcbN24MO+0V4DrgudQlk9LHklKSJEmSJNWHQcDdwGFhBnfs2JH8/Hwuuugi4vF4apM1EitXruRnP/sZM2bMIAiCMFMCYDrJC4n+ndJwUopZUkqSJEmSpL3xLWAC0CvM4JycHIYNG8YNN9xAmzahFltmnNdee43CwkJeffXVsFM2Ab8Abge2pCyYlEKWlJIkSZIkqS72J3mJy0ggK8yEvLw8SkpKOPTQQ1MarCkIgoDp06czfvx4Vq9eHXbaUpJb7R8iucpSajQsKSVJkiRJUm3kAFcDxUDbMBO6dOnCLbfcQt++fVOZq0natGkT9957L3fffTfbtm0LO20uyfMq30hZMKmeWVJKkiRJkqSwBgG/Ao4KM7h9+/b85Cc/4bLLLiMrK9RiS+3GkiVLuO2225g5c2bYKQng9yTPqwy9FFOKiiWlJEmSJEmqSVeS5eR3wwzOyclh6NChjB07lg4dOqQ2WYaZP38+48aNY9GiRWGnfEryrMq7gK0pCybtJUtKSZIkSZK0O/sCRSS3d2eHmZCbm0txcTFdu3ZNZa6Mtn37dv74xz/y85//nLVr14ad9i6QDzyeumRS3VlSSpIkSZKkqrKBy4FbLine truncated
|
||||
}
|
||||
},
|
||||
"cell_type": "markdown",
|
||||
"id": "a0fb6974",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"<div>\n",
|
||||
"<img src=\"attachment:g7043.png\" align=\"left\" width=\"300\"/>\n",
|
||||
"</div>"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": null,
|
||||
"id": "8c725aa9",
|
||||
"metadata": {},
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"code4 = quote\n",
|
||||
" function collect_result!(C,myC,comm)\n",
|
||||
" L = size(myC,1)\n",
|
||||
" rank = MPI.Comm_rank(comm)\n",
|
||||
" P = MPI.Comm_size(comm)\n",
|
||||
" if rank == 0\n",
|
||||
" lb = L*rank+1\n",
|
||||
" ub = L*(rank+1)\n",
|
||||
" C[lb:ub,:] = myC\n",
|
||||
" for source in 1:(P-1)\n",
|
||||
" lb = L*source+1\n",
|
||||
" ub = L*(source+1)\n",
|
||||
" MPI.Recv!(view(C,lb:ub,:),comm;source)\n",
|
||||
" end\n",
|
||||
" else\n",
|
||||
" dest = 0\n",
|
||||
" MPI.Send(myC,comm;dest)\n",
|
||||
" end\n",
|
||||
" C\n",
|
||||
" end\n",
|
||||
"end;"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"id": "d573108e",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"### Running and testing the code\n",
|
||||
"\n",
|
||||
"In the cell below, we run the parallel code and compare it against the sequential. Note that we can only compare both results on rank 0 since this is the only one that contains the result of the parallel code. We have also included a function that generates random distance tables of a given size $N$."
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": null,
|
||||
"id": "c8e9cea7",
|
||||
"metadata": {},
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"code = quote\n",
|
||||
" using MPI\n",
|
||||
" MPI.Init()\n",
|
||||
" $code1\n",
|
||||
" $code2\n",
|
||||
" $code3\n",
|
||||
" $code4\n",
|
||||
" function input_distance_table(n)\n",
|
||||
" threshold = 0.1\n",
|
||||
" mincost = 3\n",
|
||||
" maxcost = 9\n",
|
||||
" inf = 10000\n",
|
||||
" C = fill(inf,n,n)\n",
|
||||
" for j in 1:n\n",
|
||||
" for i in 1:n\n",
|
||||
" if rand() > threshold\n",
|
||||
" C[i,j] = rand(mincost:maxcost)\n",
|
||||
" end\n",
|
||||
" end\n",
|
||||
" C[j,j] = 0\n",
|
||||
" end\n",
|
||||
" C\n",
|
||||
" end\n",
|
||||
" function floyd!(C)\n",
|
||||
" n = size(C,1)\n",
|
||||
" @assert size(C,2) == n\n",
|
||||
" for k in 1:n\n",
|
||||
" for j in 1:n\n",
|
||||
" for i in 1:n\n",
|
||||
" @inbounds C[i,j] = min(C[i,j],C[i,k]+C[k,j])\n",
|
||||
" end\n",
|
||||
" end\n",
|
||||
" end\n",
|
||||
" C\n",
|
||||
" end\n",
|
||||
" comm = MPI.Comm_dup(MPI.COMM_WORLD)\n",
|
||||
" rank = MPI.Comm_rank(comm)\n",
|
||||
" if rank == 0\n",
|
||||
" N = 24\n",
|
||||
" else\n",
|
||||
" N = 0\n",
|
||||
" end\n",
|
||||
" C = input_distance_table(N)\n",
|
||||
" C_par = copy(C)\n",
|
||||
" floyd_mpi!(C_par,comm)\n",
|
||||
" if rank == 0\n",
|
||||
" C_seq = copy(C)\n",
|
||||
" floyd!(C_seq)\n",
|
||||
" if C_seq == C_par\n",
|
||||
" println(\"Test passed 🥳\")\n",
|
||||
" else\n",
|
||||
" println(\"Test failed\")\n",
|
||||
" end\n",
|
||||
" end\n",
|
||||
"end\n",
|
||||
"run(`$(mpiexec()) -np 3 julia --project=. -e $code`);"
|
||||
]
|
||||
},
|
||||
{
|
||||
@@ -909,10 +1022,12 @@
|
||||
"id": "91a772df",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"### Is this implementation correct?\n",
|
||||
"## Is this implementation correct?\n",
|
||||
"\n",
|
||||
"Point-to-point messages are *non-overtaking* (i.e. FIFO order) between the specified sender and receiver according to section 3.5 of the MPI standard 4.0.\n",
|
||||
"Unfortunately this is not enough in this case. The messages can still arrive in the wrong order if messages from different processes overtake each other."
|
||||
"In the cell above, the result of the parallel code was probably the same as for the sequential code. However, is this sufficient to assert that the code is correct? Unfortunately, it is not. In fact, **the parallel code we implemented is not correct!** There is no guarantee that this code computes the correct result. This is why:\n",
|
||||
"\n",
|
||||
"In MPI, point-to-point messages are *non-overtaking* between a given sender and receiver. Say that process 1 sends several messages to process 3. All these will arrive in FIFO order. This is according to section 3.5 of the MPI standard 4.0.\n",
|
||||
"Unfortunately, this is not enough in our case. The messages could arrive in the wrong order *from different senders*. If process 1 sends messages to process 3, and then process 2 sends other messages to process 3, it is not granted that process 3 will receive first the messages from process 1 and then from process 2 (see figure below)."
|
||||
]
|
||||
},
|
||||
{
|
||||
@@ -920,7 +1035,7 @@
|
||||
"id": "1277772d",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"If we are lucky all messages will arrive in order and we will process all rows in the right order in all processors."
|
||||
"If we are lucky all messages will arrive in order. In our parallel code, all processors would receive first row one, then row 2, then row 3, etc. The computed result will be correct."
|
||||
]
|
||||
},
|
||||
{
|
||||
@@ -943,7 +1058,7 @@
|
||||
"id": "df60e4e7",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"However, FIFO ordering is not enough. In the next figure, communication between process 1 and process 3 is particularly slow. Note that process 3 receives messages from process 1 after it receives the messages from 2 even though FIFO ordering is satisfied between any two processors."
|
||||
"However, FIFO ordering between pairs of processors is not enough to guarantee that rows arrive in consecutive order. The next figure shows a counter example. In this case, communication between process 1 and process 3 is particularly slow for some unknown reason. As result processor 3 receives first messages from processor 2, even though processor 1 sent the messages first. In our parallel code, the received rows would be first row 3, then row 4, then row 1, then row 2, which is not correct. Note however that process 3 received all messages from process 1 in the correct order (guaranteed by MPI). But this is not enough in our algorithm."
|
||||
]
|
||||
},
|
||||
{
|
||||
@@ -968,10 +1083,25 @@
|
||||
"source": [
|
||||
"### Possible solutions\n",
|
||||
"\n",
|
||||
"There are several solution to this synchronization problem:\n",
|
||||
"\n",
|
||||
"1. **Synchronous sends**: Use synchronous send MPI_SSEND. This is less efficient because we spend time waiting until each message is received. Note that the blocking send MPI_SEND used above does not guarantee that the message was received. \n",
|
||||
"2. **MPI.Barrier**: Use a barrier at the end of each iteration over $k$. This is easy to implement, but we get a synchronization overhead.\n",
|
||||
"3. **Order incoming messages**: The receiver orders the incoming messages, e.g. according to MPI.Status or the sender rank. This requires buffering and extra user code.\n",
|
||||
"4. **MPI.Bcast!**: Communicate row k using `MPI.Bcast!`. One needs to know which are the rows owned by the other ranks."
|
||||
"3. **Order incoming messages**: The receiver orders the incoming messages, e.g. using to MPI.Status to get the sender rank. This requires buffering and extra user code.\n",
|
||||
"4. **MPI.Bcast!**: Communicate row k using `MPI.Bcast!`. One needs to know which are the rows owned by the other ranks since we cannot use `MPI.ANY_SOURCE` in `MPI.Bcast!`. This is trivial however if the number of rows is multiple of the number of ranks."
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"id": "fb21edc4",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"## Summary\n",
|
||||
"\n",
|
||||
"- We learned how to parallelize Floyd's algorithm.\n",
|
||||
"- The considered strategy based on a row-wise data partition has little communication overhead if the problem size is large enough.\n",
|
||||
"- One needs to be careful in which order the messages are received to have a correct algorithm.\n",
|
||||
"- There are several strategies to solve this synchronization problem, each one with pros and cons.\n"
|
||||
]
|
||||
},
|
||||
{
|
||||
@@ -979,61 +1109,11 @@
|
||||
"id": "de96ad1b",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"## Exercise \n",
|
||||
"Rewrite the worker code of the parallel ASP algorithm so it runs correctly. Use the `MPI.Bcast!` to solve the problem of overtaking messages. Note: Only use `MPI.Bcast!`, do not use other MPI directives in addition. You can test your function with the following code cell. "
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": null,
|
||||
"id": "31194529",
|
||||
"metadata": {},
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"function floyd_par!(C,N)\n",
|
||||
" comm = MPI.Comm_dup(MPI.COMM_WORLD)\n",
|
||||
" nranks = MPI.Comm_size(comm)\n",
|
||||
" rank = MPI.Comm_rank(comm)\n",
|
||||
" T = eltype(C)\n",
|
||||
" if rank == 0\n",
|
||||
" buffer_root = Vector{T}(undef,N*N)\n",
|
||||
" buffer_root[:] = transpose(C)[:]\n",
|
||||
" else\n",
|
||||
" buffer_root = Vector{T}(undef,0)\n",
|
||||
" end \n",
|
||||
" Nw = div(N,nranks)\n",
|
||||
" buffer = Vector{T}(undef,Nw*N)\n",
|
||||
" MPI.Scatter!(buffer_root,buffer,comm;root=0)\n",
|
||||
" Cw = Matrix{T}(undef,Nw,N)\n",
|
||||
" transpose(Cw)[:] = buffer\n",
|
||||
" MPI.Barrier(comm)\n",
|
||||
" floyd_worker_bcast!(Cw,comm)\n",
|
||||
" buffer[:] = transpose(Cw)[:]\n",
|
||||
" MPI.Gather!(buffer,buffer_root,comm;root=0)\n",
|
||||
" if rank == 0\n",
|
||||
" transpose(C)[:] = buffer_root[:]\n",
|
||||
" end\n",
|
||||
" C\n",
|
||||
"end\n",
|
||||
"## Exercise\n",
|
||||
"\n",
|
||||
"@everywhere function floyd_worker_bcast!(Cw,comm)\n",
|
||||
" # Your implementation here\n",
|
||||
"end\n"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": null,
|
||||
"id": "1b7eb4c2",
|
||||
"metadata": {},
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"load = 10\n",
|
||||
"n = nworkers()*load\n",
|
||||
"C = rand_distance_table(n)\n",
|
||||
"C_seq = floyd!(copy(C))\n",
|
||||
"C_par = floyd_par!(copy(C),n)\n",
|
||||
"@test C_seq == C_par"
|
||||
"### Exercise 1\n",
|
||||
"\n",
|
||||
"Modify the `floyd_iterations!` function so that it is guaranteed that the result is computed correctly. Use `MPI.Bcast!` to solve the synchronization problem. Note: only use `MPI.Bcast!`in `floyd_iterations!`, do not use other MPI directives. You can assume that the number of rows is a multiple of the number of processes."
|
||||
]
|
||||
},
|
||||
{
|
||||
|
||||
+6773
-1212
File diff suppressed because it is too large.
Load diff
|
Before Width: | Height: | Size: 6.6 MiB After Width: | Height: | Size: 6.8 MiB |
File diff suppressed because it is too large.
Load diff
|
Before Width: | Height: | Size: 248 KiB After Width: | Height: | Size: 445 KiB |
Reference in new issue
Block a user