Fork me on GitHub

StringToDouble源码解析

由于工作需要更快速的将String转换成double,不满足于现有的API方式。故我需要去查看常用的String转Value方式以及期望能写出一个更高效的实现。

常用方法及耗时

常用的StringToDouble方法主要有以下四种:

1
2
3
4
5
6
7
8
9
10
11
@Test
public void StringToDoubleValue(){
Double.parseDouble("1.0");
Double.valueOf("1.0");
new BigDecimal("1.0").doubleValue();
try {
new DecimalFormat().parse("1.0").doubleValue();
} catch (ParseException e) {
e.printStackTrace();
}
}

其中上述四种方法的耗时统计如下(100万次循环结果,时间单位ms):

根据上图的结果可知Double类自带的parseDoublevalueOf方法速度较快。经过查看API可知:
这两个方法的具体实现都是一样的。


都是实现了FloatingDecimal类的readJavaFormatString(String in)方法
那么问题来了:
为什么方式二耗时明显低于方式一?
摆上我的代码:

1
2
3
4
5
6
7
8
9
10
11
long start1 = System.currentTimeMillis();
for (int i = 0; i < 1000000; i++) {
double d = Double.parseDouble("1.0");
}
System.out.println("方式一:" + (System.currentTimeMillis()-start1));
long start2 = System.currentTimeMillis();
for (int i = 0; i < 1000000; i++) {
double d = Double.valueOf("1.0");
}
System.out.println("方式二:" + (System.currentTimeMillis()-start2));

这是因为方式一我是在方式二上面执行的。如果换个位置,方式二先执行,那么方式一的耗时会低很多。
个人观点是:这里调用了FloatingDecimal.readJavaFormatString("1.0").doubleValue()方法。类的加载花费了一定的时间。

API中的StringToDouble方法其实一共调用了2个方法:
readJavaFormatString(String in)以及doubleValue()
这两个方法的源码我会在最后贴出

看不下去了,就两个简单的方法源码就800+行了

附件(源码)

readJavaFormatString(String in)以及doubleValue()以及相应的辅助方法和变量被我提取到一个类中

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
529
530
531
532
533
534
535
536
537
538
539
540
541
542
543
544
545
546
547
548
549
550
551
552
553
554
555
556
557
558
559
560
561
562
563
564
565
566
567
568
569
570
571
572
573
574
575
576
577
578
579
580
581
582
583
584
585
586
587
588
589
590
591
592
593
594
595
596
597
598
599
600
601
602
603
604
605
606
607
608
609
610
611
612
613
614
615
616
617
618
619
620
621
622
623
624
625
626
627
628
629
630
631
632
633
634
635
636
637
638
639
640
641
642
643
644
645
646
647
648
649
650
651
652
653
654
655
656
657
658
659
660
661
662
663
664
665
666
667
668
669
670
671
672
673
674
675
676
677
678
679
680
681
682
683
684
685
686
687
688
689
690
691
692
693
694
695
696
697
698
699
700
701
702
703
704
705
706
707
708
709
710
711
712
713
714
715
716
717
718
719
720
721
722
723
724
725
726
727
728
729
730
731
732
733
734
735
736
737
738
739
740
741
742
743
744
745
746
747
748
749
750
751
752
753
754
755
756
757
758
759
760
761
762
763
764
765
766
767
768
769
770
771
772
773
774
775
776
777
778
779
780
781
782
783
784
785
786
787
788
789
790
791
792
793
794
795
796
797
798
799
800
801
802
803
804
805
806
807
808
809
810
811
812
813
814
815
816
817
818
819
820
package MyTest;
public class MyFloatingDecimal {
boolean isExceptional;
boolean isNegative;
int decExponent;
char digits[];
int nDigits;
int bigIntExp;
int bigIntNBits;
boolean mustSetRoundDir = false;
int roundDir; // set by doubleValue
static final long signMask = 0x8000000000000000L;
static final long expMask = 0x7ff0000000000000L;
static final long fractMask= ~(signMask|expMask);
static final int expShift = 52;
static final int expBias = 1023;
static final long fractHOB = ( 1L<<expShift ); // assumed High-Order bit
static final int maxDecimalDigits = 15;
static final int maxDecimalExponent = 308;
static final int minDecimalExponent = -324;
static final int bigDecimalExponent = 324; // i.e. abs(minDecimalExponent)
static final long highbyte = 0xff00000000000000L;
static final long lowbytes = ~highbyte;
static final int singleSignMask = 0x80000000;
static final int singleExpMask = 0x7f800000;
static final int singleExpShift = 23;
static final int intDecimalDigits = 9;
private static final double small10pow[] = {
1.0e0,
1.0e1, 1.0e2, 1.0e3, 1.0e4, 1.0e5,
1.0e6, 1.0e7, 1.0e8, 1.0e9, 1.0e10,
1.0e11, 1.0e12, 1.0e13, 1.0e14, 1.0e15,
1.0e16, 1.0e17, 1.0e18, 1.0e19, 1.0e20,
1.0e21, 1.0e22
};
private static final double big10pow[] = {
1e16, 1e32, 1e64, 1e128, 1e256 };
private static final double tiny10pow[] = {
1e-16, 1e-32, 1e-64, 1e-128, 1e-256 };
private static final int maxSmallTen = small10pow.length-1;
private static final int small5pow[] = {
1,
5,
5*5,
5*5*5,
5*5*5*5,
5*5*5*5*5,
5*5*5*5*5*5,
5*5*5*5*5*5*5,
5*5*5*5*5*5*5*5,
5*5*5*5*5*5*5*5*5,
5*5*5*5*5*5*5*5*5*5,
5*5*5*5*5*5*5*5*5*5*5,
5*5*5*5*5*5*5*5*5*5*5*5,
5*5*5*5*5*5*5*5*5*5*5*5*5
};
private static final long long5pow[] = {
1L,
5L,
5L*5,
5L*5*5,
5L*5*5*5,
5L*5*5*5*5,
5L*5*5*5*5*5,
5L*5*5*5*5*5*5,
5L*5*5*5*5*5*5*5,
5L*5*5*5*5*5*5*5*5,
5L*5*5*5*5*5*5*5*5*5,
5L*5*5*5*5*5*5*5*5*5*5,
5L*5*5*5*5*5*5*5*5*5*5*5,
5L*5*5*5*5*5*5*5*5*5*5*5*5,
5L*5*5*5*5*5*5*5*5*5*5*5*5*5,
5L*5*5*5*5*5*5*5*5*5*5*5*5*5*5,
5L*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5,
5L*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5,
5L*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5,
5L*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5,
5L*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5,
5L*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5,
5L*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5,
5L*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5,
5L*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5,
5L*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5,
5L*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5*5,
};
private static final char zero[] = { '0', '0', '0', '0', '0', '0', '0', '0' };
private MyFloatingDecimal( boolean negSign, int decExponent, char []digits, int n, boolean e )
{
isNegative = negSign;
isExceptional = e;
this.decExponent = decExponent;
this.digits = digits;
this.nDigits = n;
}
public static MyFloatingDecimal
readJavaFormatString( String in ) throws NumberFormatException {
boolean isNegative = false;
boolean signSeen = false;
int decExp;
char c;
parseNumber:
try{
in = in.trim(); // don't fool around with white space.
// throws NullPointerException if null
int l = in.length();
if ( l == 0 ) throw new NumberFormatException("empty String");
int i = 0;
switch ( c = in.charAt( i ) ){
case '-':
isNegative = true;
//FALLTHROUGH
case '+':
i++;
signSeen = true;
}
// Would handle NaN and Infinity here, but it isn't
// part of the spec!
//
char[] digits = new char[ l ];
int nDigits= 0;
boolean decSeen = false;
int decPt = 0;
int nLeadZero = 0;
int nTrailZero= 0;
digitLoop:
while ( i < l ){
switch ( c = in.charAt( i ) ){
case '0':
if ( nDigits > 0 ){
nTrailZero += 1;
} else {
nLeadZero += 1;
}
break; // out of switch.
case '1':
case '2':
case '3':
case '4':
case '5':
case '6':
case '7':
case '8':
case '9':
while ( nTrailZero > 0 ){
digits[nDigits++] = '0';
nTrailZero -= 1;
}
digits[nDigits++] = c;
break; // out of switch.
case '.':
if ( decSeen ){
// already saw one ., this is the 2nd.
throw new NumberFormatException("multiple points");
}
decPt = i;
if ( signSeen ){
decPt -= 1;
}
decSeen = true;
break; // out of switch.
default:
break digitLoop;
}
i++;
}
/*
* At this point, we've scanned all the digits and decimal
* point we're going to see. Trim off leading and trailing
* zeros, which will just confuse us later, and adjust
* our initial decimal exponent accordingly.
* To review:
* we have seen i total characters.
* nLeadZero of them were zeros before any other digits.
* nTrailZero of them were zeros after any other digits.
* if ( decSeen ), then a . was seen after decPt characters
* ( including leading zeros which have been discarded )
* nDigits characters were neither lead nor trailing
* zeros, nor point
*/
/*
* special hack: if we saw no non-zero digits, then the
* answer is zero!
* Unfortunately, we feel honor-bound to keep parsing!
*/
if ( nDigits == 0 ){
digits = zero;
nDigits = 1;
if ( nLeadZero == 0 ){
// we saw NO DIGITS AT ALL,
// not even a crummy 0!
// this is not allowed.
break parseNumber; // go throw exception
}
}
/* Our initial exponent is decPt, adjusted by the number of
* discarded zeros. Or, if there was no decPt,
* then its just nDigits adjusted by discarded trailing zeros.
*/
if ( decSeen ){
decExp = decPt - nLeadZero;
} else {
decExp = nDigits+nTrailZero;
}
/*
* Look for 'e' or 'E' and an optionally signed integer.
*/
if ( (i < l) && ((c = in.charAt(i) )=='e') || (c == 'E') ){
int expSign = 1;
int expVal = 0;
int reallyBig = Integer.MAX_VALUE / 10;
boolean expOverflow = false;
switch( in.charAt(++i) ){
case '-':
expSign = -1;
//FALLTHROUGH
case '+':
i++;
}
int expAt = i;
expLoop:
while ( i < l ){
if ( expVal >= reallyBig ){
// the next character will cause integer
// overflow.
expOverflow = true;
}
switch ( c = in.charAt(i++) ){
case '0':
case '1':
case '2':
case '3':
case '4':
case '5':
case '6':
case '7':
case '8':
case '9':
expVal = expVal*10 + ( (int)c - (int)'0' );
continue;
default:
i--; // back up.
break expLoop; // stop parsing exponent.
}
}
int expLimit = bigDecimalExponent+nDigits+nTrailZero;
if ( expOverflow || ( expVal > expLimit ) ){
//
// The intent here is to end up with
// infinity or zero, as appropriate.
// The reason for yielding such a small decExponent,
// rather than something intuitive such as
// expSign*Integer.MAX_VALUE, is that this value
// is subject to further manipulation in
// doubleValue() and floatValue(), and I don't want
// it to be able to cause overflow there!
// (The only way we can get into trouble here is for
// really outrageous nDigits+nTrailZero, such as 2 billion. )
//
decExp = expSign*expLimit;
} else {
// this should not overflow, since we tested
// for expVal > (MAX+N), where N >= abs(decExp)
decExp = decExp + expSign*expVal;
}
// if we saw something not a digit ( or end of string )
// after the [Ee][+-], without seeing any digits at all
// this is certainly an error. If we saw some digits,
// but then some trailing garbage, that might be ok.
// so we just fall through in that case.
// HUMBUG
if ( i == expAt )
break parseNumber; // certainly bad
}
/*
* We parsed everything we could.
* If there are leftovers, then this is not good input!
*/
if ( i < l &&
((i != l - 1) ||
(in.charAt(i) != 'f' &&
in.charAt(i) != 'F' &&
in.charAt(i) != 'd' &&
in.charAt(i) != 'D'))) {
break parseNumber; // go throw exception
}
return new MyFloatingDecimal( isNegative, decExp, digits, nDigits, false );
} catch ( StringIndexOutOfBoundsException e ){ }
throw new NumberFormatException( in );
}
public double
doubleValue(){
int kDigits = Math.min( nDigits, maxDecimalDigits+1 );
long lValue;
double dValue;
double rValue, tValue;
roundDir = 0;
/*
* convert the lead kDigits to a long integer.
*/
// (special performance hack: start to do it using int)
int iValue = (int)digits[0]-(int)'0';
int iDigits = Math.min( kDigits, intDecimalDigits );
for ( int i=1; i < iDigits; i++ ){
iValue = iValue*10 + (int)digits[i]-(int)'0';
}
lValue = (long)iValue;
for ( int i=iDigits; i < kDigits; i++ ){
lValue = lValue*10L + (long)((int)digits[i]-(int)'0');
}
dValue = (double)lValue;
int exp = decExponent-kDigits;
/*
* lValue now contains a long integer with the value of
* the first kDigits digits of the number.
* dValue contains the (double) of the same.
*/
if ( nDigits <= maxDecimalDigits ){
/*
* possibly an easy case.
* We know that the digits can be represented
* exactly. And if the exponent isn't too outrageous,
* the whole thing can be done with one operation,
* thus one rounding error.
* Note that all our constructors trim all leading and
* trailing zeros, so simple values (including zero)
* will always end up here
*/
if (exp == 0 || dValue == 0.0)
return (isNegative)? -dValue : dValue; // small floating integer
else if ( exp >= 0 ){
if ( exp <= maxSmallTen ){
/*
* Can get the answer with one operation,
* thus one roundoff.
*/
rValue = dValue * small10pow[exp];
if ( mustSetRoundDir ){
tValue = rValue / small10pow[exp];
roundDir = ( tValue == dValue ) ? 0
:( tValue < dValue ) ? 1
: -1;
}
return (isNegative)? -rValue : rValue;
}
int slop = maxDecimalDigits - kDigits;
if ( exp <= maxSmallTen+slop ){
/*
* We can multiply dValue by 10^(slop)
* and it is still "small" and exact.
* Then we can multiply by 10^(exp-slop)
* with one rounding.
*/
dValue *= small10pow[slop];
rValue = dValue * small10pow[exp-slop];
if ( mustSetRoundDir ){
tValue = rValue / small10pow[exp-slop];
roundDir = ( tValue == dValue ) ? 0
:( tValue < dValue ) ? 1
: -1;
}
return (isNegative)? -rValue : rValue;
}
/*
* Else we have a hard case with a positive exp.
*/
} else {
if ( exp >= -maxSmallTen ){
/*
* Can get the answer in one division.
*/
rValue = dValue / small10pow[-exp];
tValue = rValue * small10pow[-exp];
if ( mustSetRoundDir ){
roundDir = ( tValue == dValue ) ? 0
:( tValue < dValue ) ? 1
: -1;
}
return (isNegative)? -rValue : rValue;
}
/*
* Else we have a hard case with a negative exp.
*/
}
}
/*
* Harder cases:
* The sum of digits plus exponent is greater than
* what we think we can do with one error.
*
* Start by approximating the right answer by,
* naively, scaling by powers of 10.
*/
if ( exp > 0 ){
if ( decExponent > maxDecimalExponent+1 ){
/*
* Lets face it. This is going to be
* Infinity. Cut to the chase.
*/
return (isNegative)? Double.NEGATIVE_INFINITY : Double.POSITIVE_INFINITY;
}
if ( (exp&15) != 0 ){
dValue *= small10pow[exp&15];
}
if ( (exp>>=4) != 0 ){
int j;
for( j = 0; exp > 1; j++, exp>>=1 ){
if ( (exp&1)!=0)
dValue *= big10pow[j];
}
/*
* The reason for the weird exp > 1 condition
* in the above loop was so that the last multiply
* would get unrolled. We handle it here.
* It could overflow.
*/
double t = dValue * big10pow[j];
if ( Double.isInfinite( t ) ){
/*
* It did overflow.
* Look more closely at the result.
* If the exponent is just one too large,
* then use the maximum finite as our estimate
* value. Else call the result infinity
* and punt it.
* ( I presume this could happen because
* rounding forces the result here to be
* an ULP or two larger than
* Double.MAX_VALUE ).
*/
t = dValue / 2.0;
t *= big10pow[j];
if ( Double.isInfinite( t ) ){
return (isNegative)? Double.NEGATIVE_INFINITY : Double.POSITIVE_INFINITY;
}
t = Double.MAX_VALUE;
}
dValue = t;
}
} else if ( exp < 0 ){
exp = -exp;
if ( decExponent < minDecimalExponent-1 ){
/*
* Lets face it. This is going to be
* zero. Cut to the chase.
*/
return (isNegative)? -0.0 : 0.0;
}
if ( (exp&15) != 0 ){
dValue /= small10pow[exp&15];
}
if ( (exp>>=4) != 0 ){
int j;
for( j = 0; exp > 1; j++, exp>>=1 ){
if ( (exp&1)!=0)
dValue *= tiny10pow[j];
}
/*
* The reason for the weird exp > 1 condition
* in the above loop was so that the last multiply
* would get unrolled. We handle it here.
* It could underflow.
*/
double t = dValue * tiny10pow[j];
if ( t == 0.0 ){
/*
* It did underflow.
* Look more closely at the result.
* If the exponent is just one too small,
* then use the minimum finite as our estimate
* value. Else call the result 0.0
* and punt it.
* ( I presume this could happen because
* rounding forces the result here to be
* an ULP or two less than
* Double.MIN_VALUE ).
*/
t = dValue * 2.0;
t *= tiny10pow[j];
if ( t == 0.0 ){
return (isNegative)? -0.0 : 0.0;
}
t = Double.MIN_VALUE;
}
dValue = t;
}
}
/*
* dValue is now approximately the result.
* The hard part is adjusting it, by comparison
* with FDBigInt arithmetic.
* Formulate the EXACT big-number result as
* bigD0 * 10^exp
*/
FDBigInt bigD0 = new FDBigInt( lValue, digits, kDigits, nDigits );
exp = decExponent - nDigits;
correctionLoop:
while(true){
/* AS A SIDE EFFECT, THIS METHOD WILL SET THE INSTANCE VARIABLES
* bigIntExp and bigIntNBits
*/
FDBigInt bigB = doubleToBigInt( dValue );
/*
* Scale bigD, bigB appropriately for
* big-integer operations.
* Naively, we multiply by powers of ten
* and powers of two. What we actually do
* is keep track of the powers of 5 and
* powers of 2 we would use, then factor out
* common divisors before doing the work.
*/
int B2, B5; // powers of 2, 5 in bigB
int D2, D5; // powers of 2, 5 in bigD
int Ulp2; // powers of 2 in halfUlp.
if ( exp >= 0 ){
B2 = B5 = 0;
D2 = D5 = exp;
} else {
B2 = B5 = -exp;
D2 = D5 = 0;
}
if ( bigIntExp >= 0 ){
B2 += bigIntExp;
} else {
D2 -= bigIntExp;
}
Ulp2 = B2;
// shift bigB and bigD left by a number s. t.
// halfUlp is still an integer.
int hulpbias;
if ( bigIntExp+bigIntNBits <= -expBias+1 ){
// This is going to be a denormalized number
// (if not actually zero).
// half an ULP is at 2^-(expBias+expShift+1)
hulpbias = bigIntExp+ expBias + expShift;
} else {
hulpbias = expShift + 2 - bigIntNBits;
}
B2 += hulpbias;
D2 += hulpbias;
// if there are common factors of 2, we might just as well
// factor them out, as they add nothing useful.
int common2 = Math.min( B2, Math.min( D2, Ulp2 ) );
B2 -= common2;
D2 -= common2;
Ulp2 -= common2;
// do multiplications by powers of 5 and 2
bigB = multPow52( bigB, B5, B2 );
FDBigInt bigD = multPow52( new FDBigInt( bigD0 ), D5, D2 );
//
// to recap:
// bigB is the scaled-big-int version of our floating-point
// candidate.
// bigD is the scaled-big-int version of the exact value
// as we understand it.
// halfUlp is 1/2 an ulp of bigB, except for special cases
// of exact powers of 2
//
// the plan is to compare bigB with bigD, and if the difference
// is less than halfUlp, then we're satisfied. Otherwise,
// use the ratio of difference to halfUlp to calculate a fudge
// factor to add to the floating value, then go 'round again.
//
FDBigInt diff;
int cmpResult;
boolean overvalue;
if ( (cmpResult = bigB.cmp( bigD ) ) > 0 ){
overvalue = true; // our candidate is too big.
diff = bigB.sub( bigD );
if ( (bigIntNBits == 1) && (bigIntExp > -expBias) ){
// candidate is a normalized exact power of 2 and
// is too big. We will be subtracting.
// For our purposes, ulp is the ulp of the
// next smaller range.
Ulp2 -= 1;
if ( Ulp2 < 0 ){
// rats. Cannot de-scale ulp this far.
// must scale diff in other direction.
Ulp2 = 0;
diff.lshiftMe( 1 );
}
}
} else if ( cmpResult < 0 ){
overvalue = false; // our candidate is too small.
diff = bigD.sub( bigB );
} else {
// the candidate is exactly right!
// this happens with surprising frequency
break correctionLoop;
}
FDBigInt halfUlp = constructPow52( B5, Ulp2 );
if ( (cmpResult = diff.cmp( halfUlp ) ) < 0 ){
// difference is small.
// this is close enough
roundDir = overvalue ? -1 : 1;
break correctionLoop;
} else if ( cmpResult == 0 ){
// difference is exactly half an ULP
// round to some other value maybe, then finish
dValue += 0.5*ulp( dValue, overvalue );
// should check for bigIntNBits == 1 here??
roundDir = overvalue ? -1 : 1;
break correctionLoop;
} else {
// difference is non-trivial.
// could scale addend by ratio of difference to
// halfUlp here, if we bothered to compute that difference.
// Most of the time ( I hope ) it is about 1 anyway.
dValue += ulp( dValue, overvalue );
if ( dValue == 0.0 || dValue == Double.POSITIVE_INFINITY )
break correctionLoop; // oops. Fell off end of range.
continue; // try again.
}
}
return (isNegative)? -dValue : dValue;
}
/*
* Make a floating double into a FDBigInt.
* This could also be structured as a FDBigInt
* constructor, but we'd have to build a lot of knowledge
* about floating-point representation into it, and we don't want to.
*
* AS A SIDE EFFECT, THIS METHOD WILL SET THE INSTANCE VARIABLES
* bigIntExp and bigIntNBits
*
*/
private FDBigInt
doubleToBigInt( double dval ){
long lbits = Double.doubleToLongBits( dval ) & ~signMask;
int binexp = (int)(lbits >>> expShift);
lbits &= fractMask;
if ( binexp > 0 ){
lbits |= fractHOB;
} else {
if ( lbits == 0L )
throw new RuntimeException("Assertion botch: doubleToBigInt(0.0)");
binexp +=1;
while ( (lbits & fractHOB ) == 0L){
lbits <<= 1;
binexp -= 1;
}
}
binexp -= expBias;
int nbits = countBits( lbits );
/*
* We now know where the high-order 1 bit is,
* and we know how many there are.
*/
int lowOrderZeros = expShift+1-nbits;
lbits >>>= lowOrderZeros;
bigIntExp = binexp+1-nbits;
bigIntNBits = nbits;
return new FDBigInt( lbits );
}
//
// a common operation
//
private static FDBigInt
multPow52( FDBigInt v, int p5, int p2 ){
if ( p5 != 0 ){
if ( p5 < small5pow.length ){
v = v.mult( small5pow[p5] );
} else {
v = v.mult( big5pow( p5 ) );
}
}
if ( p2 != 0 ){
v.lshiftMe( p2 );
}
return v;
}
//
// another common operation
//
private static FDBigInt
constructPow52( int p5, int p2 ){
FDBigInt v = new FDBigInt( big5pow( p5 ) );
if ( p2 != 0 ){
v.lshiftMe( p2 );
}
return v;
}
/*
* Compute a number that is the ULP of the given value,
* for purposes of addition/subtraction. Generally easy.
* More difficult if subtracting and the argument
* is a normalized a power of 2, as the ULP changes at these points.
*/
private static double
ulp( double dval, boolean subtracting ){
long lbits = Double.doubleToLongBits( dval ) & ~signMask;
int binexp = (int)(lbits >>> expShift);
double ulpval;
if ( subtracting && ( binexp >= expShift ) && ((lbits&fractMask) == 0L) ){
// for subtraction from normalized, powers of 2,
// use next-smaller exponent
binexp -= 1;
}
if ( binexp > expShift ){
ulpval = Double.longBitsToDouble( ((long)(binexp-expShift))<<expShift );
} else if ( binexp == 0 ){
ulpval = Double.MIN_VALUE;
} else {
ulpval = Double.longBitsToDouble( 1L<<(binexp-1) );
}
if ( subtracting ) ulpval = - ulpval;
return ulpval;
}
/*
* count number of bits from high-order 1 bit to low-order 1 bit,
* inclusive.
*/
private static int
countBits( long v ){
//
// the strategy is to shift until we get a non-zero sign bit
// then shift until we have no bits left, counting the difference.
// we do byte shifting as a hack. Hope it helps.
//
if ( v == 0L ) return 0;
while ( ( v & highbyte ) == 0L ){
v <<= 8;
}
while ( v > 0L ) { // i.e. while ((v&highbit) == 0L )
v <<= 1;
}
int n = 0;
while (( v & lowbytes ) != 0L ){
v <<= 8;
n += 8;
}
while ( v != 0L ){
v <<= 1;
n += 1;
}
return n;
}
/*
* Keep big powers of 5 handy for future reference.
*/
private static FDBigInt b5p[];
private static synchronized FDBigInt
big5pow( int p ){
if ( p < 0 )
throw new RuntimeException( "Assertion botch: negative power of 5");
if ( b5p == null ){
b5p = new FDBigInt[ p+1 ];
}else if (b5p.length <= p ){
FDBigInt t[] = new FDBigInt[ p+1 ];
System.arraycopy( b5p, 0, t, 0, b5p.length );
b5p = t;
}
if ( b5p[p] != null )
return b5p[p];
else if ( p < small5pow.length )
return b5p[p] = new FDBigInt( small5pow[p] );
else if ( p < long5pow.length )
return b5p[p] = new FDBigInt( long5pow[p] );
else {
// construct the value.
// recursively.
int q, r;
// in order to compute 5^p,
// compute its square root, 5^(p/2) and square.
// or, let q = p / 2, r = p -q, then
// 5^p = 5^(q+r) = 5^q * 5^r
q = p >> 1;
r = p - q;
FDBigInt bigq = b5p[q];
if ( bigq == null )
bigq = big5pow ( q );
if ( r < small5pow.length ){
return (b5p[p] = bigq.mult( small5pow[r] ) );
}else{
FDBigInt bigr = b5p[ r ];
if ( bigr == null )
bigr = big5pow( r );
return (b5p[p] = bigq.mult( bigr ) );
}
}
}
}

「真诚赞赏,手留余香」