xref: /reactos/base/applications/calc/fun_mpfr.c (revision 682f85ad)
1 /*
2  * ReactOS Calc (Math functions, GMP/MPFR engine)
3  *
4  * Copyright 2007-2017, Carlo Bramini
5  *
6  * This program is free software; you can redistribute it and/or
7  * modify it under the terms of the GNU Lesser General Public
8  * License as published by the Free Software Foundation; either
9  * version 2 of the License, or (at your option) any later version.
10  *
11  * This library is distributed in the hope that it will be useful,
12  * but WITHOUT ANY WARRANTY; without even the implied warranty of
13  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
14  * Lesser General Public License for more details.
15  *
16  * You should have received a copy of the GNU Lesser General Public
17  * License along with this library; if not, write to the Free Software
18  * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301  USA
19  */
20 
21 #include "calc.h"
22 
23 static void validate_rad2angle(calc_number_t *c);
24 static void validate_angle2rad(calc_number_t *c);
25 
26 void apply_int_mask(calc_number_t *r)
27 {
28     mpz_t a, mask;
29 
30     switch (calc.size) {
31     case IDC_RADIO_QWORD:
32         mpz_init_set_str(mask, "FFFFFFFFFFFFFFFF", 16);
33         break;
34     case IDC_RADIO_DWORD:
35         mpz_init_set_str(mask, "00000000FFFFFFFF", 16);
36         break;
37     case IDC_RADIO_WORD:
38         mpz_init_set_str(mask, "000000000000FFFF", 16);
39         break;
40     case IDC_RADIO_BYTE:
41         mpz_init_set_str(mask, "00000000000000FF", 16);
42         break;
43     default:
44         mpz_init_set_si(mask, -1);
45     }
46     mpz_init(a);
47     mpfr_get_z(a, r->mf, MPFR_DEFAULT_RND);
48     mpz_and(a, a, mask);
49     mpfr_set_z(r->mf, a, MPFR_DEFAULT_RND);
50     mpz_clear(a);
51     mpz_clear(mask);
52 }
53 
54 static void validate_rad2angle(calc_number_t *r)
55 {
56     mpfr_t mult, divs;
57 
58     mpfr_init(mult);
59     mpfr_init(divs);
60     switch (calc.degr) {
61     case IDC_RADIO_DEG:
62         mpfr_set_ui(mult, 180, MPFR_DEFAULT_RND);
63         mpfr_const_pi(divs, MPFR_DEFAULT_RND);
64         break;
65     case IDC_RADIO_RAD:
66         mpfr_set_ui(mult, 1, MPFR_DEFAULT_RND);
67         mpfr_set_ui(divs, 1, MPFR_DEFAULT_RND);
68         break;
69     case IDC_RADIO_GRAD:
70         mpfr_set_ui(mult, 200, MPFR_DEFAULT_RND);
71         mpfr_const_pi(divs, MPFR_DEFAULT_RND);
72         break;
73     }
74     mpfr_mul(r->mf, r->mf, mult, MPFR_DEFAULT_RND);
75     mpfr_div(r->mf, r->mf, divs, MPFR_DEFAULT_RND);
76 
77     mpfr_clear(mult);
78     mpfr_clear(divs);
79 }
80 
81 static void validate_angle2rad(calc_number_t *r)
82 {
83     mpfr_t mult, divs;
84 
85     if (!mpfr_number_p(r->mf)) {
86         calc.is_nan = TRUE;
87         return;
88     }
89     mpfr_init(mult);
90     mpfr_init(divs);
91     switch (calc.degr) {
92     case IDC_RADIO_DEG:
93         mpfr_const_pi(mult, MPFR_DEFAULT_RND);
94         mpfr_set_ui(divs, 180, MPFR_DEFAULT_RND);
95         break;
96     case IDC_RADIO_RAD:
97         mpfr_set_ui(mult, 1, MPFR_DEFAULT_RND);
98         mpfr_set_ui(divs, 1, MPFR_DEFAULT_RND);
99         break;
100     case IDC_RADIO_GRAD:
101         mpfr_const_pi(mult, MPFR_DEFAULT_RND);
102         mpfr_set_ui(divs, 200, MPFR_DEFAULT_RND);
103         break;
104     }
105     mpfr_mul(r->mf, r->mf, mult, MPFR_DEFAULT_RND);
106     mpfr_div(r->mf, r->mf, divs, MPFR_DEFAULT_RND);
107 
108     mpfr_clear(mult);
109     mpfr_clear(divs);
110 }
111 
112 static void build_rad_const(
113     mpfr_t *mp_pi,
114     mpfr_t *mp_pi_2,
115     mpfr_t *mp_3_pi_2,
116     mpfr_t *mp_2_pi)
117 {
118     mpfr_init(*mp_pi);
119     mpfr_init(*mp_pi_2);
120     mpfr_init(*mp_3_pi_2);
121     mpfr_init(*mp_2_pi);
122     mpfr_const_pi(*mp_pi, MPFR_DEFAULT_RND);
123     mpfr_div_ui(*mp_pi_2, *mp_pi, 2, MPFR_DEFAULT_RND);
124     mpfr_mul_ui(*mp_3_pi_2, *mp_pi, 3, MPFR_DEFAULT_RND);
125     mpfr_div_ui(*mp_3_pi_2, *mp_3_pi_2, 2, MPFR_DEFAULT_RND);
126     mpfr_mul_ui(*mp_2_pi, *mp_pi, 2, MPFR_DEFAULT_RND);
127 }
128 
129 void rpn_sin(calc_number_t *c)
130 {
131     mpfr_t mp_pi, mp_pi_2, mp_3_pi_2, mp_2_pi;
132 
133     validate_angle2rad(c);
134     build_rad_const(&mp_pi, &mp_pi_2, &mp_3_pi_2, &mp_2_pi);
135 
136     if (rpn_is_zero(c) || !mpfr_cmp(c->mf, mp_pi) || !mpfr_cmp(c->mf, mp_2_pi))
137         rpn_zero(c);
138     else
139     if (!mpfr_cmp(c->mf, mp_3_pi_2))
140         mpfr_set_si(c->mf, -1, MPFR_DEFAULT_RND);
141     else
142     if (!mpfr_cmp(c->mf, mp_pi_2))
143         mpfr_set_si(c->mf, 1, MPFR_DEFAULT_RND);
144     else {
145         mpfr_sin(c->mf, c->mf, MPFR_DEFAULT_RND);
146         if (!mpfr_number_p(c->mf)) calc.is_nan = TRUE;
147     }
148     mpfr_clear(mp_pi);
149     mpfr_clear(mp_pi_2);
150     mpfr_clear(mp_3_pi_2);
151     mpfr_clear(mp_2_pi);
152 }
153 void rpn_cos(calc_number_t *c)
154 {
155     mpfr_t mp_pi, mp_pi_2, mp_3_pi_2, mp_2_pi;
156 
157     validate_angle2rad(c);
158     build_rad_const(&mp_pi, &mp_pi_2, &mp_3_pi_2, &mp_2_pi);
159 
160     if (!mpfr_cmp(c->mf, mp_pi_2) || !mpfr_cmp(c->mf, mp_3_pi_2))
161         rpn_zero(c);
162     else
163     if (!mpfr_cmp(c->mf, mp_pi))
164         mpfr_set_si(c->mf, -1, MPFR_DEFAULT_RND);
165     else
166     if (!mpfr_cmp(c->mf, mp_2_pi))
167         mpfr_set_si(c->mf, 1, MPFR_DEFAULT_RND);
168     else {
169         mpfr_cos(c->mf, c->mf, MPFR_DEFAULT_RND);
170         if (!mpfr_number_p(c->mf)) calc.is_nan = TRUE;
171     }
172     mpfr_clear(mp_pi);
173     mpfr_clear(mp_pi_2);
174     mpfr_clear(mp_3_pi_2);
175     mpfr_clear(mp_2_pi);
176 }
177 void rpn_tan(calc_number_t *c)
178 {
179     mpfr_t mp_pi, mp_pi_2, mp_3_pi_2, mp_2_pi;
180 
181     validate_angle2rad(c);
182     build_rad_const(&mp_pi, &mp_pi_2, &mp_3_pi_2, &mp_2_pi);
183 
184     if (!mpfr_cmp(c->mf, mp_pi_2) || !mpfr_cmp(c->mf, mp_3_pi_2))
185         calc.is_nan = TRUE;
186     else
187     if (!mpfr_cmp(c->mf, mp_pi) || !mpfr_cmp(c->mf, mp_2_pi))
188         rpn_zero(c);
189     else {
190         mpfr_tan(c->mf, c->mf, MPFR_DEFAULT_RND);
191         if (!mpfr_number_p(c->mf)) calc.is_nan = TRUE;
192     }
193     mpfr_clear(mp_pi);
194     mpfr_clear(mp_pi_2);
195     mpfr_clear(mp_3_pi_2);
196     mpfr_clear(mp_2_pi);
197 }
198 
199 void rpn_asin(calc_number_t *c)
200 {
201     mpfr_asin(c->mf, c->mf, MPFR_DEFAULT_RND);
202     validate_rad2angle(c);
203 }
204 void rpn_acos(calc_number_t *c)
205 {
206     mpfr_acos(c->mf, c->mf, MPFR_DEFAULT_RND);
207     validate_rad2angle(c);
208 }
209 void rpn_atan(calc_number_t *c)
210 {
211     mpfr_atan(c->mf, c->mf, MPFR_DEFAULT_RND);
212     validate_rad2angle(c);
213 }
214 
215 void rpn_sinh(calc_number_t *c)
216 {
217     mpfr_sinh(c->mf, c->mf, MPFR_DEFAULT_RND);
218     if (!mpfr_number_p(c->mf)) calc.is_nan = TRUE;
219 }
220 void rpn_cosh(calc_number_t *c)
221 {
222     mpfr_cosh(c->mf, c->mf, MPFR_DEFAULT_RND);
223     if (!mpfr_number_p(c->mf)) calc.is_nan = TRUE;
224 }
225 void rpn_tanh(calc_number_t *c)
226 {
227     mpfr_tanh(c->mf, c->mf, MPFR_DEFAULT_RND);
228     if (!mpfr_number_p(c->mf)) calc.is_nan = TRUE;
229 }
230 
231 void rpn_asinh(calc_number_t *c)
232 {
233     mpfr_asinh(c->mf, c->mf, MPFR_DEFAULT_RND);
234     if (!mpfr_number_p(c->mf)) calc.is_nan = TRUE;
235 }
236 void rpn_acosh(calc_number_t *c)
237 {
238     mpfr_acosh(c->mf, c->mf, MPFR_DEFAULT_RND);
239     if (!mpfr_number_p(c->mf)) calc.is_nan = TRUE;
240 }
241 void rpn_atanh(calc_number_t *c)
242 {
243     mpfr_atanh(c->mf, c->mf, MPFR_DEFAULT_RND);
244     if (!mpfr_number_p(c->mf)) calc.is_nan = TRUE;
245 }
246 
247 void rpn_int(calc_number_t *c)
248 {
249     mpfr_trunc(c->mf, c->mf);
250 }
251 
252 void rpn_frac(calc_number_t *c)
253 {
254     mpfr_frac(c->mf, c->mf, MPFR_DEFAULT_RND);
255 }
256 
257 void rpn_reci(calc_number_t *c)
258 {
259     if (mpfr_sgn(c->mf) == 0)
260         calc.is_nan = TRUE;
261     else
262         mpfr_ui_div(c->mf, 1, c->mf, MPFR_DEFAULT_RND);
263 }
264 
265 void rpn_fact(calc_number_t *c)
266 {
267     if (mpfr_sgn(c->mf) < 0) {
268         calc.is_nan = TRUE;
269         return;
270     }
271 
272     mpfr_trunc(c->mf, c->mf);
273     if (mpfr_fits_ulong_p(c->mf, MPFR_DEFAULT_RND) == 0)
274         calc.is_nan = TRUE;
275     else {
276         mpfr_fac_ui(c->mf, mpfr_get_ui(c->mf, MPFR_DEFAULT_RND), MPFR_DEFAULT_RND);
277         if (!mpfr_number_p(c->mf)) calc.is_nan = TRUE;
278     }
279 }
280 
281 void rpn_not(calc_number_t *c)
282 {
283     mpz_t a;
284 
285     mpz_init(a);
286     mpfr_get_z(a, c->mf, MPFR_DEFAULT_RND);
287     mpz_com(a, a);
288     mpfr_set_z(c->mf, a, MPFR_DEFAULT_RND);
289     mpz_clear(a);
290 }
291 
292 void rpn_pi(calc_number_t *c)
293 {
294     mpfr_const_pi(c->mf, MPFR_DEFAULT_RND);
295 }
296 
297 void rpn_2pi(calc_number_t *c)
298 {
299     mpfr_const_pi(c->mf, MPFR_DEFAULT_RND);
300     mpfr_mul_ui(c->mf, c->mf, 2, MPFR_DEFAULT_RND);
301 }
302 
303 void rpn_sign(calc_number_t *c)
304 {
305     mpfr_mul_si(c->mf, c->mf, -1, MPFR_DEFAULT_RND);
306 }
307 
308 void rpn_exp2(calc_number_t *c)
309 {
310     mpfr_sqr(c->mf, c->mf, MPFR_DEFAULT_RND);
311     if (!mpfr_number_p(c->mf)) calc.is_nan = TRUE;
312 }
313 
314 void rpn_exp3(calc_number_t *c)
315 {
316     mpfr_pow_ui(c->mf, c->mf, 3, MPFR_DEFAULT_RND);
317     if (!mpfr_number_p(c->mf)) calc.is_nan = TRUE;
318 }
319 
320 void rpn_sqrt(calc_number_t *c)
321 {
322     mpfr_sqrt(c->mf, c->mf, MPFR_DEFAULT_RND);
323     if (!mpfr_number_p(c->mf)) calc.is_nan = TRUE;
324 }
325 
326 void rpn_cbrt(calc_number_t *c)
327 {
328     mpfr_cbrt(c->mf, c->mf, MPFR_DEFAULT_RND);
329     if (!mpfr_number_p(c->mf)) calc.is_nan = TRUE;
330 }
331 
332 void rpn_exp(calc_number_t *c)
333 {
334     mpfr_exp(c->mf, c->mf, MPFR_DEFAULT_RND);
335     if (!mpfr_number_p(c->mf)) calc.is_nan = TRUE;
336 }
337 
338 void rpn_exp10(calc_number_t *c)
339 {
340     mpfr_exp10(c->mf, c->mf, MPFR_DEFAULT_RND);
341     if (!mpfr_number_p(c->mf)) calc.is_nan = TRUE;
342 }
343 
344 void rpn_ln(calc_number_t *c)
345 {
346     mpfr_log(c->mf, c->mf, MPFR_DEFAULT_RND);
347     if (!mpfr_number_p(c->mf)) calc.is_nan = TRUE;
348 }
349 
350 void rpn_log(calc_number_t *c)
351 {
352     mpfr_log10(c->mf, c->mf, MPFR_DEFAULT_RND);
353     if (!mpfr_number_p(c->mf)) calc.is_nan = TRUE;
354 }
355 
356 static void stat_sum(mpfr_t sum)
357 {
358     statistic_t *p = calc.stat;
359 
360     mpfr_set_ui(sum, 0, MPFR_DEFAULT_RND);
361     while (p != NULL) {
362         mpfr_add(sum, sum, p->num.mf, MPFR_DEFAULT_RND);
363         p = (statistic_t *)(p->next);
364     }
365 }
366 
367 static void stat_sum2(mpfr_t sum)
368 {
369     statistic_t *p = calc.stat;
370     mpfr_t       sqr;
371 
372     mpfr_init(sqr);
373     mpfr_set_ui(sum, 0, MPFR_DEFAULT_RND);
374     while (p != NULL) {
375         mpfr_mul(sqr, p->num.mf, p->num.mf, MPFR_DEFAULT_RND);
376         mpfr_add(sum, sum, sqr, MPFR_DEFAULT_RND);
377         p = (statistic_t *)(p->next);
378     }
379     mpfr_clear(sqr);
380 }
381 
382 void rpn_ave(calc_number_t *c)
383 {
384     int     n;
385 
386     stat_sum(c->mf);
387     n = SendDlgItemMessage(calc.hStatWnd, IDC_LIST_STAT, LB_GETCOUNT, 0, 0);
388 
389     if (n)
390         mpfr_div_ui(c->mf, c->mf, n, MPFR_DEFAULT_RND);
391 
392     if (calc.base != IDC_RADIO_DEC)
393         mpfr_trunc(c->mf, c->mf);
394 }
395 
396 void rpn_ave2(calc_number_t *c)
397 {
398     int     n;
399 
400     stat_sum2(c->mf);
401     n = SendDlgItemMessage(calc.hStatWnd, IDC_LIST_STAT, LB_GETCOUNT, 0, 0);
402 
403     if (n)
404         mpfr_div_ui(c->mf, c->mf, n, MPFR_DEFAULT_RND);
405 
406     if (calc.base != IDC_RADIO_DEC)
407         mpfr_trunc(c->mf, c->mf);
408 }
409 
410 void rpn_sum(calc_number_t *c)
411 {
412     stat_sum(c->mf);
413 
414     if (calc.base != IDC_RADIO_DEC)
415         mpfr_trunc(c->mf, c->mf);
416 }
417 
418 void rpn_sum2(calc_number_t *c)
419 {
420     stat_sum2(c->mf);
421 
422     if (calc.base != IDC_RADIO_DEC)
423         mpfr_trunc(c->mf, c->mf);
424 }
425 
426 static void rpn_s_ex(calc_number_t *c, int pop_type)
427 {
428     mpfr_t        dev;
429     mpfr_t        num;
430     unsigned long n = 0;
431     statistic_t  *p = calc.stat;
432 
433     n = SendDlgItemMessage(calc.hStatWnd, IDC_LIST_STAT, LB_GETCOUNT, 0, 0);
434     if (n < 2) {
435         mpfr_set_ui(c->mf, 0, MPFR_DEFAULT_RND);
436         return;
437     }
438 
439     stat_sum(c->mf);
440     mpfr_div_ui(c->mf, c->mf, n, MPFR_DEFAULT_RND);
441 
442     mpfr_init(dev);
443     mpfr_init(num);
444 
445     mpfr_set_ui(dev, 0, MPFR_DEFAULT_RND);
446     p = calc.stat;
447     while (p != NULL) {
448         mpfr_sub(num, p->num.mf, c->mf, MPFR_DEFAULT_RND);
449         mpfr_sqr(num, num, MPFR_DEFAULT_RND);
450         mpfr_add(dev, dev, num, MPFR_DEFAULT_RND);
451         p = (statistic_t *)(p->next);
452     }
453     mpfr_div_ui(c->mf, dev, pop_type ? n-1 : n, MPFR_DEFAULT_RND);
454     mpfr_sqrt(c->mf, c->mf, MPFR_DEFAULT_RND);
455 
456     if (calc.base != IDC_RADIO_DEC)
457         mpfr_trunc(c->mf, c->mf);
458 
459     mpfr_clear(dev);
460     mpfr_clear(num);
461 }
462 
463 void rpn_s(calc_number_t *c)
464 {
465     rpn_s_ex(c, 0);
466 }
467 
468 void rpn_s_m1(calc_number_t *c)
469 {
470     rpn_s_ex(c, 1);
471 }
472 
473 void rpn_dms2dec(calc_number_t *c)
474 {
475     mpfr_t d, m, s;
476 
477     mpfr_init(d);
478     mpfr_init(m);
479     mpfr_init(s);
480 
481     mpfr_trunc(d, c->mf);
482     mpfr_frac(m, c->mf, MPFR_DEFAULT_RND);
483     mpfr_mul_ui(m, m, 100, MPFR_DEFAULT_RND);
484 
485     mpfr_frac(s, m, MPFR_DEFAULT_RND);
486     mpfr_trunc(m, m);
487     mpfr_mul_ui(s, s, 100, MPFR_DEFAULT_RND);
488     mpfr_ceil(s, s);
489 
490     mpfr_div_ui(m, m, 60, MPFR_DEFAULT_RND);
491     mpfr_div_ui(s, s, 3600, MPFR_DEFAULT_RND);
492     mpfr_add(c->mf, d, m, MPFR_DEFAULT_RND);
493     mpfr_add(c->mf, c->mf, s, MPFR_DEFAULT_RND);
494 
495     mpfr_clear(d);
496     mpfr_clear(m);
497     mpfr_clear(s);
498 }
499 
500 void rpn_dec2dms(calc_number_t *c)
501 {
502     mpfr_t d, m, s;
503 
504     mpfr_init(d);
505     mpfr_init(m);
506     mpfr_init(s);
507 
508     mpfr_trunc(d, c->mf);
509     mpfr_frac(m, c->mf, MPFR_DEFAULT_RND);
510     mpfr_mul_ui(m, m, 60, MPFR_DEFAULT_RND);
511 
512     mpfr_frac(s, m, MPFR_DEFAULT_RND);
513     mpfr_trunc(m, m);
514     mpfr_mul_ui(s, s, 60, MPFR_DEFAULT_RND);
515     mpfr_ceil(s, s);
516 
517     mpfr_div_ui(m, m, 100, MPFR_DEFAULT_RND);
518     mpfr_div_ui(s, s, 10000, MPFR_DEFAULT_RND);
519     mpfr_add(c->mf, d, m, MPFR_DEFAULT_RND);
520     mpfr_add(c->mf, c->mf, s, MPFR_DEFAULT_RND);
521 
522     mpfr_clear(d);
523     mpfr_clear(m);
524     mpfr_clear(s);
525 }
526 
527 void rpn_zero(calc_number_t *c)
528 {
529     mpfr_set_ui(c->mf, 0, MPFR_DEFAULT_RND);
530 }
531 
532 void rpn_copy(calc_number_t *dst, calc_number_t *src)
533 {
534     mpfr_set(dst->mf, src->mf, MPFR_DEFAULT_RND);
535 }
536 
537 int rpn_is_zero(calc_number_t *c)
538 {
539     return (mpfr_sgn(c->mf) == 0);
540 }
541 
542 void rpn_alloc(calc_number_t *c)
543 {
544     mpfr_init(c->mf);
545 }
546 
547 void rpn_free(calc_number_t *c)
548 {
549     mpfr_clear(c->mf);
550 }
551