Merge branch 'master' into delete-indices-dev

general/array.hpp

@@ -190,6 +190,12 @@ public:
    /// Prepend an 'el' to the array, resize if necessary.
    inline int Prepend(const T &el);
 
+   /// Insert @a els into the array at index @a i
+   inline int Insert(int i, const Array<T> &els);
+
+   /// Insert @a el into the array at index @a i
+   inline int Insert(int i, const T &el) { return Insert(i, Array<T>({el})); }
+
    /// Return the last element in the array.
    inline T &Last();
 
@@ -211,6 +217,9 @@ public:
    /// Delete the first entry with value == 'el'.
    inline void DeleteFirst(const T &el);
 
+   /// Delete entries at @a indices, and resize
+   inline void DeleteAt(const Array<int> &indices);
+
    /// Delete the whole array.
    inline void DeleteAll();
 
@@ -249,6 +258,9 @@ public:
    /// Copy sub array starting from @a offset out to the provided @a sa.
    inline void GetSubArray(int offset, int sa_size, Array<T> &sa) const;
 
+   /// Set from sub array @sa at @a offset
+   inline void SetSubArray(int offset, const Array<T> &sa);
+
    /// Prints array to stream with width elements per row.
    void Print(std::ostream &out = mfem::out, int width = 4) const;
 
@@ -873,6 +885,22 @@ inline int Array<T>::Prepend(const T &el)
    return size;
 }
 
+template<class T>
+inline int Array<T>::Insert(int i, const Array<T> &els)
+{
+   MFEM_ASSERT(i < size, "Insert index is out-of-bounds.");
+
+   const int old_size = size;
+
+   SetSize(size + els.Size());
+   for (int j = old_size-1; j >= i; j--)
+   {
+      data[j+els.Size()] = data[j];
+   }
+   SetSubArray(i, els);
+   return size;
+}
+
 template <class T>
 inline T &Array<T>::Last()
 {
@@ -935,6 +963,30 @@ inline void Array<T>::DeleteFirst(const T &el)
    }
 }
 
+template <class T>
+inline void Array<T>::DeleteAt(const Array<int> &indices)
+{
+   // Make a copy of the indices, sorted.
+   Array<int> sorted_indices(indices);
+   sorted_indices.Sort();
+
+   int rm_count = 0;
+   for (int i = 0; i < size; i++)
+   {
+      if (rm_count < sorted_indices.Size() && i == sorted_indices[rm_count])
+      {
+         rm_count++;
+      }
+      else
+      {
+         data[i-rm_count] = data[i]; // shift data rm_count
+      }
+   }
+
+   // Resize to remove tail
+   SetSize(size - rm_count);
+}
+
 template <class T>
 inline void Array<T>::DeleteAll()
 {
@@ -987,6 +1039,18 @@ inline void Array<T>::GetSubArray(int offset, int sa_size, Array<T> &sa) const
    }
 }
 
+template<class T>
+inline void Array<T>::SetSubArray(int offset, const Array<T> &sa)
+{
+   MFEM_ASSERT(offset + sa.Size() < size,
+               "Sub-array with size " << sa.Size() << " is too large to set at offset " <<
+               offset << ", given array size " << size);
+   for (int i = 0; i < sa.Size(); i++)
+   {
+      data[offset + i] = sa[i];
+   }
+}
+
 template <class T>
 inline void Array<T>::operator=(const T &a)
 {

linalg/vector.cpp

@@ -668,6 +668,28 @@ void Vector::median(const Vector &lo, const Vector &hi)
    });
 }
 
+void Vector::Insert(int offset, const Vector &sv)
+{
+   const int old_size = size;
+
+   if (sv.Size() + old_size > Capacity())
+   {
+      Vector copy = *this;
+      SetSize(size + sv.Size());
+      SetVector(copy, 0);
+   }
+   else
+   {
+      SetSize(size + sv.Size());
+   }
+
+   for (int j = old_size-1; j >= offset; j--)
+   {
+      data[j+sv.Size()] = data[j];
+   }
+   SetVector(sv, offset);
+}
+
 void Vector::GetSubVector(const Array<int> &dofs, Vector &elemvect) const
 {
    const int n = dofs.Size();

linalg/vector.hpp

@@ -171,6 +171,9 @@ public:
    /// Resize the vector to size @a s using the MemoryType of @a v.
    void SetSize(int s, const Vector &v) { SetSize(s, v.GetMemory().GetMemoryType()); }
 
+   /// Delete elements at @a indices and resize vector accordingly
+   void DeleteAt(const Array<int> &indices);
+
    /// Set the Vector data.
    /// @warning This method should be called only when OwnsData() is false.
    void SetData(real_t *d) { data.Wrap(d, data.Capacity(), false); }
@@ -396,6 +399,12 @@ public:
    /// v = median(v,lo,hi) entrywise.  Implementation assumes lo <= hi.
    void median(const Vector &lo, const Vector &hi);
 
+   /// Insert sub Vector @a sv at @a offset and resize
+   void Insert(int offset, const Vector &sv);
+
+   /// Insert @a value at @a offset and resize
+   void Insert(int offset, const real_t value) { Insert(offset, Vector({value})); }
+
    /// Extract entries listed in @a dofs to the output Vector @a elemvect.
    /** Negative dof values cause the -dof-1 position in @a elemvect to receive
        the -val in from this Vector. */
@@ -621,6 +630,29 @@ inline void Vector::SetSize(int s, MemoryType mt)
    data.UseDevice(use_dev);
 }
 
+inline void Vector::DeleteAt(const Array<int> &indices)
+{
+   // Make copy of the indices, sorted.
+   Array<int> sorted_indices(indices);
+   sorted_indices.Sort();
+
+   int rm_count = 0;
+   for (int i = 0; i < size; i++)
+   {
+      if (rm_count < sorted_indices.Size() && i == sorted_indices[rm_count])
+      {
+         rm_count++;
+      }
+      else
+      {
+         data[i-rm_count] = data[i]; // shift data rm_count
+      }
+   }
+
+   // Resize to remove tail
+   SetSize(size - rm_count);
+}
+
 inline void Vector::NewMemoryAndSize(const Memory<real_t> &mem, int s,
                                      bool own_mem)
 {

tests/unit/general/test_array.cpp

@@ -124,3 +124,51 @@ TEST_CASE("Array stl-interactions", "[Array]")
       CHECK(x[i] == y[i]);
    }
 }
+
+TEST_CASE("Array delete at indices", "[Array]")
+{
+   Array<int>       test({0,1,2,3,4,5,6,7,8});
+   Array<int> rm_indices({0,    3,4,  6,  8});
+   Array<int>     result({  1,2,    5,  7  });
+
+   test.DeleteAt(rm_indices);
+
+   REQUIRE(test.Size() == result.Size());
+
+   for (int i = 0; i < test.Size(); i++)
+   {
+      CHECK(test[i] == result[i]);
+   }
+}
+
+TEST_CASE("Array set sub array", "[Array]")
+{
+   Array<int>   test({0,1,2,3,4,5,6,7,8});
+   Array<int>     sa({      9,9,9      });
+   Array<int> result({0,1,2,9,9,9,6,7,8});
+
+   test.SetSubArray(3, sa);
+
+   for (int i = 0; i < test.Size(); i++)
+   {
+      CHECK(test[i] == result[i]);
+   }
+}
+
+TEST_CASE("Array insert", "[Array]")
+{
+   Array<int>   test({0,1,2,3,      7,8});
+   Array<int>     sa({        4,5,     });
+   Array<int> result({0,1,2,3,4,5,6,7,8});
+
+   test.Insert(4, sa);
+   test.Insert(6, 6);
+
+   REQUIRE(test.Size() == result.Size());
+
+   for (int i = 0; i < test.Size(); i++)
+   {
+      CHECK(test[i] == result[i]);
+   }
+
+}

tests/unit/linalg/test_vector.cpp

@@ -247,3 +247,36 @@ TEST_CASE("Vector Sum", "[Vector],[GPU]")
 
    REQUIRE(sum_1 == MFEM_Approx(sum_2));
 }
+
+TEST_CASE("Vector delete at indices", "[Vector]")
+{
+   Vector           test({0,1,2,3,4,5,6,7,8});
+   Array<int> rm_indices({0,    3,4,  6,  8});
+   Vector         result({  1,2,    5,  7  });
+
+   test.DeleteAt(rm_indices);
+
+   REQUIRE(test.Size() == result.Size());
+
+   for (int i = 0; i < test.Size(); i++)
+   {
+      CHECK(test[i] == result[i]);
+   }
+}
+
+TEST_CASE("Vector insert", "[Vector]")
+{
+   Vector   test({0,1,2,3,      7,8});
+   Vector     sa({        4,5,     });
+   Vector result({0,1,2,3,4,5,6,7,8});
+
+   test.Insert(4, sa);
+   test.Insert(6, 6);
+
+   REQUIRE(test.Size() == result.Size());
+
+   for (int i = 0; i < test.Size(); i++)
+   {
+      CHECK(test[i] == result[i]);
+   }
+}