{"id":11120,"date":"2024-03-28T16:00:00","date_gmt":"2024-03-28T20:00:00","guid":{"rendered":"https:\/\/blogs.mathworks.com\/cleve\/?p=11120"},"modified":"2024-03-28T17:33:53","modified_gmt":"2024-03-28T21:33:53","slug":"closest-pair-of-points-problem","status":"publish","type":"post","link":"https:\/\/blogs.mathworks.com\/cleve\/2024\/03\/28\/closest-pair-of-points-problem\/","title":{"rendered":"Closest Pair of Points Problem"},"content":{"rendered":"<div class=\"content\"><!--introduction--><p>The Closest Pair of Points problem is a standard topic in an algorithms course today, but when I taught such a course fifty years ago, the algorithm was not yet known.<\/p><!--\/introduction--><h3>Contents<\/h3><div><ul><li><a href=\"#0f62a942-c8a3-4395-8f83-9509ed1dabd0\">California Dreaming<\/a><\/li><li><a href=\"#825b1d10-1532-4f6e-8927-6d9e45dbdb49\">Closest Pair of Points<\/a><\/li><li><a href=\"#92c53c07-0f81-494a-91a7-c364eb3132ea\"><tt>Pairs<\/tt><\/a><\/li><li><a href=\"#05e723c3-e7a0-4091-bb97-463ac017b7e0\"><tt>DivCon<\/tt><\/a><\/li><li><a href=\"#387c9a5e-eee6-47d0-be8a-6ea3d8d95f8f\"><tt>Center<\/tt><\/a><\/li><li><a href=\"#53330f5c-f9a0-4060-ba34-e5e98e1c9ed4\">Complexity<\/a><\/li><li><a href=\"#5d10f4df-79d5-46b8-9375-01260abdba40\">Timing<\/a><\/li><li><a href=\"#4de823af-8e67-4218-a471-d8def3204d38\">Software<\/a><\/li><li><a href=\"#0dc6f3e8-7d3a-4381-8d34-13514fb2cad1\">References<\/a><\/li><\/ul><\/div><h4>California Dreaming<a name=\"0f62a942-c8a3-4395-8f83-9509ed1dabd0\"><\/a><\/h4><p>Imagine you are driving a car on the Harbor Freeway in southern California with typical Los Angeles traffic conditions. Among the many things you might want to know is which pair of vehicles is nearest each other.<\/p><p>This is an instance of the Closest Pair of Points problem:<\/p><div><ul><li>Given the location of <tt>n<\/tt> points in the plane, which pair of points is closest to each other?<\/li><\/ul><\/div><p><img decoding=\"async\" vspace=\"5\" hspace=\"5\" src=\"https:\/\/blogs.mathworks.com\/cleve\/files\/traffic.png\" alt=\"\"> <\/p><h4>Closest Pair of Points<a name=\"825b1d10-1532-4f6e-8927-6d9e45dbdb49\"><\/a><\/h4><p>It is convenient to represent the points by a vector of complex values. The distance between points <tt>z(k)<\/tt> and <tt>z(j)<\/tt> is then<\/p><pre>  d = abs(z(k) - z(j))<\/pre><p>Here are a few points in the unit square.  The closest pair is highlighted.<\/p><p><img decoding=\"async\" vspace=\"5\" hspace=\"5\" src=\"https:\/\/blogs.mathworks.com\/cleve\/files\/pairs.png\" alt=\"\"> <\/p><h4><tt>Pairs<\/tt><a name=\"92c53c07-0f81-494a-91a7-c364eb3132ea\"><\/a><\/h4><p>The first algorithm you might think of computes the distances between all possible pairs of points and finds the minimum.  This is a brute force approach that requires only a few lines of code.<\/p><pre>function d = Pairs(z)\r\n    % Pairs.\r\n    % d = Pairs(z) is the minimum distance between any two elements\r\n    % of the complex vector z.<\/pre><pre>    n = length(z);\r\n    d = Inf;\r\n    for k = 1:n\r\n        for j = k+1:n\r\n            if abs(z(k) - z(j)) &lt; d\r\n                d = abs(z(k) - z(j));\r\n            end\r\n        end\r\n    end\r\nend<\/pre><h4><tt>DivCon<\/tt><a name=\"05e723c3-e7a0-4091-bb97-463ac017b7e0\"><\/a><\/h4><p>DivCon stands for Divide and Conquer.  In outline, the steps are:<\/p><div><ul><li>Divide the set of points into two halves.<\/li><\/ul><\/div><div><ul><li>Recursively, find the closest pair in each half.<\/li><\/ul><\/div><div><ul><li>Consider the case when the closest pair has one point in each half.<\/li><\/ul><\/div><div><ul><li>Terminate the recursion with sets of two or three points.<\/li><\/ul><\/div><pre>function d = DivCon(z,sorted)\r\n    % DivCon.\r\n    % d = DivCon(z) is the minimum distance between any two elements\r\n    % of the complex vector z.\r\n    %\r\n    % d = DivCon(z,true) is a recursive call with ascending real(z).<\/pre><pre>    n = length(z);\r\n    if n &lt;= 3\r\n       d = Pairs(z);\r\n    else\r\n       if nargin &lt; 2 || ~sorted\r\n           [~,p] = sort(real(z));\r\n           z = z(p);\r\n       end\r\n       m = floor(n\/2);<\/pre><pre>       % Left half\r\n       dl = DivCon(z(1:m),true)<\/pre><pre>       % Right half\r\n       dr = DivCon(z(m+1:end),true);<\/pre><pre>       % Choose\r\n       d = min(dl,dr);<\/pre><pre>       % Center strip\r\n       ds = Center(z,d);\r\n       d = min(ds,d);\r\n    end\r\nend<\/pre><h4><tt>Center<\/tt><a name=\"387c9a5e-eee6-47d0-be8a-6ea3d8d95f8f\"><\/a><\/h4><p>The delicate case involves the strip of points near the center dividing line.  The width of the strip is the closest distance found in the recursion.  Any closer pair with one point in each half must be in this strip.<\/p><pre>function d = Center(z,d)\r\n    % Center(z,d) is used by DivCon to examine the\r\n    % strip of half-width d about the center point.<\/pre><pre>    n = length(z)\r\n    m = floor(n\/2);\r\n    xh = real(z(m));\r\n    [~,p] = sort(imag(z));\r\n    z = z(p);\r\n    s = [];\r\n    for i = 1:n\r\n        if abs(real(z(i)) - xh) &lt;  d\r\n            s = [s; z(i)];\r\n        end\r\n    end<\/pre><pre>    ns = length(s);\r\n    for k = 1:ns\r\n        for j = k+1:ns\r\n            if (imag(s(j)) - imag(s(k))) &lt; d &amp;&amp; abs(s(k) - s(j)) &lt; d\r\n                d = abs(s(k) - s(j));\r\n            end\r\n        end\r\n    end\r\nend<\/pre><h4>Complexity<a name=\"53330f5c-f9a0-4060-ba34-e5e98e1c9ed4\"><\/a><\/h4><p>Let <tt>n<\/tt> be the number of points.  An asymptotic execution-time complexity analysis involves <tt>n<\/tt> approaching infinity.<\/p><p>It is not hard to see that the complexity of the brute force algorithm implemented in <tt>Pairs<\/tt> is O(<tt>n<\/tt>^2).<\/p><p>There are <a href=\"https:\/\/www.google.com\/search?q=closest+pair+of+points+(divide+and+conquer)\">dozens of pages on the web<\/a> devoted to showing that the complexity of the divide and conquer algorithm implemented in <tt>DivCon<\/tt> and <tt>Center<\/tt> is O(<tt>n<\/tt>*log(<tt>n<\/tt>)).  The best page that I have seen is the <a href=\"https:\/\/www.youtube.com\/watch?v=6u_hWxbOc7E\">YouTube video by Ling Qi<\/a>. The key to the analysis is showing that the inner loop in <tt>Center<\/tt> is executed at most 7 times for any <tt>n<\/tt>.<\/p><h4>Timing<a name=\"5d10f4df-79d5-46b8-9375-01260abdba40\"><\/a><\/h4><p>We measured the execution time of <tt>Pairs(z)<\/tt> and <tt>DivCon(z)<\/tt> for <tt>n<\/tt> from 1,000 to 40,000 and computed the ratios of the two times. The complexity analysis predicts that this ratio is asymptotically<\/p><pre>  O(n\/log(n))<\/pre><p>Here are the timing results and a least square fit by <tt>n<\/tt>\/log(<tt>n<\/tt>).<\/p><p><img decoding=\"async\" vspace=\"5\" hspace=\"5\" src=\"https:\/\/blogs.mathworks.com\/cleve\/files\/fit.png\" alt=\"\"> <\/p><h4>Software<a name=\"4de823af-8e67-4218-a471-d8def3204d38\"><\/a><\/h4><p>A self-extracting MATLAB archive is available at <a href=\"https:\/\/blogs.mathworks.com\/cleve\/files\/TestDivCon_mzip.m\">https:\/\/blogs.mathworks.com\/cleve\/files\/TestDivCon_mzip.m<\/a><\/p><h4>References<a name=\"0dc6f3e8-7d3a-4381-8d34-13514fb2cad1\"><\/a><\/h4><p>Ling Qi, IDeer7, <i>Closest Pair of Points (Divide and Conquer) Explained<\/i>. <a href=\"https:\/\/www.youtube.com\/watch?v=6u_hWxbOc7E\">https:\/\/www.youtube.com\/watch?v=6u_hWxbOc7E<\/a>.<\/p><p>Cormen, Thomas H.; Leiserson, Charles E.; Rivest, Ronald L.; Stein, Clifford. <i>Introduction to Algorithms (4th ed.)<\/i>. MIT Press and McGraw-Hill. ISBN 0-262-04630-X. 1312 pp.<\/p><script language=\"JavaScript\"> <!-- \r\n    function grabCode_fc0430dede624952a83439140fe547ff() {\r\n        \/\/ Remember the title so we can use it in the new page\r\n        title = document.title;\r\n\r\n        \/\/ Break up these strings so that their presence\r\n        \/\/ in the Javascript doesn't mess up the search for\r\n        \/\/ the MATLAB code.\r\n        t1='fc0430dede624952a83439140fe547ff ' + '##### ' + 'SOURCE BEGIN' + ' #####';\r\n        t2='##### ' + 'SOURCE END' + ' #####' + ' fc0430dede624952a83439140fe547ff';\r\n    \r\n        b=document.getElementsByTagName('body')[0];\r\n        i1=b.innerHTML.indexOf(t1)+t1.length;\r\n        i2=b.innerHTML.indexOf(t2);\r\n \r\n        code_string = b.innerHTML.substring(i1, i2);\r\n        code_string = code_string.replace(\/REPLACE_WITH_DASH_DASH\/g,'--');\r\n\r\n        \/\/ Use \/x3C\/g instead of the less-than character to avoid errors \r\n        \/\/ in the XML parser.\r\n        \/\/ Use '\\x26#60;' instead of '<' so that the XML parser\r\n        \/\/ doesn't go ahead and substitute the less-than character. \r\n        code_string = code_string.replace(\/\\x3C\/g, '\\x26#60;');\r\n\r\n        copyright = 'Copyright 2024 The MathWorks, Inc.';\r\n\r\n        w = window.open();\r\n        d = w.document;\r\n        d.write('<pre>\\n');\r\n        d.write(code_string);\r\n\r\n        \/\/ Add copyright line at the bottom if specified.\r\n        if (copyright.length > 0) {\r\n            d.writeln('');\r\n            d.writeln('%%');\r\n            if (copyright.length > 0) {\r\n                d.writeln('% _' + copyright + '_');\r\n            }\r\n        }\r\n\r\n        d.write('<\/pre>\\n');\r\n\r\n        d.title = title + ' (MATLAB code)';\r\n        d.close();\r\n    }   \r\n     --> <\/script><p style=\"text-align: right; font-size: xx-small; font-weight:lighter;   font-style: italic; color: gray\"><br><a href=\"javascript:grabCode_fc0430dede624952a83439140fe547ff()\"><span style=\"font-size: x-small;        font-style: italic;\">Get \r\n      the MATLAB code <noscript>(requires JavaScript)<\/noscript><\/span><\/a><br><br>\r\n      Published with MATLAB&reg; R2023a<br><\/p><\/div><!--\r\nfc0430dede624952a83439140fe547ff ##### SOURCE BEGIN #####\r\n%% The Closest Pair of Points Problem\r\n% The Closest Pair of Points problem is a standard topic in an\r\n% algorithms course today, but when I taught such \r\n% a course fifty years ago, the algorithm was not yet known.\r\n\r\n%% California Dreaming \r\n% Imagine you are driving a car on the Harbor Freeway in southern \r\n% California with typical Los Angeles traffic conditions.\r\n% Among the many things you might want to know is which pair of \r\n% vehicles is nearest each other.  \r\n%\r\n% This is an instance of the Closest Pair of Points problem:\r\n%\r\n% * Given the location of |n| points in the plane, which pair of points\r\n% is closest to each other?\r\n%\r\n% <<traffic.png>>\r\n%\r\n%% Closest Pair of Points\r\n% It is convenient to represent the points by a vector of complex values. \r\n% The distance between points |z(k)| and |z(j)| is then\r\n%\r\n%    d = abs(z(k) - z(j))\r\n%\r\n% Here are a few points in the unit square.  The closest pair is\r\n% highlighted.\r\n%\r\n% <<pairs.png>>\r\n%\r\n\r\n%% |Pairs|\r\n% The first algorithm you might think of computes the distance\r\n% between all possible pairs of points and finds the minimum.  This is\r\n% a brute force approach that requires only a few lines of code.\r\n%\r\n%  function d = Pairs(z)\r\n%      % Pairs. \r\n%      % d = Pairs(z) is the minimum distance between any two elements\r\n%      % of the complex vector z.  \r\n%\r\n%      n = length(z);\r\n%      d = Inf; \r\n%      for k = 1:n\r\n%          for j = k+1:n\r\n%              if abs(z(k) - z(j)) < d\r\n%                  d = abs(z(k) - z(j));\r\n%              end    \r\n%          end\r\n%      end\r\n%  end\r\n\r\n%% |DivCon|\r\n% DivCon stands for Divide and Conquer.  In outline, the steps are:\r\n%\r\n% * Divide the set of points into two halves.\r\n%\r\n% * Recursively, find the closest pair in each half.\r\n%\r\n% * Consider the case when the closest pair has one point in each half. \r\n%\r\n% * Terminate the recursion with sets of two or three points.\r\n%\r\n%  function d = DivCon(z,sorted)\r\n%      % DivCon. \r\n%      % d = DivCon(z) is the minimum distance between any two elements\r\n%      % of the complex vector z. \r\n%      %\r\n%      % d = DivCon(z,true) is a recursive call with ascending real(z).\r\n%\r\n%      n = length(z);\r\n%      if n <= 3\r\n%         d = Pairs(z);\r\n%      else\r\n%         if nargin < 2 || ~sorted\r\n%             [~,p] = sort(real(z));\r\n%             z = z(p);\r\n%         end\r\n%         m = floor(n\/2);\r\n%   \r\n%         % Left half\r\n%         dl = DivCon(z(1:m),true)\r\n%   \r\n%         % Right half\r\n%         dr = DivCon(z(m+1:end),true);\r\n% \r\n%         % Choose\r\n%         d = min(dl,dr);\r\n%   \r\n%         % Center strip\r\n%         ds = Center(z,d);\r\n%         d = min(ds,d);\r\n%      end\r\n%  end\r\n\r\n%% |Center|\r\n% The delicate case involves the strip of points near the \r\n% center dividing line.  The width of the strip is the closest\r\n% distance found in the recursion.  Any closer pair with one point\r\n% in each half must be in this strip.\r\n%\r\n%  function d = Center(z,d)    \r\n%      % Center(z,d) is used by DivCon to examine the\r\n%      % strip of half-width d about the center point.\r\n%\r\n%      n = length(z) \r\n%      m = floor(n\/2);\r\n%      xh = real(z(m));\r\n%      [~,p] = sort(imag(z));\r\n%      z = z(p);\r\n%      s = [];\r\n%      for i = 1:n\r\n%          if abs(real(z(i)) - xh) <  d\r\n%              s = [s; z(i)];\r\n%          end\r\n%      end\r\n%\r\n%      ns = length(s);\r\n%      for k = 1:ns\r\n%          for j = k+1:ns\r\n%              if (imag(s(j)) - imag(s(k))) < d && abs(s(k) - s(j)) < d\r\n%                  d = abs(s(k) - s(j));\r\n%              end\r\n%          end\r\n%      end\r\n%  end\r\n\r\n%% Complexity\r\n% Let |n| be the number of points.  An asymptotic execution-time \r\n% complexity analysis involves |n| approaching infinity.\r\n%\r\n% It is not hard to see that the complexity of the brute force algorithm\r\n% implemented in |Pairs| is O(|n|^2).  \r\n%\r\n% There are \r\n% <https:\/\/www.google.com\/search?q=closest+pair+of+points+(divide+and+conquer)\r\n% dozens of pages on the web> devoted to showing that the complexity\r\n% of the divide and conquer algorithm implemented in |DivCon| and |Center|\r\n% is O(|n|*log(|n|)).  The best page that I have seen is the \r\n% <https:\/\/www.youtube.com\/watch?v=6u_hWxbOc7E\r\n% YouTube video by Ling Qi>.\r\n% The key to the analysis is showing that the inner\r\n% loop in |Center| is executed at most 7 times for any |n|.  \r\n% \r\n%% Timing\r\n% We measured the execution time of |Pairs(z)| and |DivCon(z)| for\r\n% |n| from 1,000 to 40,000 and computed the ratios of the two times.\r\n% The complexity analysis predicts that this ratio is asymptotically\r\n%\r\n%    O(n\/log(n))\r\n%\r\n% Here are the timing results and a least square fit by |n|\/log(|n|).\r\n%\r\n% <<fit.png>>\r\n\r\n%% Software\r\n% A self-extracting MATLAB archive is available at\r\n% <https:\/\/blogs.mathworks.com\/cleve\/files\/TestDivCon_mzip.m>\r\n\r\n%% References\r\n%\r\n% Ling Qi, IDeer7, _Closest Pair of Points (Divide and Conquer) Explained_.\r\n% <https:\/\/www.youtube.com\/watch?v=6u_hWxbOc7E>.\r\n%\r\n% Cormen, Thomas H.; Leiserson, Charles E.; Rivest, Ronald L.; \r\n% Stein, Clifford. _Introduction to Algorithms (4th ed.)_.\r\n% MIT Press and McGraw-Hill. ISBN 0-262-04630-X. 1312 pp.\r\n##### SOURCE END ##### fc0430dede624952a83439140fe547ff\r\n-->","protected":false},"excerpt":{"rendered":"<div class=\"overview-image\"><img src=\"https:\/\/blogs.mathworks.com\/cleve\/files\/pairs.png\" class=\"img-responsive attachment-post-thumbnail size-post-thumbnail wp-post-image\" alt=\"\" decoding=\"async\" loading=\"lazy\" \/><\/div><!--introduction--><p>The Closest Pair of Points problem is a standard topic in an algorithms course today, but when I taught such a course fifty years ago, the algorithm was not yet known.... <a class=\"read-more\" href=\"https:\/\/blogs.mathworks.com\/cleve\/2024\/03\/28\/closest-pair-of-points-problem\/\">read more >><\/a><\/p>","protected":false},"author":78,"featured_media":11090,"comment_status":"open","ping_status":"closed","sticky":false,"template":"","format":"standard","meta":[],"categories":[11,5,4,14],"tags":[],"_links":{"self":[{"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/posts\/11120"}],"collection":[{"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/posts"}],"about":[{"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/types\/post"}],"author":[{"embeddable":true,"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/users\/78"}],"replies":[{"embeddable":true,"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/comments?post=11120"}],"version-history":[{"count":3,"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/posts\/11120\/revisions"}],"predecessor-version":[{"id":11129,"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/posts\/11120\/revisions\/11129"}],"wp:featuredmedia":[{"embeddable":true,"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/media\/11090"}],"wp:attachment":[{"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/media?parent=11120"}],"wp:term":[{"taxonomy":"category","embeddable":true,"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/categories?post=11120"},{"taxonomy":"post_tag","embeddable":true,"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/tags?post=11120"}],"curies":[{"name":"wp","href":"https:\/\/api.w.org\/{rel}","templated":true}]}}