В алгоритме, который я хочу реализовать, мне нужно многократно работать с матрицей, которая образована набором строк, взятых из другой существующей матрицы. Для ясности я назову «A» исходную матрицу, а «B» - «рабочую матрицу» (состоящую из нескольких строк A).
Алгоритм работает следующим образом:
- Выберите строку
r из A - Вставьте или удалите такую строку в
B - Выполните некоторые операции с
B (если это не пусто) - Повтор
Некоторая соответствующая информация:
- Если
r необходимо вставить в B, я уверен, что его еще нет (нет риска создания дубликатов) - Если
r необходимо удалить из B, я уверен, что он уже существует (нет риска удалить несуществующие строки) - Порядок строк внутри
B не имеет значения
Моя цель - эффективно реализовать эту часть алгоритма, используя Eigen. С эффективностью я в основном обращаюсь ко времени выполнения.
Я предложил 3 разных решения, я хотел бы знать, что вы думаете о них и какие улучшения можно добавить.
Решение 1 : изменить размер и скопировать
Идея состоит в том, чтобы сохранить B в MatrixXd экземпляре, который многократно расширяется / сокращается. Поскольку порядок строк не имеет значения, я складываю новые строки внизу B. Когда мне нужно удалить строку из B, я просто копирую нижние строки вверх.
class Algorithm1 {
public:
bool show_b;
Algorithm1(const Eigen::MatrixXd& M) : A(M), B(0,M.cols()), show_b(false) {
rows.reserve(A.rows());
}
void run(const std::vector<std::pair<bool,int>>& ops) {
for(const auto& op : ops) {
if(rows.size() > 0) {
// If B is not empty, do something with it (print is just an example)
if(show_b)
std::cout << "-------------------" << std::endl << "B is:" << std::endl << B << std::endl;
}
// check if we have to insert or remove the row
if(op.first) {
// Row 'op.second' should be added to B.
B.conservativeResize(rows.size()+1, A.cols()); // add space for a column at the bottom
B.row(rows.size()) = A.row(op.second); // copy the desired row from A to B
rows.push_back(op.second); // keep track of which rows are inside B
}
else {
// Row 'op.second' should be removed from B.
auto it = std::find(rows.begin(), rows.end(), op.second); // find where the row is stored in B
int r = std::distance(rows.begin(), it); // get the actual index
for(int i=r; i<rows.size()-1; i++)
B.row(i) = B.row(i+1); // copy each row to the one above
B.conservativeResize(rows.size()-1, A.cols()); // remove the last row of B
rows.erase(it); // remove the row also from this list
}
}
}
private:
Eigen::MatrixXd A; // Original matrix
Eigen::MatrixXd B; // Working matrix
std::vector<int> rows; // list of rows in B, such that A.row(rows[i]) equals B.row(i)
};
Обратите внимание, что в принципе строки, которые будут добавлены / удалены, выбираются алгоритмом внутри. Однако, чтобы упростить его, я вручную передаю «список операций» в качестве параметра алгоритму. Каждый элемент std::vector - это pair, первый элемент которого сообщает, нужно ли добавить строку (pair.first равно true) или удалить (pair.first равно false). Индекс строки, которую нужно выбрать, указан в pair.second (индекс относится к строкам A).
Решение 2: Eigen::Map
Другое решение, о котором я мог подумать, - это сохранить данные B в «сыром буфере», который можно использовать в основном алгоритме через Eigen::Map. Чтобы иметь возможность быстро добавлять / удалять целую строку из буфера, мне пришлось использовать упорядочение основных строк для A и B.
class Algorithm2 {
public:
// Row-major matrix
typedef Eigen::Matrix<double,Eigen::Dynamic,Eigen::Dynamic,Eigen::RowMajor> RMatrixXd;
bool show_b;
Algorithm2(const Eigen::MatrixXd& M) : A(M), show_b(false) {
rows.reserve(A.rows());
data.reserve(A.size());
}
void run(const std::vector<std::pair<bool,int>>& ops) {
for(const auto& op : ops) {
if(rows.size() > 0) {
// Create the sub-matrix B using a Map
Eigen::Map<RMatrixXd> B(data.data(), rows.size(), A.cols());
// Do something with it (print is just an example)
if(show_b)
std::cout << "-------------------" << std::endl << "B is:" << std::endl << B << std::endl;
}
// check if we have to insert or remove the row
if(op.first) {
// Row 'op.second' should be added to B.
// append the row content inside the data
data.insert(
data.end(), // append at the end
A.data() + op.second*A.cols(), // from the beginning of the row
A.data() + (op.second+1)*A.cols() // until the beginning of the next one
);
rows.push_back(op.second); // keep track of which rows are inside B
}
else {
// Row 'op.second' should be removed from B.
auto it = std::find(rows.begin(), rows.end(), op.second); // find where the row is stored in B
int r = std::distance(rows.begin(), it); // get the actual index
data.erase(data.begin()+r*A.cols(), data.begin()+(r+1)*A.cols()); // remove the corresponding values from the data
rows.erase(it); // remove the row also from this list
}
}
}
private:
RMatrixXd A; // The original matrix
std::vector<int> rows; // list of rows in B, such that A.row(rows[i]) equals B.row(i)
std::vector<double> data; // content of B (row-major ordering)
};
Решение 3: новая перегрузка operator()
Используя последнюю версию Eigen, можно выбрать подматрицу напрямую через operator(). На практике матрицу B можно получить как A( rows-to-select , Eigen::all ).
class Algorithm3 {
public:
bool show_b;
Algorithm3(const Eigen::MatrixXd& M) : A(M), show_b(false) {
rows.reserve(A.rows());
}
void run(const std::vector<std::pair<bool,int>>& ops) {
for(const auto& op : ops) {
if(rows.size() > 0) {
// Create the sub-matrix B
auto B = A(rows, Eigen::all);
// Do something with it (print is just an example)
if(show_b)
std::cout << "-------------------" << std::endl << "B is:" << std::endl << B << std::endl;
}
// check if we have to insert or remove the row
if(op.first) {
// Row 'op.second' should be added to B.
rows.push_back(op.second); // keep track of which rows are inside B
}
else {
// Row 'op.second' should be removed from B.
rows.erase(std::find(rows.begin(), rows.end(), op.second)); // remove the row from the list
}
}
}
private:
Eigen::MatrixXd A; // The original matrix
std::vector<int> rows; // list of rows in B, such that A.row(rows[i]) equals B.row(i)
};
Быстрое сравнение решений
Ниже приведена простая программа для записи времени выполнения:
int main() {
Eigen::MatrixXd A(3,3);
A << 0,0,0,1,1,1,2,2,2;
std::vector<std::pair<bool,int>> ops {
// {add=true/remove=false, row-index} // which rows will now be inside B
{true, 1}, // 1
{true, 0}, // 0 1
{true, 2}, // 0 1 2
{false,0}, // 1 2
{false,2}, // 1
{true, 2}, // 1 2
{false,1}, // 2
{true, 0}, // 0 2
{true, 1}, // 0 1 2
{false,0}, // 1 2
{false,1}, // 2
{false,2}, // empty
};
// Instanciate
Algorithm1 alg1(A);
Algorithm2 alg2(A);
Algorithm3 alg3(A);
// warm-up
for(int i=0; i<1000; i++) {
alg1.run(ops);
alg2.run(ops);
alg3.run(ops);
}
// run the 3 algorithms
auto start1 = std::chrono::high_resolution_clock::now();
for(int i=0; i<1000; i++)
alg1.run(ops);
auto end1 = std::chrono::high_resolution_clock::now();
auto start2 = std::chrono::high_resolution_clock::now();
for(int i=0; i<1000; i++)
alg2.run(ops);
auto end2 = std::chrono::high_resolution_clock::now();
auto start3 = std::chrono::high_resolution_clock::now();
for(int i=0; i<1000; i++)
alg3.run(ops);
auto end3 = std::chrono::high_resolution_clock::now();
// show the results
std::cout << "Alg1 took " << std::chrono::duration_cast<std::chrono::nanoseconds>(end1-start1).count() << "ns" << std::endl;
std::cout << "Alg2 took " << std::chrono::duration_cast<std::chrono::nanoseconds>(end2-start2).count() << "ns" << std::endl;
std::cout << "Alg3 took " << std::chrono::duration_cast<std::chrono::nanoseconds>(end3-start3).count() << "ns" << std::endl;
return 0;
}
На моем ноутбуке лучшие результаты дают 2-я стратегия, затем 3-я (в 3,4 раза дольше) и 1-я (в 4,3 раза медленнее).
Что вы думаете? Имеют ли результаты смысл по сравнению с вашим опытом? Можете ли вы предложить лучшую стратегию?