[Rcpp-commits] r3147 - in pkg/RcppArmadillo: . inst inst/include/armadillo_bits
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Sat Jul 23 02:03:18 CEST 2011
Author: edd
Date: 2011-07-23 02:03:17 +0200 (Sat, 23 Jul 2011)
New Revision: 3147
Modified:
pkg/RcppArmadillo/ChangeLog
pkg/RcppArmadillo/DESCRIPTION
pkg/RcppArmadillo/inst/NEWS
pkg/RcppArmadillo/inst/include/armadillo_bits/Col_bones.hpp
pkg/RcppArmadillo/inst/include/armadillo_bits/Col_meat.hpp
pkg/RcppArmadillo/inst/include/armadillo_bits/Mat_bones.hpp
pkg/RcppArmadillo/inst/include/armadillo_bits/Mat_meat.hpp
pkg/RcppArmadillo/inst/include/armadillo_bits/Row_bones.hpp
pkg/RcppArmadillo/inst/include/armadillo_bits/Row_meat.hpp
pkg/RcppArmadillo/inst/include/armadillo_bits/arma_static_assert.hpp
pkg/RcppArmadillo/inst/include/armadillo_bits/arma_version.hpp
pkg/RcppArmadillo/inst/include/armadillo_bits/auxlib_bones.hpp
pkg/RcppArmadillo/inst/include/armadillo_bits/auxlib_meat.hpp
pkg/RcppArmadillo/inst/include/armadillo_bits/config.hpp
pkg/RcppArmadillo/inst/include/armadillo_bits/diskio_meat.hpp
pkg/RcppArmadillo/inst/include/armadillo_bits/fn_as_scalar.hpp
pkg/RcppArmadillo/inst/include/armadillo_bits/fn_norm.hpp
pkg/RcppArmadillo/inst/include/armadillo_bits/fn_svd.hpp
pkg/RcppArmadillo/inst/include/armadillo_bits/fn_toeplitz.hpp
pkg/RcppArmadillo/inst/include/armadillo_bits/glue_relational_meat.hpp
pkg/RcppArmadillo/inst/include/armadillo_bits/glue_times_meat.hpp
pkg/RcppArmadillo/inst/include/armadillo_bits/glue_toeplitz_bones.hpp
pkg/RcppArmadillo/inst/include/armadillo_bits/glue_toeplitz_meat.hpp
pkg/RcppArmadillo/inst/include/armadillo_bits/op_strans_meat.hpp
pkg/RcppArmadillo/inst/include/armadillo_bits/promote_type.hpp
pkg/RcppArmadillo/inst/include/armadillo_bits/span.hpp
pkg/RcppArmadillo/inst/include/armadillo_bits/strip.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
pkg/RcppArmadillo/inst/include/armadillo_bits/unwrap.hpp
Log:
Release 0.2.27 with Armadillo 2.1.91
Modified: pkg/RcppArmadillo/ChangeLog
===================================================================
--- pkg/RcppArmadillo/ChangeLog 2011-07-20 15:27:17 UTC (rev 3146)
+++ pkg/RcppArmadillo/ChangeLog 2011-07-23 00:03:17 UTC (rev 3147)
@@ -1,3 +1,9 @@
+2011-07-22 Dirk Eddelbuettel <edd at debian.org>
+
+ * DESCRIPTION: Release 0.2.27
+
+ * inst/include/*: Updated to release 2.1.91 of Armadillo
+
2011-07-17 Dirk Eddelbuettel <edd at debian.org>
* DESCRIPTION: Release 0.2.26
Modified: pkg/RcppArmadillo/DESCRIPTION
===================================================================
--- pkg/RcppArmadillo/DESCRIPTION 2011-07-20 15:27:17 UTC (rev 3146)
+++ pkg/RcppArmadillo/DESCRIPTION 2011-07-23 00:03:17 UTC (rev 3147)
@@ -1,7 +1,7 @@
Package: RcppArmadillo
Type: Package
Title: Rcpp integration for Armadillo templated linear algebra library
-Version: 0.2.26
+Version: 0.2.27
Date: $Date$
Author: Romain Francois, Dirk Eddelbuettel and Doug Bates
Maintainer: Romain Francois, Dirk Eddelbuettel and Doug Bates <RcppArmadillo-authors at r-enthusiasts.com>
@@ -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 2.0.2). Thus users do not need to
+ Armadillo library (currently version 2.1.91). 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
===================================================================
--- pkg/RcppArmadillo/inst/NEWS 2011-07-20 15:27:17 UTC (rev 3146)
+++ pkg/RcppArmadillo/inst/NEWS 2011-07-23 00:03:17 UTC (rev 3147)
@@ -1,3 +1,14 @@
+0.2.27 2011-07-22
+
+ o Upgraded to Armadillo release 2.1.91 "v2.2 beta 1"
+
+ * faster multiplication of small matrices
+ * faster trans()
+ * faster handling of submatrices by norm()
+ * added economical singular value decomposition: svd_thin()
+ * added circ_toeplitz()
+ * added .is_colvec() & .is_rowvec()
+
0.2.26 2011-07-17
o Upgraded to Armadillo release 2.0.2
Modified: pkg/RcppArmadillo/inst/include/armadillo_bits/Col_bones.hpp
===================================================================
--- pkg/RcppArmadillo/inst/include/armadillo_bits/Col_bones.hpp 2011-07-20 15:27:17 UTC (rev 3146)
+++ pkg/RcppArmadillo/inst/include/armadillo_bits/Col_bones.hpp 2011-07-23 00:03:17 UTC (rev 3147)
@@ -62,7 +62,13 @@
arma_inline subview_col<eT> subvec(const u32 in_row1, const u32 in_row2);
arma_inline const subview_col<eT> subvec(const u32 in_row1, const u32 in_row2) const;
+ arma_inline subview_col<eT> subvec(const span& row_span);
+ arma_inline const subview_col<eT> subvec(const span& row_span) const;
+ //arma_inline subview_col<eT> operator()(const span& row_span);
+ //arma_inline const subview_col<eT> operator()(const span& row_span) const;
+
+
inline void shed_row (const u32 row_num);
inline void shed_rows(const u32 in_row1, const u32 in_row2);
Modified: pkg/RcppArmadillo/inst/include/armadillo_bits/Col_meat.hpp
===================================================================
--- pkg/RcppArmadillo/inst/include/armadillo_bits/Col_meat.hpp 2011-07-20 15:27:17 UTC (rev 3146)
+++ pkg/RcppArmadillo/inst/include/armadillo_bits/Col_meat.hpp 2011-07-23 00:03:17 UTC (rev 3147)
@@ -374,6 +374,74 @@
+template<typename eT>
+arma_inline
+subview_col<eT>
+Col<eT>::subvec(const span& row_span)
+ {
+ arma_extra_debug_sigprint();
+
+ const bool row_all = row_span.whole;
+
+ const u32 local_n_rows = Mat<eT>::n_rows;
+
+ const u32 in_row1 = row_all ? 0 : row_span.a;
+ const u32 in_row2 = row_span.b;
+ const u32 subvec_n_rows = row_all ? local_n_rows : in_row2 - in_row1 + 1;
+
+ arma_debug_check( ( row_all ? false : ((in_row1 > in_row2) || (in_row2 >= local_n_rows)) ), "Col::subvec(): indices out of bounds or incorrectly used");
+
+ return subview_col<eT>(*this, 0, in_row1, subvec_n_rows);
+ }
+
+
+
+template<typename eT>
+arma_inline
+const subview_col<eT>
+Col<eT>::subvec(const span& row_span) const
+ {
+ arma_extra_debug_sigprint();
+
+ const bool row_all = row_span.whole;
+
+ const u32 local_n_rows = Mat<eT>::n_rows;
+
+ const u32 in_row1 = row_all ? 0 : row_span.a;
+ const u32 in_row2 = row_span.b;
+ const u32 subvec_n_rows = row_all ? local_n_rows : in_row2 - in_row1 + 1;
+
+ arma_debug_check( ( row_all ? false : ((in_row1 > in_row2) || (in_row2 >= local_n_rows)) ), "Col::subvec(): indices out of bounds or incorrectly used");
+
+ return subview_col<eT>(*this, 0, in_row1, subvec_n_rows);
+ }
+
+
+
+// template<typename eT>
+// arma_inline
+// subview_col<eT>
+// Col<eT>::operator()(const span& row_span)
+// {
+// arma_extra_debug_sigprint();
+//
+// return subvec(row_span);
+// }
+//
+//
+//
+// template<typename eT>
+// arma_inline
+// const subview_col<eT>
+// Col<eT>::operator()(const span& row_span) const
+// {
+// arma_extra_debug_sigprint();
+//
+// return subvec(row_span);
+// }
+
+
+
//! remove specified row
template<typename eT>
inline
@@ -421,7 +489,7 @@
arrayops::copy( &(X_mem[n_keep_front]), &(t_mem[in_row2+1]), n_keep_back);
}
- steal_mem(X);
+ Mat<eT>::steal_mem(X);
}
@@ -465,7 +533,7 @@
arrayops::inplace_set( &(out_mem[row_num]), eT(0), N );
}
- steal_mem(out);
+ Mat<eT>::steal_mem(out);
}
}
Modified: pkg/RcppArmadillo/inst/include/armadillo_bits/Mat_bones.hpp
===================================================================
--- pkg/RcppArmadillo/inst/include/armadillo_bits/Mat_bones.hpp 2011-07-20 15:27:17 UTC (rev 3146)
+++ pkg/RcppArmadillo/inst/include/armadillo_bits/Mat_bones.hpp 2011-07-23 00:03:17 UTC (rev 3147)
@@ -420,10 +420,9 @@
inline row_iterator end_row (const u32 row_num);
inline const_row_iterator end_row (const u32 row_num) const;
- // inline void clear();
- // inline void resize(const u32 in_n_elem);
- // arma_inline bool empty() const;
- // arma_inline u32 size() const;
+ inline void clear();
+ inline bool empty() const;
+ inline u32 size() const;
template<u32 fixed_n_rows, u32 fixed_n_cols>
class fixed : public Mat<eT>
Modified: pkg/RcppArmadillo/inst/include/armadillo_bits/Mat_meat.hpp
===================================================================
--- pkg/RcppArmadillo/inst/include/armadillo_bits/Mat_meat.hpp 2011-07-20 15:27:17 UTC (rev 3146)
+++ pkg/RcppArmadillo/inst/include/armadillo_bits/Mat_meat.hpp 2011-07-23 00:03:17 UTC (rev 3147)
@@ -568,6 +568,7 @@
if( (x_mem_state == 0) && (x_n_elem > arma_config::mat_prealloc) && (layout_ok == true) )
{
reset();
+ // note: calling reset() also prevents fixed size matrices from changing size or using non-local memory
access::rw(n_rows) = x_n_rows;
access::rw(n_cols) = x_n_cols;
@@ -5022,7 +5023,40 @@
+//! resets this matrix to an empty matrix
template<typename eT>
+inline
+void
+Mat<eT>::clear()
+ {
+ reset();
+ }
+
+
+
+//! returns true if the matrix has no elements
+template<typename eT>
+inline
+bool
+Mat<eT>::empty() const
+ {
+ return (n_elem == 0);
+ }
+
+
+
+//! returns the number of elements in this matrix
+template<typename eT>
+inline
+u32
+Mat<eT>::size() const
+ {
+ return n_elem;
+ }
+
+
+
+template<typename eT>
template<u32 fixed_n_rows, u32 fixed_n_cols>
arma_inline
void
Modified: pkg/RcppArmadillo/inst/include/armadillo_bits/Row_bones.hpp
===================================================================
--- pkg/RcppArmadillo/inst/include/armadillo_bits/Row_bones.hpp 2011-07-20 15:27:17 UTC (rev 3146)
+++ pkg/RcppArmadillo/inst/include/armadillo_bits/Row_bones.hpp 2011-07-23 00:03:17 UTC (rev 3147)
@@ -62,7 +62,13 @@
arma_inline subview_row<eT> subvec(const u32 in_col1, const u32 in_col2);
arma_inline const subview_row<eT> subvec(const u32 in_col1, const u32 in_col2) const;
+ arma_inline subview_row<eT> subvec(const span& col_span);
+ arma_inline const subview_row<eT> subvec(const span& col_span) const;
+ // arma_inline subview_row<eT> operator()(const span& col_span);
+ // arma_inline const subview_row<eT> operator()(const span& col_span) const;
+
+
inline void shed_col (const u32 col_num);
inline void shed_cols(const u32 in_col1, const u32 in_col2);
Modified: pkg/RcppArmadillo/inst/include/armadillo_bits/Row_meat.hpp
===================================================================
--- pkg/RcppArmadillo/inst/include/armadillo_bits/Row_meat.hpp 2011-07-20 15:27:17 UTC (rev 3146)
+++ pkg/RcppArmadillo/inst/include/armadillo_bits/Row_meat.hpp 2011-07-23 00:03:17 UTC (rev 3147)
@@ -349,6 +349,74 @@
+template<typename eT>
+arma_inline
+subview_row<eT>
+Row<eT>::subvec(const span& col_span)
+ {
+ arma_extra_debug_sigprint();
+
+ const bool col_all = col_span.whole;
+
+ const u32 local_n_cols = Mat<eT>::n_cols;
+
+ const u32 in_col1 = col_all ? 0 : col_span.a;
+ const u32 in_col2 = col_span.b;
+ const u32 subvec_n_cols = col_all ? local_n_cols : in_col2 - in_col1 + 1;
+
+ arma_debug_check( ( col_all ? false : ((in_col1 > in_col2) || (in_col2 >= local_n_cols)) ), "Row::subvec(): indices out of bounds or incorrectly used");
+
+ return subview_row<eT>(*this, 0, in_col1, subvec_n_cols);
+ }
+
+
+
+template<typename eT>
+arma_inline
+const subview_row<eT>
+Row<eT>::subvec(const span& col_span) const
+ {
+ arma_extra_debug_sigprint();
+
+ const bool col_all = col_span.whole;
+
+ const u32 local_n_cols = Mat<eT>::n_cols;
+
+ const u32 in_col1 = col_all ? 0 : col_span.a;
+ const u32 in_col2 = col_span.b;
+ const u32 subvec_n_cols = col_all ? local_n_cols : in_col2 - in_col1 + 1;
+
+ arma_debug_check( ( col_all ? false : ((in_col1 > in_col2) || (in_col2 >= local_n_cols)) ), "Row::subvec(): indices out of bounds or incorrectly used");
+
+ return subview_row<eT>(*this, 0, in_col1, subvec_n_cols);
+ }
+
+
+
+// template<typename eT>
+// arma_inline
+// subview_row<eT>
+// Row<eT>::operator()(const span& col_span)
+// {
+// arma_extra_debug_sigprint();
+//
+// return subvec(col_span);
+// }
+//
+//
+//
+// template<typename eT>
+// arma_inline
+// const subview_row<eT>
+// Row<eT>::operator()(const span& col_span) const
+// {
+// arma_extra_debug_sigprint();
+//
+// return subvec(col_span);
+// }
+
+
+
//! remove specified columns
template<typename eT>
inline
@@ -396,7 +464,7 @@
arrayops::copy( &(X_mem[n_keep_front]), &(t_mem[in_col2+1]), n_keep_back);
}
- steal_mem(X);
+ Mat<eT>::steal_mem(X);
}
@@ -440,7 +508,7 @@
arrayops::inplace_set( &(out_mem[col_num]), eT(0), N );
}
- steal_mem(out);
+ Mat<eT>::steal_mem(out);
}
}
Modified: pkg/RcppArmadillo/inst/include/armadillo_bits/arma_static_assert.hpp
===================================================================
--- pkg/RcppArmadillo/inst/include/armadillo_bits/arma_static_assert.hpp 2011-07-20 15:27:17 UTC (rev 3146)
+++ pkg/RcppArmadillo/inst/include/armadillo_bits/arma_static_assert.hpp 2011-07-23 00:03:17 UTC (rev 3147)
@@ -17,13 +17,21 @@
//! Classes for primitive compile time assertions and checks (until the next version of C++)
-template<bool> struct arma_static_assert;
-template<> struct arma_static_assert<true> {};
+template<bool x>
+struct arma_static_assert
+ {
+ static const char
+ static_error[ x ? +1 : -1 ];
+ };
-template<bool> struct arma_static_check;
-template<> struct arma_static_check<false> {};
+template<bool x>
+struct arma_static_check
+ {
+ static const char
+ static_error[ x ? -1 : +1 ];
+ };
Modified: pkg/RcppArmadillo/inst/include/armadillo_bits/arma_version.hpp
===================================================================
--- pkg/RcppArmadillo/inst/include/armadillo_bits/arma_version.hpp 2011-07-20 15:27:17 UTC (rev 3146)
+++ pkg/RcppArmadillo/inst/include/armadillo_bits/arma_version.hpp 2011-07-23 00:03:17 UTC (rev 3147)
@@ -17,12 +17,12 @@
#define ARMA_VERSION_MAJOR 2
-#define ARMA_VERSION_MINOR 0
-#define ARMA_VERSION_PATCH 2
-#define ARMA_VERSION_NAME "Carnivorous Sugar Glider"
-// http://en.wikipedia.org/wiki/Sugar_glider
+#define ARMA_VERSION_MINOR 1
+#define ARMA_VERSION_PATCH 91
+#define ARMA_VERSION_NAME "2.2 beta 1"
+
struct arma_version
{
static const unsigned int major = ARMA_VERSION_MAJOR;
Modified: pkg/RcppArmadillo/inst/include/armadillo_bits/auxlib_bones.hpp
===================================================================
--- pkg/RcppArmadillo/inst/include/armadillo_bits/auxlib_bones.hpp 2011-07-20 15:27:17 UTC (rev 3146)
+++ pkg/RcppArmadillo/inst/include/armadillo_bits/auxlib_bones.hpp 2011-07-23 00:03:17 UTC (rev 3147)
@@ -162,7 +162,13 @@
template<typename T, typename T1>
inline static bool svd(Mat< std::complex<T> >& U, Col<T>& S, Mat< std::complex<T> >& V, const Base< std::complex<T>, T1>& X);
+ template<typename eT, typename T1>
+ inline static bool svd_thin(Mat<eT>& U, Col<eT>& S, Mat<eT>& V, const Base<eT,T1>& X, const char mode);
+ template<typename T, typename T1>
+ inline static bool svd_thin(Mat< std::complex<T> >& U, Col<T>& S, Mat< std::complex<T> >& V, const Base< std::complex<T>, T1>& X, const char mode);
+
+
//
// solve
Modified: pkg/RcppArmadillo/inst/include/armadillo_bits/auxlib_meat.hpp
===================================================================
--- pkg/RcppArmadillo/inst/include/armadillo_bits/auxlib_meat.hpp 2011-07-20 15:27:17 UTC (rev 3146)
+++ pkg/RcppArmadillo/inst/include/armadillo_bits/auxlib_meat.hpp 2011-07-23 00:03:17 UTC (rev 3147)
@@ -2114,6 +2114,297 @@
+template<typename eT, typename T1>
+inline
+bool
+auxlib::svd_thin(Mat<eT>& U, Col<eT>& S, Mat<eT>& V, const Base<eT,T1>& X, const char mode)
+ {
+ arma_extra_debug_sigprint();
+
+ #if defined(ARMA_USE_LAPACK)
+ {
+ Mat<eT> A(X.get_ref());
+
+ blas_int m = A.n_rows;
+ blas_int n = A.n_cols;
+ blas_int lda = A.n_rows;
+
+ S.set_size( static_cast<u32>((std::min)(m,n)) );
+
+ blas_int ldu = 0;
+ blas_int ldvt = 0;
+
+ char jobu;
+ char jobvt;
+
+ switch(mode)
+ {
+ case 'l':
+ jobu = 'S';
+ jobvt = 'N';
+
+ ldu = m;
+ ldvt = 1;
+
+ U.set_size( static_cast<u32>(ldu), static_cast<u32>((std::min)(m,n)) );
+ V.reset();
+
+ break;
+
+
+ case 'r':
+ jobu = 'N';
+ jobvt = 'S';
+
+ ldu = 1;
+ ldvt = (std::min)(m,n);
+
+ U.reset();
+ V.set_size( static_cast<u32>(ldvt), static_cast<u32>(n) );
+
+ break;
+
+
+ case 'b':
+ jobu = 'S';
+ jobvt = 'S';
+
+ ldu = m;
+ ldvt = (std::min)(m,n);
+
+ U.set_size( static_cast<u32>(ldu), static_cast<u32>((std::min)(m,n)) );
+ V.set_size( static_cast<u32>(ldvt), static_cast<u32>(n) );
+
+ break;
+
+
+ default:
+ U.reset();
+ S.reset();
+ V.reset();
+ return false;
+ }
+
+
+ if(A.is_empty())
+ {
+ U.eye();
+ S.reset();
+ V.eye();
+ return true;
+ }
+
+
+ blas_int lwork = 2 * (std::max)(blas_int(1), (std::max)( (3*(std::min)(m,n) + (std::max)(m,n)), 5*(std::min)(m,n) ) );
+ blas_int info = 0;
+
+
+ podarray<eT> work( static_cast<u32>(lwork) );
+
+ // let gesvd_() calculate the optimum size of the workspace
+ blas_int lwork_tmp = -1;
+
+ lapack::gesvd<eT>
+ (
+ &jobu, &jobvt,
+ &m, &n,
+ A.memptr(), &lda,
+ S.memptr(),
+ U.memptr(), &ldu,
+ V.memptr(), &ldvt,
+ work.memptr(), &lwork_tmp,
+ &info
+ );
+
+ if(info == 0)
+ {
+ blas_int proposed_lwork = static_cast<blas_int>(work[0]);
+ if(proposed_lwork > lwork)
+ {
+ lwork = proposed_lwork;
+ work.set_size( static_cast<u32>(lwork) );
+ }
+
+ lapack::gesvd<eT>
+ (
+ &jobu, &jobvt,
+ &m, &n,
+ A.memptr(), &lda,
+ S.memptr(),
+ U.memptr(), &ldu,
+ V.memptr(), &ldvt,
+ work.memptr(), &lwork,
+ &info
+ );
+
+ op_strans::apply(V,V); // op_strans will work out that an in-place transpose can be done
+ }
+
+ return (info == 0);
+ }
+ #else
+ {
+ arma_ignore(U);
+ arma_ignore(S);
+ arma_ignore(V);
+ arma_ignore(X);
+ arma_ignore(mode);
+ arma_stop("svd(): use of LAPACK needs to be enabled");
+ return false;
+ }
+ #endif
+ }
+
+
+
+template<typename T, typename T1>
+inline
+bool
+auxlib::svd_thin(Mat< std::complex<T> >& U, Col<T>& S, Mat< std::complex<T> >& V, const Base< std::complex<T>, T1>& X, const char mode)
+ {
+ arma_extra_debug_sigprint();
+
+ typedef std::complex<T> eT;
+
+ #if defined(ARMA_USE_LAPACK)
+ {
+ Mat<eT> A(X.get_ref());
+
+ blas_int m = A.n_rows;
+ blas_int n = A.n_cols;
+ blas_int lda = A.n_rows;
+
+ S.set_size( static_cast<u32>((std::min)(m,n)) );
+
+ blas_int ldu = 0;
+ blas_int ldvt = 0;
+
+ char jobu;
+ char jobvt;
+
+ switch(mode)
+ {
+ case 'l':
+ jobu = 'S';
+ jobvt = 'N';
+
+ ldu = m;
+ ldvt = 1;
+
+ U.set_size( static_cast<u32>(ldu), static_cast<u32>((std::min)(m,n)) );
+ V.reset();
+
+ break;
+
+
+ case 'r':
+ jobu = 'N';
+ jobvt = 'S';
+
+ ldu = 1;
+ ldvt = (std::min)(m,n);
+
+ U.reset();
+ V.set_size( static_cast<u32>(ldvt), static_cast<u32>(n) );
+
+ break;
+
+
+ case 'b':
+ jobu = 'S';
+ jobvt = 'S';
+
+ ldu = m;
+ ldvt = (std::min)(m,n);
+
+ U.set_size( static_cast<u32>(ldu), static_cast<u32>((std::min)(m,n)) );
+ V.set_size( static_cast<u32>(ldvt), static_cast<u32>(n) );
+
+ break;
+
+
+ default:
+ U.reset();
+ S.reset();
+ V.reset();
+ return false;
+ }
+
+
+ if(A.is_empty())
+ {
+ U.eye();
+ S.reset();
+ V.eye();
+ return true;
+ }
+
+
+ blas_int lwork = 2 * (std::max)(blas_int(1), (std::max)( (3*(std::min)(m,n) + (std::max)(m,n)), 5*(std::min)(m,n) ) );
+ blas_int info = 0;
+
+
+ podarray<eT> work( static_cast<u32>(lwork) );
+ podarray<T> rwork( static_cast<u32>(5*(std::min)(m,n)) );
+
+ // let gesvd_() calculate the optimum size of the workspace
+ blas_int lwork_tmp = -1;
+
+ lapack::cx_gesvd<T>
+ (
+ &jobu, &jobvt,
+ &m, &n,
+ A.memptr(), &lda,
+ S.memptr(),
+ U.memptr(), &ldu,
+ V.memptr(), &ldvt,
+ work.memptr(), &lwork_tmp,
+ rwork.memptr(),
+ &info
+ );
+
+ if(info == 0)
+ {
+ blas_int proposed_lwork = static_cast<blas_int>(real(work[0]));
+ if(proposed_lwork > lwork)
+ {
+ lwork = proposed_lwork;
+ work.set_size( static_cast<u32>(lwork) );
+ }
+
+ lapack::cx_gesvd<T>
+ (
+ &jobu, &jobvt,
+ &m, &n,
+ A.memptr(), &lda,
+ S.memptr(),
+ U.memptr(), &ldu,
+ V.memptr(), &ldvt,
+ work.memptr(), &lwork,
+ rwork.memptr(),
+ &info
+ );
+
+ op_htrans::apply(V,V); // op_strans will work out that an in-place transpose can be done
+ }
+
+ return (info == 0);
+ }
+ #else
+ {
+ arma_ignore(U);
+ arma_ignore(S);
+ arma_ignore(V);
+ arma_ignore(X);
+ arma_ignore(mode);
+ arma_stop("svd(): use of LAPACK needs to be enabled");
+ return false;
+ }
+ #endif
+ }
+
+
+
//! Solve a system of linear equations.
//! Assumes that A.n_rows = A.n_cols and B.n_rows = A.n_rows
template<typename eT>
Modified: pkg/RcppArmadillo/inst/include/armadillo_bits/config.hpp
===================================================================
--- pkg/RcppArmadillo/inst/include/armadillo_bits/config.hpp 2011-07-20 15:27:17 UTC (rev 3146)
+++ pkg/RcppArmadillo/inst/include/armadillo_bits/config.hpp 2011-07-23 00:03:17 UTC (rev 3147)
@@ -75,3 +75,19 @@
#undef ARMA_USE_ATLAS
#undef ARMA_ATLAS_INCLUDE_DIR
#endif
+
+#if defined(ARMA_DONT_USE_LAPACK)
+ #undef ARMA_USE_LAPACK
+#endif
+
+#if defined(ARMA_DONT_USE_BLAS)
+ #undef ARMA_USE_BLAS
+#endif
+
+#if defined(ARMA_DONT_PRINT_LOGIC_ERRORS)
+ #undef ARMA_PRINT_LOGIC_ERRORS
+#endif
+
+#if defined(ARMA_DONT_PRINT_RUNTIME_ERRORS)
+ #undef ARMA_PRINT_RUNTIME_ERRORS
+#endif
Modified: pkg/RcppArmadillo/inst/include/armadillo_bits/diskio_meat.hpp
===================================================================
--- pkg/RcppArmadillo/inst/include/armadillo_bits/diskio_meat.hpp 2011-07-20 15:27:17 UTC (rev 3146)
+++ pkg/RcppArmadillo/inst/include/armadillo_bits/diskio_meat.hpp 2011-07-23 00:03:17 UTC (rev 3147)
@@ -1680,6 +1680,8 @@
return false;
}
}
+
+ return false;
}
@@ -2292,6 +2294,8 @@
return false;
}
}
+
+ return false;
}
Modified: pkg/RcppArmadillo/inst/include/armadillo_bits/fn_as_scalar.hpp
===================================================================
--- pkg/RcppArmadillo/inst/include/armadillo_bits/fn_as_scalar.hpp 2011-07-20 15:27:17 UTC (rev 3146)
+++ pkg/RcppArmadillo/inst/include/armadillo_bits/fn_as_scalar.hpp 2011-07-23 00:03:17 UTC (rev 3147)
@@ -86,7 +86,7 @@
const u32 B_n_rows = (tmp2.do_trans == false) ? B.n_rows : B.n_cols;
const u32 B_n_cols = (tmp2.do_trans == false) ? B.n_cols : B.n_rows;
- const eT val = tmp1.val * tmp2.val;
+ const eT val = tmp1.get_val() * tmp2.get_val();
arma_debug_check( (A_n_rows != 1) || (B_n_cols != 1) || (A_n_cols != B_n_rows), "as_scalar(): incompatible dimensions" );
@@ -136,7 +136,7 @@
const u32 C_n_rows = (tmp3.do_trans == false) ? C.n_rows : C.n_cols;
const u32 C_n_cols = (tmp3.do_trans == false) ? C.n_cols : C.n_rows;
- const eT val = tmp1.val * tmp2.val * tmp3.val;
+ const eT val = tmp1.get_val() * tmp2.get_val() * tmp3.get_val();
arma_debug_check
(
@@ -192,7 +192,7 @@
const u32 C_n_rows = (tmp3.do_trans == false) ? C.n_rows : C.n_cols;
const u32 C_n_cols = (tmp3.do_trans == false) ? C.n_cols : C.n_rows;
- const eT val = tmp1.val * tmp2.val * tmp3.val;
+ const eT val = tmp1.get_val() * tmp2.get_val() * tmp3.get_val();
arma_debug_check
(
@@ -287,7 +287,7 @@
const u32 C_n_rows = (tmp3.do_trans == false) ? C.n_rows : C.n_cols;
const u32 C_n_cols = (tmp3.do_trans == false) ? C.n_cols : C.n_rows;
- const eT val = tmp1.val * tmp2.val * tmp3.val;
+ const eT val = tmp1.get_val() * tmp2.get_val() * tmp3.get_val();
arma_debug_check
(
Modified: pkg/RcppArmadillo/inst/include/armadillo_bits/fn_norm.hpp
===================================================================
--- pkg/RcppArmadillo/inst/include/armadillo_bits/fn_norm.hpp 2011-07-20 15:27:17 UTC (rev 3146)
+++ pkg/RcppArmadillo/inst/include/armadillo_bits/fn_norm.hpp 2011-07-23 00:03:17 UTC (rev 3147)
@@ -25,26 +25,50 @@
arma_extra_debug_sigprint();
typedef typename T1::pod_type T;
- typedef typename Proxy<T1>::ea_type ea_type;
-
+
T acc = T(0);
-
- ea_type P = A.get_ea();
- const u32 N = A.get_n_elem();
-
- u32 i,j;
-
- for(i=0, j=1; j<N; i+=2, j+=2)
+
+ if(Proxy<T1>::prefer_at_accessor == false)
{
- acc += std::abs(P[i]);
- acc += std::abs(P[j]);
+ typename Proxy<T1>::ea_type P = A.get_ea();
+
+ const u32 N = A.get_n_elem();
+
+ u32 i,j;
+
+ for(i=0, j=1; j<N; i+=2, j+=2)
+ {
+ acc += std::abs(P[i]);
+ acc += std::abs(P[j]);
+ }
+
+ if(i < N)
+ {
+ acc += std::abs(P[i]);
+ }
}
-
- if(i < N)
+ else
{
- acc += std::abs(P[i]);
+ const u32 n_rows = A.get_n_rows();
+ const u32 n_cols = A.get_n_cols();
+
+ for(u32 col=0; col<n_cols; ++col)
+ {
+ u32 i,j;
+
+ for(i=0, j=1; j<n_rows; i+=2, j+=2)
+ {
+ acc += std::abs(A.at(i,col));
+ acc += std::abs(A.at(j,col));
+ }
+
+ if(i < n_rows)
+ {
+ acc += std::abs(A.at(i,col));
+ }
+ }
}
-
+
return acc;
}
@@ -59,30 +83,59 @@
arma_extra_debug_sigprint();
arma_ignore(junk);
- typedef typename T1::pod_type T;
- typedef typename Proxy<T1>::ea_type ea_type;
+ typedef typename T1::pod_type T;
T acc = T(0);
- ea_type P = A.get_ea();
- const u32 N = A.get_n_elem();
-
- u32 i,j;
-
- for(i=0, j=1; j<N; i+=2, j+=2)
+ if(Proxy<T1>::prefer_at_accessor == false)
{
- const T tmp_i = P[i];
- const T tmp_j = P[j];
+ typename Proxy<T1>::ea_type P = A.get_ea();
- acc += tmp_i * tmp_i;
- acc += tmp_j * tmp_j;
+ const u32 N = A.get_n_elem();
+
+ u32 i,j;
+
+ for(i=0, j=1; j<N; i+=2, j+=2)
+ {
+ const T tmp_i = P[i];
+ const T tmp_j = P[j];
+
+ acc += tmp_i * tmp_i;
+ acc += tmp_j * tmp_j;
+ }
+
+ if(i < N)
+ {
+ const T tmp_i = P[i];
+
+ acc += tmp_i * tmp_i;
+ }
}
-
- if(i < N)
+ else
{
- const T tmp_i = P[i];
+ const u32 n_rows = A.get_n_rows();
+ const u32 n_cols = A.get_n_cols();
- acc += tmp_i * tmp_i;
+ for(u32 col=0; col<n_cols; ++col)
+ {
+ u32 i,j;
+
+ for(i=0, j=1; j<n_rows; i+=2, j+=2)
+ {
+ const T tmp_i = A.at(i,col);
+ const T tmp_j = A.at(j,col);
+
+ acc += tmp_i * tmp_i;
+ acc += tmp_j * tmp_j;
+ }
+
+ if(i < n_rows)
+ {
+ const T tmp_i = A.at(i,col);
+
+ acc += tmp_i * tmp_i;
+ }
+ }
}
return std::sqrt(acc);
@@ -99,19 +152,34 @@
arma_extra_debug_sigprint();
arma_ignore(junk);
- typedef typename T1::pod_type T;
- typedef typename Proxy<T1>::ea_type ea_type;
+ typedef typename T1::pod_type T;
T acc = T(0);
- ea_type P = A.get_ea();
- const u32 N = A.get_n_elem();
-
- for(u32 i=0; i<N; ++i)
+ if(Proxy<T1>::prefer_at_accessor == false)
{
- const T tmp = std::abs(P[i]);
- acc += tmp*tmp;
+ typename Proxy<T1>::ea_type P = A.get_ea();
+
+ const u32 N = A.get_n_elem();
+
+ for(u32 i=0; i<N; ++i)
+ {
+ const T tmp = std::abs(P[i]);
+ acc += tmp*tmp;
+ }
}
+ else
+ {
+ const u32 n_rows = A.get_n_rows();
+ const u32 n_cols = A.get_n_cols();
+
+ for(u32 col=0; col<n_cols; ++col)
+ for(u32 row=0; row<n_rows; ++row)
+ {
+ const T tmp = std::abs(A.at(row,col));
+ acc += tmp*tmp;
+ }
+ }
return std::sqrt(acc);
}
@@ -126,25 +194,39 @@
{
arma_extra_debug_sigprint();
- typedef typename T1::pod_type T;
- typedef typename Proxy<T1>::ea_type ea_type;
+ typedef typename T1::pod_type T;
T acc = T(0);
- ea_type P = A.get_ea();
- const u32 N = A.get_n_elem();
-
- u32 i,j;
-
- for(i=0, j=1; j<N; i+=2, j+=2)
+ if(Proxy<T1>::prefer_at_accessor == false)
{
- acc += std::pow(std::abs(P[i]), k);
- acc += std::pow(std::abs(P[j]), k);
+ typename Proxy<T1>::ea_type P = A.get_ea();
+
+ const u32 N = A.get_n_elem();
+
+ u32 i,j;
+
+ for(i=0, j=1; j<N; i+=2, j+=2)
+ {
+ acc += std::pow(std::abs(P[i]), k);
+ acc += std::pow(std::abs(P[j]), k);
+ }
+
+ if(i < N)
+ {
+ acc += std::pow(std::abs(P[i]), k);
+ }
}
-
- if(i < N)
+ else
{
- acc += std::pow(std::abs(P[i]), k);
+ const u32 n_rows = A.get_n_rows();
+ const u32 n_cols = A.get_n_cols();
+
+ for(u32 col=0; col<n_cols; ++col)
+ for(u32 row=0; row<n_rows; ++row)
+ {
+ acc += std::pow(std::abs(A.at(row,col)), k);
+ }
}
return std::pow(acc, T(1)/T(k));
Modified: pkg/RcppArmadillo/inst/include/armadillo_bits/fn_svd.hpp
===================================================================
--- pkg/RcppArmadillo/inst/include/armadillo_bits/fn_svd.hpp 2011-07-20 15:27:17 UTC (rev 3146)
+++ pkg/RcppArmadillo/inst/include/armadillo_bits/fn_svd.hpp 2011-07-23 00:03:17 UTC (rev 3147)
@@ -110,4 +110,51 @@
+template<typename T1>
+inline
+bool
+svd_thin
+ (
+ Mat<typename T1::elem_type>& U,
+ Col<typename T1::pod_type >& S,
+ Mat<typename T1::elem_type>& V,
+ const Base<typename T1::elem_type,T1>& X,
+ const char mode = 'b',
+ const typename arma_blas_type_only<typename T1::elem_type>::result* junk = 0
+ )
+ {
+ arma_extra_debug_sigprint();
+ arma_ignore(junk);
+
+ typedef typename T1::elem_type eT;
+
+ arma_debug_check
+ (
+ ( ((void*)(&U) == (void*)(&S)) || (&U == &V) || ((void*)(&S) == (void*)(&V)) ),
+ "svd_thin(): two or more output objects are the same object"
+ );
+
+ arma_debug_check
+ (
+ ( (mode != 'l') && (mode != 'r') && (mode != 'b') ),
+ "svd_thin(): parameter 'mode' is incorrect"
+ );
+
+
+ // auxlib::svd_thin() makes an internal copy of X
+ const bool status = auxlib::svd_thin(U, S, V, X, mode);
+
+ if(status == false)
+ {
+ U.reset();
+ S.reset();
+ V.reset();
+ arma_bad("svd_thin(): failed to converge", false);
+ }
+
+ return status;
+ }
+
+
+
//! @}
Modified: pkg/RcppArmadillo/inst/include/armadillo_bits/fn_toeplitz.hpp
===================================================================
--- pkg/RcppArmadillo/inst/include/armadillo_bits/fn_toeplitz.hpp 2011-07-20 15:27:17 UTC (rev 3146)
+++ pkg/RcppArmadillo/inst/include/armadillo_bits/fn_toeplitz.hpp 2011-07-23 00:03:17 UTC (rev 3147)
@@ -1,5 +1,7 @@
-// Copyright (C) 2010 NICTA (www.nicta.com.au)
-// Copyright (C) 2010 Conrad Sanderson
+// Copyright (C) 2010-2011 NICTA (www.nicta.com.au)
+// Copyright (C) 2010-2011 Conrad Sanderson
+// Copyright (C) 2011 Alcatel Lucent
+// Copyright (C) 2011 Gerhard Schreiber
//
// This file is part of the Armadillo C++ library.
// It is provided without any warranty of fitness
@@ -40,4 +42,16 @@
+template<typename T1>
+inline
+Glue<T1, T1, glue_toeplitz_circ>
+circ_toeplitz(const Base<typename T1::elem_type,T1>& X)
+ {
+ arma_extra_debug_sigprint();
+
+ return Glue<T1, T1, glue_toeplitz_circ>( X.get_ref(), X.get_ref() );
+ }
+
+
+
//! @}
Modified: pkg/RcppArmadillo/inst/include/armadillo_bits/glue_relational_meat.hpp
===================================================================
--- pkg/RcppArmadillo/inst/include/armadillo_bits/glue_relational_meat.hpp 2011-07-20 15:27:17 UTC (rev 3146)
+++ pkg/RcppArmadillo/inst/include/armadillo_bits/glue_relational_meat.hpp 2011-07-23 00:03:17 UTC (rev 3147)
@@ -1,5 +1,5 @@
-// Copyright (C) 2009-2010 NICTA (www.nicta.com.au)
-// Copyright (C) 2009-2010 Conrad Sanderson
+// Copyright (C) 2009-2011 NICTA (www.nicta.com.au)
+// Copyright (C) 2009-2011 Conrad Sanderson
//
// This file is part of the Armadillo C++ library.
// It is provided without any warranty of fitness
@@ -16,6 +16,102 @@
+#undef operator_rel
+#undef operator_str
+
+#undef arma_applier_mat
+#undef arma_applier_cube
+
+
+#define arma_applier_mat(operator_rel, operator_str) \
+ {\
+ const Proxy<T1> P1(X.A);\
+ const Proxy<T2> P2(X.B);\
+ \
+ arma_debug_assert_same_size(P1, P2, operator_str);\
+ \
+ const u32 n_rows = P1.get_n_rows();\
+ const u32 n_cols = P1.get_n_cols();\
+ \
+ out.set_size(n_rows, n_cols);\
+ \
+ u32* out_mem = out.memptr();\
+ \
+ const bool prefer_at_accessor = (Proxy<T1>::prefer_at_accessor || Proxy<T2>::prefer_at_accessor);\
+ \
+ if(prefer_at_accessor == false)\
+ {\
+ typename Proxy<T1>::ea_type A = P1.get_ea();\
+ typename Proxy<T2>::ea_type B = P2.get_ea();\
+ \
+ const u32 n_elem = out.n_elem;\
+ \
+ for(u32 i=0; i<n_elem; ++i)\
+ {\
+ out_mem[i] = (A[i] operator_rel B[i]) ? u32(1) : u32(0);\
+ }\
+ }\
+ else\
+ {\
+ u32 count = 0;\
+ \
+ for(u32 col=0; col<n_cols; ++col)\
+ {\
+ for(u32 row=0; row<n_rows; ++row, ++count)\
+ {\
+ out_mem[count] = (P1.at(row,col) operator_rel P2.at(row,col)) ? u32(1) : u32(0);\
+ }\
+ }\
+ }\
+ }
+
+
+
+
+#define arma_applier_cube(operator_rel, operator_str) \
+ {\
+ const ProxyCube<T1> P1(X.A);\
+ const ProxyCube<T2> P2(X.B);\
+ \
+ arma_debug_assert_same_size(P1, P2, operator_str);\
+ \
+ const u32 n_rows = P1.get_n_rows();\
+ const u32 n_cols = P1.get_n_cols();\
+ const u32 n_slices = P1.get_n_slices();\
+ \
+ out.set_size(n_rows, n_cols, n_slices);\
+ \
+ u32* out_mem = out.memptr();\
+ \
+ const bool prefer_at_accessor = (ProxyCube<T1>::prefer_at_accessor || ProxyCube<T2>::prefer_at_accessor);\
+ \
+ if(prefer_at_accessor == false)\
+ {\
+ typename ProxyCube<T1>::ea_type A = P1.get_ea();\
+ typename ProxyCube<T2>::ea_type B = P2.get_ea();\
+ \
+ const u32 n_elem = out.n_elem;\
+ \
+ for(u32 i=0; i<n_elem; ++i)\
+ {\
+ out_mem[i] = (A[i] operator_rel B[i]) ? u32(1) : u32(0);\
+ }\
+ }\
+ else\
+ {\
+ u32 count = 0;\
+ \
[TRUNCATED]
To get the complete diff run:
svnlook diff /svnroot/rcpp -r 3147
More information about the Rcpp-commits
mailing list