val_ptr);
const size_type n_values = col_ptr - cols.data();
if (n_values > 0)
- global_matrix.add(
- row, n_values, cols.data(), vals.data(), false, true);
+ global_matrix.add(row,
+ n_values,
+ cols.data(),
+ vals.data(),
+ /* elide zero additions */ false,
+ /* sorted by column index */ true);
}
else
internal::AffineConstraints::resolve_matrix_row(
const bool elide_zero_values = true);
/**
- * Set several elements in the specified row of the matrix with column
- * indices as given by <tt>col_indices</tt> to the respective value.
+ * Add the provided values to several elements in the specified row of the
+ * matrix with column indices as given by <tt>col_indices</tt>.
*
* The optional parameter <tt>elide_zero_values</tt> can be used to specify
* whether zero values should be added anyway or these should be filtered
++post_diag;
}
- // add indices before diagonal
+ // Add indices before diagonal. Because the input array
+ // is sorted, and because the entries in this matrix row
+ // are sorted, we can just linearly walk the colnums array
+ // and the input array in parallel, stopping whenever the
+ // former matches the column index of the next index in
+ // the input array:
size_type counter = 1;
for (size_type i = 0; i < diag; ++i)
{
val_ptr[counter] += values[i];
}
- // add indices after diagonal
+ // Then do the same to add indices after the diagonal:
for (size_type i = post_diag; i < n_cols; ++i)
{
while (this_cols[counter] < col_indices[i] &&
}
else
{
+ // Use the same algorithm as above, but because the matrix is
+ // not square, we can now do without the split for diagonal/
+ // entries before the diagional/entries are the diagonal.
size_type counter = 0;
for (size_type i = 0; i < n_cols; ++i)
{