[Rcpp-commits] r4311 - in pkg/RcppArmadillo: . debian inst inst/doc/unitTests inst/include inst/include/armadillo_bits
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Sat Apr 20 03:15:30 CEST 2013
Author: edd
Date: 2013-04-20 03:15:28 +0200 (Sat, 20 Apr 2013)
New Revision: 4311
Added:
pkg/RcppArmadillo/inst/include/armadillo_bits/fft_engine.hpp
pkg/RcppArmadillo/inst/include/armadillo_bits/fn_fft.hpp
pkg/RcppArmadillo/inst/include/armadillo_bits/op_fft_bones.hpp
pkg/RcppArmadillo/inst/include/armadillo_bits/op_fft_meat.hpp
Modified:
pkg/RcppArmadillo/ChangeLog
pkg/RcppArmadillo/DESCRIPTION
pkg/RcppArmadillo/debian/changelog
pkg/RcppArmadillo/debian/control
pkg/RcppArmadillo/inst/NEWS.Rd
pkg/RcppArmadillo/inst/doc/unitTests/RcppArmadillo-unitTests.R
pkg/RcppArmadillo/inst/include/armadillo
pkg/RcppArmadillo/inst/include/armadillo_bits/Mat_bones.hpp
pkg/RcppArmadillo/inst/include/armadillo_bits/Mat_meat.hpp
pkg/RcppArmadillo/inst/include/armadillo_bits/SpMat_bones.hpp
pkg/RcppArmadillo/inst/include/armadillo_bits/SpMat_meat.hpp
pkg/RcppArmadillo/inst/include/armadillo_bits/arma_version.hpp
pkg/RcppArmadillo/inst/include/armadillo_bits/config.hpp
pkg/RcppArmadillo/inst/include/armadillo_bits/fn_misc.hpp
pkg/RcppArmadillo/inst/include/armadillo_bits/spglue_times_meat.hpp
pkg/RcppArmadillo/inst/include/armadillo_bits/subview_bones.hpp
pkg/RcppArmadillo/inst/include/armadillo_bits/subview_cube_bones.hpp
pkg/RcppArmadillo/inst/include/armadillo_bits/subview_cube_meat.hpp
pkg/RcppArmadillo/inst/include/armadillo_bits/subview_meat.hpp
pkg/RcppArmadillo/inst/include/armadillo_bits/traits.hpp
Log:
RcppArmadillo release 0.3.810.0 containing Armadillo 3.810.0
Modified: pkg/RcppArmadillo/ChangeLog
===================================================================
--- pkg/RcppArmadillo/ChangeLog 2013-04-16 18:38:26 UTC (rev 4310)
+++ pkg/RcppArmadillo/ChangeLog 2013-04-20 01:15:28 UTC (rev 4311)
@@ -1,3 +1,10 @@
+2013-04-19 Dirk Eddelbuettel <edd at debian.org>
+
+ * DESCRIPTION: Release 0.3.810.0
+ * inst/NEWS: Release 0.3.810.0
+ * inst/include/*: Upgraded to new release 3.810.0 of Armadillo
+ * debian/*: Similarly updated for new release to Debian
+
2013-04-05 Dirk Eddelbuettel <edd at debian.org>
* debian/: Package will be added to Debian as well
Modified: pkg/RcppArmadillo/DESCRIPTION
===================================================================
--- pkg/RcppArmadillo/DESCRIPTION 2013-04-16 18:38:26 UTC (rev 4310)
+++ pkg/RcppArmadillo/DESCRIPTION 2013-04-20 01:15:28 UTC (rev 4311)
@@ -1,7 +1,7 @@
Package: RcppArmadillo
Type: Package
Title: Rcpp integration for Armadillo templated linear algebra library
-Version: 0.3.800.1
+Version: 0.3.810.0
Date: $Date$
Author: Romain Francois, Dirk Eddelbuettel and Doug Bates
Maintainer: Dirk Eddelbuettel <edd at debian.org>
@@ -21,7 +21,7 @@
(due to speed and/or integration capabilities), rather than another language.
.
The RcppArmadillo package includes the header files from the templated
- Armadillo library (currently version 3.800.1). Thus users do not need to
+ Armadillo library (currently version 3.810.0). Thus users do not need to
install Armadillo itself in order to use RcppArmadillo.
.
This Armadillo integration provides a nice illustration of the
Modified: pkg/RcppArmadillo/debian/changelog
===================================================================
--- pkg/RcppArmadillo/debian/changelog 2013-04-16 18:38:26 UTC (rev 4310)
+++ pkg/RcppArmadillo/debian/changelog 2013-04-20 01:15:28 UTC (rev 4311)
@@ -1,3 +1,9 @@
+r-cran-rcpparmadillo (0.3.810.0-1) unstable; urgency=low
+
+ * New upstream release
+
+ -- Dirk Eddelbuettel <edd at debian.org> Fri, 19 Apr 2013 19:56:04 -0500
+
r-cran-rcpparmadillo (0.3.800.1-1) unstable; urgency=low
* Initial Debian release (Closes: #704746)
Modified: pkg/RcppArmadillo/debian/control
===================================================================
--- pkg/RcppArmadillo/debian/control 2013-04-16 18:38:26 UTC (rev 4310)
+++ pkg/RcppArmadillo/debian/control 2013-04-20 01:15:28 UTC (rev 4311)
@@ -25,7 +25,7 @@
(due to speed and/or integration capabilities), rather than another language.
.
The RcppArmadillo package includes the header files from the templated
- Armadillo library (currently version 3.800.1). Thus users do not need to
+ Armadillo library (currently version 3.810.0). Thus users do not need to
install Armadillo itself in order to use RcppArmadillo.
.
This Armadillo integration provides a nice illustration of the
Modified: pkg/RcppArmadillo/inst/NEWS.Rd
===================================================================
--- pkg/RcppArmadillo/inst/NEWS.Rd 2013-04-16 18:38:26 UTC (rev 4310)
+++ pkg/RcppArmadillo/inst/NEWS.Rd 2013-04-20 01:15:28 UTC (rev 4311)
@@ -2,6 +2,21 @@
\title{News for Package 'RcppArmadillo'}
\newcommand{\cpkg}{\href{http://CRAN.R-project.org/package=#1}{\pkg{#1}}}
+\section{Changes in RcppArmadillo version 0.3.810.0 (2013-04-19)}{
+ \itemize{
+ \item Upgraded to Armadillo release Version 3.810.0 (Newell Highway)
+ \itemize{
+ \item added fast Fourier transform: \code{fft()}
+ \item added handling of \code{.imbue()} and \code{.transform()} by
+ submatrices and subcubes
+ \item added batch insertion constructors for sparse matrices
+ \item minor fix for multiplication of complex sparse matrices
+ }
+ \item Updated sample() function and test again contributed by
+ Christian Gunning
+ }
+}
+
\section{Changes in RcppArmadillo version 0.3.800.1 (2013-03-12)}{
\itemize{
\item Upgraded to Armadillo release Version 3.800.1 (Miami Beach)
Modified: pkg/RcppArmadillo/inst/doc/unitTests/RcppArmadillo-unitTests.R
===================================================================
--- pkg/RcppArmadillo/inst/doc/unitTests/RcppArmadillo-unitTests.R 2013-04-16 18:38:26 UTC (rev 4310)
+++ pkg/RcppArmadillo/inst/doc/unitTests/RcppArmadillo-unitTests.R 2013-04-20 01:15:28 UTC (rev 4311)
@@ -1,20 +1,22 @@
pkg <- "RcppArmadillo"
# load this package
-require( pkg, character.only = TRUE )
+require(pkg, character.only = TRUE)
#load RUnit
-runit <- "RUnit" ; require( runit, character.only = TRUE )
-if( file.exists( "unitTests-results" ) ){ unlink("unitTests-results", recursive = TRUE ) }
-dir.create( "unitTests-results" )
+runit <- "RUnit"
+require(runit, character.only = TRUE)
+if (file.exists("unitTests-results")) { unlink("unitTests-results", recursive = TRUE ) }
+dir.create("unitTests-results")
+pathRcppArmadilloTests <<- system.file("unitTests", package = pkg)
path <- system.file("unitTests", package = pkg)
testSuite <- defineTestSuite(name=paste(pkg, "unit testing"), dirs = path)
tests <- runTestSuite(testSuite)
printHTMLProtocol(tests, fileName= sprintf( "unitTests-results/%s-unitTests.html" , pkg ) )
printTextProtocol(tests, fileName= sprintf( "unitTests-results/%s-unitTests.txt" , pkg ) )
-if( file.exists( "/tmp" ) ){
- file.copy( sprintf( "unitTests-results/%s-unitTests.txt" , pkg ) , "/tmp", overwrite = TRUE )
- file.copy( sprintf( "unitTests-results/%s-unitTests.html", pkg ) , "/tmp", overwrite = TRUE )
+if (file.exists( "/tmp")) {
+ file.copy(sprintf("unitTests-results/%s-unitTests.txt" , pkg), "/tmp", overwrite = TRUE)
+ file.copy(sprintf("unitTests-results/%s-unitTests.html", pkg), "/tmp", overwrite = TRUE)
}
Modified: pkg/RcppArmadillo/inst/include/armadillo
===================================================================
--- pkg/RcppArmadillo/inst/include/armadillo 2013-04-16 18:38:26 UTC (rev 4310)
+++ pkg/RcppArmadillo/inst/include/armadillo 2013-04-20 01:15:28 UTC (rev 4311)
@@ -225,6 +225,7 @@
#include "armadillo_bits/op_hist_bones.hpp"
#include "armadillo_bits/op_unique_bones.hpp"
#include "armadillo_bits/op_toeplitz_bones.hpp"
+ #include "armadillo_bits/op_fft_bones.hpp"
#include "armadillo_bits/glue_times_bones.hpp"
#include "armadillo_bits/glue_mixed_bones.hpp"
@@ -408,6 +409,7 @@
#include "armadillo_bits/fn_hist.hpp"
#include "armadillo_bits/fn_histc.hpp"
#include "armadillo_bits/fn_unique.hpp"
+ #include "armadillo_bits/fn_fft.hpp"
#include "armadillo_bits/fn_speye.hpp"
#include "armadillo_bits/fn_spones.hpp"
@@ -418,8 +420,9 @@
// misc stuff
#include "armadillo_bits/hdf5_misc.hpp"
-
-
+ #include "armadillo_bits/fft_engine.hpp"
+
+
//
// class meat
@@ -499,6 +502,7 @@
#include "armadillo_bits/op_hist_meat.hpp"
#include "armadillo_bits/op_unique_meat.hpp"
#include "armadillo_bits/op_toeplitz_meat.hpp"
+ #include "armadillo_bits/op_fft_meat.hpp"
#include "armadillo_bits/glue_times_meat.hpp"
#include "armadillo_bits/glue_mixed_meat.hpp"
Modified: pkg/RcppArmadillo/inst/include/armadillo_bits/Mat_bones.hpp
===================================================================
--- pkg/RcppArmadillo/inst/include/armadillo_bits/Mat_bones.hpp 2013-04-16 18:38:26 UTC (rev 4310)
+++ pkg/RcppArmadillo/inst/include/armadillo_bits/Mat_bones.hpp 2013-04-20 01:15:28 UTC (rev 4311)
@@ -184,7 +184,10 @@
template<typename T1> arma_inline subview_elem1<eT,T1> elem(const Base<uword,T1>& a);
template<typename T1> arma_inline const subview_elem1<eT,T1> elem(const Base<uword,T1>& a) const;
+ template<typename T1> arma_inline subview_elem1<eT,T1> operator()(const Base<uword,T1>& a);
+ template<typename T1> arma_inline const subview_elem1<eT,T1> operator()(const Base<uword,T1>& a) const;
+
template<typename T1, typename T2> arma_inline subview_elem2<eT,T1,T2> elem(const Base<uword,T1>& ri, const Base<uword,T2>& ci);
template<typename T1, typename T2> arma_inline const subview_elem2<eT,T1,T2> elem(const Base<uword,T1>& ri, const Base<uword,T2>& ci) const;
@@ -194,6 +197,7 @@
template<typename T1, typename T2> arma_inline subview_elem2<eT,T1,T2> operator()(const Base<uword,T1>& ri, const Base<uword,T2>& ci);
template<typename T1, typename T2> arma_inline const subview_elem2<eT,T1,T2> operator()(const Base<uword,T1>& ri, const Base<uword,T2>& ci) const;
+
template<typename T1> arma_inline subview_elem2<eT,T1,T1> rows(const Base<uword,T1>& ri);
template<typename T1> arma_inline const subview_elem2<eT,T1,T1> rows(const Base<uword,T1>& ri) const;
Modified: pkg/RcppArmadillo/inst/include/armadillo_bits/Mat_meat.hpp
===================================================================
--- pkg/RcppArmadillo/inst/include/armadillo_bits/Mat_meat.hpp 2013-04-16 18:38:26 UTC (rev 4310)
+++ pkg/RcppArmadillo/inst/include/armadillo_bits/Mat_meat.hpp 2013-04-20 01:15:28 UTC (rev 4311)
@@ -2727,6 +2727,32 @@
template<typename eT>
+template<typename T1>
+arma_inline
+subview_elem1<eT,T1>
+Mat<eT>::operator()(const Base<uword,T1>& a)
+ {
+ arma_extra_debug_sigprint();
+
+ return subview_elem1<eT,T1>(*this, a);
+ }
+
+
+
+template<typename eT>
+template<typename T1>
+arma_inline
+const subview_elem1<eT,T1>
+Mat<eT>::operator()(const Base<uword,T1>& a) const
+ {
+ arma_extra_debug_sigprint();
+
+ return subview_elem1<eT,T1>(*this, a);
+ }
+
+
+
+template<typename eT>
template<typename T1, typename T2>
arma_inline
subview_elem2<eT,T1,T2>
Modified: pkg/RcppArmadillo/inst/include/armadillo_bits/SpMat_bones.hpp
===================================================================
--- pkg/RcppArmadillo/inst/include/armadillo_bits/SpMat_bones.hpp 2013-04-16 18:38:26 UTC (rev 4310)
+++ pkg/RcppArmadillo/inst/include/armadillo_bits/SpMat_bones.hpp 2013-04-20 01:15:28 UTC (rev 4311)
@@ -64,8 +64,21 @@
inline SpMat(const std::string& text);
inline const SpMat& operator=(const std::string& text);
inline SpMat(const SpMat<eT>& x);
+
+ /**
+ * Construct a sparse matrix and insert a large number of points at once.
+ * locations.row[0] should be row indices, locations.row[1] should be column
+ * indices, and values should be the corresponding values. If sorted is true,
+ * then it is assumed that the locations and values are already sorted in
+ * column-major ordering.
+ *
+ * The size is automatically determined to be the minimum size to fit all
+ * points in the locations list.
+ */
+ template<typename T1, typename T2> inline SpMat(const Base<uword,T1>& locations, const Base<eT,T2>& values, const bool sort_locations = true);
+ //! Same as the other batch insertion constructor, but with manual size specification.
+ template<typename T1, typename T2> inline SpMat(const Base<uword,T1>& locations, const Base<eT,T2>& values, const uword n_rows, const uword n_cols, const bool sort_locations = true);
-
inline const SpMat& operator=(const eT val); //! Sets size to 1x1.
inline const SpMat& operator*=(const eT val);
inline const SpMat& operator/=(const eT val);
Modified: pkg/RcppArmadillo/inst/include/armadillo_bits/SpMat_meat.hpp
===================================================================
--- pkg/RcppArmadillo/inst/include/armadillo_bits/SpMat_meat.hpp 2013-04-16 18:38:26 UTC (rev 4310)
+++ pkg/RcppArmadillo/inst/include/armadillo_bits/SpMat_meat.hpp 2013-04-20 01:15:28 UTC (rev 4311)
@@ -164,6 +164,250 @@
+//! Insert a large number of values at once.
+//! locations.row[0] should be row indices, locations.row[1] should be column indices,
+//! and values should be the corresponding values.
+//! If sort_locations is false, then it is assumed that the locations and values
+//! are already sorted in column-major ordering.
+template<typename eT>
+template<typename T1, typename T2>
+inline
+SpMat<eT>::SpMat(const Base<uword,T1>& locations_expr, const Base<eT,T2>& vals_expr, const bool sort_locations)
+ : n_rows(0)
+ , n_cols(0)
+ , n_elem(0)
+ , n_nonzero(0)
+ , vec_state(0)
+ , values(NULL)
+ , row_indices(NULL)
+ , col_ptrs(NULL)
+ {
+ arma_extra_debug_sigprint_this(this);
+
+ const unwrap<T1> locs_tmp( locations_expr.get_ref() );
+ const Mat<uword>& locs = locs_tmp.M;
+
+ const unwrap<T2> vals_tmp( vals_expr.get_ref() );
+ const Mat<eT>& vals = vals_tmp.M;
+
+ arma_debug_check( (vals.is_vec() == false), "SpMat::SpMat(): given 'values' object is not a vector" );
+
+ arma_debug_check((locs.n_cols != vals.n_elem), "SpMat::SpMat(): number of locations is different than number of values");
+
+ // If there are no elements in the list, max() will fail.
+ if (locs.n_cols == 0)
+ {
+ init(0, 0);
+ return;
+ }
+
+ arma_debug_check((locs.n_rows != 2), "SpMat::SpMat(): locations matrix must have two rows");
+
+ // Automatically determine size (and check if it's sorted).
+ uvec bounds = arma::max(locs, 1);
+ init(bounds[0] + 1, bounds[1] + 1);
+
+ // Resize to correct number of elements.
+ mem_resize(vals.n_elem);
+
+ // Reset column pointers to zero.
+ arrayops::inplace_set(access::rwp(col_ptrs), uword(0), n_cols + 1);
+
+ bool actually_sorted = true;
+ if(sort_locations == true)
+ {
+ // sort_index() uses std::sort() which may use quicksort... so we better
+ // make sure it's not already sorted before taking an O(N^2) sort penalty.
+ for (uword i = 1; i < locs.n_cols; ++i)
+ {
+ if ((locs.at(1, i) < locs.at(1, i - 1)) || (locs.at(1, i) == locs.at(1, i - 1) && locs.at(0, i) <= locs.at(0, i - 1)))
+ {
+ actually_sorted = false;
+ break;
+ }
+ }
+
+ if(actually_sorted == false)
+ {
+ // This may not be the fastest possible implementation but it maximizes code reuse.
+ Col<uword> abslocs(locs.n_cols);
+
+ for (uword i = 0; i < locs.n_cols; ++i)
+ {
+ abslocs[i] = locs.at(1, i) * n_rows + locs.at(0, i);
+ }
+
+ // Now we will sort with sort_index().
+ uvec sorted_indices = sort_index(abslocs); // Ascending sort.
+
+ // Now we add the elements in this sorted order.
+ for (uword i = 0; i < sorted_indices.n_elem; ++i)
+ {
+ arma_debug_check((locs.at(0, sorted_indices[i]) >= n_rows), "SpMat::SpMat(): invalid row index");
+ arma_debug_check((locs.at(1, sorted_indices[i]) >= n_cols), "SpMat::SpMat(): invalid column index");
+
+ access::rw(values[i]) = vals[sorted_indices[i]];
+ access::rw(row_indices[i]) = locs.at(0, sorted_indices[i]);
+
+ access::rw(col_ptrs[locs.at(1, sorted_indices[i]) + 1])++;
+ }
+ }
+ }
+
+ if( (sort_locations == false) || (actually_sorted == true) )
+ {
+ // Now set the values and row indices correctly.
+ // Increment the column pointers in each column (so they are column "counts").
+ for (uword i = 0; i < vals.n_elem; ++i)
+ {
+ arma_debug_check((locs.at(0, i) >= n_rows), "SpMat::SpMat(): invalid row index");
+ arma_debug_check((locs.at(1, i) >= n_cols), "SpMat::SpMat(): invalid column index");
+
+ // Check ordering in debug mode.
+ if(i > 0)
+ {
+ arma_debug_check
+ (
+ ( (locs.at(1, i) < locs.at(1, i - 1)) || (locs.at(1, i) == locs.at(1, i - 1) && locs.at(0, i) < locs.at(0, i - 1)) ),
+ "SpMat::SpMat(): out of order points; either pass sort_locations = true, or sort points in column-major ordering"
+ );
+ arma_debug_check((locs.at(1, i) == locs.at(1, i - 1) && locs.at(0, i) == locs.at(0, i - 1)), "SpMat::SpMat(): two identical point locations in list");
+ }
+
+ access::rw(values[i]) = vals[i];
+ access::rw(row_indices[i]) = locs.at(0, i);
+
+ access::rw(col_ptrs[locs.at(1, i) + 1])++;
+ }
+ }
+
+ // Now fix the column pointers.
+ for (uword i = 0; i <= n_cols; ++i)
+ {
+ access::rw(col_ptrs[i + 1]) += col_ptrs[i];
+ }
+ }
+
+
+
+//! Insert a large number of values at once.
+//! locations.row[0] should be row indices, locations.row[1] should be column indices,
+//! and values should be the corresponding values.
+//! If sort_locations is false, then it is assumed that the locations and values
+//! are already sorted in column-major ordering.
+//! In this constructor the size is explicitly given.
+template<typename eT>
+template<typename T1, typename T2>
+inline
+SpMat<eT>::SpMat(const Base<uword,T1>& locations_expr, const Base<eT,T2>& vals_expr, const uword in_n_rows, const uword in_n_cols, const bool sort_locations)
+ : n_rows(0)
+ , n_cols(0)
+ , n_elem(0)
+ , n_nonzero(0)
+ , vec_state(0)
+ , values(NULL)
+ , row_indices(NULL)
+ , col_ptrs(NULL)
+ {
+ arma_extra_debug_sigprint_this(this);
+
+ init(in_n_rows, in_n_cols);
+
+ const unwrap<T1> locs_tmp( locations_expr.get_ref() );
+ const Mat<uword>& locs = locs_tmp.M;
+
+ const unwrap<T2> vals_tmp( vals_expr.get_ref() );
+ const Mat<eT>& vals = vals_tmp.M;
+
+ arma_debug_check( (vals.is_vec() == false), "SpMat::SpMat(): given 'values' object is not a vector" );
+
+ arma_debug_check((locs.n_rows != 2), "SpMat::SpMat(): locations matrix must have two rows");
+
+ arma_debug_check((locs.n_cols != vals.n_elem), "SpMat::SpMat(): number of locations is different than number of values");
+
+ // Resize to correct number of elements.
+ mem_resize(vals.n_elem);
+
+ // Reset column pointers to zero.
+ arrayops::inplace_set(access::rwp(col_ptrs), uword(0), n_cols + 1);
+
+ bool actually_sorted = true;
+ if(sort_locations == true)
+ {
+ // sort_index() uses std::sort() which may use quicksort... so we better
+ // make sure it's not already sorted before taking an O(N^2) sort penalty.
+ for (uword i = 1; i < locs.n_cols; ++i)
+ {
+ if ((locs.at(1, i) < locs.at(1, i - 1)) || (locs.at(1, i) == locs.at(1, i - 1) && locs.at(0, i) <= locs.at(0, i - 1)))
+ {
+ actually_sorted = false;
+ break;
+ }
+ }
+
+ if(actually_sorted == false)
+ {
+ // This may not be the fastest possible implementation but it maximizes code reuse.
+ Col<uword> abslocs(locs.n_cols);
+
+ for (uword i = 0; i < locs.n_cols; ++i)
+ {
+ abslocs[i] = locs.at(1, i) * n_rows + locs.at(0, i);
+ }
+
+ // Now we will sort with sort_index().
+ uvec sorted_indices = sort_index(abslocs); // Ascending sort.
+
+ // Now we add the elements in this sorted order.
+ for (uword i = 0; i < sorted_indices.n_elem; ++i)
+ {
+ arma_debug_check((locs.at(0, sorted_indices[i]) >= n_rows), "SpMat::SpMat(): invalid row index");
+ arma_debug_check((locs.at(1, sorted_indices[i]) >= n_cols), "SpMat::SpMat(): invalid column index");
+
+ access::rw(values[i]) = vals[sorted_indices[i]];
+ access::rw(row_indices[i]) = locs.at(0, sorted_indices[i]);
+
+ access::rw(col_ptrs[locs.at(1, sorted_indices[i]) + 1])++;
+ }
+ }
+ }
+
+ if( (sort_locations == false) || (actually_sorted == true) )
+ {
+ // Now set the values and row indices correctly.
+ // Increment the column pointers in each column (so they are column "counts").
+ for (uword i = 0; i < vals.n_elem; ++i)
+ {
+ arma_debug_check((locs.at(0, i) >= n_rows), "SpMat::SpMat(): invalid row index");
+ arma_debug_check((locs.at(1, i) >= n_cols), "SpMat::SpMat(): invalid column index");
+
+ // Check ordering in debug mode.
+ if(i > 0)
+ {
+ arma_debug_check
+ (
+ ( (locs.at(1, i) < locs.at(1, i - 1)) || (locs.at(1, i) == locs.at(1, i - 1) && locs.at(0, i) < locs.at(0, i - 1)) ),
+ "SpMat::SpMat(): out of order points; either pass sort_locations = true or sort points in column-major ordering"
+ );
+ arma_debug_check((locs.at(1, i) == locs.at(1, i - 1) && locs.at(0, i) == locs.at(0, i - 1)), "SpMat::SpMat(): two identical point locations in list");
+ }
+
+ access::rw(values[i]) = vals[i];
+ access::rw(row_indices[i]) = locs.at(0, i);
+
+ access::rw(col_ptrs[locs.at(1, i) + 1])++;
+ }
+ }
+
+ // Now fix the column pointers.
+ for (uword i = 0; i <= n_cols; ++i)
+ {
+ access::rw(col_ptrs[i + 1]) += col_ptrs[i];
+ }
+ }
+
+
+
/**
* Simple operators with plain values. These operate on every value in the
* matrix, so a sparse matrix += 1 will turn all those zeroes into ones. Be
Modified: pkg/RcppArmadillo/inst/include/armadillo_bits/arma_version.hpp
===================================================================
--- pkg/RcppArmadillo/inst/include/armadillo_bits/arma_version.hpp 2013-04-16 18:38:26 UTC (rev 4310)
+++ pkg/RcppArmadillo/inst/include/armadillo_bits/arma_version.hpp 2013-04-20 01:15:28 UTC (rev 4311)
@@ -12,9 +12,9 @@
#define ARMA_VERSION_MAJOR 3
-#define ARMA_VERSION_MINOR 800
-#define ARMA_VERSION_PATCH 1
-#define ARMA_VERSION_NAME "Miami Beach"
+#define ARMA_VERSION_MINOR 810
+#define ARMA_VERSION_PATCH 0
+#define ARMA_VERSION_NAME "Newell Highway"
Modified: pkg/RcppArmadillo/inst/include/armadillo_bits/config.hpp
===================================================================
--- pkg/RcppArmadillo/inst/include/armadillo_bits/config.hpp 2013-04-16 18:38:26 UTC (rev 4310)
+++ pkg/RcppArmadillo/inst/include/armadillo_bits/config.hpp 2013-04-20 01:15:28 UTC (rev 4311)
@@ -17,7 +17,7 @@
#if !defined(ARMA_USE_BLAS)
// #define ARMA_USE_BLAS
//// Uncomment the above line if you have BLAS or a high-speed replacement for BLAS,
-//// such as GotoBLAS, Intel's MKL, AMD's ACML, or the Accelerate framework.
+//// such as OpenBLAS, GotoBLAS, Intel's MKL, AMD's ACML, or the Accelerate framework.
//// BLAS is used for matrix multiplication.
//// Without BLAS, matrix multiplication will still work, but might be slower.
#endif
Added: pkg/RcppArmadillo/inst/include/armadillo_bits/fft_engine.hpp
===================================================================
--- pkg/RcppArmadillo/inst/include/armadillo_bits/fft_engine.hpp (rev 0)
+++ pkg/RcppArmadillo/inst/include/armadillo_bits/fft_engine.hpp 2013-04-20 01:15:28 UTC (rev 4311)
@@ -0,0 +1,409 @@
+// This Source Code Form is a compilation of:
+// (1) source code written by Conrad Sanderson, and
+// (2) a modified form of source code referred to as "kissfft.hh".
+//
+// This compilation is Copyright (C) 2013 Conrad Sanderson
+// and is subject to the terms of the Mozilla Public License, v. 2.0.
+//
+// The source code that is distinct and separate from "kissfft.hh"
+// is Copyright (C) 2013 Conrad Sanderson and is subject to the
+// terms of the Mozilla Public License, v. 2.0.
+//
+// If a copy of the MPL was not distributed with this file,
+// You can obtain one at http://mozilla.org/MPL/2.0/.
+//
+// The original "kissfft.hh" source code is licensed under a 3-clause BSD license,
+// as follows:
+//
+// Copyright (c) 2003-2010 Mark Borgerding
+//
+// All rights reserved.
+//
+// Redistribution and use in source and binary forms, with or without modification,
+// are permitted provided that the following conditions are met:
+//
+// * Redistributions of source code must retain the above copyright notice,
+// this list of conditions and the following disclaimer.
+//
+// * Redistributions in binary form must reproduce the above copyright notice,
+// this list of conditions and the following disclaimer in the documentation
+// and/or other materials provided with the distribution.
+//
+// * Neither the author nor the names of any contributors may be used to endorse or promote
+// products derived from this software without specific prior written permission.
+//
+// THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY EXPRESS
+// OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY
+// AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER
+// OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY,
+// OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS
+// OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON
+// ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE
+// OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
+
+
+
+template<typename cx_type, uword fixed_N, bool> struct store {};
+
+template<typename cx_type, uword fixed_N>
+struct store<cx_type, fixed_N, true>
+ {
+ static const uword N = fixed_N;
+
+ arma_aligned cx_type coeffs_array[fixed_N];
+
+ inline store() {}
+ inline store(uword) {}
+
+ arma_inline cx_type* coeffs_ptr() { return &coeffs_array[0]; }
+ arma_inline const cx_type* coeffs_ptr() const { return &coeffs_array[0]; }
+ };
+
+
+
+template<typename cx_type, uword fixed_N>
+struct store<cx_type, fixed_N, false>
+ {
+ const uword N;
+
+ podarray<cx_type> coeffs_array;
+
+ inline store() : N(0) {}
+ inline store(uword in_N) : N(in_N) { coeffs_array.set_size(N); }
+
+ arma_inline cx_type* coeffs_ptr() { return coeffs_array.memptr(); }
+ arma_inline const cx_type* coeffs_ptr() const { return coeffs_array.memptr(); }
+ };
+
+
+
+template<typename cx_type, bool inverse, uword fixed_N = 0>
+class fft_engine : public store<cx_type, fixed_N, (fixed_N > 0)>
+ {
+ public:
+
+ typedef typename get_pod_type<cx_type>::result T;
+
+ using store<cx_type, fixed_N, (fixed_N > 0)>::N;
+ using store<cx_type, fixed_N, (fixed_N > 0)>::coeffs_ptr;
+
+ podarray<uword> residue;
+ podarray<uword> radix;
+
+ podarray<cx_type> tmp_array;
+
+
+ template<bool fill>
+ inline
+ uword
+ calc_radix()
+ {
+ uword i = 0;
+
+ for(uword n = N, r=4; n >= 2; ++i)
+ {
+ while( (n % r) > 0 )
+ {
+ switch(r)
+ {
+ case 2: r = 3; break;
+ case 4: r = 2; break;
+ default: r += 2; break;
+ }
+
+ if(r*r > n) { r = n; }
+ }
+
+ n /= r;
+
+ if(fill)
+ {
+ residue[i] = n;
+ radix[i] = r;
+ }
+ }
+
+ return i;
+ }
+
+
+
+ inline
+ fft_engine(const uword in_N)
+ : store< cx_type, fixed_N, (fixed_N > 0) >(in_N)
+ {
+ arma_extra_debug_sigprint();
+
+ const uword len = calc_radix<false>();
+
+ residue.set_size(len);
+ radix.set_size(len);
+
+ calc_radix<true>();
+
+
+ // calculate the constant coefficients
+
+ cx_type* coeffs = coeffs_ptr();
+
+ const T k = T( (inverse) ? +2 : -2 ) * std::acos( T(-1) ) / T(N);
+
+ for(uword i=0; i < N; ++i) { coeffs[i] = std::exp( cx_type(T(0), i*k) ); }
+ }
+
+
+
+ arma_hot
+ inline
+ void
+ butterfly_2(cx_type* Y, const uword stride, const uword m)
+ {
+ arma_extra_debug_sigprint();
+
+ const cx_type* coeffs = coeffs_ptr();
+
+ for(uword i=0; i < m; ++i)
+ {
+ const cx_type t = Y[i+m] * coeffs[i*stride];
+
+ Y[i+m] = Y[i] - t;
+ Y[i ] += t;
+ }
+ }
+
+
+
+ arma_hot
+ inline
+ void
+ butterfly_3(cx_type* Y, const uword stride, const uword m)
+ {
+ arma_extra_debug_sigprint();
+
+ arma_aligned cx_type tmp[5];
+
+ cx_type* coeffs1 = coeffs_ptr();
+ cx_type* coeffs2 = coeffs1;
+
+ const T coeff_sm_imag = coeffs1[stride*m].imag();
+
+ const uword n = m*2;
+
+ // TODO: rearrange the indices within tmp[] into a more sane order
+
+ for(uword i = m; i > 0; --i)
+ {
+ tmp[1] = Y[m] * (*coeffs1);
+ tmp[2] = Y[n] * (*coeffs2);
+
+ tmp[0] = tmp[1] - tmp[2];
+ tmp[0] *= coeff_sm_imag;
+
+ tmp[3] = tmp[1] + tmp[2];
+
+ Y[m] = cx_type( (Y[0].real() - (0.5*tmp[3].real())), (Y[0].imag() - (0.5*tmp[3].imag())) );
+
+ Y[0] += tmp[3];
+
+
+ Y[n] = cx_type( (Y[m].real() + tmp[0].imag()), (Y[m].imag() - tmp[0].real()) );
+
+ Y[m] += cx_type( -tmp[0].imag(), tmp[0].real() );
+
+ Y++;
+
+ coeffs1 += stride;
+ coeffs2 += stride*2;
+ }
+ }
+
+
+
+ arma_hot
+ inline
+ void
+ butterfly_4(cx_type* Y, const uword stride, const uword m)
+ {
+ arma_extra_debug_sigprint();
+
+ arma_aligned cx_type tmp[7];
+
+ const cx_type* coeffs = coeffs_ptr();
+
+ const uword m2 = m*2;
+ const uword m3 = m*3;
+
+ // TODO: rearrange the indices within tmp[] into a more sane order
+
+ for(uword i=0; i < m; ++i)
+ {
+ tmp[0] = Y[i + m ] * coeffs[i*stride ];
+ tmp[2] = Y[i + m3] * coeffs[i*stride*3];
+ tmp[3] = tmp[0] + tmp[2];
+
+ //tmp[4] = tmp[0] - tmp[2];
+ //tmp[4] = (inverse) ? cx_type( -(tmp[4].imag()), tmp[4].real() ) : cx_type( tmp[4].imag(), -tmp[4].real() );
+
+ tmp[4] = (inverse)
+ ? cx_type( (tmp[2].imag() - tmp[0].imag()), (tmp[0].real() - tmp[2].real()) )
+ : cx_type( (tmp[0].imag() - tmp[2].imag()), (tmp[2].real() - tmp[0].real()) );
+
+ tmp[1] = Y[i + m2] * coeffs[i*stride*2];
+ tmp[5] = Y[i] - tmp[1];
+
+
+ Y[i ] += tmp[1];
+ Y[i + m2] = Y[i] - tmp[3];
+ Y[i ] += tmp[3];
+ Y[i + m ] = tmp[5] + tmp[4];
+ Y[i + m3] = tmp[5] - tmp[4];
+ }
+ }
+
+
+
+ inline
+ arma_hot
+ void
+ butterfly_5(cx_type* Y, const uword stride, const uword m)
+ {
+ arma_extra_debug_sigprint();
+
+ arma_aligned cx_type tmp[13];
+
+ const cx_type* coeffs = coeffs_ptr();
+
+ const T a_real = coeffs[stride*1*m].real();
+ const T a_imag = coeffs[stride*1*m].imag();
+
+ const T b_real = coeffs[stride*2*m].real();
+ const T b_imag = coeffs[stride*2*m].imag();
+
+ cx_type* Y0 = Y;
+ cx_type* Y1 = Y + 1*m;
+ cx_type* Y2 = Y + 2*m;
+ cx_type* Y3 = Y + 3*m;
+ cx_type* Y4 = Y + 4*m;
+
+ for(uword i=0; i < m; ++i)
+ {
+ tmp[0] = (*Y0);
+
+ tmp[1] = (*Y1) * coeffs[stride*1*i];
+ tmp[2] = (*Y2) * coeffs[stride*2*i];
+ tmp[3] = (*Y3) * coeffs[stride*3*i];
+ tmp[4] = (*Y4) * coeffs[stride*4*i];
+
+ tmp[7] = tmp[1] + tmp[4];
+ tmp[8] = tmp[2] + tmp[3];
+ tmp[9] = tmp[2] - tmp[3];
+ tmp[10] = tmp[1] - tmp[4];
+
+ (*Y0) += tmp[7];
+ (*Y0) += tmp[8];
+
+ tmp[5] = tmp[0] + cx_type( ( (tmp[7].real() * a_real) + (tmp[8].real() * b_real) ), ( (tmp[7].imag() * a_real) + (tmp[8].imag() * b_real) ) );
+
+ tmp[6] = cx_type( ( (tmp[10].imag() * a_imag) + (tmp[9].imag() * b_imag) ), ( -(tmp[10].real() * a_imag) - (tmp[9].real() * b_imag) ) );
+
+ (*Y1) = tmp[5] - tmp[6];
+ (*Y4) = tmp[5] + tmp[6];
+
+ tmp[11] = tmp[0] + cx_type( ( (tmp[7].real() * b_real) + (tmp[8].real() * a_real) ), ( (tmp[7].imag() * b_real) + (tmp[8].imag() * a_real) ) );
+
+ tmp[12] = cx_type( ( -(tmp[10].imag() * b_imag) + (tmp[9].imag() * a_imag) ), ( (tmp[10].real() * b_imag) - (tmp[9].real() * a_imag) ) );
+
+ (*Y2) = tmp[11] + tmp[12];
+ (*Y3) = tmp[11] - tmp[12];
+
+ Y0++;
+ Y1++;
+ Y2++;
+ Y3++;
+ Y4++;
+ }
+ }
+
+
+
+ arma_hot
+ inline
+ void
+ butterfly_N(cx_type* Y, const uword stride, const uword m, const uword r)
+ {
+ arma_extra_debug_sigprint();
+
+ const cx_type* coeffs = coeffs_ptr();
+
+ tmp_array.set_min_size(r);
+ cx_type* tmp = tmp_array.memptr();
+
+ for(uword u=0; u < m; ++u)
+ {
+ uword k = u;
+
+ for(uword v=0; v < r; ++v)
+ {
+ tmp[v] = Y[k];
+ k += m;
+ }
+
+ k = u;
+
+ for(uword v=0; v < r; ++v)
+ {
+ Y[k] = tmp[0];
+
+ uword j = 0;
+
+ for(uword w=1; w < r; ++w)
+ {
+ j += stride * k;
+
+ if(j >= N) { j -= N; }
+
+ Y[k] += tmp[w] * coeffs[j];
+ }
+
+ k += m;
+ }
+ }
+ }
+
+
+
+ inline
+ void
+ run(cx_type* Y, const cx_type* X, const uword stage = 0, const uword stride = 1)
+ {
+ arma_extra_debug_sigprint();
+
+ const uword m = residue[stage];
+ const uword r = radix[stage];
+
+ const cx_type *Y_end = Y + r*m;
+
+ if(m == 1)
+ {
+ for(cx_type* Yi = Y; Yi != Y_end; Yi++, X += stride) { (*Yi) = (*X); }
+ }
+ else
+ {
+ const uword next_stage = stage + 1;
+ const uword next_stride = stride * r;
+
+ for(cx_type* Yi = Y; Yi != Y_end; Yi += m, X += stride) { run(Yi, X, next_stage, next_stride); }
+ }
+
+ switch(r)
+ {
+ case 2: butterfly_2(Y, stride, m ); break;
+ case 3: butterfly_3(Y, stride, m ); break;
+ case 4: butterfly_4(Y, stride, m ); break;
+ case 5: butterfly_5(Y, stride, m ); break;
+ default: butterfly_N(Y, stride, m, r); break;
+ }
+ }
+
+
+ };
Added: pkg/RcppArmadillo/inst/include/armadillo_bits/fn_fft.hpp
===================================================================
--- pkg/RcppArmadillo/inst/include/armadillo_bits/fn_fft.hpp (rev 0)
[TRUNCATED]
To get the complete diff run:
svnlook diff /svnroot/rcpp -r 4311
More information about the Rcpp-commits
mailing list