Welcome to mirror list, hosted at ThFree Co, Russian Federation.

git.blender.org/blender.git - Unnamed repository; edit this file 'description' to name the repository.
summaryrefslogtreecommitdiff
diff options
context:
space:
mode:
Diffstat (limited to 'extern/ceres/internal/ceres/schur_eliminator_impl.h')
-rw-r--r--extern/ceres/internal/ceres/schur_eliminator_impl.h336
1 files changed, 177 insertions, 159 deletions
diff --git a/extern/ceres/internal/ceres/schur_eliminator_impl.h b/extern/ceres/internal/ceres/schur_eliminator_impl.h
index f2535880f15..bd0881eec1e 100644
--- a/extern/ceres/internal/ceres/schur_eliminator_impl.h
+++ b/extern/ceres/internal/ceres/schur_eliminator_impl.h
@@ -48,23 +48,23 @@
// This include must come before any #ifndef check on Ceres compile options.
#include "ceres/internal/port.h"
-#ifdef CERES_USE_OPENMP
-#include <omp.h>
-#endif
-
#include <algorithm>
#include <map>
+
+#include "Eigen/Dense"
#include "ceres/block_random_access_matrix.h"
#include "ceres/block_sparse_matrix.h"
#include "ceres/block_structure.h"
#include "ceres/internal/eigen.h"
#include "ceres/internal/fixed_array.h"
-#include "ceres/internal/scoped_ptr.h"
+#include "ceres/invert_psd_matrix.h"
#include "ceres/map_util.h"
+#include "ceres/parallel_for.h"
#include "ceres/schur_eliminator.h"
+#include "ceres/scoped_thread_token.h"
#include "ceres/small_blas.h"
#include "ceres/stl_util.h"
-#include "Eigen/Dense"
+#include "ceres/thread_token_provider.h"
#include "glog/logging.h"
namespace ceres {
@@ -76,14 +76,16 @@ SchurEliminator<kRowBlockSize, kEBlockSize, kFBlockSize>::~SchurEliminator() {
}
template <int kRowBlockSize, int kEBlockSize, int kFBlockSize>
-void
-SchurEliminator<kRowBlockSize, kEBlockSize, kFBlockSize>::
-Init(int num_eliminate_blocks, const CompressedRowBlockStructure* bs) {
+void SchurEliminator<kRowBlockSize, kEBlockSize, kFBlockSize>::Init(
+ int num_eliminate_blocks,
+ bool assume_full_rank_ete,
+ const CompressedRowBlockStructure* bs) {
CHECK_GT(num_eliminate_blocks, 0)
<< "SchurComplementSolver cannot be initialized with "
<< "num_eliminate_blocks = 0.";
num_eliminate_blocks_ = num_eliminate_blocks;
+ assume_full_rank_ete_ = assume_full_rank_ete;
const int num_col_blocks = bs->cols.size();
const int num_row_blocks = bs->rows.size();
@@ -102,6 +104,13 @@ Init(int num_eliminate_blocks, const CompressedRowBlockStructure* bs) {
lhs_num_rows += bs->cols[i].size;
}
+ // TODO(sameeragarwal): Now that we may have subset block structure,
+ // we need to make sure that we account for the fact that somep
+ // point blocks only have a "diagonal" row and nothing more.
+ //
+ // This likely requires a slightly different algorithm, which works
+ // off of the number of elimination blocks.
+
int r = 0;
// Iterate over the row blocks of A, and detect the chunks. The
// matrix should already have been ordered so that all rows
@@ -143,15 +152,12 @@ Init(int num_eliminate_blocks, const CompressedRowBlockStructure* bs) {
++chunk.size;
}
- CHECK_GT(chunk.size, 0);
+ CHECK_GT(chunk.size, 0); // This check will need to be resolved.
r += chunk.size;
}
const Chunk& chunk = chunks_.back();
uneliminated_row_begins_ = chunk.start + chunk.size;
- if (num_threads_ > 1) {
- random_shuffle(chunks_.begin(), chunks_.end());
- }
buffer_.reset(new double[buffer_size_ * num_threads_]);
@@ -163,46 +169,51 @@ Init(int num_eliminate_blocks, const CompressedRowBlockStructure* bs) {
STLDeleteElements(&rhs_locks_);
rhs_locks_.resize(num_col_blocks - num_eliminate_blocks_);
for (int i = 0; i < num_col_blocks - num_eliminate_blocks_; ++i) {
- rhs_locks_[i] = new Mutex;
+ rhs_locks_[i] = new std::mutex;
}
}
template <int kRowBlockSize, int kEBlockSize, int kFBlockSize>
void
SchurEliminator<kRowBlockSize, kEBlockSize, kFBlockSize>::
-Eliminate(const BlockSparseMatrix* A,
+Eliminate(const BlockSparseMatrixData& A,
const double* b,
const double* D,
BlockRandomAccessMatrix* lhs,
double* rhs) {
if (lhs->num_rows() > 0) {
lhs->SetZero();
- VectorRef(rhs, lhs->num_rows()).setZero();
+ if (rhs) {
+ VectorRef(rhs, lhs->num_rows()).setZero();
+ }
}
- const CompressedRowBlockStructure* bs = A->block_structure();
+ const CompressedRowBlockStructure* bs = A.block_structure();
const int num_col_blocks = bs->cols.size();
// Add the diagonal to the schur complement.
if (D != NULL) {
-#pragma omp parallel for num_threads(num_threads_) schedule(dynamic)
- for (int i = num_eliminate_blocks_; i < num_col_blocks; ++i) {
- const int block_id = i - num_eliminate_blocks_;
- int r, c, row_stride, col_stride;
- CellInfo* cell_info = lhs->GetCell(block_id, block_id,
- &r, &c,
- &row_stride, &col_stride);
- if (cell_info != NULL) {
- const int block_size = bs->cols[i].size;
- typename EigenTypes<Eigen::Dynamic>::ConstVectorRef
- diag(D + bs->cols[i].position, block_size);
-
- CeresMutexLock l(&cell_info->m);
- MatrixRef m(cell_info->values, row_stride, col_stride);
- m.block(r, c, block_size, block_size).diagonal()
- += diag.array().square().matrix();
- }
- }
+ ParallelFor(
+ context_,
+ num_eliminate_blocks_,
+ num_col_blocks,
+ num_threads_,
+ [&](int i) {
+ const int block_id = i - num_eliminate_blocks_;
+ int r, c, row_stride, col_stride;
+ CellInfo* cell_info = lhs->GetCell(block_id, block_id, &r, &c,
+ &row_stride, &col_stride);
+ if (cell_info != NULL) {
+ const int block_size = bs->cols[i].size;
+ typename EigenTypes<Eigen::Dynamic>::ConstVectorRef diag(
+ D + bs->cols[i].position, block_size);
+
+ std::lock_guard<std::mutex> l(cell_info->m);
+ MatrixRef m(cell_info->values, row_stride, col_stride);
+ m.block(r, c, block_size, block_size).diagonal() +=
+ diag.array().square().matrix();
+ }
+ });
}
// Eliminate y blocks one chunk at a time. For each chunk, compute
@@ -218,79 +229,78 @@ Eliminate(const BlockSparseMatrix* A,
// z blocks that share a row block/residual term with the y
// block. EliminateRowOuterProduct does the corresponding operation
// for the lhs of the reduced linear system.
-#pragma omp parallel for num_threads(num_threads_) schedule(dynamic)
- for (int i = 0; i < chunks_.size(); ++i) {
-#ifdef CERES_USE_OPENMP
- int thread_id = omp_get_thread_num();
-#else
- int thread_id = 0;
-#endif
- double* buffer = buffer_.get() + thread_id * buffer_size_;
- const Chunk& chunk = chunks_[i];
- const int e_block_id = bs->rows[chunk.start].cells.front().block_id;
- const int e_block_size = bs->cols[e_block_id].size;
-
- VectorRef(buffer, buffer_size_).setZero();
-
- typename EigenTypes<kEBlockSize, kEBlockSize>::Matrix
- ete(e_block_size, e_block_size);
+ ParallelFor(
+ context_,
+ 0,
+ int(chunks_.size()),
+ num_threads_,
+ [&](int thread_id, int i) {
+ double* buffer = buffer_.get() + thread_id * buffer_size_;
+ const Chunk& chunk = chunks_[i];
+ const int e_block_id = bs->rows[chunk.start].cells.front().block_id;
+ const int e_block_size = bs->cols[e_block_id].size;
+
+ VectorRef(buffer, buffer_size_).setZero();
+
+ typename EigenTypes<kEBlockSize, kEBlockSize>::Matrix
+ ete(e_block_size, e_block_size);
+
+ if (D != NULL) {
+ const typename EigenTypes<kEBlockSize>::ConstVectorRef
+ diag(D + bs->cols[e_block_id].position, e_block_size);
+ ete = diag.array().square().matrix().asDiagonal();
+ } else {
+ ete.setZero();
+ }
- if (D != NULL) {
- const typename EigenTypes<kEBlockSize>::ConstVectorRef
- diag(D + bs->cols[e_block_id].position, e_block_size);
- ete = diag.array().square().matrix().asDiagonal();
- } else {
- ete.setZero();
- }
+ FixedArray<double, 8> g(e_block_size);
+ typename EigenTypes<kEBlockSize>::VectorRef gref(g.data(),
+ e_block_size);
+ gref.setZero();
+
+ // We are going to be computing
+ //
+ // S += F'F - F'E(E'E)^{-1}E'F
+ //
+ // for each Chunk. The computation is broken down into a number of
+ // function calls as below.
+
+ // Compute the outer product of the e_blocks with themselves (ete
+ // = E'E). Compute the product of the e_blocks with the
+ // corresponding f_blocks (buffer = E'F), the gradient of the terms
+ // in this chunk (g) and add the outer product of the f_blocks to
+ // Schur complement (S += F'F).
+ ChunkDiagonalBlockAndGradient(
+ chunk, A, b, chunk.start, &ete, g.data(), buffer, lhs);
+
+ // Normally one wouldn't compute the inverse explicitly, but
+ // e_block_size will typically be a small number like 3, in
+ // which case its much faster to compute the inverse once and
+ // use it to multiply other matrices/vectors instead of doing a
+ // Solve call over and over again.
+ typename EigenTypes<kEBlockSize, kEBlockSize>::Matrix inverse_ete =
+ InvertPSDMatrix<kEBlockSize>(assume_full_rank_ete_, ete);
+
+ // For the current chunk compute and update the rhs of the reduced
+ // linear system.
+ //
+ // rhs = F'b - F'E(E'E)^(-1) E'b
+
+ if (rhs) {
+ FixedArray<double, 8> inverse_ete_g(e_block_size);
+ MatrixVectorMultiply<kEBlockSize, kEBlockSize, 0>(
+ inverse_ete.data(),
+ e_block_size,
+ e_block_size,
+ g.data(),
+ inverse_ete_g.data());
+ UpdateRhs(chunk, A, b, chunk.start, inverse_ete_g.data(), rhs);
+ }
- FixedArray<double, 8> g(e_block_size);
- typename EigenTypes<kEBlockSize>::VectorRef gref(g.get(), e_block_size);
- gref.setZero();
-
- // We are going to be computing
- //
- // S += F'F - F'E(E'E)^{-1}E'F
- //
- // for each Chunk. The computation is broken down into a number of
- // function calls as below.
-
- // Compute the outer product of the e_blocks with themselves (ete
- // = E'E). Compute the product of the e_blocks with the
- // corresonding f_blocks (buffer = E'F), the gradient of the terms
- // in this chunk (g) and add the outer product of the f_blocks to
- // Schur complement (S += F'F).
- ChunkDiagonalBlockAndGradient(
- chunk, A, b, chunk.start, &ete, g.get(), buffer, lhs);
-
- // Normally one wouldn't compute the inverse explicitly, but
- // e_block_size will typically be a small number like 3, in
- // which case its much faster to compute the inverse once and
- // use it to multiply other matrices/vectors instead of doing a
- // Solve call over and over again.
- typename EigenTypes<kEBlockSize, kEBlockSize>::Matrix inverse_ete =
- ete
- .template selfadjointView<Eigen::Upper>()
- .llt()
- .solve(Matrix::Identity(e_block_size, e_block_size));
-
- // For the current chunk compute and update the rhs of the reduced
- // linear system.
- //
- // rhs = F'b - F'E(E'E)^(-1) E'b
-
- FixedArray<double, 8> inverse_ete_g(e_block_size);
- MatrixVectorMultiply<kEBlockSize, kEBlockSize, 0>(
- inverse_ete.data(),
- e_block_size,
- e_block_size,
- g.get(),
- inverse_ete_g.get());
-
- UpdateRhs(chunk, A, b, chunk.start, inverse_ete_g.get(), rhs);
-
- // S -= F'E(E'E)^{-1}E'F
- ChunkOuterProduct(bs, inverse_ete, buffer, chunk.buffer_layout, lhs);
- }
+ // S -= F'E(E'E)^{-1}E'F
+ ChunkOuterProduct(
+ thread_id, bs, inverse_ete, buffer, chunk.buffer_layout, lhs);
+ });
// For rows with no e_blocks, the schur complement update reduces to
// S += F'F.
@@ -300,19 +310,25 @@ Eliminate(const BlockSparseMatrix* A,
template <int kRowBlockSize, int kEBlockSize, int kFBlockSize>
void
SchurEliminator<kRowBlockSize, kEBlockSize, kFBlockSize>::
-BackSubstitute(const BlockSparseMatrix* A,
+BackSubstitute(const BlockSparseMatrixData& A,
const double* b,
const double* D,
const double* z,
double* y) {
- const CompressedRowBlockStructure* bs = A->block_structure();
-#pragma omp parallel for num_threads(num_threads_) schedule(dynamic)
- for (int i = 0; i < chunks_.size(); ++i) {
+ const CompressedRowBlockStructure* bs = A.block_structure();
+ const double* values = A.values();
+
+ ParallelFor(
+ context_,
+ 0,
+ int(chunks_.size()),
+ num_threads_,
+ [&](int i) {
const Chunk& chunk = chunks_[i];
const int e_block_id = bs->rows[chunk.start].cells.front().block_id;
const int e_block_size = bs->cols[e_block_id].size;
- double* y_ptr = y + bs->cols[e_block_id].position;
+ double* y_ptr = y + bs->cols[e_block_id].position;
typename EigenTypes<kEBlockSize>::VectorRef y_block(y_ptr, e_block_size);
typename EigenTypes<kEBlockSize, kEBlockSize>::Matrix
@@ -325,7 +341,6 @@ BackSubstitute(const BlockSparseMatrix* A,
ete.setZero();
}
- const double* values = A->values();
for (int j = 0; j < chunk.size; ++j) {
const CompressedRow& row = bs->rows[chunk.start + j];
const Cell& e_cell = row.cells.front();
@@ -333,9 +348,9 @@ BackSubstitute(const BlockSparseMatrix* A,
FixedArray<double, 8> sj(row.block.size);
- typename EigenTypes<kRowBlockSize>::VectorRef(sj.get(), row.block.size) =
- typename EigenTypes<kRowBlockSize>::ConstVectorRef
- (b + bs->rows[chunk.start + j].block.position, row.block.size);
+ typename EigenTypes<kRowBlockSize>::VectorRef(sj.data(), row.block.size) =
+ typename EigenTypes<kRowBlockSize>::ConstVectorRef(
+ b + bs->rows[chunk.start + j].block.position, row.block.size);
for (int c = 1; c < row.cells.size(); ++c) {
const int f_block_id = row.cells[c].block_id;
@@ -345,23 +360,24 @@ BackSubstitute(const BlockSparseMatrix* A,
MatrixVectorMultiply<kRowBlockSize, kFBlockSize, -1>(
values + row.cells[c].position, row.block.size, f_block_size,
z + lhs_row_layout_[r_block],
- sj.get());
+ sj.data());
}
MatrixTransposeVectorMultiply<kRowBlockSize, kEBlockSize, 1>(
values + e_cell.position, row.block.size, e_block_size,
- sj.get(),
+ sj.data(),
y_ptr);
MatrixTransposeMatrixMultiply
<kRowBlockSize, kEBlockSize, kRowBlockSize, kEBlockSize, 1>(
- values + e_cell.position, row.block.size, e_block_size,
- values + e_cell.position, row.block.size, e_block_size,
- ete.data(), 0, 0, e_block_size, e_block_size);
+ values + e_cell.position, row.block.size, e_block_size,
+ values + e_cell.position, row.block.size, e_block_size,
+ ete.data(), 0, 0, e_block_size, e_block_size);
}
- ete.llt().solveInPlace(y_block);
- }
+ y_block =
+ InvertPSDMatrix<kEBlockSize>(assume_full_rank_ete_, ete) * y_block;
+ });
}
// Update the rhs of the reduced linear system. Compute
@@ -371,17 +387,17 @@ template <int kRowBlockSize, int kEBlockSize, int kFBlockSize>
void
SchurEliminator<kRowBlockSize, kEBlockSize, kFBlockSize>::
UpdateRhs(const Chunk& chunk,
- const BlockSparseMatrix* A,
+ const BlockSparseMatrixData& A,
const double* b,
int row_block_counter,
const double* inverse_ete_g,
double* rhs) {
- const CompressedRowBlockStructure* bs = A->block_structure();
+ const CompressedRowBlockStructure* bs = A.block_structure();
+ const double* values = A.values();
+
const int e_block_id = bs->rows[chunk.start].cells.front().block_id;
const int e_block_size = bs->cols[e_block_id].size;
-
int b_pos = bs->rows[row_block_counter].block.position;
- const double* values = A->values();
for (int j = 0; j < chunk.size; ++j) {
const CompressedRow& row = bs->rows[row_block_counter + j];
const Cell& e_cell = row.cells.front();
@@ -398,7 +414,7 @@ UpdateRhs(const Chunk& chunk,
const int block_id = row.cells[c].block_id;
const int block_size = bs->cols[block_id].size;
const int block = block_id - num_eliminate_blocks_;
- CeresMutexLock l(rhs_locks_[block]);
+ std::lock_guard<std::mutex> l(*rhs_locks_[block]);
MatrixTransposeVectorMultiply<kRowBlockSize, kFBlockSize, 1>(
values + row.cells[c].position,
row.block.size, block_size,
@@ -432,14 +448,15 @@ void
SchurEliminator<kRowBlockSize, kEBlockSize, kFBlockSize>::
ChunkDiagonalBlockAndGradient(
const Chunk& chunk,
- const BlockSparseMatrix* A,
+ const BlockSparseMatrixData& A,
const double* b,
int row_block_counter,
typename EigenTypes<kEBlockSize, kEBlockSize>::Matrix* ete,
double* g,
double* buffer,
BlockRandomAccessMatrix* lhs) {
- const CompressedRowBlockStructure* bs = A->block_structure();
+ const CompressedRowBlockStructure* bs = A.block_structure();
+ const double* values = A.values();
int b_pos = bs->rows[row_block_counter].block.position;
const int e_block_size = ete->rows();
@@ -448,7 +465,6 @@ ChunkDiagonalBlockAndGradient(
// contribution of its F blocks to the Schur complement, the
// contribution of its E block to the matrix EE' (ete), and the
// corresponding block in the gradient vector.
- const double* values = A->values();
for (int j = 0; j < chunk.size; ++j) {
const CompressedRow& row = bs->rows[row_block_counter + j];
@@ -464,12 +480,13 @@ ChunkDiagonalBlockAndGradient(
values + e_cell.position, row.block.size, e_block_size,
ete->data(), 0, 0, e_block_size, e_block_size);
- // g += E_i' b_i
- MatrixTransposeVectorMultiply<kRowBlockSize, kEBlockSize, 1>(
- values + e_cell.position, row.block.size, e_block_size,
- b + b_pos,
- g);
-
+ if (b) {
+ // g += E_i' b_i
+ MatrixTransposeVectorMultiply<kRowBlockSize, kEBlockSize, 1>(
+ values + e_cell.position, row.block.size, e_block_size,
+ b + b_pos,
+ g);
+ }
// buffer = E'F. This computation is done by iterating over the
// f_blocks for each row in the chunk.
@@ -495,7 +512,8 @@ ChunkDiagonalBlockAndGradient(
template <int kRowBlockSize, int kEBlockSize, int kFBlockSize>
void
SchurEliminator<kRowBlockSize, kEBlockSize, kFBlockSize>::
-ChunkOuterProduct(const CompressedRowBlockStructure* bs,
+ChunkOuterProduct(int thread_id,
+ const CompressedRowBlockStructure* bs,
const Matrix& inverse_ete,
const double* buffer,
const BufferLayoutType& buffer_layout,
@@ -507,11 +525,6 @@ ChunkOuterProduct(const CompressedRowBlockStructure* bs,
const int e_block_size = inverse_ete.rows();
BufferLayoutType::const_iterator it1 = buffer_layout.begin();
-#ifdef CERES_USE_OPENMP
- int thread_id = omp_get_thread_num();
-#else
- int thread_id = 0;
-#endif
double* b1_transpose_inverse_ete =
chunk_outer_product_buffer_.get() + thread_id * buffer_size_;
@@ -535,7 +548,7 @@ ChunkOuterProduct(const CompressedRowBlockStructure* bs,
&row_stride, &col_stride);
if (cell_info != NULL) {
const int block2_size = bs->cols[it2->first].size;
- CeresMutexLock l(&cell_info->m);
+ std::lock_guard<std::mutex> l(cell_info->m);
MatrixMatrixMultiply
<kFBlockSize, kEBlockSize, kEBlockSize, kFBlockSize, -1>(
b1_transpose_inverse_ete, block1_size, e_block_size,
@@ -552,14 +565,18 @@ ChunkOuterProduct(const CompressedRowBlockStructure* bs,
template <int kRowBlockSize, int kEBlockSize, int kFBlockSize>
void
SchurEliminator<kRowBlockSize, kEBlockSize, kFBlockSize>::
-NoEBlockRowsUpdate(const BlockSparseMatrix* A,
+NoEBlockRowsUpdate(const BlockSparseMatrixData& A,
const double* b,
int row_block_counter,
BlockRandomAccessMatrix* lhs,
double* rhs) {
- const CompressedRowBlockStructure* bs = A->block_structure();
- const double* values = A->values();
+ const CompressedRowBlockStructure* bs = A.block_structure();
+ const double* values = A.values();
for (; row_block_counter < bs->rows.size(); ++row_block_counter) {
+ NoEBlockRowOuterProduct(A, row_block_counter, lhs);
+ if (!rhs) {
+ continue;
+ }
const CompressedRow& row = bs->rows[row_block_counter];
for (int c = 0; c < row.cells.size(); ++c) {
const int block_id = row.cells[c].block_id;
@@ -570,7 +587,6 @@ NoEBlockRowsUpdate(const BlockSparseMatrix* A,
b + row.block.position,
rhs + lhs_row_layout_[block]);
}
- NoEBlockRowOuterProduct(A, row_block_counter, lhs);
}
}
@@ -582,7 +598,7 @@ NoEBlockRowsUpdate(const BlockSparseMatrix* A,
// one difference. It does not use any of the template
// parameters. This is because the algorithm used for detecting the
// static structure of the matrix A only pays attention to rows with
-// e_blocks. This is becase rows without e_blocks are rare and
+// e_blocks. This is because rows without e_blocks are rare and
// typically arise from regularization terms in the original
// optimization problem, and have a very different structure than the
// rows with e_blocks. Including them in the static structure
@@ -592,12 +608,13 @@ NoEBlockRowsUpdate(const BlockSparseMatrix* A,
template <int kRowBlockSize, int kEBlockSize, int kFBlockSize>
void
SchurEliminator<kRowBlockSize, kEBlockSize, kFBlockSize>::
-NoEBlockRowOuterProduct(const BlockSparseMatrix* A,
+NoEBlockRowOuterProduct(const BlockSparseMatrixData& A,
int row_block_index,
BlockRandomAccessMatrix* lhs) {
- const CompressedRowBlockStructure* bs = A->block_structure();
+ const CompressedRowBlockStructure* bs = A.block_structure();
+ const double* values = A.values();
+
const CompressedRow& row = bs->rows[row_block_index];
- const double* values = A->values();
for (int i = 0; i < row.cells.size(); ++i) {
const int block1 = row.cells[i].block_id - num_eliminate_blocks_;
DCHECK_GE(block1, 0);
@@ -608,7 +625,7 @@ NoEBlockRowOuterProduct(const BlockSparseMatrix* A,
&r, &c,
&row_stride, &col_stride);
if (cell_info != NULL) {
- CeresMutexLock l(&cell_info->m);
+ std::lock_guard<std::mutex> l(cell_info->m);
// This multiply currently ignores the fact that this is a
// symmetric outer product.
MatrixTransposeMatrixMultiply
@@ -628,7 +645,7 @@ NoEBlockRowOuterProduct(const BlockSparseMatrix* A,
&row_stride, &col_stride);
if (cell_info != NULL) {
const int block2_size = bs->cols[row.cells[j].block_id].size;
- CeresMutexLock l(&cell_info->m);
+ std::lock_guard<std::mutex> l(cell_info->m);
MatrixTransposeMatrixMultiply
<Eigen::Dynamic, Eigen::Dynamic, Eigen::Dynamic, Eigen::Dynamic, 1>(
values + row.cells[i].position, row.block.size, block1_size,
@@ -639,18 +656,19 @@ NoEBlockRowOuterProduct(const BlockSparseMatrix* A,
}
}
-// For a row with an e_block, compute the contribition S += F'F. This
+// For a row with an e_block, compute the contribution S += F'F. This
// function has the same structure as NoEBlockRowOuterProduct, except
// that this function uses the template parameters.
template <int kRowBlockSize, int kEBlockSize, int kFBlockSize>
void
SchurEliminator<kRowBlockSize, kEBlockSize, kFBlockSize>::
-EBlockRowOuterProduct(const BlockSparseMatrix* A,
+EBlockRowOuterProduct(const BlockSparseMatrixData& A,
int row_block_index,
BlockRandomAccessMatrix* lhs) {
- const CompressedRowBlockStructure* bs = A->block_structure();
+ const CompressedRowBlockStructure* bs = A.block_structure();
+ const double* values = A.values();
+
const CompressedRow& row = bs->rows[row_block_index];
- const double* values = A->values();
for (int i = 1; i < row.cells.size(); ++i) {
const int block1 = row.cells[i].block_id - num_eliminate_blocks_;
DCHECK_GE(block1, 0);
@@ -661,7 +679,7 @@ EBlockRowOuterProduct(const BlockSparseMatrix* A,
&r, &c,
&row_stride, &col_stride);
if (cell_info != NULL) {
- CeresMutexLock l(&cell_info->m);
+ std::lock_guard<std::mutex> l(cell_info->m);
// block += b1.transpose() * b1;
MatrixTransposeMatrixMultiply
<kRowBlockSize, kFBlockSize, kRowBlockSize, kFBlockSize, 1>(
@@ -681,7 +699,7 @@ EBlockRowOuterProduct(const BlockSparseMatrix* A,
&row_stride, &col_stride);
if (cell_info != NULL) {
// block += b1.transpose() * b2;
- CeresMutexLock l(&cell_info->m);
+ std::lock_guard<std::mutex> l(cell_info->m);
MatrixTransposeMatrixMultiply
<kRowBlockSize, kFBlockSize, kRowBlockSize, kFBlockSize, 1>(
values + row.cells[i].position, row.block.size, block1_size,