Skip to content

Commit

Permalink
Merge pull request #3 from nlmixr2/3-strict
Browse files Browse the repository at this point in the history
STRICT_R_HEADERS
  • Loading branch information
mattfidler authored Aug 25, 2024
2 parents 2b52707 + 481d51d commit a23f045
Show file tree
Hide file tree
Showing 7 changed files with 57 additions and 45 deletions.
2 changes: 2 additions & 0 deletions NEWS.md
Original file line number Diff line number Diff line change
Expand Up @@ -5,6 +5,8 @@
the `rxode2` internals use the function pointers, while the models
built by `rxode2` use the binary linkage.

* As requested by `CRAN` use strict headers.

# PreciseSums 0.6

* Updated linkage for C.
Expand Down
3 changes: 2 additions & 1 deletion src/Makevars
Original file line number Diff line number Diff line change
@@ -1,2 +1,3 @@
PKG_CFLAGS = -DSTRICT_R_HEADERS=1
SOURCES_C = sum.c prod.c init.c
OBJECTS = $(SOURCES_C:.c=.o)
OBJECTS = $(SOURCES_C:.c=.o)
1 change: 1 addition & 0 deletions src/PreciseSums.h
Original file line number Diff line number Diff line change
Expand Up @@ -11,6 +11,7 @@
#define PS_Prod 2
#define PS_LogifyProd 3

#define R_NO_REMAP
#include <R.h>
#include <Rinternals.h>
#include <Rdefines.h>
Expand Down
1 change: 1 addition & 0 deletions src/PreciseSumsPtr.h
Original file line number Diff line number Diff line change
Expand Up @@ -11,6 +11,7 @@
#define PS_Prod 2
#define PS_LogifyProd 3

#define R_NO_REMAP
#include <R.h>
#include <Rinternals.h>
#include <Rdefines.h>
Expand Down
19 changes: 10 additions & 9 deletions src/init.c
Original file line number Diff line number Diff line change
@@ -1,5 +1,6 @@
#include <stdio.h>
#include <stdlib.h>
#define R_NO_REMAP
#include <R.h>
#include <Rinternals.h>
#include <Rmath.h> //Rmath includes math.
Expand Down Expand Up @@ -40,23 +41,23 @@ SEXP _PreciseSumsPtr(void) {
SEXP PreciseSums_sum_get_ptr = PROTECT(R_MakeExternalPtrFn((DL_FUNC)&PreciseSums_sum_get, R_NilValue, R_NilValue)); pro++;
SEXP PreciseSums_prod_get_ptr = PROTECT(R_MakeExternalPtrFn((DL_FUNC)&PreciseSums_prod_get, R_NilValue, R_NilValue)); pro++;

SEXP ret = PROTECT(allocVector(VECSXP, 6)); pro++;
SEXP ret = PROTECT(Rf_allocVector(VECSXP, 6)); pro++;
SET_VECTOR_ELT(ret, 0, PreciseSums_sum_ptr);
SET_VECTOR_ELT(ret, 1, PreciseSums_prod_ptr);
SET_VECTOR_ELT(ret, 2, PreciseSums_sum_r_ptr);
SET_VECTOR_ELT(ret, 3, PreciseSums_prod_r_ptr);
SET_VECTOR_ELT(ret, 4, PreciseSums_sum_get_ptr);
SET_VECTOR_ELT(ret, 5, PreciseSums_prod_get_ptr);

SEXP retn = PROTECT(allocVector(STRSXP, 6)); pro++;
SET_STRING_ELT(retn, 0, mkChar("PreciseSums_sum"));
SET_STRING_ELT(retn, 1, mkChar("PreciseSums_prod"));
SET_STRING_ELT(retn, 2, mkChar("PreciseSums_sum_r"));
SET_STRING_ELT(retn, 3, mkChar("PreciseSums_prod_r"));
SET_STRING_ELT(retn, 4, mkChar("PreciseSums_sum_get"));
SET_STRING_ELT(retn, 5, mkChar("PreciseSums_prod_get"));
SEXP retn = PROTECT(Rf_allocVector(STRSXP, 6)); pro++;
SET_STRING_ELT(retn, 0, Rf_mkChar("PreciseSums_sum"));
SET_STRING_ELT(retn, 1, Rf_mkChar("PreciseSums_prod"));
SET_STRING_ELT(retn, 2, Rf_mkChar("PreciseSums_sum_r"));
SET_STRING_ELT(retn, 3, Rf_mkChar("PreciseSums_prod_r"));
SET_STRING_ELT(retn, 4, Rf_mkChar("PreciseSums_sum_get"));
SET_STRING_ELT(retn, 5, Rf_mkChar("PreciseSums_prod_get"));

setAttrib(ret, R_NamesSymbol, retn);
Rf_setAttrib(ret, R_NamesSymbol, retn);

UNPROTECT(pro);
return ret;
Expand Down
23 changes: 14 additions & 9 deletions src/prod.c
Original file line number Diff line number Diff line change
@@ -1,21 +1,26 @@
#include <stdio.h>
#include <stdlib.h>
#include <stdarg.h>
#define R_NO_REMAP
#include <R.h>
#include <Rinternals.h>
#include <Rmath.h> //Rmath includes math.

// https://github.com/wch/r-source/blob/41ca75a4ffa1a4ec66f7697820fcf229b6a64f8e/src/nmath/qgamma.c#L122
#define LN_EPS -36.043653389117156

extern double PreciseSums_sum (double *input, int n);

extern double PreciseSums_pairwise_add_DOUBLE(double *a, int n);

double PreciseSums_safe_log(double x){
if (x <= 0){
// Warning?
return log(DBL_EPSILON);
return LN_EPS;
} else {
return log(x);
}
return LN_EPS;
}

int PreciseSums_prod_type = 3;
Expand All @@ -38,7 +43,7 @@ extern double PreciseSums_prod_ld(double *input, int n){
long double p = 1;
for (int i = 0; i < n; i++){
if (input[i] == 0){
return 0.0;
return 0.0;
}
p *= input[i];
}
Expand All @@ -49,7 +54,7 @@ extern double PreciseSums_prod_d(double *input, int n){
double p = 1;
for (int i = 0; i < n; i++){
if (input[i] == 0){
return 0.0;
return 0.0;
}
p *= input[i];
}
Expand All @@ -70,9 +75,9 @@ extern double PreciseSums_prod_logify_r(double *input, double *p, int n){
}

extern double PreciseSums_prod_logify(double *input, int n){
double *p = Calloc(n,double);
double *p = R_Calloc(n,double);
double ret= PreciseSums_prod_logify_r(input, p, n);
Free(p);
R_Free(p);
return ret;
}

Expand Down Expand Up @@ -110,9 +115,9 @@ extern double PreciseSums_prod(double *input, int n){


SEXP _psProd(SEXP input){
int len = length(input);
int len = Rf_length(input);
double *dinput = REAL(input);
SEXP rets = PROTECT(allocVector(REALSXP,1));
SEXP rets = PROTECT(Rf_allocVector(REALSXP,1));
REAL(rets)[0] = PreciseSums_prod(dinput, len);
UNPROTECT(1);
return rets;
Expand All @@ -121,13 +126,13 @@ SEXP _psProd(SEXP input){
extern double PreciseSums_prodV(int n, ...){
va_list valist;
va_start(valist, n);
double *p = Calloc(n, double);
double *p = R_Calloc(n, double);
for (int i = 0; i < n; i++){
p[i] = va_arg(valist, double);
}
va_end(valist);
double s = PreciseSums_prod(p, n);
Free(p);
R_Free(p);
return s;
}

Expand Down
53 changes: 27 additions & 26 deletions src/sum.c
Original file line number Diff line number Diff line change
@@ -1,6 +1,7 @@
#include <stdio.h>
#include <stdlib.h>
#include <stdarg.h>
#define R_NO_REMAP
#include <R.h>
#include <Rinternals.h>
#include <Rmath.h> //Rmath includes math.
Expand Down Expand Up @@ -87,9 +88,9 @@ PreciseSums_pairwise_add_DOUBLE(double *a, int n)
}

SEXP _psPairwiseSum(SEXP input){
int len = length(input);
int len = Rf_length(input);
double *dinput = REAL(input);
SEXP rets = PROTECT(allocVector(REALSXP,1));
SEXP rets = PROTECT(Rf_allocVector(REALSXP,1));
REAL(rets)[0] = PreciseSums_pairwise_add_DOUBLE(dinput, len);
UNPROTECT(1);
return rets;
Expand All @@ -102,7 +103,7 @@ extern double PreciseSums_KahanSum(double *input, int len){
int i;

for (i = 0; i < len; i++){
y = (double)input[i] - c;
y = (double)input[i] - c;
t = sum + y;
c = (t - sum) - y; // (t - sum) cancels the high-order part of y; subtracting y recovers negative (low part of y)
sum = t; // Algebraically, c should always be zero. Beware overly-aggressive optimizing compilers!
Expand All @@ -111,9 +112,9 @@ extern double PreciseSums_KahanSum(double *input, int len){
}

SEXP _psKahanSum(SEXP input){
int len = length(input);
int len = Rf_length(input);
double *dinput = REAL(input);
SEXP rets = PROTECT(allocVector(REALSXP,1));
SEXP rets = PROTECT(Rf_allocVector(REALSXP,1));
REAL(rets)[0] = PreciseSums_KahanSum(dinput, len);
UNPROTECT(1);
return rets;
Expand Down Expand Up @@ -159,9 +160,9 @@ extern double PreciseSums_KleinSum(double *input, int len){
}

SEXP _psNeumaierSum(SEXP input){
int len = length(input);
int len = Rf_length(input);
double *dinput = REAL(input);
SEXP rets = PROTECT(allocVector(REALSXP,1));
SEXP rets = PROTECT(Rf_allocVector(REALSXP,1));
REAL(rets)[0] = PreciseSums_NeumaierSum(dinput, len);
UNPROTECT(1);
return rets;
Expand Down Expand Up @@ -194,16 +195,16 @@ extern double PreciseSums_Python_fsum_r(double *iterable, int iterable_len, doub
p[i++] = lo;
x = hi;
}
n = i;

n = i;
if (x != 0.0) {
if (!R_FINITE(x)) {
/* a nonfinite x could arise either as
a result of intermediate overflow, or
as a result of a nan or inf in the
summands */
if (R_FINITE(xsave) || ISNAN(xsave)) {
if (m > 0) Free(p);
if (m > 0) R_Free(p);
return PreciseSums_KleinSum(iterable, iterable_len);
/* error("intermediate overflow in fsum"); */
} else {
Expand All @@ -217,9 +218,9 @@ extern double PreciseSums_Python_fsum_r(double *iterable, int iterable_len, doub
//&& _fsum_realloc(&p, n, ps, &m)
// Doubles the size of array.
m += m;
p = Realloc(p, m, double);
p = R_Realloc(p, m, double);
} else if (m < 0 && n >= -m){
/* if (m > 0) Free(p); */
/* if (m > 0) R_Free(p); */
return PreciseSums_KleinSum(iterable, iterable_len);
/* error("The size of the saved partials is too small to calculate the sum."); */
}
Expand All @@ -229,8 +230,8 @@ extern double PreciseSums_Python_fsum_r(double *iterable, int iterable_len, doub
}
if (special_sum != 0.0) {
if (ISNAN(inf_sum)){
if (m > 0) Free(p);
error("-inf + inf in fsum");
if (m > 0) R_Free(p);
Rf_error("-inf + inf in fsum");
}
sum = special_sum;
return sum;
Expand All @@ -246,15 +247,15 @@ extern double PreciseSums_Python_fsum_r(double *iterable, int iterable_len, doub
x = hi;
y = p[--n];
if (fabs(y) >= fabs(x)){
if (m > 0) Free(p);
if (m > 0) R_Free(p);
return PreciseSums_KleinSum(iterable, iterable_len);
/* Rprintf("Partial Sums:\n"); */
/* for (i = 0; i < j; i++){ */
/* Rprintf("p[%d] = %f\n",i,p[i]); */
/* } */
/* Rprintf("Assertion Error:\n"); */
/* Rprintf("fabs(y) >= fabs(x) or %f >= %f\n",fabs(y),fabs(x)); */
/* if (m > 0) Free(p); */
/* if (m > 0) R_Free(p); */
/* error("Error in parital sums."); */
}
hi = x + y;
Expand Down Expand Up @@ -282,17 +283,17 @@ extern double PreciseSums_Python_fsum_r(double *iterable, int iterable_len, doub
}

extern double PreciseSums_Python_fsum(double *iterable, int iterable_len){
double *p = Calloc(NUM_PARTIALS, double);
double *p = R_Calloc(NUM_PARTIALS, double);
int m = NUM_PARTIALS;
double ret = PreciseSums_Python_fsum_r(&iterable[0], iterable_len, p, m);
Free(p);
R_Free(p);
return ret;
}

SEXP _psPythonSum(SEXP input){
int len = length(input);
int len = Rf_length(input);
double *dinput = REAL(input);
SEXP rets = PROTECT(allocVector(REALSXP,1));
SEXP rets = PROTECT(Rf_allocVector(REALSXP,1));
REAL(rets)[0] = PreciseSums_Python_fsum(dinput, len);
UNPROTECT(1);
return rets;
Expand Down Expand Up @@ -321,7 +322,7 @@ extern double PreciseSums_sum (double *input, int n){
}
//PreciseSums_KahanSum;
// PreciseSums_NeumaierSum
error("Unknown sum type.");
Rf_error("Unknown sum type.");
return 0;
}

Expand All @@ -347,7 +348,7 @@ extern double PreciseSums_sum_r(double *input, int n, double *p, int m, int type
}
//PreciseSums_KahanSum;
// PreciseSums_NeumaierSum
error("Unknown sum type.");
Rf_error("Unknown sum type.");
return 0;
}

Expand All @@ -366,9 +367,9 @@ SEXP _psSetSum(SEXP input){
}

SEXP _psSum(SEXP input){
int len = length(input);
int len = Rf_length(input);
double *dinput = REAL(input);
SEXP rets = PROTECT(allocVector(REALSXP,1));
SEXP rets = PROTECT(Rf_allocVector(REALSXP,1));
REAL(rets)[0] = PreciseSums_sum(dinput, len);
UNPROTECT(1);
return rets;
Expand All @@ -377,13 +378,13 @@ SEXP _psSum(SEXP input){
extern double PreciseSums_sumV(int n, ...){
va_list valist;
va_start(valist, n);
double *p = Calloc(n, double);
double *p = R_Calloc(n, double);
for (unsigned int i = 0; i < n; i++){
p[i] = va_arg(valist, double);
}
va_end(valist);
double s = PreciseSums_sum(p, n);
Free(p);
R_Free(p);
return s;
}

Expand Down

0 comments on commit a23f045

Please sign in to comment.