1/*
2 *  Licensed to the Apache Software Foundation (ASF) under one or more
3 *  contributor license agreements.  See the NOTICE file distributed with
4 *  this work for additional information regarding copyright ownership.
5 *  The ASF licenses this file to You under the Apache License, Version 2.0
6 *  (the "License"); you may not use this file except in compliance with
7 *  the License.  You may obtain a copy of the License at
8 *
9 *     http://www.apache.org/licenses/LICENSE-2.0
10 *
11 *  Unless required by applicable law or agreed to in writing, software
12 *  distributed under the License is distributed on an "AS IS" BASIS,
13 *  WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
14 *  See the License for the specific language governing permissions and
15 *  limitations under the License.
16 */
17
18#include <stdlib.h>
19#include <string.h>
20#include <math.h>
21#include "JNIHelp.h"
22#include "commonDblParce.h"
23#include "cbigint.h"
24
25#if defined(LINUX) || defined(FREEBSD) || defined(ZOS)
26#define USE_LL
27#endif
28
29#define LOW_I32_FROM_VAR(u64)     LOW_I32_FROM_LONG64(u64)
30#ifdef HY_LITTLE_ENDIAN
31#define LOW_I32_FROM_PTR(ptr64) (*(I_32 *) (ptr64))
32#else
33#define LOW_I32_FROM_PTR(ptr64) (*(((I_32 *) (ptr64)) + 1))
34#endif
35#define HIGH_I32_FROM_VAR(u64)    HIGH_I32_FROM_LONG64(u64)
36#define HIGH_I32_FROM_PTR(u64ptr) HIGH_I32_FROM_LONG64_PTR(u64ptr)
37
38#define MAX_ACCURACY_WIDTH 8
39
40#define DEFAULT_WIDTH MAX_ACCURACY_WIDTH
41
42JNIEXPORT jfloat JNICALL
43Java_org_apache_harmony_luni_util_FloatingPointParser_parseFltImpl (JNIEnv * env,
44                                                        jclass clazz,
45                                                        jstring s, jint e);
46JNIEXPORT jdouble JNICALL
47Java_org_apache_harmony_luni_util_FloatingPointParser_parseDblImpl (JNIEnv * env,
48                                                        jclass clazz,
49                                                        jstring s, jint e);
50
51jfloat createFloat1 (JNIEnv * env, U_64 * f, IDATA length, jint e);
52jfloat floatAlgorithm (JNIEnv * env, U_64 * f, IDATA length, jint e,
53                       jfloat z);
54jfloat createFloat (JNIEnv * env, const char *s, jint e);
55
56static const U_32 tens[] = {
57  0x3f800000,
58  0x41200000,
59  0x42c80000,
60  0x447a0000,
61  0x461c4000,
62  0x47c35000,
63  0x49742400,
64  0x4b189680,
65  0x4cbebc20,
66  0x4e6e6b28,
67  0x501502f9                    /* 10 ^ 10 in float */
68};
69
70#define tenToTheE(e) (*((jfloat *) (tens + (e))))
71#define LOG5_OF_TWO_TO_THE_N 11
72
73#define sizeOfTenToTheE(e) (((e) / 19) + 1)
74
75#define INFINITE_INTBITS (0x7F800000)
76#define MINIMUM_INTBITS (1)
77
78#define MANTISSA_MASK (0x007FFFFF)
79#define EXPONENT_MASK (0x7F800000)
80#define NORMAL_MASK   (0x00800000)
81#define FLOAT_TO_INTBITS(flt) (*((U_32 *)(&flt)))
82
83/* Keep a count of the number of times we decrement and increment to
84 * approximate the double, and attempt to detect the case where we
85 * could potentially toggle back and forth between decrementing and
86 * incrementing. It is possible for us to be stuck in the loop when
87 * incrementing by one or decrementing by one may exceed or stay below
88 * the value that we are looking for. In this case, just break out of
89 * the loop if we toggle between incrementing and decrementing for more
90 * than twice.
91 */
92#define INCREMENT_FLOAT(_x, _decCount, _incCount) \
93    { \
94        ++FLOAT_TO_INTBITS(_x); \
95        _incCount++; \
96        if( (_incCount > 2) && (_decCount > 2) ) { \
97            if( _decCount > _incCount ) { \
98                FLOAT_TO_INTBITS(_x) += _decCount - _incCount; \
99            } else if( _incCount > _decCount ) { \
100                FLOAT_TO_INTBITS(_x) -= _incCount - _decCount; \
101            } \
102            break; \
103        } \
104    }
105#define DECREMENT_FLOAT(_x, _decCount, _incCount) \
106    { \
107        --FLOAT_TO_INTBITS(_x); \
108        _decCount++; \
109        if( (_incCount > 2) && (_decCount > 2) ) { \
110            if( _decCount > _incCount ) { \
111                FLOAT_TO_INTBITS(_x) += _decCount - _incCount; \
112            } else if( _incCount > _decCount ) { \
113                FLOAT_TO_INTBITS(_x) -= _incCount - _decCount; \
114            } \
115            break; \
116        } \
117    }
118#define ERROR_OCCURED(x) (HIGH_I32_FROM_VAR(x) < 0)
119
120#define allocateU64(x, n) if (!((x) = (U_64*) malloc((n) * sizeof(U_64)))) goto OutOfMemory;
121#define release(r) if ((r)) free((r));
122
123jfloat
124createFloat (JNIEnv * env, const char *s, jint e)
125{
126  /* assumes s is a null terminated string with at least one
127   * character in it */
128  U_64 def[DEFAULT_WIDTH];
129  U_64 defBackup[DEFAULT_WIDTH];
130  U_64 *f, *fNoOverflow, *g, *tempBackup;
131  U_32 overflow;
132  jfloat result;
133  IDATA index = 1;
134  int unprocessedDigits = 0;
135
136  f = def;
137  fNoOverflow = defBackup;
138  *f = 0;
139  tempBackup = g = 0;
140  do
141    {
142      if (*s >= '0' && *s <= '9')
143        {
144          /* Make a back up of f before appending, so that we can
145           * back out of it if there is no more room, i.e. index >
146           * MAX_ACCURACY_WIDTH.
147           */
148          memcpy (fNoOverflow, f, sizeof (U_64) * index);
149          overflow =
150            simpleAppendDecimalDigitHighPrecision (f, index, *s - '0');
151          if (overflow)
152            {
153
154              f[index++] = overflow;
155              /* There is an overflow, but there is no more room
156               * to store the result. We really only need the top 52
157               * bits anyway, so we must back out of the overflow,
158               * and ignore the rest of the string.
159               */
160              if (index >= MAX_ACCURACY_WIDTH)
161                {
162                  index--;
163                  memcpy (f, fNoOverflow, sizeof (U_64) * index);
164                  break;
165                }
166              if (tempBackup)
167                {
168                  fNoOverflow = tempBackup;
169                }
170            }
171        }
172      else
173        index = -1;
174    }
175  while (index > 0 && *(++s) != '\0');
176
177  /* We've broken out of the parse loop either because we've reached
178   * the end of the string or we've overflowed the maximum accuracy
179   * limit of a double. If we still have unprocessed digits in the
180   * given string, then there are three possible results:
181   *   1. (unprocessed digits + e) == 0, in which case we simply
182   *      convert the existing bits that are already parsed
183   *   2. (unprocessed digits + e) < 0, in which case we simply
184   *      convert the existing bits that are already parsed along
185   *      with the given e
186   *   3. (unprocessed digits + e) > 0 indicates that the value is
187   *      simply too big to be stored as a double, so return Infinity
188   */
189  if ((unprocessedDigits = strlen (s)) > 0)
190    {
191      e += unprocessedDigits;
192      if (index > -1)
193        {
194          if (e <= 0)
195            {
196              result = createFloat1 (env, f, index, e);
197            }
198          else
199            {
200              FLOAT_TO_INTBITS (result) = INFINITE_INTBITS;
201            }
202        }
203      else
204        {
205          result = *(jfloat *) & index;
206        }
207    }
208  else
209    {
210      if (index > -1)
211        {
212          result = createFloat1 (env, f, index, e);
213        }
214      else
215        {
216        result = *(jfloat *) & index;
217        }
218    }
219
220  return result;
221
222}
223
224jfloat
225createFloat1 (JNIEnv * env, U_64 * f, IDATA length, jint e)
226{
227  IDATA numBits;
228  jdouble dresult;
229  jfloat result;
230
231  numBits = highestSetBitHighPrecision (f, length) + 1;
232  if (numBits < 25 && e >= 0 && e < LOG5_OF_TWO_TO_THE_N)
233    {
234      return ((jfloat) LOW_I32_FROM_PTR (f)) * tenToTheE (e);
235    }
236  else if (numBits < 25 && e < 0 && (-e) < LOG5_OF_TWO_TO_THE_N)
237    {
238      return ((jfloat) LOW_I32_FROM_PTR (f)) / tenToTheE (-e);
239    }
240  else if (e >= 0 && e < 39)
241    {
242      result = (jfloat) (toDoubleHighPrecision (f, length) * pow (10.0, (double) e));
243    }
244  else if (e >= 39)
245    {
246      /* Convert the partial result to make sure that the
247       * non-exponential part is not zero. This check fixes the case
248       * where the user enters 0.0e309! */
249      result = (jfloat) toDoubleHighPrecision (f, length);
250
251      if (result == 0.0)
252
253        FLOAT_TO_INTBITS (result) = MINIMUM_INTBITS;
254      else
255        FLOAT_TO_INTBITS (result) = INFINITE_INTBITS;
256    }
257  else if (e > -309)
258    {
259      int dexp;
260      U_32 fmant, fovfl;
261      U_64 dmant;
262      dresult = toDoubleHighPrecision (f, length) / pow (10.0, (double) -e);
263      if (IS_DENORMAL_DBL (dresult))
264        {
265          FLOAT_TO_INTBITS (result) = 0;
266          return result;
267        }
268      dexp = doubleExponent (dresult) + 51;
269      dmant = doubleMantissa (dresult);
270      /* Is it too small to be represented by a single-precision
271       * float? */
272      if (dexp <= -155)
273        {
274          FLOAT_TO_INTBITS (result) = 0;
275          return result;
276        }
277      /* Is it a denormalized single-precision float? */
278      if ((dexp <= -127) && (dexp > -155))
279        {
280          /* Only interested in 24 msb bits of the 53-bit double mantissa */
281          fmant = (U_32) (dmant >> 29);
282          fovfl = ((U_32) (dmant & 0x1FFFFFFF)) << 3;
283          while ((dexp < -127) && ((fmant | fovfl) != 0))
284            {
285              if ((fmant & 1) != 0)
286                {
287                  fovfl |= 0x80000000;
288                }
289              fovfl >>= 1;
290              fmant >>= 1;
291              dexp++;
292            }
293          if ((fovfl & 0x80000000) != 0)
294            {
295              if ((fovfl & 0x7FFFFFFC) != 0)
296                {
297                  fmant++;
298                }
299              else if ((fmant & 1) != 0)
300                {
301                  fmant++;
302                }
303            }
304          else if ((fovfl & 0x40000000) != 0)
305            {
306              if ((fovfl & 0x3FFFFFFC) != 0)
307                {
308                  fmant++;
309                }
310            }
311          FLOAT_TO_INTBITS (result) = fmant;
312        }
313      else
314        {
315          result = (jfloat) dresult;
316        }
317    }
318
319  /* Don't go straight to zero as the fact that x*0 = 0 independent
320   * of x might cause the algorithm to produce an incorrect result.
321   * Instead try the min  value first and let it fall to zero if need
322   * be.
323   */
324  if (e <= -309 || FLOAT_TO_INTBITS (result) == 0)
325    FLOAT_TO_INTBITS (result) = MINIMUM_INTBITS;
326
327  return floatAlgorithm (env, f, length, e, (jfloat) result);
328}
329
330#if defined(WIN32)
331/* disable global optimizations on the microsoft compiler for the
332 * floatAlgorithm function otherwise it won't properly compile */
333#pragma optimize("g",off)
334#endif
335
336/* The algorithm for the function floatAlgorithm() below can be found
337 * in:
338 *
339 *      "How to Read Floating-Point Numbers Accurately", William D.
340 *      Clinger, Proceedings of the ACM SIGPLAN '90 Conference on
341 *      Programming Language Design and Implementation, June 20-22,
342 *      1990, pp. 92-101.
343 *
344 * There is a possibility that the function will end up in an endless
345 * loop if the given approximating floating-point number (a very small
346 * floating-point whose value is very close to zero) straddles between
347 * two approximating integer values. We modified the algorithm slightly
348 * to detect the case where it oscillates back and forth between
349 * incrementing and decrementing the floating-point approximation. It
350 * is currently set such that if the oscillation occurs more than twice
351 * then return the original approximation.
352 */
353jfloat
354floatAlgorithm (JNIEnv * env, U_64 * f, IDATA length, jint e, jfloat z)
355{
356  U_64 m;
357  IDATA k, comparison, comparison2;
358  U_64 *x, *y, *D, *D2;
359  IDATA xLength, yLength, DLength, D2Length;
360  IDATA decApproxCount, incApproxCount;
361
362  x = y = D = D2 = 0;
363  xLength = yLength = DLength = D2Length = 0;
364  decApproxCount = incApproxCount = 0;
365
366  do
367    {
368      m = floatMantissa (z);
369      k = floatExponent (z);
370
371      if (x && x != f)
372          free(x);
373
374      release (y);
375      release (D);
376      release (D2);
377
378      if (e >= 0 && k >= 0)
379        {
380          xLength = sizeOfTenToTheE (e) + length;
381          allocateU64 (x, xLength);
382          memset (x + length, 0, sizeof (U_64) * (xLength - length));
383          memcpy (x, f, sizeof (U_64) * length);
384          timesTenToTheEHighPrecision (x, xLength, e);
385
386          yLength = (k >> 6) + 2;
387          allocateU64 (y, yLength);
388          memset (y + 1, 0, sizeof (U_64) * (yLength - 1));
389          *y = m;
390          simpleShiftLeftHighPrecision (y, yLength, k);
391        }
392      else if (e >= 0)
393        {
394          xLength = sizeOfTenToTheE (e) + length + ((-k) >> 6) + 1;
395          allocateU64 (x, xLength);
396          memset (x + length, 0, sizeof (U_64) * (xLength - length));
397          memcpy (x, f, sizeof (U_64) * length);
398          timesTenToTheEHighPrecision (x, xLength, e);
399          simpleShiftLeftHighPrecision (x, xLength, -k);
400
401          yLength = 1;
402          allocateU64 (y, 1);
403          *y = m;
404        }
405      else if (k >= 0)
406        {
407          xLength = length;
408          x = f;
409
410          yLength = sizeOfTenToTheE (-e) + 2 + (k >> 6);
411          allocateU64 (y, yLength);
412          memset (y + 1, 0, sizeof (U_64) * (yLength - 1));
413          *y = m;
414          timesTenToTheEHighPrecision (y, yLength, -e);
415          simpleShiftLeftHighPrecision (y, yLength, k);
416        }
417      else
418        {
419          xLength = length + ((-k) >> 6) + 1;
420          allocateU64 (x, xLength);
421          memset (x + length, 0, sizeof (U_64) * (xLength - length));
422          memcpy (x, f, sizeof (U_64) * length);
423          simpleShiftLeftHighPrecision (x, xLength, -k);
424
425          yLength = sizeOfTenToTheE (-e) + 1;
426          allocateU64 (y, yLength);
427          memset (y + 1, 0, sizeof (U_64) * (yLength - 1));
428          *y = m;
429          timesTenToTheEHighPrecision (y, yLength, -e);
430        }
431
432      comparison = compareHighPrecision (x, xLength, y, yLength);
433      if (comparison > 0)
434        {                       /* x > y */
435          DLength = xLength;
436          allocateU64 (D, DLength);
437          memcpy (D, x, DLength * sizeof (U_64));
438          subtractHighPrecision (D, DLength, y, yLength);
439        }
440      else if (comparison)
441        {                       /* y > x */
442          DLength = yLength;
443          allocateU64 (D, DLength);
444          memcpy (D, y, DLength * sizeof (U_64));
445          subtractHighPrecision (D, DLength, x, xLength);
446        }
447      else
448        {                       /* y == x */
449          DLength = 1;
450          allocateU64 (D, 1);
451          *D = 0;
452        }
453
454      D2Length = DLength + 1;
455      allocateU64 (D2, D2Length);
456      m <<= 1;
457      multiplyHighPrecision (D, DLength, &m, 1, D2, D2Length);
458      m >>= 1;
459
460      comparison2 = compareHighPrecision (D2, D2Length, y, yLength);
461      if (comparison2 < 0)
462        {
463          if (comparison < 0 && m == NORMAL_MASK)
464            {
465              simpleShiftLeftHighPrecision (D2, D2Length, 1);
466              if (compareHighPrecision (D2, D2Length, y, yLength) > 0)
467                {
468                  DECREMENT_FLOAT (z, decApproxCount, incApproxCount);
469                }
470              else
471                {
472                  break;
473                }
474            }
475          else
476            {
477              break;
478            }
479        }
480      else if (comparison2 == 0)
481        {
482          if ((m & 1) == 0)
483            {
484              if (comparison < 0 && m == NORMAL_MASK)
485                {
486                  DECREMENT_FLOAT (z, decApproxCount, incApproxCount);
487                }
488              else
489                {
490                  break;
491                }
492            }
493          else if (comparison < 0)
494            {
495              DECREMENT_FLOAT (z, decApproxCount, incApproxCount);
496              break;
497            }
498          else
499            {
500              INCREMENT_FLOAT (z, decApproxCount, incApproxCount);
501              break;
502            }
503        }
504      else if (comparison < 0)
505        {
506          DECREMENT_FLOAT (z, decApproxCount, incApproxCount);
507        }
508      else
509        {
510          if (FLOAT_TO_INTBITS (z) == EXPONENT_MASK)
511            break;
512          INCREMENT_FLOAT (z, decApproxCount, incApproxCount);
513        }
514    }
515  while (1);
516
517  if (x && x != f)
518      free(x);
519  release (y);
520  release (D);
521  release (D2);
522  return z;
523
524OutOfMemory:
525  if (x && x != f)
526      free(x);
527  release (y);
528  release (D);
529  release (D2);
530
531  FLOAT_TO_INTBITS (z) = -2;
532
533  return z;
534}
535
536#if defined(WIN32)
537#pragma optimize("",on)         /*restore optimizations */
538#endif
539
540JNIEXPORT jfloat JNICALL
541Java_org_apache_harmony_luni_util_FloatingPointParser_parseFltImpl (JNIEnv * env,
542                                                        jclass clazz,
543                                                        jstring s, jint e)
544{
545  jfloat flt;
546  const char *str = (*env)->GetStringUTFChars (env, s, 0);
547  flt = createFloat (env, str, e);
548  (*env)->ReleaseStringUTFChars (env, s, str);
549
550  if (((I_32) FLOAT_TO_INTBITS (flt)) >= 0)
551    {
552      return flt;
553    }
554  else if (((I_32) FLOAT_TO_INTBITS (flt)) == (I_32) - 1)
555    {                           /* NumberFormatException */
556      jniThrowException(env, "java/lang/NumberFormatException", "");
557    }
558  else
559    {                           /* OutOfMemoryError */
560      jniThrowException(env, "java/lang/OutOfMemoryError", "");
561    }
562
563  return 0.0;
564}
565
566JNIEXPORT jdouble JNICALL
567Java_org_apache_harmony_luni_util_FloatingPointParser_parseDblImpl (JNIEnv * env,
568                                                        jclass clazz,
569                                                        jstring s, jint e)
570{
571  jdouble dbl;
572  const char *str = (*env)->GetStringUTFChars (env, s, 0);
573  dbl = createDouble (env, str, e);
574  (*env)->ReleaseStringUTFChars (env, s, str);
575
576  if (!ERROR_OCCURED (dbl))
577    {
578      return dbl;
579    }
580  else if (LOW_I32_FROM_VAR (dbl) == (I_32) - 1)
581    {                           /* NumberFormatException */
582      jniThrowException(env, "java/lang/NumberFormatException", "");
583    }
584  else
585    {                           /* OutOfMemoryError */
586      jniThrowException(env, "java/lang/OutOfMemoryError", "");
587    }
588
589  return 0.0;
590}
591
592/*
593 * JNI registration
594 */
595static JNINativeMethod gMethods[] = {
596    /* NAME,          SIGNATURE,                FUNCPTR */
597    { "parseFltImpl", "(Ljava/lang/String;I)F",
598        Java_org_apache_harmony_luni_util_FloatingPointParser_parseFltImpl },
599    { "parseDblImpl", "(Ljava/lang/String;I)D",
600        Java_org_apache_harmony_luni_util_FloatingPointParser_parseDblImpl },
601};
602int register_org_apache_harmony_luni_util_fltparse(JNIEnv *env)
603{
604    return jniRegisterNativeMethods(env,
605               "org/apache/harmony/luni/util/FloatingPointParser",
606                gMethods, NELEM(gMethods));
607}
608