Custom exp and log (#44)

* custom exp and log

* uncomment debug stuff in Makevars

* replace a few more calls

* use exp from openlibm

* use log from openlibm

* use log

* use it everywhere

* cleanup

* remove duplicated function

* remove newlines

Author: vandenman

.Rbuildignore

@@ -1,3 +1,5 @@
+^renv$
+^renv\.lock$
 ^.*\.Rproj$
 ^\.Rproj\.user$
 

R/RcppExports.R

@@ -1,6 +1,18 @@
 # Generated by using Rcpp::compileAttributes() -> do not edit by hand
 # Generator token: 10BE3573-1514-4C36-9D1C-5A225CD40393
 
+get_explog_switch <- function() {
+    .Call(`_bgms_get_explog_switch`)
+}
+
+rcpp_ieee754_exp <- function(x) {
+    .Call(`_bgms_rcpp_ieee754_exp`, x)
+}
+
+rcpp_ieee754_log <- function(x) {
+    .Call(`_bgms_rcpp_ieee754_log`, x)
+}
+
 sample_omrf_gibbs <- function(no_states, no_variables, no_categories, interactions, thresholds, iter) {
     .Call(`_bgms_sample_omrf_gibbs`, no_states, no_variables, no_categories, interactions, thresholds, iter)
 }

bgms.Rproj

@@ -16,5 +16,6 @@ StripTrailingWhitespace: Yes
 
 BuildType: Package
 PackageUseDevtools: Yes
+PackageCleanBeforeInstall: No
 PackageInstallArgs: --no-multiarch --with-keep.source
 PackageCheckArgs: --as-cran

src/Makevars.win

@@ -1 +1 @@
-PKG_LIBS = $(LAPACK_LIBS) $(BLAS_LIBS) $(FLIBS)
+PKG_LIBS = $(LAPACK_LIBS) $(BLAS_LIBS) $(FLIBS)

src/RcppExports.cpp

@@ -11,6 +11,38 @@ Rcpp::Rostream<true>&  Rcpp::Rcout = Rcpp::Rcpp_cout_get();
 Rcpp::Rostream<false>& Rcpp::Rcerr = Rcpp::Rcpp_cerr_get();
 #endif
 
+// get_explog_switch
+Rcpp::String get_explog_switch();
+RcppExport SEXP _bgms_get_explog_switch() {
+BEGIN_RCPP
+    Rcpp::RObject rcpp_result_gen;
+    Rcpp::RNGScope rcpp_rngScope_gen;
+    rcpp_result_gen = Rcpp::wrap(get_explog_switch());
+    return rcpp_result_gen;
+END_RCPP
+}
+// rcpp_ieee754_exp
+Rcpp::NumericVector rcpp_ieee754_exp(Rcpp::NumericVector x);
+RcppExport SEXP _bgms_rcpp_ieee754_exp(SEXP xSEXP) {
+BEGIN_RCPP
+    Rcpp::RObject rcpp_result_gen;
+    Rcpp::RNGScope rcpp_rngScope_gen;
+    Rcpp::traits::input_parameter< Rcpp::NumericVector >::type x(xSEXP);
+    rcpp_result_gen = Rcpp::wrap(rcpp_ieee754_exp(x));
+    return rcpp_result_gen;
+END_RCPP
+}
+// rcpp_ieee754_log
+Rcpp::NumericVector rcpp_ieee754_log(Rcpp::NumericVector x);
+RcppExport SEXP _bgms_rcpp_ieee754_log(SEXP xSEXP) {
+BEGIN_RCPP
+    Rcpp::RObject rcpp_result_gen;
+    Rcpp::RNGScope rcpp_rngScope_gen;
+    Rcpp::traits::input_parameter< Rcpp::NumericVector >::type x(xSEXP);
+    rcpp_result_gen = Rcpp::wrap(rcpp_ieee754_log(x));
+    return rcpp_result_gen;
+END_RCPP
+}
 // sample_omrf_gibbs
 IntegerMatrix sample_omrf_gibbs(int no_states, int no_variables, IntegerVector no_categories, NumericMatrix interactions, NumericMatrix thresholds, int iter);
 RcppExport SEXP _bgms_sample_omrf_gibbs(SEXP no_statesSEXP, SEXP no_variablesSEXP, SEXP no_categoriesSEXP, SEXP interactionsSEXP, SEXP thresholdsSEXP, SEXP iterSEXP) {
@@ -143,6 +175,9 @@ END_RCPP
 }
 
 static const R_CallMethodDef CallEntries[] = {
+    {"_bgms_get_explog_switch", (DL_FUNC) &_bgms_get_explog_switch, 0},
+    {"_bgms_rcpp_ieee754_exp", (DL_FUNC) &_bgms_rcpp_ieee754_exp, 1},
+    {"_bgms_rcpp_ieee754_log", (DL_FUNC) &_bgms_rcpp_ieee754_log, 1},
     {"_bgms_sample_omrf_gibbs", (DL_FUNC) &_bgms_sample_omrf_gibbs, 6},
     {"_bgms_sample_bcomrf_gibbs", (DL_FUNC) &_bgms_sample_bcomrf_gibbs, 8},
     {"_bgms_gibbs_sampler", (DL_FUNC) &_bgms_gibbs_sampler, 28},

src/custom_exp.cpp

@@ -0,0 +1,29 @@
+#include "Rcpp.h"
+#include "explog_switch.h"
+
+// [[Rcpp::export]]
+Rcpp::String get_explog_switch() {
+#if USE_CUSTOM_LOG
+      return "custom";
+  #else
+      return "standard";
+  #endif
+}
+
+// [[Rcpp::export]]
+Rcpp::NumericVector rcpp_ieee754_exp(Rcpp::NumericVector x) {
+  Rcpp::NumericVector y(x.size());
+  for (int i = 0; i < x.size(); i++) {
+    y[i] = MY_EXP(x[i]);
+  }
+  return y;
+}
+
+// [[Rcpp::export]]
+Rcpp::NumericVector rcpp_ieee754_log(Rcpp::NumericVector x) {
+  Rcpp::NumericVector y(x.size());
+  for (int i = 0; i < x.size(); i++) {
+    y[i] = MY_LOG(x[i]);
+  }
+  return y;
+}

src/data_simulation.cpp

@@ -1,4 +1,5 @@
 #include <Rcpp.h>
+#include "explog_switch.h"
 using namespace Rcpp;
 
 // [[Rcpp::export]]
@@ -53,7 +54,7 @@ IntegerMatrix sample_omrf_gibbs(int no_states,
         for(int category = 0; category < no_categories[variable]; category++) {
           exponent = thresholds(variable, category);
           exponent += (category + 1) * rest_score;
-          cumsum += std::exp(exponent);
+          cumsum += MY_EXP(exponent);
           probabilities[category + 1] = cumsum;
         }
 
@@ -132,7 +133,7 @@ IntegerMatrix sample_bcomrf_gibbs(int no_states,
               (category - reference_category[variable]);
             //The pairwise interactions
             exponent += category * rest_score;
-            cumsum += std::exp(exponent);
+            cumsum += MY_EXP(exponent);
             probabilities[category] = cumsum;
           }
         } else {
@@ -141,7 +142,7 @@ IntegerMatrix sample_bcomrf_gibbs(int no_states,
           for(int category = 0; category < no_categories[variable]; category++) {
             exponent = thresholds(variable, category);
             exponent += (category + 1) * rest_score;
-            cumsum += std::exp(exponent);
+            cumsum += MY_EXP(exponent);
             probabilities[category + 1] = cumsum;
           }
         }