l4/pkg/ocaml/contrib/byterun/floats.c

   1 /***********************************************************************/
   2 /*                                                                     */
   3 /*                           Objective Caml                            */
   4 /*                                                                     */
   5 /*            Xavier Leroy, projet Cristal, INRIA Rocquencourt         */
   6 /*                                                                     */
   7 /*  Copyright 1996 Institut National de Recherche en Informatique et   */
   8 /*  en Automatique.  All rights reserved.  This file is distributed    */
   9 /*  under the terms of the GNU Library General Public License, with    */
  10 /*  the special exception on linking described in file ../LICENSE.     */
  11 /*                                                                     */
  12 /***********************************************************************/
  13
  14 /* $Id: floats.c 8976 2008-08-02 11:02:28Z xleroy $ */
  15
  16 /* The interface of this file is in "mlvalues.h" and "alloc.h" */
  17
  18 #include <math.h>
  19 #include <stdio.h>
  20 #include <stdlib.h>
  21 #include <string.h>
  22
  23 #include "alloc.h"
  24 #include "fail.h"
  25 #include "memory.h"
  26 #include "mlvalues.h"
  27 #include "misc.h"
  28 #include "reverse.h"
  29 #include "stacks.h"
  30
  31 #ifdef ARCH_ALIGN_DOUBLE
  32
  33 CAMLexport double caml_Double_val(value val)
  34 {
  35   union { value v[2]; double d; } buffer;
  36
  37   Assert(sizeof(double) == 2 * sizeof(value));
  38   buffer.v[0] = Field(val, 0);
  39   buffer.v[1] = Field(val, 1);
  40   return buffer.d;
  41 }
  42
  43 CAMLexport void caml_Store_double_val(value val, double dbl)
  44 {
  45   union { value v[2]; double d; } buffer;
  46
  47   Assert(sizeof(double) == 2 * sizeof(value));
  48   buffer.d = dbl;
  49   Field(val, 0) = buffer.v[0];
  50   Field(val, 1) = buffer.v[1];
  51 }
  52
  53 #endif
  54
  55 CAMLexport value caml_copy_double(double d)
  56 {
  57   value res;
  58
  59 #define Setup_for_gc
  60 #define Restore_after_gc
  61   Alloc_small(res, Double_wosize, Double_tag);
  62 #undef Setup_for_gc
  63 #undef Restore_after_gc
  64   Store_double_val(res, d);
  65   return res;
  66 }
  67
  68 CAMLprim value caml_format_float(value fmt, value arg)
  69 {
  70 #define MAX_DIGITS 350
  71 /* Max number of decimal digits in a "natural" (not artificially padded)
  72    representation of a float. Can be quite big for %f format.
  73    Max exponent for IEEE format is 308 decimal digits.
  74    Rounded up for good measure. */
  75   char format_buffer[MAX_DIGITS + 20];
  76   int prec, i;
  77   char * p;
  78   char * dest;
  79   value res;
  80
  81   prec = MAX_DIGITS;
  82   for (p = String_val(fmt); *p != 0; p++) {
  83     if (*p >= '0' && *p <= '9') {
  84       i = atoi(p) + MAX_DIGITS;
  85       if (i > prec) prec = i;
  86       break;
  87     }
  88   }
  89   for( ; *p != 0; p++) {
  90     if (*p == '.') {
  91       i = atoi(p+1) + MAX_DIGITS;
  92       if (i > prec) prec = i;
  93       break;
  94     }
  95   }
  96   if (prec < sizeof(format_buffer)) {
  97     dest = format_buffer;
  98   } else {
  99     dest = caml_stat_alloc(prec);
 100   }
 101   sprintf(dest, String_val(fmt), Double_val(arg));
 102   res = caml_copy_string(dest);
 103   if (dest != format_buffer) {
 104     caml_stat_free(dest);
 105   }
 106   return res;
 107 }
 108
 109 /*CAMLprim*/ value caml_float_of_substring(value vs, value idx, value l)
 110 {
 111   char parse_buffer[64];
 112   char * buf, * src, * dst, * end;
 113   mlsize_t len, lenvs;
 114   double d;
 115   intnat flen = Long_val(l);
 116   intnat fidx = Long_val(idx);
 117
 118   lenvs = caml_string_length(vs);
 119   len =
 120     fidx >= 0 && fidx < lenvs && flen > 0 && flen <= lenvs - fidx
 121     ? flen : 0;
 122   buf = len < sizeof(parse_buffer) ? parse_buffer : caml_stat_alloc(len + 1);
 123   src = String_val(vs) + fidx;
 124   dst = buf;
 125   while (len--) {
 126     char c = *src++;
 127     if (c != '_') *dst++ = c;
 128   }
 129   *dst = 0;
 130   if (dst == buf) goto error;
 131   d = strtod((const char *) buf, &end);
 132   if (end != dst) goto error;
 133   if (buf != parse_buffer) caml_stat_free(buf);
 134   return caml_copy_double(d);
 135  error:
 136   if (buf != parse_buffer) caml_stat_free(buf);
 137   caml_failwith("float_of_string");
 138 }
 139
 140 CAMLprim value caml_float_of_string(value vs)
 141 {
 142   char parse_buffer[64];
 143   char * buf, * src, * dst, * end;
 144   mlsize_t len;
 145   double d;
 146
 147   len = caml_string_length(vs);
 148   buf = len < sizeof(parse_buffer) ? parse_buffer : caml_stat_alloc(len + 1);
 149   src = String_val(vs);
 150   dst = buf;
 151   while (len--) {
 152     char c = *src++;
 153     if (c != '_') *dst++ = c;
 154   }
 155   *dst = 0;
 156   if (dst == buf) goto error;
 157   d = strtod((const char *) buf, &end);
 158   if (end != dst) goto error;
 159   if (buf != parse_buffer) caml_stat_free(buf);
 160   return caml_copy_double(d);
 161  error:
 162   if (buf != parse_buffer) caml_stat_free(buf);
 163   caml_failwith("float_of_string");
 164 }
 165
 166 CAMLprim value caml_int_of_float(value f)
 167 {
 168   return Val_long((intnat) Double_val(f));
 169 }
 170
 171 CAMLprim value caml_float_of_int(value n)
 172 {
 173   return caml_copy_double((double) Long_val(n));
 174 }
 175
 176 CAMLprim value caml_neg_float(value f)
 177 {
 178   return caml_copy_double(- Double_val(f));
 179 }
 180
 181 CAMLprim value caml_abs_float(value f)
 182 {
 183   return caml_copy_double(fabs(Double_val(f)));
 184 }
 185
 186 CAMLprim value caml_add_float(value f, value g)
 187 {
 188   return caml_copy_double(Double_val(f) + Double_val(g));
 189 }
 190
 191 CAMLprim value caml_sub_float(value f, value g)
 192 {
 193   return caml_copy_double(Double_val(f) - Double_val(g));
 194 }
 195
 196 CAMLprim value caml_mul_float(value f, value g)
 197 {
 198   return caml_copy_double(Double_val(f) * Double_val(g));
 199 }
 200
 201 CAMLprim value caml_div_float(value f, value g)
 202 {
 203   return caml_copy_double(Double_val(f) / Double_val(g));
 204 }
 205
 206 CAMLprim value caml_exp_float(value f)
 207 {
 208   return caml_copy_double(exp(Double_val(f)));
 209 }
 210
 211 CAMLprim value caml_floor_float(value f)
 212 {
 213   return caml_copy_double(floor(Double_val(f)));
 214 }
 215
 216 CAMLprim value caml_fmod_float(value f1, value f2)
 217 {
 218   return caml_copy_double(fmod(Double_val(f1), Double_val(f2)));
 219 }
 220
 221 CAMLprim value caml_frexp_float(value f)
 222 {
 223   CAMLparam1 (f);
 224   CAMLlocal2 (res, mantissa);
 225   int exponent;
 226
 227   mantissa = caml_copy_double(frexp (Double_val(f), &exponent));
 228   res = caml_alloc_tuple(2);
 229   Field(res, 0) = mantissa;
 230   Field(res, 1) = Val_int(exponent);
 231   CAMLreturn (res);
 232 }
 233
 234 CAMLprim value caml_ldexp_float(value f, value i)
 235 {
 236   return caml_copy_double(ldexp(Double_val(f), Int_val(i)));
 237 }
 238
 239 CAMLprim value caml_log_float(value f)
 240 {
 241   return caml_copy_double(log(Double_val(f)));
 242 }
 243
 244 CAMLprim value caml_log10_float(value f)
 245 {
 246   return caml_copy_double(log10(Double_val(f)));
 247 }
 248
 249 CAMLprim value caml_modf_float(value f)
 250 {
 251 #if __SC__
 252   _float_eval frem;       /* Problem with Apple's <math.h> */
 253 #else
 254   double frem;
 255 #endif
 256   CAMLparam1 (f);
 257   CAMLlocal3 (res, quo, rem);
 258
 259   quo = caml_copy_double(modf (Double_val(f), &frem));
 260   rem = caml_copy_double(frem);
 261   res = caml_alloc_tuple(2);
 262   Field(res, 0) = quo;
 263   Field(res, 1) = rem;
 264   CAMLreturn (res);
 265 }
 266
 267 CAMLprim value caml_sqrt_float(value f)
 268 {
 269   return caml_copy_double(sqrt(Double_val(f)));
 270 }
 271
 272 CAMLprim value caml_power_float(value f, value g)
 273 {
 274   return caml_copy_double(pow(Double_val(f), Double_val(g)));
 275 }
 276
 277 CAMLprim value caml_sin_float(value f)
 278 {
 279   return caml_copy_double(sin(Double_val(f)));
 280 }
 281
 282 CAMLprim value caml_sinh_float(value f)
 283 {
 284   return caml_copy_double(sinh(Double_val(f)));
 285 }
 286
 287 CAMLprim value caml_cos_float(value f)
 288 {
 289   return caml_copy_double(cos(Double_val(f)));
 290 }
 291
 292 CAMLprim value caml_cosh_float(value f)
 293 {
 294   return caml_copy_double(cosh(Double_val(f)));
 295 }
 296
 297 CAMLprim value caml_tan_float(value f)
 298 {
 299   return caml_copy_double(tan(Double_val(f)));
 300 }
 301
 302 CAMLprim value caml_tanh_float(value f)
 303 {
 304   return caml_copy_double(tanh(Double_val(f)));
 305 }
 306
 307 CAMLprim value caml_asin_float(value f)
 308 {
 309   return caml_copy_double(asin(Double_val(f)));
 310 }
 311
 312 CAMLprim value caml_acos_float(value f)
 313 {
 314   return caml_copy_double(acos(Double_val(f)));
 315 }
 316
 317 CAMLprim value caml_atan_float(value f)
 318 {
 319   return caml_copy_double(atan(Double_val(f)));
 320 }
 321
 322 CAMLprim value caml_atan2_float(value f, value g)
 323 {
 324   return caml_copy_double(atan2(Double_val(f), Double_val(g)));
 325 }
 326
 327 CAMLprim value caml_ceil_float(value f)
 328 {
 329   return caml_copy_double(ceil(Double_val(f)));
 330 }
 331
 332 CAMLprim value caml_eq_float(value f, value g)
 333 {
 334   return Val_bool(Double_val(f) == Double_val(g));
 335 }
 336
 337 CAMLprim value caml_neq_float(value f, value g)
 338 {
 339   return Val_bool(Double_val(f) != Double_val(g));
 340 }
 341
 342 CAMLprim value caml_le_float(value f, value g)
 343 {
 344   return Val_bool(Double_val(f) <= Double_val(g));
 345 }
 346
 347 CAMLprim value caml_lt_float(value f, value g)
 348 {
 349   return Val_bool(Double_val(f) < Double_val(g));
 350 }
 351
 352 CAMLprim value caml_ge_float(value f, value g)
 353 {
 354   return Val_bool(Double_val(f) >= Double_val(g));
 355 }
 356
 357 CAMLprim value caml_gt_float(value f, value g)
 358 {
 359   return Val_bool(Double_val(f) > Double_val(g));
 360 }
 361
 362 CAMLprim value caml_float_compare(value vf, value vg)
 363 {
 364   double f = Double_val(vf);
 365   double g = Double_val(vg);
 366   if (f == g) return Val_int(0);
 367   if (f < g) return Val_int(-1);
 368   if (f > g) return Val_int(1);
 369   /* One or both of f and g is NaN.  Order according to the
 370      convention NaN = NaN and NaN < x for all other floats x. */
 371   if (f == f) return Val_int(1);  /* f is not NaN, g is NaN */
 372   if (g == g) return Val_int(-1); /* g is not NaN, f is NaN */
 373   return Val_int(0);              /* both f and g are NaN */
 374 }
 375
 376 enum { FP_normal, FP_subnormal, FP_zero, FP_infinite, FP_nan };
 377
 378 CAMLprim value caml_classify_float(value vd)
 379 {
 380   /* Cygwin 1.3 has problems with fpclassify (PR#1293), so don't use it */
 381 #if defined(fpclassify) && !defined(__CYGWIN32__) && !defined(__MINGW32__)
 382   switch (fpclassify(Double_val(vd))) {
 383   case FP_NAN:
 384     return Val_int(FP_nan);
 385   case FP_INFINITE:
 386     return Val_int(FP_infinite);
 387   case FP_ZERO:
 388     return Val_int(FP_zero);
 389   case FP_SUBNORMAL:
 390     return Val_int(FP_subnormal);
 391   default: /* case FP_NORMAL */
 392     return Val_int(FP_normal);
 393   }
 394 #else
 395   union {
 396     double d;
 397 #if defined(ARCH_BIG_ENDIAN) || (defined(__arm__) && !defined(__ARM_EABI__))
 398     struct { uint32 h; uint32 l; } i;
 399 #else
 400     struct { uint32 l; uint32 h; } i;
 401 #endif
 402   } u;
 403   uint32 h, l;
 404
 405   u.d = Double_val(vd);
 406   h = u.i.h;  l = u.i.l;
 407   l = l | (h & 0xFFFFF);
 408   h = h & 0x7FF00000;
 409   if ((h | l) == 0)
 410     return Val_int(FP_zero);
 411   if (h == 0)
 412     return Val_int(FP_subnormal);
 413   if (h == 0x7FF00000) {
 414     if (l == 0)
 415       return Val_int(FP_infinite);
 416     else
 417       return Val_int(FP_nan);
 418   }
 419   return Val_int(FP_normal);
 420 #endif
 421 }
 422
 423 /* The [caml_init_ieee_float] function should initialize floating-point hardware
 424    so that it behaves as much as possible like the IEEE standard.
 425    In particular, return special numbers like Infinity and NaN instead
 426    of signalling exceptions.  Currently, everyone is in IEEE mode
 427    at program startup, except FreeBSD prior to 4.0R. */
 428
 429 #ifdef __FreeBSD__
 430 #include <osreldate.h>
 431 #if (__FreeBSD_version < 400017)
 432 #include <floatingpoint.h>
 433 #endif
 434 #endif
 435
 436 void caml_init_ieee_floats(void)
 437 {
 438 #if defined(__FreeBSD__) && (__FreeBSD_version < 400017)
 439   fpsetmask(0);
 440 #endif
 441 }