Skip to content

Instantly share code, notes, and snippets.

@mclements
Last active February 21, 2023 15:15
Show Gist options
  • Save mclements/68aa486f6ef4cfc9634076dad357e483 to your computer and use it in GitHub Desktop.
Save mclements/68aa486f6ef4cfc9634076dad357e483 to your computer and use it in GitHub Desktop.
Error with linking Mercury with a static file
minimal_fails: minimal_fails.m
mmc --make -v -l:libRmath.a -L `pkg-config --variable=libdir libRmath` minimal_fails && ./minimal_fails
minimal_ok: minimal_ok.m
mmc --make -v -l:libRmath.a -L `pkg-config --variable=libdir libRmath` minimal_ok && ./minimal_ok
minimal_fixed: minimal_fixed.m
mmc --make -v -l:libRmath.a -L `pkg-config --variable=libdir libRmath` minimal_fixed && ./minimal_fixed
:- module minimal_fails.
:- interface.
:- import_module float, io.
:- func dwilcox(float,float,float,int) = float.
:- func pnorm(float,float,float,int,int) = float.
:- pred main(io::di, io::uo) is det.
:- implementation.
:- pragma foreign_proc("C",
dwilcox(P::in,M::in,N::in,Log::in) = (Dwilcox_::out),
[promise_pure, will_not_call_mercury],
"Dwilcox_ = dwilcox(P,M,N,Log);").
:- pragma foreign_proc("C",
pnorm(Q::in,Mean::in,Sd::in,Lower::in,Log::in) = (Pnorm_::out),
[promise_pure, will_not_call_mercury],
"Pnorm_ = pnorm(Q,Mean,Sd,Lower,Log);").
:- pragma foreign_decl("C", "
#define MATHLIB_STANDALONE
#include \"Rmath.h\"
typedef unsigned int Int32;
").
:- pragma foreign_code("C",
"void set_seed(Int32 inseed, Int32 ignored) { }
double unif_rand(void) { return 0.5; }
").
main(!IO) :-
P = pnorm(1.96, 0.0, 1.0, 1, 0),
print_line(P, !IO).
:- module minimal_fixed.
:- interface.
:- import_module float, io.
:- func dwilcox(float,float,float,int) = float.
:- func pnorm(float,float,float,int,int) = float.
:- pred main(io::di, io::uo) is det.
:- implementation.
:- pragma foreign_proc("C",
dwilcox(P::in,M::in,N::in,Log::in) = (Dwilcox_::out),
[promise_pure, will_not_call_mercury],
"Dwilcox_ = dwilcox(P,M,N,Log);").
:- pragma foreign_proc("C",
pnorm(Q::in,Mean::in,Sd::in,Lower::in,Log::in) = (Pnorm_::out),
[promise_pure, will_not_call_mercury],
"Pnorm_ = pnorm(Q,Mean,Sd,Lower,Log);").
:- pragma foreign_decl("C", "
#define MATHLIB_STANDALONE
#include \"Rmath.h\"
#include <math.h>
#include <stdint.h>
typedef unsigned int Int32;
").
:- pragma foreign_code("C",
"void set_seed(Int32 inseed, Int32 ignored) { }
double unif_rand(void) { return 0.5; }
//copied from src/nmath/standalone/sunif.c:
//generate a random non-negative integer < 2 ^ bits in 16 bit chunks
static double rbits(int bits)
{
int_least64_t v = 0;
for (int n = 0; n <= bits; n += 16) {
int v1 = (int) floor(unif_rand() * 65536);
v = 65536 * v + v1;
}
// mask out the bits in the result that are not needed
return (double) (v & ((1L << bits) - 1));
}
double R_unif_index(double dn)
{
// rejection sampling from integers below the next larger power of two
if (dn <= 0)
return 0.0;
int bits = (int) ceil(log2(dn));
double dv;
do { dv = rbits(bits); } while (dn <= dv);
return dv;
}
").
main(!IO) :-
P = pnorm(1.96, 0.0, 1.0, 1, 0),
print_line(P, !IO).
:- module minimal_ok.
:- interface.
:- import_module float, io.
:- func pnorm(float,float,float,int,int) = float.
:- pred main(io::di, io::uo) is det.
:- implementation.
:- func dwilcox(float,float,float,int) = float. %% THIS HAS MOVED FROM THE INTERFACE
:- pragma foreign_proc("C",
dwilcox(P::in,M::in,N::in,Log::in) = (Dwilcox_::out),
[promise_pure, will_not_call_mercury],
"Dwilcox_ = dwilcox(P,M,N,Log);").
:- pragma foreign_proc("C",
pnorm(Q::in,Mean::in,Sd::in,Lower::in,Log::in) = (Pnorm_::out),
[promise_pure, will_not_call_mercury],
"Pnorm_ = pnorm(Q,Mean,Sd,Lower,Log);").
:- pragma foreign_decl("C", "
#define MATHLIB_STANDALONE
#include \"Rmath.h\"
typedef unsigned int Int32;
").
:- pragma foreign_code("C",
"void set_seed(Int32 inseed, Int32 ignored) { }
double unif_rand(void) { return 0.5; }
").
main(!IO) :-
P = pnorm(1.96, 0.0, 1.0, 1, 0), %% expects approximately 1-0.05/2.0=0.975
print_line(P, !IO).
@mclements
Copy link
Author

For the rmath library in Mercury (https://github.com/mclements/mercury-rmath), I am trying to use C code for the Mersenne Twister random number generator to replace the default random number generator. This involves replacing two C functions: double unif_rand(void) and void set_seed(unsigned int, unsigned int); for details, see https://rstudio.github.io/r-manuals/r-admin/The-standalone-Rmath-library.html.

I have got this working -- except if any of four functions (dwilcox, pwilcox, qwilcox, rwilcox) are included in the interface section. I then get a linking error (make minimal_fails):

/usr/bin/ld: /usr/lib/libRmath.a(std_unif.o): in function `set_seed':
(.text+0x0): multiple definition of `set_seed'; Mercury/os/minimal_fails.o:minimal_fails.c:(.text+0xb0): first defined here
/usr/bin/ld: /usr/lib/libRmath.a(std_unif.o): in function `unif_rand':
(.text+0x40): multiple definition of `unif_rand'; Mercury/os/minimal_fails.o:minimal_fails.c:(.text+0xc0): first defined here

If I move those functions to the implementation section, then the linking is okay (make minimal_ok).

Can anyone suggest why this fails and then works?

@mclements
Copy link
Author

mclements commented Feb 13, 2023

I have now added a call to the library using pnorm. This has the same behaviour as previously:(.

For make minimal_ok, the output is 0.9750021048517796.

@mclements
Copy link
Author

I have now added minimal_fixed.m, which fixes the issue with [dpqr]wilcox. In outline, rwilcox() calls R_rand_index() -- which needs to be re-defined for a user-defined random number generator. I have passed this information back to the R-devel mailing list.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment