diff options
Diffstat (limited to 'double-error.c')
| -rw-r--r-- | double-error.c | 39 |
1 files changed, 24 insertions, 15 deletions
diff --git a/double-error.c b/double-error.c index dfbb34e..4e52b18 100644 --- a/double-error.c +++ b/double-error.c @@ -4,13 +4,18 @@ #include <float.h> #include <stdbool.h> -typedef double vtype; -typedef double etype; +typedef long double vtype; +typedef long double etype; #define EPS 2*DBL_EPSILON #define MIN 2*DBL_MIN typedef struct {vtype v; etype e;} double_error; +vtype vabsolute(vtype v){ + //return fabs(v); + return fabsl(v); +} + double_error int_to_valuetype(int i){ double_error de; de.v = (vtype) i; @@ -33,7 +38,7 @@ bool exptype_is_infinite(double_error d){ return d.v == INFINITY; } double_error sum(double_error de1, double_error de2){ double_error de; de.v = de1.v + de2.v; - de.e = de1.e + de2.e + ( fabs(de1.v) + fabs(de2.v) + de1.e + de2.e )*EPS + MIN; + de.e = de1.e + de2.e + ( vabsolute(de1.v) + vabsolute(de2.v) + de1.e + de2.e )*EPS + MIN; return de; } @@ -45,11 +50,11 @@ double_error difference(double_error de1, double_error de2){ } bool is_greater_certainly(double_error de1, double_error de2){ - return de1.v-de1.e > de2.v+de2.e + EPS*(fabs(de1.v)+fabs(de1.e)+fabs(de2.v)+fabs(de2.e)) + MIN; + return de1.v-de1.e > de2.v+de2.e + EPS*(vabsolute(de1.v)+vabsolute(de1.e)+vabsolute(de2.v)+vabsolute(de2.e)) + MIN; } bool is_greater_possibly(double_error de1, double_error de2){ - return de1.v+de1.e + EPS*(fabs(de1.v)+fabs(de1.e)+fabs(de2.v)+fabs(de2.e)) + MIN > de2.v-de2.e; + return de1.v+de1.e + EPS*(vabsolute(de1.v)+vabsolute(de1.e)+vabsolute(de2.v)+vabsolute(de2.e)) + MIN > de2.v-de2.e; } double_error maximum(double_error de1, double_error de2){ @@ -68,26 +73,26 @@ double_error maximum(double_error de1, double_error de2){ double_error product(double_error de1, double_error de2){ double_error de; de.v = de1.v * de2.v; - vtype p12 = fabs(de1.v)*de2.e*(1+EPS) + MIN; - vtype p21 = fabs(de2.v)*de1.e*(1+EPS) + MIN; + vtype p12 = vabsolute(de1.v)*de2.e*(1+EPS) + MIN; + vtype p21 = vabsolute(de2.v)*de1.e*(1+EPS) + MIN; vtype p22 = de1.e*de2.e*(1+EPS) + MIN; etype s = (p12 + p21 + p22)*(1+2*EPS); - de.e = s + ( fabs(de.v) + s )*EPS + MIN; + de.e = s + ( vabsolute(de.v) + s )*EPS + MIN; return de; } double_error ratio(double_error de1, double_error de2){ double_error de2i; de2i.v = 1/de2.v; - if(de2.e > fabs(de2.v)/2) de2i.e = INFINITY; - else de2i.e = 2*de2.e/(de2.v*de2.v) + fabs(de2i.v)*EPS + MIN; //2* is extra margin in case de2.e is near de2.v/2 + if(de2.e > vabsolute(de2.v)/2) de2i.e = INFINITY; + else de2i.e = 2*de2.e/(de2.v*de2.v) + vabsolute(de2i.v)*EPS + MIN; //2* is extra margin in case de2.e is near de2.v/2 return product(de1,de2i); } double_error absolute(double_error de){ double_error dea; - dea.v = fabs(de.v); + dea.v = vabsolute(de.v); dea.e = de.e; return dea; } @@ -134,9 +139,11 @@ int valuetype_to_string(char* s, double_error de){ if(de.e == 0.0){ /*If number really is an int, print accordingly to give cleaner output for example when printing the function itself, which often will be integer valued without errors.*/ if((vtype) ((int) de.v) == de.v) sprintf(s,"%d", (int)de.v); - else sprintf(s,"%4.3f…",de.v); + //else sprintf(s,"%4.3f…",de.v); + else sprintf(s,"%4.3Lf…",de.v); } - else sprintf(s,"%4.3f… +/- %6.1e",de.v,de.e); + //else sprintf(s,"%4.3f… +/- %6.1e",de.v,de.e); + else sprintf(s,"%4.3Lf… +/- %6.1e",de.v,de.e); return 0; } @@ -149,9 +156,11 @@ int valuetype_to_latex(char* s, double_error de){ if(de.e == 0.0){ /*If number really is an int, print accordingly to give cleaner output for example when printing the function itself, which often will be integer valued without errors.*/ if((vtype) ((int) de.v) == de.v) sprintf(s,"%d", (int)de.v); - else sprintf(s,"%4.3f\\ldots",de.v); + //else sprintf(s,"%4.3f\\ldots",de.v); + else sprintf(s,"%4.3Lf\\ldots",de.v); } - else sprintf(s,"%4.3f\\ldots\\pm\\texttt{%1.0e}",de.v,de.e); + //else sprintf(s,"%4.3f\\ldots\\pm\\texttt{%1.0e}",de.v,de.e); + else sprintf(s,"%4.3Lf\\ldots\\pm\\texttt{%1.0e}",de.v,de.e); return 0; } |