g.cpp

#include <Rcpp.h>
using namespace Rcpp;

// [[Rcpp::export]]
NumericVector g(NumericMatrix x) {
    double sigma = 0.25;
    double S0 = 1;
    double Time = 1;
    double r = 0;
    double K = 1;
    int nrow = x.nrow(), ncol = x.ncol();
    NumericVector out(ncol);
    NumericMatrix temp(nrow,ncol);
    for (int j = 0; j < ncol; j++) {
        double acc = 0;
        for(int i = 0; i < nrow; i++) {
            acc += x(i,j);
            temp(i,j) = acc;
        }
    }
    for (int i = 0; i < nrow; i++) {
        temp(i,_) = S0*exp(sigma*temp(i,_)+(r-pow(sigma,2)/2)*i/nrow*Time) ; // Works until here but seems to do not give the exact same numerical value
    }
    for (int j = 0; j < ncol; j++) { //Breaks from here
        out(j) = 0;
        for(int i = 0; i < nrow; i++) {
            out(j) += temp(i,j);
        }
        out(j) *= 1/nrow;
        if ( out(j) > K ) {out(j) -= K;} else {out(j)=0;}
    }
    return out;
}