-
Notifications
You must be signed in to change notification settings - Fork 0
/
My_exp.cpp
162 lines (128 loc) · 4.14 KB
/
My_exp.cpp
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
#include <stdio.h>
#include <inttypes.h>
#include <iostream>
//#include <gmp.h>
//#include <mpfr.h>
#include <fstream>
#include <random>
float* readDataFromFile(const char* filename, int& size) {
std::ifstream file(filename);
if (!file.is_open()) {
std::cerr << "Error opening file " << filename << std::endl;
return nullptr;
}
// Îïðåäåëåíèå ðàçìåðà ìàññèâà äàííûõ
file >> size;
// Âûäåëåíèå ïàìÿòè ïîä äèíàìè÷åñêèé ìàññèâ
float* data = new float[size];
// Ñ÷èòûâàíèå äàííûõ èç ôàéëà â ìàññèâ
for (int i = 0; i < size; ++i) {
file >> data[i];
}
file.close();
return data;
}
void update_exponent(float& num, int k) {
// Ïåðâûì äåëîì, ïîëó÷èì ïðåäñòàâëåíèå ÷èñëà â âèäå öåëîãî äëÿ äîñòóïà ê ýêñïîíåíòå
uint32_t* num_as_int = reinterpret_cast<uint32_t*>(&num);
// Èçâëå÷åíèå ýêñïîíåíòû
uint32_t exponent = (*num_as_int & 0x7F800000) >> 23;
// Îáíîâëåíèå ýêñïîíåíòû
exponent += k;
// Ïðîâåðêà íà ïåðåïîëíåíèå
/* if (exponent >= 255)
exponent = 255;
else if (exponent <= 0)
exponent = 0;
*/
// Óñòàíîâêà îáíîâëåííîé ýêñïîíåíòû
*num_as_int = (*num_as_int & 0x807FFFFFu) | (exponent << 23);
}
typedef union {
float f;
uint16_t i16[2];//bigendian xranenieÔ
uint32_t i32;
} binary32;
#define HI 1
#define LO 0
#define FLOAT2INT(_ri, _rf, x) \
{ binary32 t; \
t.f = (x + 12582912.0f); \
_rf = (t.f - 12582912.0f);\
_ri = t.i16[LO]; }
void float_to_int(float& x, float& rf, uint32_t& fint)
{
float tmp = (x + 12582912.0f);
rf = tmp - 12582912.0f;
fint = reinterpret_cast<uint32_t&>(tmp);
fint = (fint & 0xFFFF);
}
float my_exp(float x) {
float y, kf;
uint32_t ki;
float Log2 = (float)0x1.62e43p-1;
float Log2h = (float)0xb.17200p-4;
float Log2l = (float)0x1.7f7d1cf8p-20;
float InvLog2 = (float)0x1.715476p0;
// Here should be the tests for exceptional cases
float X = x * InvLog2;
float_to_int(X, kf, ki);
//y = x - kf * Log2;
y = (x - kf * Log2h) - kf * Log2l;
float result = (float)0x1p0 + y * ((float)0x1p0
+ y * ((float)0x1.fffff8p-2
+ y * ((float)0x1.55548ep-3
+ y * ((float)0x1.555b98p-5
+ y * ((float)0x1.123bccp-7
+ y * (float)0x1.6850e4p-10)))));
uint32_t resint = reinterpret_cast<float&>(result);
update_exponent(result, ki);
return result;
}
float my_exp0(float x) {
float y, kf;
int16_t k;
binary32 r;
float Log2 = (float)0x1.62e43p-1;
float Log2h = (float)0xb.17200p-4;
float Log2l = (float)0x1.7f7d1cf8p-20;
float InvLog2 = (float)0x1.715476p0;
// Here should be the tests for exceptional cases
FLOAT2INT(k, kf, x * InvLog2);
//y = x - kf * Log2;
y = (x - kf * Log2h) - kf * Log2l;
r.f = (float)0x1p0 + y * ((float)0x1p0
+ y * ((float)0x1.fffff8p-2
+ y * ((float)0x1.55548ep-3
+ y * ((float)0x1.555b98p-5
+ y * ((float)0x1.123bccp-7
+ y * (float)0x1.6850e4p-10)))));
std::cout << "r.f = " << r.f << std::endl;
r.i16[HI] += k << 7; //Exponent update
return r.f;
}
int main() {
binary32 ref, r;
float x;
int32_t diff, max_diff = 0;
r.f = (float)0x2.c5c85fdf473dep8;//0x2.c5c85fp3f; // Íà÷àëüíîå çíà÷åíèå äëÿ ãåíåðàöèè ñëó÷àéíûõ ÷èñåë
ref.f = 0x1.0p-10f; // Îãðàíè÷èâàþùåå çíà÷åíèå
std::random_device rd;
std::mt19937 gen(rd());
std::uniform_int_distribution<uint32_t> distrib(ref.i32, r.i32); // Ðàâíîìåðíîå ðàñïðåäåëåíèå
for (int i = 0; i < 10; i++) {
// r.i32 = distrib(gen);
x = r.f;
r.f = my_exp0(x);
ref.f = exp(x);
diff = r.i32 - ref.i32; // Ñðàâíåíèå ðåçóëüòàòîâ
if (abs(diff) > max_diff) {
max_diff = abs(diff);
}
// printf("x=%.12ef\nmy_exp=%.12ef\nexpf=%.12ef\ndiff=%f max_diff=%f\n",
//x, r.f, ref.f, diff, max_diff);
printf("x=%.12ef\nmy_exp=%.12ef\nexpf=%.12ef\ndiff=%i max_diff=%i\n",
x, r.f, ref.f, diff, max_diff);
}
return 0;
}