2
votes

below is the R code where functions stl_example_model and stl_example_jacobian in a dll are referenced.

library(RcppSundials)
library(microbenchmark)

stl_example_model = getNativeSymbolInfo(name = "example_model_stl",
                                        PACKAGE = "RcppSundials")$address
stl_example_jacobian = getNativeSymbolInfo(name = "example_jacobian_stl",
                                           PACKAGE = "RcppSundials")$address

simulationstl = wrap_cvodes(times = 1:5, 
                            states_ = rep(1.0,5), 
                            parameters_ = 0.1, 
                            forcings_data_ =list(cbind(1:3600,1:3600)),
            ..., 
                            model_ = example_model, jacobian_ = example_jacobian) 

I have sourced another function (test.cpp) using the Rcpp::sourceCpp. How do I a make a reference to this function (similar to stl_example_model) using Rcpp?

Test.cpp

#define ARMA_DONT_USE_CXX11
#include <RcppArmadillo.h>
#include <Rcpp.h>
#include <array>
#include <vector>
using namespace std;
using namespace Rcpp;

extern "C" {
  array<vector<double>, 2> example_model_stl(const double& t, const vector<double>& states, 
            const vector<double>& parameters, const vector<double>& forcings) {

    vector<double> derivatives(states.size());  
    for(int i = 0; i < states.size(); i++) {
      derivatives[i] = -states[i]*parameters[0];
    }
    vector<double> observed{forcings[0]};
    array<vector<double>, 2> output{derivatives, observed};
    return output;
  }

  arma::mat example_jacobian_stl(const double& t, const vector<double>& states, 
            const vector<double>& parameters, const vector<double>& forcings) {
    arma::mat output = arma::eye(states.size(), states.size());
    output = -parameters[0]*output;
    return output;
  }

  array<vector<double>, 2> example_dae_stl(const double& t, const vector<double>& states, 
            const vector<double>& derivatives, const vector<double>& parameters, 
            const vector<double>& forcings) {

    vector<double> residues(states.size());  
    residues[0] = -0.04*states[0] + 1e4*states[1]*states[2];
    residues[1] = -residues[0] - 3e7*states[1]*states[1] - derivatives[1];
    residues[0] = residues[0] - derivatives[0];
    residues[2] = states[0] + states[1] + states[2] - 1.0;
    vector<double> observed{forcings[0]};
    array<vector<double>, 2> output{residues, observed};
    return output;
  }

};
1
The RcppSundials package is dormant. I'm not sure it would be wise to count on this package for an analysis. - coatless
Besides that, you know how to access code from it -- as it is in package form. So do the same with your code: Make it into a package. - Dirk Eddelbuettel
@coatless. I know but this is really about sundials not the R package wrapper. I hope the wrapper dont break when updating the Sundials library. - Oyvind Foshaug
@DirkEddelbuettel . I finally grasped what you have been posting a few times. I created a package using RcppArmadillo, That works very good. - Oyvind Foshaug

1 Answers

1
votes

It does make a lot of sense to use a package as mentioned by Dirk in the comments. However, it is not strictly necessary.

  1. The PACKAGE argument of getNativeSymbols is actually optional. All loaded DLLs are searched if it is omitted.

  2. Compiling your example produces the warning

    No Rcpp::export attributes or RCPP_MODULE declarations found in source

    As a consequence the produced DLL is not loaded into R's memory. You can change this by adding an exported dummy function.

Example:

// [[Rcpp::depends("RcppArmadillo")]]
#define ARMA_DONT_USE_CXX11
#include <RcppArmadillo.h>
#include <array>
#include <vector>
using namespace std;
using namespace Rcpp;

extern "C" {
    array<vector<double>, 2> example_model_stl(const double& t, const vector<double>& states, 
                                               const vector<double>& parameters, const vector<double>& forcings) {

        vector<double> derivatives(states.size());  
        for(int i = 0; i < states.size(); i++) {
            derivatives[i] = -states[i]*parameters[0];
        }
        vector<double> observed{forcings[0]};
        array<vector<double>, 2> output{derivatives, observed};
        return output;
    }

    arma::mat example_jacobian_stl(const double& t, const vector<double>& states, 
                                   const vector<double>& parameters, const vector<double>& forcings) {
        arma::mat output = arma::eye(states.size(), states.size());
        output = -parameters[0]*output;
        return output;
    }

    array<vector<double>, 2> example_dae_stl(const double& t, const vector<double>& states, 
                                             const vector<double>& derivatives, const vector<double>& parameters, 
                                             const vector<double>& forcings) {

        vector<double> residues(states.size());  
        residues[0] = -0.04*states[0] + 1e4*states[1]*states[2];
        residues[1] = -residues[0] - 3e7*states[1]*states[1] - derivatives[1];
        residues[0] = residues[0] - derivatives[0];
        residues[2] = states[0] + states[1] + states[2] - 1.0;
        vector<double> observed{forcings[0]};
        array<vector<double>, 2> output{residues, observed};
        return output;
    }

};

// [[Rcpp::export]]
void dummy() {}

/***R
getNativeSymbolInfo("example_model_stl")$address
*/

Output:

> Rcpp::sourceCpp('58881676.cpp')

> getNativeSymbolInfo("example_model_stl")$address
<pointer: 0x7fd438ed5f80>
attr(,"class")
[1] "NativeSymbol"

This way it should be possible to use such functions together with RcppSundials.