# Cpp11 (R package) vendoring

**pacha.dev/blog**, and kindly contributed to R-bloggers]. (You can report issue about the content on this page here)

Want to share your content on R-bloggers? click here if you have a blog, or here if you don't.

*R and Shiny Training: If you find this blog to be interesting, please note that I offer personalized and group-based training sessions that may be reserved through Buy me a Coffee. Additionally, I provide training services in the Spanish language and am available to discuss means by which I may contribute to your Shiny project.*

*Updated 2023-05-23: I explicitly mention to add the vendoring after creating the project instead of assuming the readers know I run a commented line in the install script. Also, if we include “cpp11.hpp”, then we don’t need to call “doubles.hpp” and “matrix.hpp” because it loads that and also “list.hpp”, “external_pointer.hpp” and other headers. Therefore, I simplified and kept just the call to doubles and matrix in order to load the minimal things required to build the package.*

# Motivation

In my previous post A step by step guide to write an R package that uses C++ code (Ubuntu) I explained the steps in my journey to use C++ function from R.

Now I will explain how I figured out how to use vendoring, which means that I copy the code for the dependencies into your project’s source tree, and in this case this means I copied the C++ headers provided by the `cpp11`

package into my own R package. This ensures the dependency code is fixed and stable until it is updated.

The advantage of vendoring is that changes to the `cpp11`

package could never break your existing code. The disadvantage is that you no longer get bugfixes and new features until you manually run `cpp_vendor()`

in your project.

# Creating an R package with C++ code and vendoring

Similar to the previous post, I started with `create_package("~/github/cpp11gaussjordan")`

. I am using VSCode but all my steps also apply to RStudio.

After opening `~/github/cpp11gaussjordan`

I run `use_cpp11()`

to have a readily availably skeleton in my project.

I run `use_apache_licence()`

to have a `LICENSE`

file and indicate in `DESCRIPTION`

that my package is distributed under the Apache License.

Then I run `cpp_vendor()`

to copy the C++ headers into `inst/include`

.

You can find the final result on my GitHub profile.

## R side

I run `use_r("cpp11gaussjordan-package")`

and added the following contents.

#' @useDynLib cpp11gaussjordan, .registration = TRUE NULL #' Invert (some) square matrices #' @export #' @param A numeric matrix #' @return numeric matrix #' @examples #' A <- matrix(c(2,1,3,-1), nrow = 2, ncol = 2) #' invert_matrix(A) invert_matrix <- function(A) { invert_matrix_(A) } #' Solve (some) linear systems #' @export #' @param A numeric matrix #' @param b numeric matrix #' @return numeric matrix #' @examples #' A <- matrix(c(2,1,3,-1), nrow = 2, ncol = 2) #' b <- matrix(c(7,4), nrow = 2, ncol = 1) solve_system <- function(A,b) { solve_system_(A, b) }

Here `solve_system()`

is the “end-user” function, and `solve_system_()`

is the “development” function. The “development” functions were written in C++ and `cpp11`

handles adding the `.Call()`

when required to link the C++ code and be able to use the C++ functions from R.

## C++ side

I replaced `code.cpp`

contents with the following lines.

#include <cpp11.hpp> #include <cpp11/doubles.hpp> using namespace cpp11; [[cpp11::register]] doubles_matrix<> invert_matrix_(doubles_matrix<> A) { // Obtain (X'X)^-1 via Gauss-Jordan // Check dimensions int N = A.nrow(); int M = A.ncol(); if (N != M) { stop("X must be a square matrix"); } // Copy the matrix writable::doubles_matrix<> Acopy(N, N); for (int i = 0; i < N; i++) { for (int j = 0; j < N; j++) { Acopy(i, j) = A(i, j); } } // Create the identity matrix as a starting point for Gauss-Jordan writable::doubles_matrix<> Ainv(N, N); for (int i = 0; i < N; i++) { for (int j = 0; j < N; j++) { if (i == j) { Ainv(i, j) = 1.0; } else { Ainv(i, j) = 0.0; } } } // Overwrite Ainv by steps with the inverse of A // in other words, find the echelon form of A for (int i = 0; i < M; i++) { double a = Acopy(i, i); for (int j = 0; j < M; j++) { Acopy(i, j) /= a; Ainv(i, j) /= a; } for (int j = 0; j < M; j++) { if (i != j) { a = Acopy(j, i); for (int k = 0; k < M; k++) { Acopy(j, k) -= Acopy(i, k) * a; Ainv(j, k) -= Ainv(i, k) * a; } } } } return Ainv; } [[cpp11::register]] doubles_matrix<> solve_system_(doubles_matrix<> A, doubles_matrix<> b) { // Use Gauss-Jordan to solve the system Ax = b // Check dimensions int N1 = A.nrow(); int M1 = A.ncol(); if (N1 != M1) { stop("A must be a square matrix"); } int N2 = b.nrow(); int M2 = b.ncol(); if (N1 != N2) { stop("b must have the same number of rows as A"); } if (M2 != 1) { stop("b must be a column vector"); } // Obtain the inverse // writable::doubles_matrix<> Acopy(N1,N1); // for (int i = 0; i < N1; i++) // { // for (int j = 0; j < N1; j++) // { // Acopy(i, j) = A(i, j); // } // } // doubles_matrix<> Ainv = invert_matrix_(Acopy); // I don´t need to create a copy in this case doubles_matrix<> Ainv = invert_matrix_(A); // Multiply Ainv by b writable::doubles_matrix<> x(N1, 1); for (int i = 0; i < N1; i++) { x(i, 0) = 0.0; for (int j = 0; j < N1; j++) { x(i, 0) += Ainv(i, j) * b(j, 0); } } return x; }

I also needed a `Makevars`

file to indicate the compiler (either `clang`

or `gcc`

) to use the vendored headers. The content is the following.

PKG_CPPFLAGS = -I../inst/include

Because of the vendoring option, I also had to remove the `LinkingTo: cpp11`

line from `DESCRIPTION`

.

## Putting all together

I created a `vscode-install.r`

file in the root of the project, and I added it to `.Rbuildignore`

. It contains the following lines.

# cpp_vendor() # run only when updating C++ headers clean_dll() cpp_register() document() install() # load_all()

To test that the functions work, after running `install()`

I solved a simple system.

> library(cpp11gaussjordan) > A <- matrix(c(2,1,3,-1), nrow = 2, ncol = 2) > b <- matrix(c(7,4), nrow = 2, ncol = 1) > solve_system(A,b) [,1] [1,] 3.8 [2,] -0.2 > solve(A) %*% b [,1] [1,] 3.8 [2,] -0.2

# References

# Ask a lot

Blog posts like this are summaries of what worked for me. I asked so much on Stackoverflow, that they created the “r-lib-cpp11” for my questions. Please ask there, there is probably something I figured how to solve after many hours on Google.

**leave a comment**for the author, please follow the link and comment on their blog:

**pacha.dev/blog**.

R-bloggers.com offers

**daily e-mail updates**about R news and tutorials about learning R and many other topics. Click here if you're looking to post or find an R/data-science job.

Want to share your content on R-bloggers? click here if you have a blog, or here if you don't.