Ian Jauslin
summaryrefslogtreecommitdiff
diff options
context:
space:
mode:
authorIan Jauslin <ian.jauslin@roma1.infn.it>2015-10-07 12:51:41 +0000
committerIan Jauslin <ian.jauslin@roma1.infn.it>2015-10-07 13:00:23 +0000
commit469bdc80712dbf9c12562059dc4594620b59a076 (patch)
treec6da40a884899110d102d82a7a778f2b3afae702 /src/coefficient.c
parente7aa6859f08565d58684fa4b9c40fed716f0ba17 (diff)
Support MPFR floats in numkondov1.4
Remove '-D' option (error tolerance) in numkondo
Diffstat (limited to 'src/coefficient.c')
-rw-r--r--src/coefficient.c47
1 files changed, 46 insertions, 1 deletions
diff --git a/src/coefficient.c b/src/coefficient.c
index c26b668..d4cb6eb 100644
--- a/src/coefficient.c
+++ b/src/coefficient.c
@@ -18,6 +18,10 @@ limitations under the License.
#include <stdio.h>
#include <stdlib.h>
+#include <stdarg.h>
+// define MPFR_USE_VA_LIST to enable the use of mpfr_inits and mpfr_clears
+#define MPFR_USE_VA_LIST
+#include <mpfr.h>
#include "definitions.cpp"
#include "rational.h"
#include "istring.h"
@@ -721,7 +725,7 @@ int evalcoef(RCC rccs, Coefficient coef, long double* out){
int i,j;
long double num_factor;
- *out=0;
+ *out=0.;
// for each monomial
for(i=0;i<coef.length;i++){
@@ -737,3 +741,44 @@ int evalcoef(RCC rccs, Coefficient coef, long double* out){
}
return(0);
}
+
+// evaluate a coefficient on a vector (using mpfr floats)
+int evalcoef_mpfr(RCC_mpfr rccs, Coefficient coef, mpfr_t out){
+ int i,j;
+ mpfr_t num_factor;
+ // tmp number (do not initialize Z)
+ mpfr_t x, y, Z;
+
+ // init numbers
+ mpfr_inits(num_factor, x, y, (mpfr_ptr) NULL);
+
+ mpfr_init(out);
+ mpfr_set_zero(out, 1);
+
+ // for each monomial
+ for(i=0;i<coef.length;i++){
+ // product of factors
+ mpfr_set_flt(num_factor, 1., MPFR_RNDN);
+ for(j=0;j<coef.factors[i].length;j++){
+ mpfr_mul(x,num_factor,rccs.values[intlist_find_err(rccs.indices,rccs.length,coef.factors[i].values[j])], MPFR_RNDN);
+ mpfr_set(num_factor,x, MPFR_RNDN);
+ }
+ // denominator
+ if(coef.denoms[i].power>0){
+ mpfr_pow_si(y, rccs.values[intlist_find_err(rccs.indices,rccs.length,coef.denoms[i].index)], coef.denoms[i].power, MPFR_RNDN);
+ mpfr_div(x, num_factor, y, MPFR_RNDN);
+ mpfr_set(num_factor, x, MPFR_RNDN);
+ }
+
+ number_mpfr_val(Z, coef.nums[i]);
+ mpfr_mul(x, num_factor, Z, MPFR_RNDN);
+ mpfr_add(y, x, out, MPFR_RNDN);
+ mpfr_set(out, y, MPFR_RNDN);
+ mpfr_clear(Z);
+ }
+
+ // free numbers
+ mpfr_clears(num_factor, x, y, (mpfr_ptr)NULL);
+
+ return(0);
+}