Re: [Rcpp-devel] tracking volatile bug
Gabor Grothendieck has written at Sat, 22 Nov 2014 19:17:27 -0500
On Sat, Nov 22, 2014 at 1:43 PM, Sokol Sergueiserguei.so...@gmail.com wrote:
Let try to stick with regular rcpp code
(file: matrix_norm.cpp):
//[[Rcpp::depends(RcppArmadillo)]]
#include RcppArmadillo.h
using namespace Rcpp;
using namespace arma;
// [[Rcpp::export]]
double nmat(mat A) {
Function Matrix_norm_r_=Environment(package:Matrix)[norm];
double res=asdouble(Matrix_norm_r_(A, 1));
return res;
}
When called in R as:
library(Rcpp)
library(Matrix)
sourceCpp(matrix_norm.cpp)
gctorture(TRUE)
nmat(as.matrix(pi))
it gives an error:
Erreur : 'getCharCE' doit être appelé sur un CHARSXP
(my English translation: Error: 'getCharCE' must be called on CHARSXP)
Something was irregular on my side here?
Serguei.
Try replacing the line that sets res with:
double res=asdouble(Matrix_norm_r_(NumericMatrix(wrap(A)),
CharacterVector::create(1)));
Be sure to try it on a fresh session since errors from prior runs
could have messed up R.
Yes, this version works too.
The least change that leads to a working code that I have found (thanks
to Martin) is
double res=asdouble(Matrix_norm_r_(A, CharacterVector(1)));
but it is still unclear what is wrong with passing a plain string.
In quickref, all examples pass just numbers without any prior conversion
to NumericVector or IntegerVector. Strings are exception to this?
I have an impression that they are supposed to work too without
any explicit conversion. If you set gctorture(FALSE), it works and
there is no warning or error message that plain strings are not supported
in this place.
So, gctorture() reveals that somewhere there is an un-protected pointer
which is getting corrupted. If it is needed to be said, I don't pretend
that this unprotected pointer is in Rcpp code. May be it is in R, may be
in Matrix, it is still unclear.
___
Rcpp-devel mailing list
Rcpp-devel@lists.r-forge.r-project.org
https://lists.r-forge.r-project.org/cgi-bin/mailman/listinfo/rcpp-devel
Re: [Rcpp-devel] tracking volatile bug
In quickref, all examples pass just numbers without any prior conversion
to NumericVector or IntegerVector. Strings are exception to this?
Strings are different than numbers in that they don't occupy a fixed
amount of space, so some memory management is needed for them.
Bill Dunlap
TIBCO Software
wdunlap tibco.com
On Sun, Nov 23, 2014 at 11:07 AM, Sokol Serguei serguei.so...@gmail.com
wrote:
Gabor Grothendieck has written at Sat, 22 Nov 2014 19:17:27 -0500
On Sat, Nov 22, 2014 at 1:43 PM, Sokol Sergueiserguei.so...@gmail.com
wrote:
Let try to stick with regular rcpp code
(file: matrix_norm.cpp):
//[[Rcpp::depends(RcppArmadillo)]]
#include RcppArmadillo.h
using namespace Rcpp;
using namespace arma;
// [[Rcpp::export]]
double nmat(mat A) {
Function Matrix_norm_r_=Environment(package:Matrix)[norm];
double res=asdouble(Matrix_norm_r_(A, 1));
return res;
}
When called in R as:
library(Rcpp)
library(Matrix)
sourceCpp(matrix_norm.cpp)
gctorture(TRUE)
nmat(as.matrix(pi))
it gives an error:
Erreur : 'getCharCE' doit être appelé sur un CHARSXP
(my English translation: Error: 'getCharCE' must be called on CHARSXP)
Something was irregular on my side here?
Serguei.
Try replacing the line that sets res with:
double res=asdouble(Matrix_norm_r_(NumericMatrix(wrap(A)),
CharacterVector::create(1)));
Be sure to try it on a fresh session since errors from prior runs
could have messed up R.
Yes, this version works too.
The least change that leads to a working code that I have found (thanks to
Martin) is
double res=asdouble(Matrix_norm_r_(A, CharacterVector(1)));
but it is still unclear what is wrong with passing a plain string.
In quickref, all examples pass just numbers without any prior conversion
to NumericVector or IntegerVector. Strings are exception to this?
I have an impression that they are supposed to work too without
any explicit conversion. If you set gctorture(FALSE), it works and
there is no warning or error message that plain strings are not supported
in this place.
So, gctorture() reveals that somewhere there is an un-protected pointer
which is getting corrupted. If it is needed to be said, I don't pretend
that this unprotected pointer is in Rcpp code. May be it is in R, may be
in Matrix, it is still unclear.
___
Rcpp-devel mailing list
Rcpp-devel@lists.r-forge.r-project.org
https://lists.r-forge.r-project.org/cgi-bin/mailman/listinfo/rcpp-devel
___
Rcpp-devel mailing list
Rcpp-devel@lists.r-forge.r-project.org
https://lists.r-forge.r-project.org/cgi-bin/mailman/listinfo/rcpp-devel
Re: [Rcpp-devel] tracking volatile bug
Martin Morgan has written at Sun, 23 Nov 2014 14:09:27 -0800
On 11/23/2014 11:07 AM, Sokol Serguei wrote:
Gabor Grothendieck has written at Sat, 22 Nov 2014 19:17:27 -0500
On Sat, Nov 22, 2014 at 1:43 PM, Sokol
Sergueiserguei.so...@gmail.com wrote:
Let try to stick with regular rcpp code
(file: matrix_norm.cpp):
//[[Rcpp::depends(RcppArmadillo)]]
#include RcppArmadillo.h
using namespace Rcpp;
using namespace arma;
// [[Rcpp::export]]
double nmat(mat A) {
Function Matrix_norm_r_=Environment(package:Matrix)[norm];
double res=asdouble(Matrix_norm_r_(A, 1));
return res;
}
...
but obviously the other pairlist templates need protection, too. Maybe
it's enough, instead, to add protection to grow
template typename T
SEXP grow(const T head, SEXP tail) {
ShieldSEXP y(tail);
return internal::grow__dispatch( typename
traits::is_namedT::type(), head, y );
}
inline SEXP grow( const char* head, SEXP tail ) {
ShieldSEXP y(tail);
return grow( Rf_mkString(head), y ) ;
}
I confirm that this patch solved the problem. Nice find Martin.
Serguei.
___
Rcpp-devel mailing list
Rcpp-devel@lists.r-forge.r-project.org
https://lists.r-forge.r-project.org/cgi-bin/mailman/listinfo/rcpp-devel
Re: [Rcpp-devel] tracking volatile bug
On 21 November 2014 at 23:50, Sokol Serguei wrote: | If I did not use the protection in my Rcpp code it's because | I did not see it in any example of quick ref guide: | http://cran.r-project.org/web/packages/Rcpp/vignettes/Rcpp-quickref.pdf | I didn't see any mention of PROTECT necessity neither in | http://dirk.eddelbuettel.com/code/rcpp/Rcpp-introduction.pdf | especially in the part dedicated to the use of wrap(). I am sorry but that is somewhat uncalled for. We rarely advocate a programming style advocating SEXP manipulation (apart from rare cases, or in the guts of our code). You programmed the plain C API for R, apart from using some utility functions from Rcpp to access functions. For this programming style, you should go to r-devel rather than coming here. What you do is pretty different from standard Rcpp use, or what we recommend. Blaming us for not documenting what we consider inappropriate style is not only circular but also continues to misunderstand why we do what we do here, and how we do it. | I supposed (very naively, it is clear now) that necessary wrapping and | protection was automagicaly added behind the stage by Rcpp. It of course happes automagically if and only if you use Rcpp datatypes. From what I understand you did not use those, and hence got no magic. Dirk -- http://dirk.eddelbuettel.com | @eddelbuettel | e...@debian.org ___ Rcpp-devel mailing list Rcpp-devel@lists.r-forge.r-project.org https://lists.r-forge.r-project.org/cgi-bin/mailman/listinfo/rcpp-devel
Re: [Rcpp-devel] tracking volatile bug
Dirk Eddelbuettel has written at Sat, 22 Nov 2014 09:57:36 -0600
On 21 November 2014 at 23:50, Sokol Serguei wrote:
| If I did not use the protection in my Rcpp code it's because
| I did not see it in any example of quick ref guide:
| http://cran.r-project.org/web/packages/Rcpp/vignettes/Rcpp-quickref.pdf
| I didn't see any mention of PROTECT necessity neither in
| http://dirk.eddelbuettel.com/code/rcpp/Rcpp-introduction.pdf
| especially in the part dedicated to the use of wrap().
I am sorry but that is somewhat uncalled for.
We rarely advocate a programming style advocating SEXP manipulation (apart
from rare cases, or in the guts of our code). You programmed the plain C API
for R, apart from using some utility functions from Rcpp to access functions.
For this programming style, you should go to r-devel rather than coming here.
What you do is pretty different from standard Rcpp use, or what we recommend.
If I moved to this style it's because an Rcpp style (at least how I
understood it) did not work for me.
Obviously, if I can get a pure Rcpp code working in this case I would
be perfectly happy
with that.
Blaming us
It was not a blaming, just citing my sources.
for not documenting what we consider inappropriate style is not
only circular but also continues to misunderstand why we do what we do here,
and how we do it.
May be I still misunderstand something but from what I understood
it is perfectly regular to call an R function from rcpp. Well, you say
it is not
recommended but if there is no rcpp counterpart neither CCallable interface
what are the alternatives?
| I supposed (very naively, it is clear now) that necessary wrapping and
| protection was automagicaly added behind the stage by Rcpp.
It of course happes automagically if and only if you use Rcpp datatypes.
From what I understand you did not use those, and hence got no magic.
Let try to stick with regular rcpp code
(file: matrix_norm.cpp):
//[[Rcpp::depends(RcppArmadillo)]]
#include RcppArmadillo.h
using namespace Rcpp;
using namespace arma;
// [[Rcpp::export]]
double nmat(mat A) {
Function Matrix_norm_r_=Environment(package:Matrix)[norm];
double res=asdouble(Matrix_norm_r_(A, 1));
return res;
}
When called in R as:
library(Rcpp)
library(Matrix)
sourceCpp(matrix_norm.cpp)
gctorture(TRUE)
nmat(as.matrix(pi))
it gives an error:
Erreur : 'getCharCE' doit être appelé sur un CHARSXP
(my English translation: Error: 'getCharCE' must be called on CHARSXP)
Something was irregular on my side here?
Serguei.
___
Rcpp-devel mailing list
Rcpp-devel@lists.r-forge.r-project.org
https://lists.r-forge.r-project.org/cgi-bin/mailman/listinfo/rcpp-devel
Re: [Rcpp-devel] tracking volatile bug
On Sat, Nov 22, 2014 at 1:43 PM, Sokol Serguei serguei.so...@gmail.com wrote:
Let try to stick with regular rcpp code
(file: matrix_norm.cpp):
//[[Rcpp::depends(RcppArmadillo)]]
#include RcppArmadillo.h
using namespace Rcpp;
using namespace arma;
// [[Rcpp::export]]
double nmat(mat A) {
Function Matrix_norm_r_=Environment(package:Matrix)[norm];
double res=asdouble(Matrix_norm_r_(A, 1));
return res;
}
When called in R as:
library(Rcpp)
library(Matrix)
sourceCpp(matrix_norm.cpp)
gctorture(TRUE)
nmat(as.matrix(pi))
it gives an error:
Erreur : 'getCharCE' doit être appelé sur un CHARSXP
(my English translation: Error: 'getCharCE' must be called on CHARSXP)
Something was irregular on my side here?
Serguei.
Try replacing the line that sets res with:
double res=asdouble(Matrix_norm_r_(NumericMatrix(wrap(A)),
CharacterVector::create(1)));
Be sure to try it on a fresh session since errors from prior runs
could have messed up R.
--
Statistics Software Consulting
GKX Group, GKX Associates Inc.
tel: 1-877-GKX-GROUP
email: ggrothendieck at gmail.com
___
Rcpp-devel mailing list
Rcpp-devel@lists.r-forge.r-project.org
https://lists.r-forge.r-project.org/cgi-bin/mailman/listinfo/rcpp-devel
Re: [Rcpp-devel] tracking volatile bug
Well, I could narrow the issue but not yet resolve it.
Le 18/11/2014 20:46, William Dunlap a écrit :
I will try gctorture(TRUE) suggested by Martin.
I'll start with this easy part. Unfortunately valgrind
didn't detect any wrong access to memory.
Now, the difficult part.
The most reduced code in cpp producing an error under
gctorture(TRUE) is the following:
8file matrix_norm.cpp
//[[Rcpp::depends(RcppArmadillo)]]
#include RcppArmadillo.h
using namespace Rcpp;
using namespace arma;
// [[Rcpp::export]]
SEXP nmat(mat A) {
Function Matrix_norm_r_=Environment(package:Matrix)[norm];
SEXP s=wrap(1);
Rf_PrintValue(s);
SEXP res=Matrix_norm_r_(wrap(A), s);
return res;
}
8
It is compiled with CXXFLAGS=-g -O0 -std=c++11 in ~/.R/Makevars
The R-3.1.2 itself is also compiled with -g -O0 -std=c++11 by
gcc version 4.8.2
How to reproduce ?
$ R -d gdb
(gdb) run
library(Rcpp)
library(Matrix)
sourceCpp(matrix_norm.cpp)
[Ctrl-C]
(gdb) b nmat
Breakpoint 1 at 0x7079c253: file matrix_norm.cpp, line 8.
(gdb) c
Continuing.
nmat(as.matrix(pi))
Breakpoint 1, nmat (A=...) at matrix_norm.cpp:8
8 Function Matrix_norm_r_=Environment(package:Matrix)[norm];
(gdb) c
Continuing.
[1] 1
[1] 3.141593
This is a correct result. Now, push it to an error with gctorture
gctorture(TRUE)
nmat(as.matrix(pi))
Breakpoint 1, nmat (A=...) at matrix_norm.cpp:8
8 Function Matrix_norm_r_=Environment(package:Matrix)[norm];
(gdb) c
Continuing.
CHARSXP: NA
Erreur : erreur d'évaluation de l'argument 'type' lors de la sélection d'une
méthode pour la fonction 'norm' : Erreur : type 'char' indisponible dans 'eval'
The wrong behavior is manifested in printing 'CHARSXP: NA'
It means that the string 1 was not wrapped correctly. The error messages
(here in french but never mind) are just a consequence of the wrong string
parameter passed to a call Matrix::norm(). I must say that in a function
without evoking the Matrix environment but doing just a plain
wrap(1), the wrap() is working as expected.
At one moment I suspected that the problem was in Rf_mkString() where there
is a call to unprotected mkChar(), so in R-3.1.2/src/include/Rinlinedfuns.h:682
I replaced
SET_STRING_ELT(t, (R_xlen_t)0, mkChar(s));
by
SET_STRING_ELT(t, (R_xlen_t)0, PROTECT(mkChar(s)));
and consequent UNPROTECT(1) by UNPROTECT(2)
But ... adding this protection did not resolve the issue.
nmat(as.matrix(pi))
Breakpoint 1, nmat (A=...) at matrix_norm.cpp:8
8 Function Matrix_norm_r_=Environment(package:Matrix)[norm];
(gdb) n
9 SEXP s=wrap(1);
(gdb) s
Rcpp::wrap (v=0x707a1a5a 1)
at /home/local/src/R-3.1.2/library/Rcpp/include/Rcpp/internal/wrap.h:930
930 if (v != NULL)
(gdb)
931 return Rf_mkString(v) ;
(gdb) s
Rf_mkString (s=0x707a1a5a 1) at ../../src/include/Rinlinedfuns.h:681
681 PROTECT(t = allocVector(STRSXP, (R_xlen_t)1));
(gdb) n
682 SET_STRING_ELT(t, (R_xlen_t)0, PROTECT(mkChar(s)));
(gdb)
683 UNPROTECT(2);
(gdb)
684 return t;
(gdb) call Rf_PrintValue(t)
CHARSXP: NA
instead of a correct 1.
My problem to locate the place where a pointer corruption occurs is that
if I dig deep enough to see what is going on on the stack of protected values
and during memory allocation, it ... produce a good result, i.e. a string
vector with a single element 1.
nmat(as.matrix(pi))
Breakpoint 1, nmat (A=...) at matrix_norm.cpp:8
8 Function Matrix_norm_r_=Environment(package:Matrix)[norm];
(gdb) n
9 SEXP s=wrap(1);
(gdb) s
Rcpp::wrap (v=0x707a1a5a 1)
at /home/local/src/R-3.1.2/library/Rcpp/include/Rcpp/internal/wrap.h:930
930 if (v != NULL)
(gdb)
931 return Rf_mkString(v) ;
(gdb)
Rf_mkString (s=0x707a1a5a 1) at ../../src/include/Rinlinedfuns.h:681
681 PROTECT(t = allocVector(STRSXP, (R_xlen_t)1));
(gdb) n
682 SET_STRING_ELT(t, (R_xlen_t)0, PROTECT(mkChar(s)));
(gdb) s
Rf_mkChar (name=0x707a1a5a 1) at envir.c:3444
3444size_t len = strlen(name);
(gdb) p name
$8 = 0x707a1a5a 1
(gdb) finish
Run till exit from #0 Rf_mkChar (name=0x707a1a5a 1) at envir.c:3444
0x778ef4f6 in Rf_mkString (s=0x707a1a5a 1) at
../../src/include/Rinlinedfuns.h:682
682 SET_STRING_ELT(t, (R_xlen_t)0, PROTECT(mkChar(s)));
Value returned is $9 = (struct SEXPREC *) 0x48bf4c8
(gdb) call Rf_PrintValue(0x48bf4c8)
CHARSXP: 1
(gdb) n
683 UNPROTECT(2);
(gdb)
684 return t;
(gdb) call Rf_PrintValue(t)
[1] 1
(gdb) c
Continuing.
[1] 1
[1] 3.141593
But immediately after:
nmat(as.matrix(pi))
Breakpoint 1, nmat (A=...) at matrix_norm.cpp:8
8 Function Matrix_norm_r_=Environment(package:Matrix)[norm];
(gdb) c
Continuing.
CHARSXP: NA
Erreur : erreur d'évaluation de l'argument 'type' lors de la sélection d'une
méthode pour la fonction 'norm' : Erreur : type 'char' indisponible dans 'eval'
Any thoughts how to locate this problem with pointer
Re: [Rcpp-devel] tracking volatile bug
Thanks Martin,
Your PROTECT did the job :)
Martin Morgan has written at Fri, 21 Nov 2014 13:33:47 -0800
On 11/21/2014 09:23 AM, Serguei Sokol wrote:
Well, I could narrow the issue but not yet resolve it.
Le 18/11/2014 20:46, William Dunlap a écrit :
I will try gctorture(TRUE) suggested by Martin.
I'll start with this easy part. Unfortunately valgrind
didn't detect any wrong access to memory.
Now, the difficult part.
The most reduced code in cpp producing an error under
gctorture(TRUE) is the following:
8file matrix_norm.cpp
//[[Rcpp::depends(RcppArmadillo)]]
#include RcppArmadillo.h
using namespace Rcpp;
using namespace arma;
// [[Rcpp::export]]
SEXP nmat(mat A) {
Function Matrix_norm_r_=Environment(package:Matrix)[norm];
SEXP s=wrap(1);
the reduced example is great. s is not PROTECT()ed, Rf_PrintValue(s)
almost certainly decides that it needs memory, and re-uses the SEXP
occupied by s for its own internal purposes. So simply
SEXP s = PROTECT(wrap(1));
// ...
UNPROTECT(1)
return res;
If I did not use the protection in my Rcpp code it's because
I did not see it in any example of quick ref guide:
http://cran.r-project.org/web/packages/Rcpp/vignettes/Rcpp-quickref.pdf
I didn't see any mention of PROTECT necessity neither in
http://dirk.eddelbuettel.com/code/rcpp/Rcpp-introduction.pdf
especially in the part dedicated to the use of wrap().
You _could_ figure out where Rf_PrintValue over-writes the unprotected
s SEXP,
It turns out that the memory is getting corrupted well before
Rf_PrintValue. It
was here just to witness the devastation. The corruption occurred inside
Rf_mkString()
as shown in the debug session.
but that wouldn't be helpful in any way; R in general and
Rf_PrintValue in particular is allowed to write to unprotected SEXPs.
It's also not relevant in some ways what the error is -- invalid
'type' or anything else -- just that the error occurs. Also, while not
strictly true, it is much more likely that protection errors are in
user code (as here), then in package code, and finally in R code itself.
For what it's worth...
Rf_PrintValue(s);
SEXP res=Matrix_norm_r_(wrap(A), s);
I think wrap(A) needs (in principle) to be PROTECT()ed, too.
For the sake of safety, I'll do it.
Sometimes in simplifying the problem we introduce new errors; I think
your wrap()'s in the original code were immediately returned to the
user, so not in need of protection.
In the original code, there were not wraping at all for R function call.
I supposed (very naively, it is clear now) that necessary wrapping and
protection was automagicaly added behind the stage by Rcpp.
I would set a breakpoint on Rf_error, then trigger the error under
gctorture(TRUE), then in gdb backtrace to localize which specific line
in your code triggers the error, and simplify around that.
That was I have done to find out that Rf_mkString() call was the last part
of Rcpp where corruption occurred.
Thank again.
Serguei.
___
Rcpp-devel mailing list
Rcpp-devel@lists.r-forge.r-project.org
https://lists.r-forge.r-project.org/cgi-bin/mailman/listinfo/rcpp-devel
[Rcpp-devel] tracking volatile bug
Hi everybody,
I am facing a volatile bug which appears and disappears in
identical calls on identical data.
My RcppArmadilla code is expm_cpp() (it is obtained by rex2arma tool
which I have presented before, cf. attached file). It is compared and
benchmarked vs its R counterpart expm.higham() (also attached). Both
are computing matrix exponential by a same particular method.
The results are identical for R and Rcpp on small and medium sized
matrices nxn n=10:100 but on large matrices (n 800) various error messages
can interrupt (or not) a normal run of expm_cpp().
Sometimes message says (in French) type 'char' indisponible dans 'eval',
I suppose in English it must be unimplemented type 'char' in 'eval'
I have seen (thanks to Google) a similar error message in the following snippet:
/* call this R command: source(FileName) */
int errorOccurred;
SEXP e = lang2(install(source), mkString(FileName));
/* mkChar instead of mkString would lead to this runtime error:
* Error in source(FileName) : unimplemented type 'char' in 'eval' */
R_tryEval(e, R_GlobalEnv, errorOccurred);
which suggests that somewhere in Rcpp or RcppArmadillo there is
a mkChar() call instead of mkString().
Other times, error message can say something like
argument type[1]='x' must be one of 'M','1','O','I','F' or 'E'
or argument type[1]='character' must be a one-letter character string
This latter message is somewhat volatile per se. The part of message
just after type[1]= can be 'character' (as above) or 'method' or 'ANY' etc.
I have found these messages in the Matrix package
https://r-forge.r-project.org/scm/viewvc.php/pkg/Matrix/src/Mutils.c?view=markuproot=matrixpathrev=2614
function char La_norm_type(const char *typstr)
Seemingly, the argument *typstr is getting corrupted somewhere
on the road.
It is useless to say that debugging is of no help here.
If run in a debugger, the program stops normally with
or without error messages cited above.
I have also tried to low the level of gcc optimization
both in compiling R and the Rcpp code but it didn't help.
Anybody has an experience in tracking down similar cases?
This example can be run as:
library(expm)
library(Rcpp)
source(expm_cpp.R)
source(expm.higham.R)
n=1000
As=matrix(rnorm(n*n), n)
stopifnot(diff(range(expm.higham(As, balancing=TRUE)-expm_cpp(As,
balancing=TRUE))) 1.e-14)
The last command may be run several times before an error shows up.
Cheers,
Serguei.
cppFunction(depends='RcppArmadillo', rebuild=TRUE, includes='
template typename T
inline unsigned which_max(T v) {
unsigned i;
v.max(i);
return i+1;
}
templatetypename T
inline unsigned which_min(T v) {
unsigned i;
v.min(i);
return i+1;
}
',
using namespace arma;
using namespace Rcpp;
SEXP expm_cpp(
NumericMatrix A_in_,
bool balancing) {
// auxiliary functions
Environment base_env_r_=Environment::base_env();
Function rep_r_=base_env_r_[\rep\];
Function c_r_=base_env_r_[\c\];
// External R function declarations
Function expm_balance_r_=Environment(\package:expm\)[\balance\];
Function Matrix_norm_r_=Environment(\package:Matrix\)[\norm\];
// Input variable declarations and conversion
mat A(A_in_.begin(), A_in_.nrow(), A_in_.ncol(), false);
// Output and intermediate variable declarations
mat A2;
mat B;
mat B2;
mat B4;
mat B6;
List baP;
List baS;
mat C;
vec c_;
ivec d;
vec dd;
int i;
mat I;
int k;
int l;
int n;
double nA;
mat P;
ivec pp;
double s;
vec t;
mat tt;
mat U;
mat V;
mat X;
// Translated code starts here
d=ivec(IntegerVector::create(A.n_rows, A.n_cols));
if (d.size() != 2 || d.at(0) != d.at(1)) stop(\'A' must be a square
matrix\);
n=d.at(0);
if (n = 1) return wrap(exp(A));
if (balancing) {
baP=asList(expm_balance_r_(A, \P\));
baS=asList(expm_balance_r_(asmat(baP[\z\]), \S\));
A=asmat(baS[\z\]);
}
nA=asdouble(Matrix_norm_r_(A, \1\));
I=eyemat(n, n);
if (nA = 2.1) {
t=vec({0.015, 0.25, 0.95, 2.1});
l=which_max(nA = t);
C=join_vert(join_vert(join_vert(vec({120, 60, 12, 1, 0, 0, 0, 0, 0,
0}).st(), vec({30240, 15120, 3360, 420, 30, 1, 0, 0, 0, 0}).st()),
vec({17297280, 8648640, 1995840, 277200, 25200, 1512, 56, 1, 0, 0}).st()),
vec({17643225600, 8821612800, 2075673600, 302702400, 30270240, 2162160, 110880,
3960, 90, 1}).st());
A2=A*A;
P=I;
U=C.at((l)-1, 1) * I;
V=C.at((l)-1, 0) * I;
for (int incr_arma_=(1 = l ? 1 : -1), k=1; k != l+incr_arma_;
k+=incr_arma_) {
P=P*A2;
U=U + C((l)-1, ((2 * k) + 2)-1) * P;
V=V + C((l)-1, ((2 * k) + 1)-1) * P;
}
U=A*U;
X=solve(V - U, V + U);
} else {
s=log2(nA / 5.4);
B=A;
if (s 0) {
s=ceil(s);
B=B / (pow(2, s));
}
c_=vec({6476475253248, 3238237626624, 7771770303897600,
1187353796428800, 129060195264000, 10559470521600,
Re: [Rcpp-devel] tracking volatile bug
Serguei, First, I think we mentioned to you after your first post on rex2arma that all serious Rcpp development is generally done in packages. You should probably embrace that mode too. Many of us prefer it. Second, 'expm' does not need to be reimplemented. See my RcppKalman package (on GitHub only) to see how to access the expm code (from the package of the same name by Goulet, Dutang,Maechler et all -- all serious people) directly from RcppArmadillo. Third, your code calls back into R via Rcpp::Environment. As a first test I'd try to do without. The recommendation usually is minimally self-contained reproducible examples for bug report. I, like many other people, have found that a good number of times the bug becomes apparent while preparing such a report. Maybe someone will have time to go through your somewhat-non-standard in detail. Cheers, Dirk -- http://dirk.eddelbuettel.com | @eddelbuettel | e...@debian.org ___ Rcpp-devel mailing list Rcpp-devel@lists.r-forge.r-project.org https://lists.r-forge.r-project.org/cgi-bin/mailman/listinfo/rcpp-devel
Re: [Rcpp-devel] tracking volatile bug
Le 18 nov. 2014 à 14:42, Serguei Sokol serguei.so...@gmail.com a écrit : Hi everybody, I am facing a volatile bug which appears and disappears in identical calls on identical data. My RcppArmadilla code is expm_cpp() (it is obtained by rex2arma tool which I have presented before, cf. attached file). It is compared and benchmarked vs its R counterpart expm.higham() (also attached). Both are computing matrix exponential by a same particular method. The results are identical for R and Rcpp on small and medium sized matrices nxn n=10:100 but on large matrices (n 800) various error messages can interrupt (or not) a normal run of expm_cpp(). Sometimes message says (in French) type 'char' indisponible dans 'eval', I suppose in English it must be unimplemented type 'char' in 'eval' I have seen (thanks to Google) a similar error message in the following snippet: /* call this R command: source(FileName) */ int errorOccurred; SEXP e = lang2(install(source), mkString(FileName)); /* mkChar instead of mkString would lead to this runtime error: * Error in source(FileName) : unimplemented type 'char' in 'eval' */ R_tryEval(e, R_GlobalEnv, errorOccurred); e is not protected here. Does the problem go away if you protect it: SEXP e = PROTECT( lang2(install(source), mkString(FileName))) ; Or more R/C++ idiomatically, which is both nicer and safer. Language e( source, FileName ) ; Romain which suggests that somewhere in Rcpp or RcppArmadillo there is a mkChar() call instead of mkString(). Other times, error message can say something like argument type[1]='x' must be one of 'M','1','O','I','F' or 'E' or argument type[1]='character' must be a one-letter character string This latter message is somewhat volatile per se. The part of message just after type[1]= can be 'character' (as above) or 'method' or 'ANY' etc. That kind of problem usually is related to premature GC. I have found these messages in the Matrix package https://r-forge.r-project.org/scm/viewvc.php/pkg/Matrix/src/Mutils.c?view=markuproot=matrixpathrev=2614 function char La_norm_type(const char *typstr) Seemingly, the argument *typstr is getting corrupted somewhere on the road. It is useless to say that debugging is of no help here. If run in a debugger, the program stops normally with or without error messages cited above. I have also tried to low the level of gcc optimization both in compiling R and the Rcpp code but it didn't help. Anybody has an experience in tracking down similar cases? This example can be run as: library(expm) library(Rcpp) source(expm_cpp.R) source(expm.higham.R) n=1000 As=matrix(rnorm(n*n), n) stopifnot(diff(range(expm.higham(As, balancing=TRUE)-expm_cpp(As, balancing=TRUE))) 1.e-14) The last command may be run several times before an error shows up. Cheers, Serguei. expm_cpp.Rexpm.higham.R___ Rcpp-devel mailing list Rcpp-devel@lists.r-forge.r-project.org https://lists.r-forge.r-project.org/cgi-bin/mailman/listinfo/rcpp-devel ___ Rcpp-devel mailing list Rcpp-devel@lists.r-forge.r-project.org https://lists.r-forge.r-project.org/cgi-bin/mailman/listinfo/rcpp-devel
Re: [Rcpp-devel] tracking volatile bug
Le 18/11/2014 16:39, Romain François a écrit : Le 18 nov. 2014 à 14:42, Serguei Sokol serguei.so...@gmail.com a écrit : /* call this R command: source(FileName) */ int errorOccurred; SEXP e = lang2(install(source), mkString(FileName)); /* mkChar instead of mkString would lead to this runtime error: * Error in source(FileName) : unimplemented type 'char' in 'eval' */ R_tryEval(e, R_GlobalEnv, errorOccurred); e is not protected here. Does the problem go away if you protect it: SEXP e = PROTECT( lang2(install(source), mkString(FileName))) ; you mean SEXP e = PROTECT(lang2(install(source), mkChar(expm.higham.R))); (because mkString() had no problem)? No, it does not, the problem is still there. Or more R/C++ idiomatically, which is both nicer and safer. Language e( source, FileName ) ; This works as expected. We can suppose that behind the stage FileName is wrap()-ed here with mkString() and not with mkChar(). That's why it's working. which suggests that somewhere in Rcpp or RcppArmadillo there is a mkChar() call instead of mkString(). Other times, error message can say something like argument type[1]='x' must be one of 'M','1','O','I','F' or 'E' or argument type[1]='character' must be a one-letter character string This latter message is somewhat volatile per se. The part of message just after type[1]= can be 'character' (as above) or 'method' or 'ANY' etc. That kind of problem usually is related to premature GC. I will try gctorture(TRUE) suggested by Martin. Thanks, Serguei. ___ Rcpp-devel mailing list Rcpp-devel@lists.r-forge.r-project.org https://lists.r-forge.r-project.org/cgi-bin/mailman/listinfo/rcpp-devel
Re: [Rcpp-devel] tracking volatile bug
I will try gctorture(TRUE) suggested by Martin. Try running it under valgrind as well. E.g., % R --debugger=valgrind ==15338== Memcheck, a memory error detector ==15338== Copyright (C) 2002-2011, and GNU GPL'd, by Julian Seward et al. ==15338== Using Valgrind-3.7.0 and LibVEX; rerun with -h for copyright info ==15338== Command: /opt/sw/R/R-3.1.1.atlas1/lib/R/bin/exec/R ==15338== ==15338== Conditional jump or move depends on uninitialised value(s) ==15338==at 0x40188B6: index (strchr.S:56) ==15338==by 0x4007C92: expand_dynamic_string_token (dl-load.c:431) ==15338==by 0x40085EE: _dl_map_object (dl-load.c:2501) ==15338==by 0x400188D: map_doit (rtld.c:633) ==15338==by 0x400F175: _dl_catch_error (dl-error.c:178) ... gctorture(TRUE) # call your test code here Bill Dunlap TIBCO Software wdunlap tibco.com On Tue, Nov 18, 2014 at 10:10 AM, Serguei Sokol serguei.so...@gmail.com wrote: Le 18/11/2014 16:39, Romain François a écrit : Le 18 nov. 2014 à 14:42, Serguei Sokol serguei.so...@gmail.com a écrit : /* call this R command: source(FileName) */ int errorOccurred; SEXP e = lang2(install(source), mkString(FileName)); /* mkChar instead of mkString would lead to this runtime error: * Error in source(FileName) : unimplemented type 'char' in 'eval' */ R_tryEval(e, R_GlobalEnv, errorOccurred); e is not protected here. Does the problem go away if you protect it: SEXP e = PROTECT( lang2(install(source), mkString(FileName))) ; you mean SEXP e = PROTECT(lang2(install(source), mkChar(expm.higham.R))); (because mkString() had no problem)? No, it does not, the problem is still there. Or more R/C++ idiomatically, which is both nicer and safer. Language e( source, FileName ) ; This works as expected. We can suppose that behind the stage FileName is wrap()-ed here with mkString() and not with mkChar(). That's why it's working. which suggests that somewhere in Rcpp or RcppArmadillo there is a mkChar() call instead of mkString(). Other times, error message can say something like argument type[1]='x' must be one of 'M','1','O','I','F' or 'E' or argument type[1]='character' must be a one-letter character string This latter message is somewhat volatile per se. The part of message just after type[1]= can be 'character' (as above) or 'method' or 'ANY' etc. That kind of problem usually is related to premature GC. I will try gctorture(TRUE) suggested by Martin. Thanks, Serguei. ___ Rcpp-devel mailing list Rcpp-devel@lists.r-forge.r-project.org https://lists.r-forge.r-project.org/cgi-bin/mailman/listinfo/rcpp-devel ___ Rcpp-devel mailing list Rcpp-devel@lists.r-forge.r-project.org https://lists.r-forge.r-project.org/cgi-bin/mailman/listinfo/rcpp-devel
