2014-07-07 1 views
1

, 벡터 fv의 첫번째 소자가 매트릭스 tm 모든 컬럼 1 요소 승산 된 행렬 열 (2)에 의해 제 등RcppEigen -이 행렬 곱셈에있어 잘못된 점은 무엇입니까? I 간단한 마드 제품려고

최소 예 : 이것은 단순히 첫 번째 열을 방송보다는 tm를 반환해야

set.seed(123) 
tm <- matrix(rnorm(25,2,1),nrow=5) 
fv <- rep(1,5) 

tm 

     [,1]  [,2]  [,3]  [,4]  [,5] 
[1,] 1.439524 3.7150650 3.224082 3.78691314 0.9321763 
[2,] 1.769823 2.4609162 2.359814 2.49785048 1.7820251 
[3,] 3.558708 0.7349388 2.400771 0.03338284 0.9739956 
[4,] 2.070508 1.3131471 2.110683 2.70135590 1.2711088 
[5,] 2.129288 1.5543380 1.444159 1.52720859 1.3749607 


library(inline) 
etest <- cxxfunction(signature(tm="NumericMatrix", 
           fv="NumericVector"), 
        plugin="RcppEigen", 
        body=" 
NumericVector fvv(fv); 
NumericMatrix tmm(tm); 

const Eigen::Map<Eigen::MatrixXd> ttm(as<Eigen::Map<Eigen::MatrixXd> >(tmm)); 
const Eigen::Map<Eigen::VectorXd> ffv(as<Eigen::Map<Eigen::VectorXd> >(fvv)); 

Eigen::MatrixXd prod = ttm*ffv.transpose(); 
return(wrap(prod)); 
        ") 

etest(tm,fv) 

     [,1]  [,2]  [,3]  [,4]  [,5] 
[1,] 1.439524 1.439524 1.439524 1.439524 1.439524 
[2,] 1.769823 1.769823 1.769823 1.769823 1.769823 
[3,] 3.558708 3.558708 3.558708 3.558708 3.558708 
[4,] 2.070508 2.070508 2.070508 2.070508 2.070508 
[5,] 2.129288 2.129288 2.129288 2.129288 2.129288 

, 나는 그것이 실제로하고있는 생각하는지 모르겠어요. 나는 명백한 것을 생략하고 있는가?

편집 : etest(tm,diag(fv)) 내가 원하는 것을 제공하지만 이것이 고유 안에서 수행 할 수 있어야합니까?

답변

3

수정 해 주셔서 감사합니다. 5x1 행렬에 5x1 벡터를 곱하면 5x5를 생성 할 수 없습니다. 여기에 단위 행렬이 필요합니다. diag(rep(1,5))diag(5)과 같습니다. 다음 그냥에 sourceCpp("nameOfTheFile.cpp")을 실행

#include <RcppEigen.h> 

using namespace Eigen; 
using namespace Rcpp; 

// [[Rcpp::depends(RcppEigen)]] 

// [[Rcpp::export]] 
MatrixXd etest(const Map<MatrixXd> ttm)); 
    const MatrixXd id = MatrixXd::Identity(ttm.rows(), ttm.cols()); 

    Rcout << "ttm\n " << ttm << std::endl; 
    Rcout << "id\n " << id << std::endl; 
    MatrixXd res = ttm*id; 
    Rcout << "res\n " << res << std::endl; 

    return(res); 
} 

을하고 : 당신과 같은 tm

R> sourceCpp("/tmp/hada.cpp") 
R> etest(tm) 
ttm 
    1.43952 3.71506 3.22408 3.78691 0.932176 
    1.76982 2.46092 2.35981 2.49785 1.78203 
    3.55871 0.734939 2.40077 0.0333828 0.973996 
    2.07051 1.31315 2.11068 2.70136 1.27111 
    2.12929 1.55434 1.44416 1.52721 1.37496 
id 
    1 0 0 0 0 
0 1 0 0 0 
0 0 1 0 0 
0 0 0 1 0 
0 0 0 0 1 
res 
    1.43952 3.71506 3.22408 3.78691 0.932176 
    1.76982 2.46092 2.35981 2.49785 1.78203 
    3.55871 0.734939 2.40077 0.0333828 0.973996 
    2.07051 1.31315 2.11068 2.70136 1.27111 
    2.12929 1.55434 1.44416 1.52721 1.37496 
     [,1]  [,2] [,3]  [,4]  [,5] 
[1,] 1.43952 3.715065 3.22408 3.7869131 0.932176 
[2,] 1.76982 2.460916 2.35981 2.4978505 1.782025 
[3,] 3.55871 0.734939 2.40077 0.0333828 0.973996 
[4,] 2.07051 1.313147 2.11068 2.7013559 1.271109 
[5,] 2.12929 1.554338 1.44416 1.5272086 1.374961 
R> 

고유치의 문서에서 간단한 검색은 어떻게해야 MatrixXd::Identity()을 사용하여 다음과 같은 것을 제시 했다.

편집 1은 : 실제 마드 제품에 대해 ... 말하자면, 옥타브에서 루틴을 복사하거나 어딘가에 또 다른 하나를 찾으려면 단지 수

편집 2 : 더 간결한 인터페이스 배치 Map<MatrixXd> 직접 함수 서명.

+0

감사합니다. Dirk. 나는 더 명확하지 않다는 점에 유감스럽게 생각합니다.이 예는 좀 더 일반적인 문제는 아니지만 실례입니다. 제가 정말로하고 싶은 것은 임의의 벡터를'fv'로 전달하고 컬럼 단위의 곱셈을하는 것입니다. 필자는 선형 대수학의 표준 곱셈은 아니라고 인정하지만, 행렬 감각으로 할 수 있다면 요소 별 순수 Rcpp 솔루션보다 약간의 이득을 기대하고 있습니다. http://eigen.tuxfamily.org/dox/group__QuickRefPage.html –

+0

우리는 동시에 편집 한 것처럼 보이므로 추가 한 두 줄을 참조하십시오. 다른 작업이 필요한 경우 필요에 따라 구현하십시오. –

+0

감사. 실제로이 작업을 제대로 수행하지 못하는 것 같아서 조작원을 학대했습니다. 나는 대각선으로 더 놀아 볼 것입니다. –