Overwriting entries from a matrix in cpp11

[pachamaltese]

Hi!

I have some Rcpp code to fill entries in a matrix, here's a MWE of what it does:

#include <Rcpp.h>

using namespace Rcpp;

void invert_tri(NumericMatrix &M, int K) {
    for(int i=0 ; i<K ; ++i) {
        for(int j=i+1 ; j<K ; ++j) {
            M(j, i) = M(i, j);
        }
    }
}

// [[Rcpp::export]]
NumericMatrix multiply_entries(NumericMatrix M, int factor) {
    int K = M.ncol();

    for(int i=0 ; i<K ; ++i) {
        for(int j=0 ; j<K ; ++j) {
            M(i, j) = M(i, j) * factor;
        }
    }

    invert_tri(M, K);

    return M;
}

Then in R I would use it like this:

M = matrix(c(1,NA,2,4), nrow = 2, ncol = 2)

multiply_entries(M, 2L)

This changes

     [,1] [,2]
[1,]    1    2
[2,]   NA    4

into

     [,1] [,2]
[1,]    2    4
[2,]    4    8

Then, with cpp11 I would sketch the following code

#include <cpp11.hpp>
#include <cpp11/doubles.hpp>

using namespace cpp11;

void invert_tri(writable::doubles_matrix<> &M, int K) {
    for(int i=0 ; i<K ; ++i) {
        for(int j=i+1 ; j<K ; ++j) {
            M(j, i) = M(i, j);
        }
    }
}

[[cpp11::register]]
doubles_matrix<> multiply_entries(writable::doubles_matrix<> M, int factor) {
    int K = M.ncol();

    for(int i=0 ; i<K ; ++i) {
        for(int j=0 ; j<K ; ++j) {
            M(i, j) = M(i, j) * factor;
        }
    }

    invert_tri(M, K);

    return M;
}

But here doing M(j, i) += M(i, j) would fail and with = fails because of implicit deletion. Do you have any ideas for an alternative approach to translate the working code.

[technocrat]

@pachamaltese

alternative approach

library(Rcpp)

cppFunction('int add(int x, int y, int z) {
  int sum = x + y + z;
  return sum;
}')

add
#> function (x, y, z) 
#> .Call(<pointer: 0x1051c8d00>, x, y, z)
add(1, 2, 3)
#> [1] 6

^{Created on 2022-12-25 with reprex v2.0.2}

[pachamaltese]

thanks, the problem is that for the 3x3 case I really need to keep the array
my MWE was maybe too reduced

[technocrat]

Sorry. I understood you have working C++ code and was demonstrating that its possible to call it with Rcpp, using a different example, since I'm too lazy to track down the hpp files needed for your codr.

[pachamaltese]

A solution (maybe not efficient) was

M(j,i)+ = (-1.0 * M(j,i)) + (M(i, j) * factor);