c++ - 利用犰狳的对称性

Question

有没有办法利用犰狳中得到的矩阵是对称的这一事实？我经常有这种形式的产品A * B * A.t()（其中B是对称的和正定的），我想知道是否有某种方法可以通过使用这个事实来更快地进行乘法运算。我做了一些简单的尝试，但它们要么同样快，要么慢得多。看到提供速度改进的功能（或解释为什么无法获得）会很有趣。

#include <RcppArmadillo.h>
// [[Rcpp::depends(RcppArmadillo)]]
// [[Rcpp::export]]
arma::mat m1(const arma::mat & A, const arma::mat & B) {
  arma::mat C = (A * B) * A.t();
  return C;
}

// [[Rcpp::export]]
arma::mat m2(const arma::mat & A, const arma::mat & B) {
  arma::mat C1 = (A * B);
  arma::mat C2 = arma::mat(arma::size(C1), arma::fill::zeros);
  arma::uword n = C1.n_rows;

  for (arma::uword i = 0; i < n; i++) {
    for (arma::uword j = 0; j <= i; j++) {
      C2(i, j) = arma::as_scalar(C1.row(i) * A.row(j).t());
    }
  }
  C2 = arma::symmatl(C2);
  return C2;
}

// [[Rcpp::export]]
arma::mat m3(const arma::mat & A, const arma::mat & B) {
  arma::mat C1 = (A * B);
  arma::mat C2 = arma::mat(arma::size(C1), arma::fill::zeros);
  arma::uword n = C1.n_rows;

  for (arma::uword i = 0; i < n; i++) {
    for (arma::uword j = 0; j <= i; j++) {
      C2.submat(i, j, i, j) = C1.row(i) * A.row(j).t();
    }
  }
  C2 = arma::symmatl(C2);
  return C2;
}

/*** R
n <- 10
A <- matrix(rnorm(n*n), n, n)
B <- crossprod(matrix(rnorm(n*n), n, n))


microbenchmark::microbenchmark(m1(A, B), m2(A, B), m3(A, B))

all.equal(m1(A, B), m2(A, B))
all.equal(m1(A, B), m3(A, B))
*/

标签： c++rcpparmadillo