Approximate infinite series with a guaranteed error margin!
Installing this package in R from this GitHub page is straightforward. Make sure you have the package devtool
installed and that your Operating System has a C compiler that R can access. This is automatically true in a standard installation in Linux, requires Xcode or gcc (depends on MacOS version) in Mac, and Rtools in Windows. Then, just run:
devtools::install_github("GuidoAMoreira/sumR")
It is possible to define a series at the C level and use a Rcpp
wrapper function. For this to be possible, all that is needed is the use of the depends Rcpp attribute and including the adequate header. The following code at the R level exports a function that sums a series.
library(Rcpp)
sourceCpp(code='
#include <Rcpp.h>
// [[Rcpp::depends(sumR)]]
#include <sumR.h>
long double some_series(long n, double *p)
{
long double out = n * log1pl(-p[0]);
return out;
}
// [[Rcpp::export]]
double sum_series(double param)
{
double parameter[1];
long double r;
long n; // Number of iterations. Does not require initialization.
parameter[0] = param;
r = infiniteSum(some_series, parameter, log1p(-parameter[0]), 0, exp(-35), 100000, 0, &n);
Rcpp::Rcout << "Summation took " << n << " iterations to converge.\\n";
return (double)r;
}
')
sum_series(0.08)
The available functions and their arguments are listed below.
Package sumR
facilitates using its low-level C or C++ function in other packages. In order to use it, these steps are necessary:
- Make sure that the DESCRIPTION file in your package includes
sumR
in its LinkingTo and Imports fields. - Make sure that the NAMESPACE file in your packages includes a line with
import(sumR)
. If you are using the roxygen2 documentation package, this can be achieved by adding#' @import sumR
in one of your R files, such as mypackage-package.R, before runningroxygen2::roxygenize()
. - Include sumR's API header file, sumRAPI.h, in your C or C++ file that will use the desired sumR low-level function.
- Use sumR's functions in your code.
The following code exemplifies a C file in a package after steps 1. and 2. above were taken. See Writing R Extensions to learn about the SEXP
type and related macros and functions.
#include <Rinternals.h>
#include <Rmath.h> // Required for the log1p and log1pl functions
#include <sumR.h>
long double some_series(long n, double *p)
{
long double out = n * log1pl(-p[0]);
return out;
}
SEXP sum_series(SEXP param)
{
double parameter[1];
long double r;
long n; // Number of iterations. Doesn't require initialization.
parameter[0] = REAL(param)[0];
r = infiniteSum(some_series, parameter, log1p(-parameter[0]), 0, exp(-35), 100000, 0, &n);
Rprintf("Summation took %d iterations to converge.\n", n);
return Rf_ScalarReal((double)r);
}
Then your package can have an R wrapper function such as:
#' Wrapper function that sums a series for a given parameter
#' @export
sumSeries <- function(p) .Call("sum_series", param = p, PACKAGE = "mypackage")
The this function can be tested after the package has been installed with:
sumSeries(0.08)
The interfaced functions from sumR are:
long double infiniteSum(long double logFun(long k, double *Theta), double *params, double logL, int alternating, double eps, long maxIter, long n0, long* n);
long double infiniteSumToThreshold(long double logFun(long k, double *Theta), double *params, int alternating, double eps, long maxIter, long n0, long* n);
long double infiniteErrorBoundingPairs(long double logFun(long k, double *Theta), double *params, double logL, double eps, long maxIter, long n0, long* n);
long double infiniteBatches(long double logFun(long k, double *Theta), double *params, long batch_size, double eps, long maxIter, long n0, long* n);
long double sumNTimes(long double logFun(long k, double *Theta), double *params, long n, long n0);
Function infiniteSum
dispatches the arguments to infiniteSumToThreshold
, infiniteErrorBoundingPairs
or infiniteBatches
depending on the value of logL
and returns the result of the respectively chosen function. Namely, it is the first if logL
< log(0.5), the second if logL
< 0 and the third otherwise.
See the help documentation in the sumR package for information about the interfaced function arguments. sumNTimes
is documented under finiteSum
.
When making a wrapper function, we have found that manually typecasting the result of the low-level C function to double before passing it to R is more stable in some systems than straight up using Rf_ScalarReal on the long double result.
Since all sumR
code is in C, a C++ functor cannot be passed as the logFun argument, not even in a templated wrapper.
None of the C functions do error checking. It is the user's responsibility to pass adequate arguments. Pay particular mind to the C stack limit (the functions create an array of maxIter
size).
If defining an R function with cppFunction
, be mindful that argument k
of logFun
is a single long
. If multiple values are expected, they should be dealt with in the wrapper.