Построение SparseMatrix в Eigen

Я строю разреженную линейную систему с несколькими (мягкими) ограничениями. Я конвертирую код, который использовался для построения матрицы с помощью boost :: ublas, в Eigen. Boost: ublas имеет удобный способ создания разреженной матрицы с известным (или предполагаемым) количеством ненулевых значений и имеет достаточно быстрый оператор (int row, int col) для обновления ее элементов.

Проблема в следующем:

• Использование SparseMatrix :: setFromTriplets:
  Моя система имеет много ограничений. В качестве наивного, «слегка» преувеличенного примера допустим, что у меня есть разреженная матрица 100x100 с 500 nnz, но с 1 миллиардом избыточных ограничений (т.е. ненулевые коэффициенты изменены миллиард раз). setFromTriplets требует, чтобы я сохранил 1 миллиард коэффициентов, большинство из которых будут суммированы, чтобы сформировать мой набор из 500 ненулевых коэффициентов. Это не очень эффективно и не очень удобно для памяти. Конечно, я могу заменить свой std :: vector на std :: map и выполнить накопление ограничений вручную, но при этом каким-то образом упускается смысл наличия разреженного класса матриц, и это тоже было бы неэффективно.

• Использование SparseMatrix :: insert (i, j, val):
  Не работает, если элемент уже присутствует. Моя проблема в том, чтобы иметь возможность накапливать коэффициенты, которые уже есть.

• using SparseMatrix :: coeffRef (i, j):
  Это работает, и это будет функция, которую я ищу. Однако он на несколько порядков медленнее, чем boost :: ublas. Я удивлен, что не нашел для этого лучшей функции. Я думал, что это происходит из-за того, что количество ненулевых символов заранее не известно и приводит к многократным перераспределениям (что и происходит на практике). Однако использование SparseMatrix :: reserve () не имеет никакого эффекта, поскольку это функция, которая работает только для сжатых матриц (комментарий в источнике говорит: «Эта функция не имеет смысла в несжатом режиме» перед утверждением). .. и, как говорится в документации, «вставка нового элемента в SparseMatrix преобразует его позже в несжатый режим».

Каков наиболее эффективный способ построить разреженную матрицу в Eigen, сохранив при этом возможность обновлять ее коэффициенты несколько раз?

Спасибо

[РЕДАКТИРОВАТЬ: пример использования: матрица 10x10 с 10 ненулевыми значениями. Для простоты матрица диагональная]

SparseMatrix<double> mat(10, 10);
for (int i=0; i<10; i++) {
  for (int j=0; j<1000000; j++) {
    mat.coeffRef(i, i) += rand()%10;
  }
}

=> работает, но на порядки медленнее, чем оператор ublas () (для больших матриц и более реалистичных настроек, конечно).

std::vector<Eigen::Triplet<double> > triplets(10000000);
int k=0;
for (int i=0; i<10; i++) {
  for (int j=0; j<1000000; j++) {
    triplets[k++] = Eigen::Triplet<double>(i,i,rand()%10);
  }
}
SparseMatrix<double> mat(10, 10);
mat.setFromTriplets(triplets.begin(), triplets.end());

=> не дружественный к памяти ...