1 /*
2 * Copyright (c) 2017 Andrew Kelley
3 *
4 * This file is part of zig, which is MIT licensed.
5 * See http://opensource.org/licenses/MIT
6 */
7
8 #include "bigfloat.hpp"
9 #include "bigint.hpp"
10 #include "buffer.hpp"
11 #include "softfloat.hpp"
12 #include "softfloat_ext.hpp"
13 #include "parse_f128.h"
14 #include <stdio.h>
15 #include <math.h>
16 #include <errno.h>
17
18
bigfloat_init_128(BigFloat * dest,float128_t x)19 void bigfloat_init_128(BigFloat *dest, float128_t x) {
20 dest->value = x;
21 }
22
bigfloat_init_16(BigFloat * dest,float16_t x)23 void bigfloat_init_16(BigFloat *dest, float16_t x) {
24 f16_to_f128M(x, &dest->value);
25 }
26
bigfloat_init_32(BigFloat * dest,float x)27 void bigfloat_init_32(BigFloat *dest, float x) {
28 float32_t f32_val;
29 memcpy(&f32_val, &x, sizeof(float));
30 f32_to_f128M(f32_val, &dest->value);
31 }
32
bigfloat_init_64(BigFloat * dest,double x)33 void bigfloat_init_64(BigFloat *dest, double x) {
34 float64_t f64_val;
35 memcpy(&f64_val, &x, sizeof(double));
36 f64_to_f128M(f64_val, &dest->value);
37 }
38
bigfloat_init_bigfloat(BigFloat * dest,const BigFloat * x)39 void bigfloat_init_bigfloat(BigFloat *dest, const BigFloat *x) {
40 memcpy(&dest->value, &x->value, sizeof(float128_t));
41 }
42
bigfloat_init_bigint(BigFloat * dest,const BigInt * op)43 void bigfloat_init_bigint(BigFloat *dest, const BigInt *op) {
44 ui32_to_f128M(0, &dest->value);
45 if (op->digit_count == 0)
46 return;
47
48 float128_t base;
49 ui64_to_f128M(UINT64_MAX, &base);
50 float128_t one_f128;
51 ui32_to_f128M(1, &one_f128);
52 f128M_add(&base, &one_f128, &base);
53
54 const uint64_t *digits = bigint_ptr(op);
55
56 for (size_t i = op->digit_count - 1;;) {
57 float128_t digit_f128;
58 ui64_to_f128M(digits[i], &digit_f128);
59
60 f128M_mulAdd(&dest->value, &base, &digit_f128, &dest->value);
61
62 if (i == 0) {
63 if (op->is_negative) {
64 f128M_neg(&dest->value, &dest->value);
65 }
66 return;
67 }
68 i -= 1;
69 }
70 }
71
bigfloat_init_buf(BigFloat * dest,const uint8_t * buf_ptr)72 Error bigfloat_init_buf(BigFloat *dest, const uint8_t *buf_ptr) {
73 char *str_begin = (char *)buf_ptr;
74 char *str_end;
75
76 errno = 0;
77 dest->value = parse_f128(str_begin, &str_end);
78 if (errno) {
79 return ErrorOverflow;
80 }
81
82 return ErrorNone;
83 }
84
bigfloat_add(BigFloat * dest,const BigFloat * op1,const BigFloat * op2)85 void bigfloat_add(BigFloat *dest, const BigFloat *op1, const BigFloat *op2) {
86 f128M_add(&op1->value, &op2->value, &dest->value);
87 }
88
bigfloat_negate(BigFloat * dest,const BigFloat * op)89 void bigfloat_negate(BigFloat *dest, const BigFloat *op) {
90 f128M_neg(&op->value, &dest->value);
91 }
92
bigfloat_sub(BigFloat * dest,const BigFloat * op1,const BigFloat * op2)93 void bigfloat_sub(BigFloat *dest, const BigFloat *op1, const BigFloat *op2) {
94 f128M_sub(&op1->value, &op2->value, &dest->value);
95 }
96
bigfloat_mul(BigFloat * dest,const BigFloat * op1,const BigFloat * op2)97 void bigfloat_mul(BigFloat *dest, const BigFloat *op1, const BigFloat *op2) {
98 f128M_mul(&op1->value, &op2->value, &dest->value);
99 }
100
bigfloat_div(BigFloat * dest,const BigFloat * op1,const BigFloat * op2)101 void bigfloat_div(BigFloat *dest, const BigFloat *op1, const BigFloat *op2) {
102 f128M_div(&op1->value, &op2->value, &dest->value);
103 }
104
bigfloat_div_trunc(BigFloat * dest,const BigFloat * op1,const BigFloat * op2)105 void bigfloat_div_trunc(BigFloat *dest, const BigFloat *op1, const BigFloat *op2) {
106 f128M_div(&op1->value, &op2->value, &dest->value);
107 f128M_roundToInt(&dest->value, softfloat_round_minMag, false, &dest->value);
108 }
109
bigfloat_div_floor(BigFloat * dest,const BigFloat * op1,const BigFloat * op2)110 void bigfloat_div_floor(BigFloat *dest, const BigFloat *op1, const BigFloat *op2) {
111 f128M_div(&op1->value, &op2->value, &dest->value);
112 f128M_roundToInt(&dest->value, softfloat_round_min, false, &dest->value);
113 }
114
bigfloat_rem(BigFloat * dest,const BigFloat * op1,const BigFloat * op2)115 void bigfloat_rem(BigFloat *dest, const BigFloat *op1, const BigFloat *op2) {
116 f128M_rem(&op1->value, &op2->value, &dest->value);
117 }
118
bigfloat_mod(BigFloat * dest,const BigFloat * op1,const BigFloat * op2)119 void bigfloat_mod(BigFloat *dest, const BigFloat *op1, const BigFloat *op2) {
120 f128M_rem(&op1->value, &op2->value, &dest->value);
121 f128M_add(&dest->value, &op2->value, &dest->value);
122 f128M_rem(&dest->value, &op2->value, &dest->value);
123 }
124
bigfloat_append_buf(Buf * buf,const BigFloat * op)125 void bigfloat_append_buf(Buf *buf, const BigFloat *op) {
126 const size_t extra_len = 100;
127 size_t old_len = buf_len(buf);
128 buf_resize(buf, old_len + extra_len);
129
130 // TODO actually print f128
131 float64_t f64_value = f128M_to_f64(&op->value);
132 double double_value;
133 memcpy(&double_value, &f64_value, sizeof(double));
134
135 int len = snprintf(buf_ptr(buf) + old_len, extra_len, "%f", double_value);
136 assert(len > 0);
137 buf_resize(buf, old_len + len);
138 }
139
bigfloat_cmp(const BigFloat * op1,const BigFloat * op2)140 Cmp bigfloat_cmp(const BigFloat *op1, const BigFloat *op2) {
141 if (f128M_lt(&op1->value, &op2->value)) {
142 return CmpLT;
143 } else if (f128M_eq(&op1->value, &op2->value)) {
144 return CmpEQ;
145 } else {
146 return CmpGT;
147 }
148 }
149
bigfloat_to_f16(const BigFloat * bigfloat)150 float16_t bigfloat_to_f16(const BigFloat *bigfloat) {
151 return f128M_to_f16(&bigfloat->value);
152 }
153
bigfloat_to_f32(const BigFloat * bigfloat)154 float bigfloat_to_f32(const BigFloat *bigfloat) {
155 float32_t f32_value = f128M_to_f32(&bigfloat->value);
156 float result;
157 memcpy(&result, &f32_value, sizeof(float));
158 return result;
159 }
160
bigfloat_to_f64(const BigFloat * bigfloat)161 double bigfloat_to_f64(const BigFloat *bigfloat) {
162 float64_t f64_value = f128M_to_f64(&bigfloat->value);
163 double result;
164 memcpy(&result, &f64_value, sizeof(double));
165 return result;
166 }
167
bigfloat_to_f128(const BigFloat * bigfloat)168 float128_t bigfloat_to_f128(const BigFloat *bigfloat) {
169 return bigfloat->value;
170 }
171
bigfloat_cmp_zero(const BigFloat * bigfloat)172 Cmp bigfloat_cmp_zero(const BigFloat *bigfloat) {
173 float128_t zero_float;
174 ui32_to_f128M(0, &zero_float);
175 if (f128M_lt(&bigfloat->value, &zero_float)) {
176 return CmpLT;
177 } else if (f128M_eq(&bigfloat->value, &zero_float)) {
178 return CmpEQ;
179 } else {
180 return CmpGT;
181 }
182 }
183
bigfloat_has_fraction(const BigFloat * bigfloat)184 bool bigfloat_has_fraction(const BigFloat *bigfloat) {
185 float128_t floored;
186 f128M_roundToInt(&bigfloat->value, softfloat_round_minMag, false, &floored);
187 return !f128M_eq(&floored, &bigfloat->value);
188 }
189
bigfloat_sqrt(BigFloat * dest,const BigFloat * op)190 void bigfloat_sqrt(BigFloat *dest, const BigFloat *op) {
191 f128M_sqrt(&op->value, &dest->value);
192 }
193
bigfloat_min(BigFloat * dest,const BigFloat * op1,const BigFloat * op2)194 void bigfloat_min(BigFloat *dest, const BigFloat *op1, const BigFloat *op2) {
195 if (bigfloat_is_nan(op1)) {
196 bigfloat_init_bigfloat(dest, op2);
197 } else if (bigfloat_is_nan(op2)) {
198 bigfloat_init_bigfloat(dest, op1);
199 } else if (f128M_lt(&op1->value, &op2->value)) {
200 bigfloat_init_bigfloat(dest, op1);
201 } else {
202 bigfloat_init_bigfloat(dest, op2);
203 }
204 }
205
bigfloat_max(BigFloat * dest,const BigFloat * op1,const BigFloat * op2)206 void bigfloat_max(BigFloat *dest, const BigFloat *op1, const BigFloat *op2) {
207 if (bigfloat_is_nan(op1)) {
208 bigfloat_init_bigfloat(dest, op2);
209 } else if (bigfloat_is_nan(op2)) {
210 bigfloat_init_bigfloat(dest, op1);
211 } else if (f128M_lt(&op1->value, &op2->value)) {
212 bigfloat_init_bigfloat(dest, op2);
213 } else {
214 bigfloat_init_bigfloat(dest, op1);
215 }
216 }
217
bigfloat_is_nan(const BigFloat * op)218 bool bigfloat_is_nan(const BigFloat *op) {
219 return f128M_isSignalingNaN(&op->value);
220 }
221