[Rcpp-commits] r2564 - pkg/RcppGSL/inst/doc/RcppGSL
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Sun Nov 28 19:50:13 CET 2010
Author: edd
Date: 2010-11-28 19:50:13 +0100 (Sun, 28 Nov 2010)
New Revision: 2564
Modified:
pkg/RcppGSL/inst/doc/RcppGSL/RcppGSL.Rnw
Log:
a bunch more edits
Modified: pkg/RcppGSL/inst/doc/RcppGSL/RcppGSL.Rnw
===================================================================
--- pkg/RcppGSL/inst/doc/RcppGSL/RcppGSL.Rnw 2010-11-28 04:55:53 UTC (rev 2563)
+++ pkg/RcppGSL/inst/doc/RcppGSL/RcppGSL.Rnw 2010-11-28 18:50:13 UTC (rev 2564)
@@ -14,6 +14,7 @@
}
\usepackage{vmargin}
\setmargrb{0.75in}{0.75in}{0.75in}{0.75in}
+\usepackage{booktabs} % fancier \hrule
\newcommand{\proglang}[1]{\textsf{#1}}
\newcommand{\pkg}[1]{{\fontseries{b}\selectfont #1}}
@@ -145,7 +146,7 @@
double chisq;
RcppGSL::vector<double> coef(k); // to hold the coefficient vector
- RcppGSL::matrix<double> cov(k,k); // and the covariance matrix
+ RcppGSL::matrix<double> cov(k,k); // and the covariance matrix
// the actual fit requires working memory we allocate and free
gsl_multifit_linear_workspace *work = gsl_multifit_linear_alloc (n, k);
@@ -168,7 +169,7 @@
// we cannot take advantage of automatic memory management
coef.free(); cov.free(); y.free(); X.free();
- return res; // return the result list to R
+ return res; // return the result list to R
} catch( std::exception &ex ) {
forward_exception_to_r( ex );
@@ -179,6 +180,24 @@
}
@
+We first initialize a \textsl{RcppGSL} vector and matrix, each templated to
+the standard numeric type \texttt{double} (and the GSL supports other types
+ranging from lower precision floating point to signed and unsigned integers
+as well as complex numbers). We the reserve another vector and matrix to
+hold the resulting coefficient estimates as well as the estimate of the
+covariance matrix. Next, we allocate workspace using a GSL routine, fit the
+linear model and free the workspace. The next step involves extracting the
+diagonal element from the covariance matrix. We then employ a so-called
+iterator---a common \proglang{C++} idiom from the Standard Template Library
+(STL)---to iterate over the vector of diagonal and transforming it by
+applying the square root function to compute our standard error of the
+estimate. Finally we create a named list with the return value before we free
+temporary memory allocation (a step that has to be done because the
+underlying objects are really \proglang{C} objects conforming to the
+\pkg{GSL} interface and hence without the automatic memory management we
+could have with \proglang{C++} vector or matrix structures as used through
+the \pkg{Rcpp} package) and return the result to \proglang{R}.
+
We should note that \pkg{RcppArmadillo} \citep{CRAN:RcppArmadillo} implements
a matching \texttt{fastLm} function using the Armadillo library by
\cite{Sanderson:2010:Armadillo}, and can do so with more compact code due to
@@ -186,8 +205,13 @@
\section{Vectors}
-\subsection{GSL Vectors}
+This section details the different vector represenations, starting with their
+definition inside the \pkg{GSL}. We then discuss our layering before showing
+how the two types map. A discussion of read-only `vector view' classes
+concludes the section.
+\subsection{\pkg{GSL} Vectors}
+
\pkg{GSL} defines various vector types to manipulate one-dimensionnal
data, similar to \proglang{R} arrays. For example the \verb|gsl_vector| and \verb|gsl_vector_int|
structs are defined as:
@@ -204,8 +228,8 @@
typedef struct {
size_t size;
size_t stride;
- int *data;
- gsl_block_int *block;
+ int * data;
+ gsl_block_int * block;
int owner;
}
gsl_vector_int;
@@ -215,67 +239,57 @@
<<lang=cpp>>=
int i;
+gsl_vector * v = gsl_vector_alloc (3); // allocate a gsl_vector of size 3
-// allocate a gsl_vector of size 3
-gsl_vector * v = gsl_vector_alloc (3);
-
-// fill the vector
-for (i = 0; i < 3; i++) {
+for (i = 0; i < 3; i++) { // fill the vector
gsl_vector_set (v, i, 1.23 + i);
}
-// access elements
-double sum = 0.0 ;
+double sum = 0.0 ; // access elements
for (i = 0; i < 3; i++) {
- sum += gsl_vector_set( v, i ) ;
+ sum += gsl_vector_set( v, i ) ;
}
-// free the memory
-gsl_vector_free (v);
+gsl_vector_free (v); // free the memory
@
\subsection{RcppGSL::vector}
\pkg{RcppGSL} defines the template \texttt{RcppGSL::vector<T>} to manipulate
-\verb|gsl_vector| pointers taking advantage of C++ templates. Using the
-template the previous example becomes:
+\verb|gsl_vector| pointers taking advantage of C++ templates. Using this
+template type, the previous example now becomes:
<<lang=cpp>>=
int i;
+RcppGSL::vector<double> v(3); // allocate a gsl_vector of size 3
-// allocate a gsl_vector of size 3
-RcppGSL::vector<double> v(3);
-
-// fill the vector
-for (i = 0; i < 3; i++) {
+for (i = 0; i < 3; i++) { // fill the vector
v[i] = 1.23 + i ;
}
-// access elements
-double sum = 0.0 ;
+double sum = 0.0 ; // access elements
for (i = 0; i < 3; i++) {
- sum += v[i] ;
+ sum += v[i] ;
}
-// free the memory
-v.free() ;
+v.free() ; // free the memory
@
The class \texttt{RcppGSL::vector<double>} is a smart pointer, that can be used
anywhere where a raw pointer \verb|gsl_vector| can be used, such as the
\verb|gsl_vector_set| and \verb|gsl_vector_get| functions above.
-Beyond the convenience of a nicer syntax for allocation and release of memory,
-the \texttt{RcppGSL::vector} template faciliates interchange of \pkg{GSL}
-vectors with \pkg{Rcpp} concepts. The following example defines a \texttt{.Call}
-compatible function called \verb|sum_gsl_vector_int|
-that operates on a \verb|gsl_vector_int| through the \texttt{RcppGSL::vector<int>}
-template specialization:
+Beyond the convenience of a nicer syntax for allocation and release of
+memory, the \texttt{RcppGSL::vector} template faciliates interchange of
+\pkg{GSL} vectors with \pkg{Rcpp} objects, and hence \pkg{R} objects. The
+following example defines a \texttt{.Call} compatible function called
+\verb|sum_gsl_vector_int| that operates on a \verb|gsl_vector_int| through
+the \texttt{RcppGSL::vector<int>} template specialization:
<<lang=cpp>>=
RCPP_FUNCTION_1( int, sum_gsl_vector_int, RcppGSL::vector<int> vec){
int res = std::accumulate( vec.begin(), vec.end(), 0 ) ;
- vec.free() ;
+ vec.free() ; // we need to free vec after use
return res ;
}
@
@@ -292,20 +306,17 @@
<<lang=cpp>>=
RCPP_FUNCTION_1( double, gsl_vector_sum_2, Rcpp::List data ){
- // grab "x" as a gsl_vector through
- // the RcppGSL::vector<double> class
+ // grab "x" as a gsl_vector through the RcppGSL::vector<double> class
RcppGSL::vector<double> x = data["x"] ;
- // grab "y" as a gsl_vector through
- // the RcppGSL::vector<int> class
+ // grab "y" as a gsl_vector through the RcppGSL::vector<int> class
RcppGSL::vector<int> y = data["y"] ;
double res = 0.0 ;
for( size_t i=0; i< x->size; i++){
res += x[i] * y[i] ;
}
- // as usual with GSL, we need to explicitely free the
- // memory
+ // as usual with GSL, we need to explicitely free the memory
x.free() ;
y.free() ;
@@ -324,40 +335,42 @@
\subsection{Mapping}
-\begin{table}
-\centering
-\begin{small}
-\begin{tabular}{ll}
-\hline
-gsl vector & RcppGSL \\
-\hline
-\texttt{gsl\_vector} & \texttt{RcppGSL::vector<double>} \\
-\texttt{gsl\_vector\_int} & \texttt{RcppGSL::vector<int>} \\
-\texttt{gsl\_vector\_float} & \texttt{RcppGSL::vector<float>} \\
-\texttt{gsl\_vector\_long} & \texttt{RcppGSL::vector<long>} \\
-\texttt{gsl\_vector\_char} & \texttt{RcppGSL::vector<char>} \\
-\texttt{gsl\_vector\_complex} & \texttt{RcppGSL::vector<gsl\_complex>} \\
-\texttt{gsl\_vector\_complex\_float} & \texttt{RcppGSL::vector<gsl\_complex\_float>} \\
-\texttt{gsl\_vector\_complex\_long\_double} & \texttt{RcppGSL::vector<gsl\_complex\_long\_double>} \\
-\texttt{gsl\_vector\_long\_double} & \texttt{RcppGSL::vector<long double>} \\
-\texttt{gsl\_vector\_short} & \texttt{RcppGSL::vector<short>} \\
-\texttt{gsl\_vector\_uchar} & \texttt{RcppGSL::vector<unsigned char>} \\
-\texttt{gsl\_vector\_uint} & \texttt{RcppGSL::vector<unsigned int>} \\
-\texttt{gsl\_vector\_ushort} & \texttt{RcppGSL::vector<insigned short>} \\
-\texttt{gsl\_vector\_ulong} & \texttt{RcppGSL::vector<unsigned long>} \\
+Table~\ref{tab:mappingVectors} shows the mapping between types defined by the
+\pkg{GSL} and their corresponding types in the \pkg{RcppGSL} package.
-\hline
-\end{tabular}
-\end{small}
-\caption{Correspondance between \pkg{GSL} vector types and templates defined
-in \pkg{RcppGSL}.}
+\begin{table}[htb]
+ \centering
+ \begin{small}
+ \begin{tabular}{ll}
+ \toprule
+ gsl vector & RcppGSL \\
+ \midrule
+ \texttt{gsl\_vector} & \texttt{RcppGSL::vector<double>} \\
+ \texttt{gsl\_vector\_int} & \texttt{RcppGSL::vector<int>} \\
+ \texttt{gsl\_vector\_float} & \texttt{RcppGSL::vector<float>} \\
+ \texttt{gsl\_vector\_long} & \texttt{RcppGSL::vector<long>} \\
+ \texttt{gsl\_vector\_char} & \texttt{RcppGSL::vector<char>} \\
+ \texttt{gsl\_vector\_complex} & \texttt{RcppGSL::vector<gsl\_complex>} \\
+ \texttt{gsl\_vector\_complex\_float} & \texttt{RcppGSL::vector<gsl\_complex\_float>} \\
+ \texttt{gsl\_vector\_complex\_long\_double} & \texttt{RcppGSL::vector<gsl\_complex\_long\_double>} \\
+ \texttt{gsl\_vector\_long\_double} & \texttt{RcppGSL::vector<long double>} \\
+ \texttt{gsl\_vector\_short} & \texttt{RcppGSL::vector<short>} \\
+ \texttt{gsl\_vector\_uchar} & \texttt{RcppGSL::vector<unsigned char>} \\
+ \texttt{gsl\_vector\_uint} & \texttt{RcppGSL::vector<unsigned int>} \\
+ \texttt{gsl\_vector\_ushort} & \texttt{RcppGSL::vector<insigned short>} \\
+ \texttt{gsl\_vector\_ulong} & \texttt{RcppGSL::vector<unsigned long>} \\
+ \bottomrule
+ \end{tabular}
+ \end{small}
+ \caption{Correspondance between \pkg{GSL} vector types and templates defined in \pkg{RcppGSL}.}
+ \label{tab:mappingVectors}
\end{table}
\subsection{ Vector Views}
-Several \pkg{GSL} algorithms return \pkg{GSL} vector views as result. \pkg{RcppGSL}
-defines the template class \texttt{RcppGSL::vector\_view} to handle
-vector views using C++ syntax.
+Several \pkg{GSL} algorithms return \pkg{GSL} vector views as their result
+type. \pkg{RcppGSL} defines the template class \texttt{RcppGSL::vector\_view}
+to handle vector views using \proglang{C++} syntax.
<<lang=cpp>>=
extern "C" SEXP test_gsl_vector_view(){
@@ -366,61 +379,61 @@
for( int i=0 ; i<n; i++){
v[i] = i ;
}
- RcppGSL::vector_view<double> v_even = gsl_vector_subvector_with_stride(v, 0, 2, n/2);
- RcppGSL::vector_view<double> v_odd = gsl_vector_subvector_with_stride(v, 1, 2, n/2);
+ RcppGSL::vector_view<double> v_even = gsl_vector_subvector_with_stride(v,0,2,n/2);
+ RcppGSL::vector_view<double> v_odd = gsl_vector_subvector_with_stride(v,1,2,n/2);
List res = List::create(
_["even"] = v_even,
_["odd" ] = v_odd
) ;
- // we only need to free v, the views don t own data
- v.free() ;
-
+ v.free() ; // we only need to free v, the views do not own data
return res ;
}
@
-As with vectors, c++ objects \texttt{RcppGSL::vector\_view} can be implicitely
-converted to their associated gsl view type, as listed in table \ref{tabviews}
-so that they can be passed to compatible gsl algorithms.
+As with vectors, \proglang{C++} objects of type
+\texttt{RcppGSL::vector\_view} can be converted implicitly to their
+associated \pkg{GSL} view type. Table~\ref{tab:mappingVectorViews} displays
+the pairwise correspondance so that the \proglang{C++} objects can be passed
+to compatible \pkg{GSL} algorithms.
-\begin{table}
-\centering
-\begin{small}
-\begin{tabular}{ll}
-\hline
-gsl vector views & RcppGSL \\
-\hline
-\texttt{gsl\_vector\_view} & \texttt{RcppGSL::vector\_view<double>} \\
-\texttt{gsl\_vector\_view\_int} & \texttt{RcppGSL::vector\_view<int>} \\
-\texttt{gsl\_vector\_view\_float} & \texttt{RcppGSL::vector\_view<float>} \\
-\texttt{gsl\_vector\_view\_long} & \texttt{RcppGSL::vector\_view<long>} \\
-\texttt{gsl\_vector\_view\_char} & \texttt{RcppGSL::vector\_view<char>} \\
-\texttt{gsl\_vector\_view\_complex} & \texttt{RcppGSL::vector\_view<gsl\_complex>} \\
-\texttt{gsl\_vector\_view\_complex\_float} & \texttt{RcppGSL::vector\_view<gsl\_complex\_float>} \\
-\texttt{gsl\_vector\_view\_complex\_long\_double} & \texttt{RcppGSL::vector\_view<gsl\_complex\_long\_double>} \\
-\texttt{gsl\_vector\_view\_long\_double} & \texttt{RcppGSL::vector\_view<long double>} \\
-\texttt{gsl\_vector\_view\_short} & \texttt{RcppGSL::vector\_view<short>} \\
-\texttt{gsl\_vector\_view\_uchar} & \texttt{RcppGSL::vector\_view<unsigned char>} \\
-\texttt{gsl\_vector\_view\_uint} & \texttt{RcppGSL::vector\_view<unsigned int>} \\
-\texttt{gsl\_vector\_view\_ushort} & \texttt{RcppGSL::vector\_view<insigned short>} \\
-\texttt{gsl\_vector\_view\_ulong} & \texttt{RcppGSL::vector\_view<unsigned long>} \\
-
-\hline
-\end{tabular}
-\end{small}
-\caption{\label{tabviews}Correspondance between \pkg{GSL} vector view types and templates defined
-in \pkg{RcppGSL}.}
+\begin{table}[htb]
+ \centering
+ \begin{small}
+ \begin{tabular}{ll}
+ \toprule
+ gsl vector views & RcppGSL \\
+ \midrule
+ \texttt{gsl\_vector\_view} & \texttt{RcppGSL::vector\_view<double>} \\
+ \texttt{gsl\_vector\_view\_int} & \texttt{RcppGSL::vector\_view<int>} \\
+ \texttt{gsl\_vector\_view\_float} & \texttt{RcppGSL::vector\_view<float>} \\
+ \texttt{gsl\_vector\_view\_long} & \texttt{RcppGSL::vector\_view<long>} \\
+ \texttt{gsl\_vector\_view\_char} & \texttt{RcppGSL::vector\_view<char>} \\
+ \texttt{gsl\_vector\_view\_complex} & \texttt{RcppGSL::vector\_view<gsl\_complex>} \\
+ \texttt{gsl\_vector\_view\_complex\_float} & \texttt{RcppGSL::vector\_view<gsl\_complex\_float>} \\
+ \texttt{gsl\_vector\_view\_complex\_long\_double} & \texttt{RcppGSL::vector\_view<gsl\_complex\_long\_double>} \\
+ \texttt{gsl\_vector\_view\_long\_double} & \texttt{RcppGSL::vector\_view<long double>} \\
+ \texttt{gsl\_vector\_view\_short} & \texttt{RcppGSL::vector\_view<short>} \\
+ \texttt{gsl\_vector\_view\_uchar} & \texttt{RcppGSL::vector\_view<unsigned char>} \\
+ \texttt{gsl\_vector\_view\_uint} & \texttt{RcppGSL::vector\_view<unsigned int>} \\
+ \texttt{gsl\_vector\_view\_ushort} & \texttt{RcppGSL::vector\_view<insigned short>} \\
+ \texttt{gsl\_vector\_view\_ulong} & \texttt{RcppGSL::vector\_view<unsigned long>} \\
+ \bottomrule
+ \end{tabular}
+ \end{small}
+ \caption{Correspondance between \pkg{GSL} vector view types and templates defined
+ in \pkg{RcppGSL}.}
+ \label{tab:mappingVectorViews}
\end{table}
\section{Matrices}
-\texttt{GSL} defines a set of matrix data types : \texttt{gsl\_matrix},
-\texttt{gsl\_matrix\_int} etc ... for which \pkg{RcppGSL} also defines
+The \pkg{GSL} also defines a set of matrix data types : \texttt{gsl\_matrix},
+\texttt{gsl\_matrix\_int} etc ... for which \pkg{RcppGSL} defines a corresponding
convenience \proglang{C++} wrapper generated by the \texttt{RcppGSL::matrix}
template.
-\subsection{creating matrices}
+\subsection{Creating matrices}
The \texttt{RcppGSL::matrix} template exposes three constructors.
@@ -435,52 +448,61 @@
matrix( int nrow, int ncol)
@
-\subsection{implicit conversion}
+\subsection{Implicit conversion}
\texttt{RcppGSL::matrix} defines implicit conversion to a pointer to
-the associated GSL matrix type, as well as dereferencing operators, making
+the associated \pkg{GSL} matrix type, as well as dereferencing operators, making
the class \texttt{RcppGSL::matrix} look and feel like a pointer to a GSL
matrix type.
<<lang=cpp>>=
- gsltype* data ;
- operator gsltype*(){ return data ; }
- gsltype* operator->() const { return data; }
- gsltype& operator*() const { return *data; }
+ gsltype* data ;
+ operator gsltype*(){ return data ; }
+ gsltype* operator->() const { return data; }
+ gsltype& operator*() const { return *data; }
@
-\subsection{indexing}
+\subsection{Indexing}
-Indexing of GSL matrices is usually the task of the functions
+Indexing of \pkg{GSL} matrices is usually the task of the functions
\texttt{gsl\_matrix\_get}, \texttt{gsl\_matrix\_int\_get}, ... and
\texttt{gsl\_matrix\_set}, \texttt{gsl\_matrix\_int\_set}, ...
\pkg{RcppGSL} takes advantage of both operator overloading and templates
-to make indexing a GSL matrix much more convenient.
+to make indexing a \pkg{GSL} matrix much more convenient.
<<lang=cpp>>=
-// create a matrix of size 10x10
-RcppGSL::matrix<int> mat(10,10) ;
+RcppGSL::matrix<int> mat(10,10); // create a matrix of size 10x10
-// fill the diagonal
-for( int i=0; i<10: i++){
- mat(i,i) = i ;
+for( int i=0; i<10: i++) { // fill the diagonal
+ mat(i,i) = i ;
}
@
\subsection{Methods}
-The \texttt{RcppGSL::matrix} type also defines :
-\begin{itemize}
-\item nrow : extract the number of rows
-\item ncol : extract the number of columns
-\item size : extract the number of elements
-\item free : release the memory
-\end{itemize}
+The \texttt{RcppGSL::matrix} type also defines the following member functions:
+\begin{quote}
+ \begin{itemize}
+ \item[\texttt{nrow}] extracts the number of rows
+ \item[\texttt{ncol}] extract the number of columns
+ \item[\texttt{size}] extracts the number of elements
+ \item[\texttt{free}] releases the memory
+ \end{itemize}
+\end{quote}
-\section{References}
+\section{Summary}
+The GNU Scientific Library (GSL) by \citet{GSL} offers a very comprehensive
+collection of rigorously developed and tested functions for applied
+scientific computing under a common Open Source license. This has lead to
+widespread deployment of \pkg{GSL} among a number of disciplines.
+
+Using the automatic wrapping and converters offered by the \pkg{RcppGSL}
+package presented here, \proglang{R} users and programmers can now deploy
+algorithmns provided by the \pkg{GSL} with greater ease.
+
\bibliographystyle{plainnat}
\bibliography{RcppGSL}
More information about the Rcpp-commits
mailing list