[Lme4-commits] r1509 - in pkg/lme4Eigen: . R data src

[noreply at r-forge.r-project.org]

Author: dmbates
Date: 2012-01-20 01:00:12 +0100 (Fri, 20 Jan 2012)
New Revision: 1509

Added:
   pkg/lme4Eigen/R/sparsegrid.R
   pkg/lme4Eigen/data/GQN.rda
Modified:
   pkg/lme4Eigen/NAMESPACE
   pkg/lme4Eigen/src/external.cpp
Log:
Create the sparse grids for Gauss-Hermite quadrature in several dimensions.


Modified: pkg/lme4Eigen/NAMESPACE
===================================================================
--- pkg/lme4Eigen/NAMESPACE	2012-01-16 01:57:29 UTC (rev 1508)
+++ pkg/lme4Eigen/NAMESPACE	2012-01-20 00:00:12 UTC (rev 1509)
@@ -85,6 +85,7 @@
 
 # and the rest (S3 generics; regular functions):
 export(GHrule,
+       GQdk,
        NelderMead,
        VarCorr,
        bootMer,

Added: pkg/lme4Eigen/R/sparsegrid.R
===================================================================
--- pkg/lme4Eigen/R/sparsegrid.R	                        (rev 0)
+++ pkg/lme4Eigen/R/sparsegrid.R	2012-01-20 00:00:12 UTC (rev 1509)
@@ -0,0 +1,15 @@
+GQdk <- function(d=1L, k=1L) {
+    stopifnot(0L < (d <- as.integer(d)[1]),
+              d <= 20L,
+              0L < (k <- as.integer(k)[1]),
+              k <= length(GQNd <- GQN[[d]]))
+    tmat   <- t(GQNd[[k]])
+    dseq   <- seq_len(d)
+    rperms <- lapply(.Call(allPerm_int, dseq + 1L), function(v) c(1L, v))
+    unname(unique(t(do.call(cbind,
+                            lapply(as.data.frame(t(cbind(1,
+                                                         as.matrix(do.call(expand.grid,
+                                                                           lapply(dseq,
+                                                                                  function(i) c(-1,1))))))),
+                                   "*", e2=do.call(cbind, lapply(rperms, function(ind) tmat[ind,])))))))
+}

Added: pkg/lme4Eigen/data/GQN.rda
===================================================================
(Binary files differ)


Property changes on: pkg/lme4Eigen/data/GQN.rda
___________________________________________________________________
Added: svn:mime-type
   + application/octet-stream

Modified: pkg/lme4Eigen/src/external.cpp
===================================================================
--- pkg/lme4Eigen/src/external.cpp	2012-01-16 01:57:29 UTC (rev 1508)
+++ pkg/lme4Eigen/src/external.cpp	2012-01-20 00:00:12 UTC (rev 1509)
@@ -10,9 +10,10 @@
 
 extern "C" {
     using     Eigen::ArrayXd;
-    typedef   Eigen::Map<Eigen::MatrixXd>     MMat;
-    typedef   Eigen::Map<Eigen::VectorXd>     MVec;
-    typedef   Eigen::Map<Eigen::VectorXi>    MiVec;
+    typedef   Eigen::VectorXi               iVec;
+    typedef   Eigen::Map<Eigen::MatrixXd>   MMat;
+    typedef   Eigen::Map<Eigen::VectorXd>   MVec;
+    typedef   Eigen::Map<iVec>             MiVec;
 
     using      Rcpp::CharacterVector;
     using      Rcpp::Environment;
@@ -39,12 +40,33 @@
 
     using      std::runtime_error;
 
+    // utilities
+
+    SEXP allPerm_int(SEXP v_) {
+	BEGIN_RCPP;
+	iVec     v(as<iVec>(v_));   // forces a copy
+	int     sz(v.size());
+	std::vector<iVec> vec;
+	
+	std::sort(v.data(), v.data() + sz);
+	do {
+	    vec.push_back(iVec(v));
+	} while (std::next_permutation(v.data(), v.data() + sz));
+	
+	int  nperm(vec.size());
+	List allPerm(nperm);
+	for (int j = 0; j < nperm; ++j) allPerm[j] = wrap(vec[j]);
+	return allPerm;
+	END_RCPP;
+    }
+
     SEXP Eigen_SSE() {
 	BEGIN_RCPP;
 	return wrap(Eigen::SimdInstructionSetsInUse());
 	END_RCPP;
     }
 
+
     // generalized linear model (and generalized linear mixed model) response
 
     SEXP glm_Create(SEXP fam, SEXP y, SEXP weights, SEXP offset, SEXP mu,
@@ -766,6 +788,8 @@
 
     CALLDEF(Eigen_SSE, 0),
 
+    CALLDEF(allPerm_int, 1),
+
     CALLDEF(glm_Create, 10),	// generate external pointer
 
     CALLDEF(glm_setN, 2),	// setters