Subversion Repositories planix.SVN

Rev

Rev 33 | Details | Compare with Previous | Last modification | View Log | RSS feed

Rev Author Line No. Line
26 7u83 1 
#include "os.h"
2 
#include <mp.h>
3 
#include "dat.h"
4 
 
5 
double
6 
mptod(mpint *a)
7 
{
8 
	u64int v;
9 
	mpdigit w, r;
10 
	int sf, i, n, m, s;
11 
	FPdbleword x;
12 
 
13 
	if(a->top == 0) return 0.0;
14 
	sf = mpsignif(a);
15 
	if(sf > 1024) return Inf(a->sign);
16 
	i = a->top - 1;
17 
	v = a->p[i];
18 
	n = sf & Dbits - 1;
19 
	n |= n - 1 & Dbits;
20 
	r = 0;
21 
	if(n > 54){
22 
		s = n - 54;
23 
		r = v & (1<<s) - 1;
24 
		v >>= s;
25 
	}
26 
	while(n < 54){
27 
		if(--i < 0)
28 
			w = 0;
29 
		else
30 
			w = a->p[i];
31 
		m = 54 - n;
32 
		if(m > Dbits) m = Dbits;
33 
		s = Dbits - m & Dbits - 1;
34 
		v = v << m | w >> s;
35 
		r = w & (1<<s) - 1;
36 
		n += m;
37 
	}
38 
	if((v & 3) == 1){
39 
		while(--i >= 0)
40 
			r |= a->p[i];
41 
		if(r != 0)
42 
			v++;
43 
	}else
44 
		v++;
45 
	v >>= 1;
46 
	while((v >> 53) != 0){
47 
		v >>= 1;
48 
		if(++sf > 1024)
49 
			return Inf(a->sign);
50 
	}
51 
	x.lo = v;
52 
	x.hi = (u32int)(v >> 32) & (1<<20) - 1 | sf + 1022 << 20 | a->sign & 1<<31;
53 
	return x.x;
54 
}
55 
 
56 
mpint *
57 
dtomp(double d, mpint *a)
58 
{
59 
	FPdbleword x;
60 
	uvlong v;
61 
	int e;
62 
 
63 
	if(a == nil)
64 
		a = mpnew(0);
65 
	x.x = d;
66 
	e = x.hi >> 20 & 2047;
67 
	assert(e != 2047);
68 
	if(e < 1022){
69 
		mpassign(mpzero, a);
70 
		return a;
71 
	}
72 
	v = x.lo | (uvlong)(x.hi & (1<<20) - 1) << 32 | 1ULL<<52;
73 
	if(e < 1075){
74 
		v += (1ULL<<1074 - e) - (~v >> 1075 - e & 1);
75 
		v >>= 1075 - e;
76 
	}
77 
	uvtomp(v, a);
78 
	if(e > 1075)
79 
		mpleft(a, e - 1075, a);
80 
	if((int)x.hi < 0)
81 
		a->sign = -1;
82 
	return a;
83 
}