diff --git a/gtsam/linear/HessianFactor.cpp b/gtsam/linear/HessianFactor.cpp index f2130f068..0a75cc336 100644 --- a/gtsam/linear/HessianFactor.cpp +++ b/gtsam/linear/HessianFactor.cpp @@ -348,28 +348,26 @@ double HessianFactor::error(const VectorValues& c) const { /* ************************************************************************* */ void HessianFactor::updateHessian(const KeyVector& infoKeys, SymmetricBlockMatrix* info) const { - gttic(updateHessian_HessianFactor); assert(info); - // Apply updates to the upper triangle - DenseIndex nrVariablesInThisFactor = size(), nrBlocksInInfo = info->nBlocks() - 1; + gttic(updateHessian_HessianFactor); + const DenseIndex nrVariablesInThisFactor = size(); + vector slots(nrVariablesInThisFactor + 1); + for (DenseIndex j = 0; j < nrVariablesInThisFactor; ++j) + slots[j] = Slot(infoKeys, keys_[j]); + slots[nrVariablesInThisFactor] = info->nBlocks() - 1; + + // Apply updates to the upper triangle // Loop over this factor's blocks with indices (i,j) // For every block (i,j), we determine the block (I,J) in info. for (DenseIndex j = 0; j <= nrVariablesInThisFactor; ++j) { - const bool rhs = (j == nrVariablesInThisFactor); - const DenseIndex J = rhs ? nrBlocksInInfo : Slot(infoKeys, keys_[j]); - slots[j] = J; - for (DenseIndex i = 0; i <= j; ++i) { - const DenseIndex I = slots[i]; // because i<=j, slots[i] is valid. - - if (i == j) { - assert(I == J); - info->updateDiagonalBlock(I, info_.diagonalBlock(i)); - } else { - assert(i < j); - assert(I != J); - info->updateOffDiagonalBlock(I, J, info_.aboveDiagonalBlock(i, j)); - } + const DenseIndex J = slots[j]; + info->updateDiagonalBlock(J, info_.diagonalBlock(j)); + for (DenseIndex i = 0; i < j; ++i) { + const DenseIndex I = slots[i]; + assert(i < j); + assert(I != J); + info->updateOffDiagonalBlock(I, J, info_.aboveDiagonalBlock(i, j)); } } } diff --git a/gtsam/linear/JacobianFactor.cpp b/gtsam/linear/JacobianFactor.cpp index 9052cc460..1802475eb 100644 --- a/gtsam/linear/JacobianFactor.cpp +++ b/gtsam/linear/JacobianFactor.cpp @@ -602,16 +602,19 @@ void JacobianFactor::updateHessian(const KeyVector& infoKeys, // Ab_ is the augmented Jacobian matrix A, and we perform I += A'*A below DenseIndex n = Ab_.nBlocks() - 1, N = info->nBlocks() - 1; + // Pre-calculate slots + vector slots(n + 1); + for (DenseIndex j = 0; j < n; ++j) slots[j] = Slot(infoKeys, keys_[j]); + slots[n] = N; + // Apply updates to the upper triangle // Loop over blocks of A, including RHS with j==n - vector slots(n+1); for (DenseIndex j = 0; j <= n; ++j) { Eigen::Block Ab_j = Ab_(j); - const DenseIndex J = (j == n) ? N : Slot(infoKeys, keys_[j]); - slots[j] = J; + const DenseIndex J = slots[j]; // Fill off-diagonal blocks with Ai'*Aj for (DenseIndex i = 0; i < j; ++i) { - const DenseIndex I = slots[i]; // because iupdateOffDiagonalBlock(I, J, Ab_(i).transpose() * Ab_j); } // Fill diagonal block with Aj'*Aj