From: Martin Kronbichler <martin.kronbichler@it.uu.se>
Date: Mon, 13 Dec 2021 13:58:55 +0000 (+0100)
Subject: Avoid quadratic complexity in IndexSet::subtract_set
X-Git-Tag: v9.4.0-rc1~746^2
X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=e8229cb9b620ac56922a1d640ba3a12968f98bdd;p=dealii.git

Avoid quadratic complexity in IndexSet::subtract_set
---

diff --git a/doc/news/changes/minor/20211213MartinKronbichler b/doc/news/changes/minor/20211213MartinKronbichler
new file mode 100644
index 0000000000..d141a400ff
--- /dev/null
+++ b/doc/news/changes/minor/20211213MartinKronbichler
@@ -0,0 +1,4 @@
+Fixed: IndexSet::subtract_set would previously run in quadratic complexity in
+the number of ranges. This is now fixed.
+<br>
+(Martin Kronbichler, 2021/12/13)
diff --git a/source/base/index_set.cc b/source/base/index_set.cc
index 1383dd9495..774ddfacb6 100644
--- a/source/base/index_set.cc
+++ b/source/base/index_set.cc
@@ -270,7 +270,7 @@ IndexSet::subtract_set(const IndexSet &other)
   is_compressed = false;
 
 
-  // we save new ranges to be added to our IndexSet in an temporary vector and
+  // we save all new ranges to our IndexSet in an temporary vector and
   // add all of them in one go at the end.
   std::vector<Range> new_ranges;
 
@@ -282,6 +282,7 @@ IndexSet::subtract_set(const IndexSet &other)
       // advance own iterator until we get an overlap
       if (own_it->end <= other_it->begin)
         {
+          new_ranges.push_back(*own_it);
           ++own_it;
           continue;
         }
@@ -314,15 +315,11 @@ IndexSet::subtract_set(const IndexSet &other)
       // next.
     }
 
-  // Now delete all empty ranges we might
-  // have created.
-  for (std::vector<Range>::iterator it = ranges.begin(); it != ranges.end();)
-    {
-      if (it->begin >= it->end)
-        it = ranges.erase(it);
-      else
-        ++it;
-    }
+  // make sure to take over the remaining ranges
+  for (; own_it != ranges.end(); ++own_it)
+    new_ranges.push_back(*own_it);
+
+  ranges.clear();
 
   // done, now add the temporary ranges
   const std::vector<Range>::iterator end = new_ranges.end();