Skip to content
Newer
Older
100644 641 lines (545 sloc) 12.4 KB
511dc44 initial import
Laurent Sansonetti authored
1 /**********************************************************************
2
3 math.c -
4
5 $Author: akr $
6 created at: Tue Jan 25 14:12:56 JST 1994
7
8 Copyright (C) 1993-2007 Yukihiro Matsumoto
9
10 **********************************************************************/
11
12 #include "ruby/ruby.h"
13 #include <math.h>
14 #include <errno.h>
15
16 VALUE rb_mMath;
17
18 #define Need_Float(x) (x) = rb_Float(x)
19 #define Need_Float2(x,y) do {\
20 Need_Float(x);\
21 Need_Float(y);\
22 } while (0)
23
24 static void
25 domain_check(double x, char *msg)
26 {
27 while(1) {
28 if (errno) {
29 rb_sys_fail(msg);
30 }
31 if (isnan(x)) {
32 #if defined(EDOM)
33 errno = EDOM;
34 #elif defined(ERANGE)
35 errno = ERANGE;
36 #endif
37 continue;
38 }
39 break;
40 }
41 }
42
43
44 /*
45 * call-seq:
46 * Math.atan2(y, x) => float
47 *
48 * Computes the arc tangent given <i>y</i> and <i>x</i>. Returns
49 * -PI..PI.
50 *
51 */
52
53 static VALUE
54 math_atan2(VALUE obj, VALUE y, VALUE x)
55 {
56 Need_Float2(y, x);
57 return DOUBLE2NUM(atan2(RFLOAT_VALUE(y), RFLOAT_VALUE(x)));
58 }
59
60
61 /*
62 * call-seq:
63 * Math.cos(x) => float
64 *
65 * Computes the cosine of <i>x</i> (expressed in radians). Returns
66 * -1..1.
67 */
68
69 static VALUE
70 math_cos(VALUE obj, VALUE x)
71 {
72 Need_Float(x);
73 return DOUBLE2NUM(cos(RFLOAT_VALUE(x)));
74 }
75
76 /*
77 * call-seq:
78 * Math.sin(x) => float
79 *
80 * Computes the sine of <i>x</i> (expressed in radians). Returns
81 * -1..1.
82 */
83
84 static VALUE
85 math_sin(VALUE obj, VALUE x)
86 {
87 Need_Float(x);
88
89 return DOUBLE2NUM(sin(RFLOAT_VALUE(x)));
90 }
91
92
93 /*
94 * call-seq:
95 * Math.tan(x) => float
96 *
97 * Returns the tangent of <i>x</i> (expressed in radians).
98 */
99
100 static VALUE
101 math_tan(VALUE obj, VALUE x)
102 {
103 Need_Float(x);
104
105 return DOUBLE2NUM(tan(RFLOAT_VALUE(x)));
106 }
107
108 /*
109 * call-seq:
110 * Math.acos(x) => float
111 *
112 * Computes the arc cosine of <i>x</i>. Returns 0..PI.
113 */
114
115 static VALUE
116 math_acos(VALUE obj, VALUE x)
117 {
118 double d;
119
120 Need_Float(x);
121 errno = 0;
122 d = acos(RFLOAT_VALUE(x));
123 domain_check(d, "acos");
124 return DOUBLE2NUM(d);
125 }
126
127 /*
128 * call-seq:
129 * Math.asin(x) => float
130 *
131 * Computes the arc sine of <i>x</i>. Returns -{PI/2} .. {PI/2}.
132 */
133
134 static VALUE
135 math_asin(VALUE obj, VALUE x)
136 {
137 double d;
138
139 Need_Float(x);
140 errno = 0;
141 d = asin(RFLOAT_VALUE(x));
142 domain_check(d, "asin");
143 return DOUBLE2NUM(d);
144 }
145
146 /*
147 * call-seq:
148 * Math.atan(x) => float
149 *
150 * Computes the arc tangent of <i>x</i>. Returns -{PI/2} .. {PI/2}.
151 */
152
153 static VALUE
154 math_atan(VALUE obj, VALUE x)
155 {
156 Need_Float(x);
157 return DOUBLE2NUM(atan(RFLOAT_VALUE(x)));
158 }
159
160 #ifndef HAVE_COSH
161 double
162 cosh(double x)
163 {
164 return (exp(x) + exp(-x)) / 2;
165 }
166 #endif
167
168 /*
169 * call-seq:
170 * Math.cosh(x) => float
171 *
172 * Computes the hyperbolic cosine of <i>x</i> (expressed in radians).
173 */
174
175 static VALUE
176 math_cosh(VALUE obj, VALUE x)
177 {
178 Need_Float(x);
179
180 return DOUBLE2NUM(cosh(RFLOAT_VALUE(x)));
181 }
182
183 #ifndef HAVE_SINH
184 double
185 sinh(double x)
186 {
187 return (exp(x) - exp(-x)) / 2;
188 }
189 #endif
190
191 /*
192 * call-seq:
193 * Math.sinh(x) => float
194 *
195 * Computes the hyperbolic sine of <i>x</i> (expressed in
196 * radians).
197 */
198
199 static VALUE
200 math_sinh(VALUE obj, VALUE x)
201 {
202 Need_Float(x);
203 return DOUBLE2NUM(sinh(RFLOAT_VALUE(x)));
204 }
205
206 #ifndef HAVE_TANH
207 double
208 tanh(double x)
209 {
210 return sinh(x) / cosh(x);
211 }
212 #endif
213
214 /*
215 * call-seq:
216 * Math.tanh() => float
217 *
218 * Computes the hyperbolic tangent of <i>x</i> (expressed in
219 * radians).
220 */
221
222 static VALUE
223 math_tanh(VALUE obj, VALUE x)
224 {
225 Need_Float(x);
226 return DOUBLE2NUM(tanh(RFLOAT_VALUE(x)));
227 }
228
229 /*
230 * call-seq:
231 * Math.acosh(x) => float
232 *
233 * Computes the inverse hyperbolic cosine of <i>x</i>.
234 */
235
236 static VALUE
237 math_acosh(VALUE obj, VALUE x)
238 {
239 double d;
240
241 Need_Float(x);
242 errno = 0;
243 d = acosh(RFLOAT_VALUE(x));
244 domain_check(d, "acosh");
245 return DOUBLE2NUM(d);
246 }
247
248 /*
249 * call-seq:
250 * Math.asinh(x) => float
251 *
252 * Computes the inverse hyperbolic sine of <i>x</i>.
253 */
254
255 static VALUE
256 math_asinh(VALUE obj, VALUE x)
257 {
258 Need_Float(x);
259 return DOUBLE2NUM(asinh(RFLOAT_VALUE(x)));
260 }
261
262 /*
263 * call-seq:
264 * Math.atanh(x) => float
265 *
266 * Computes the inverse hyperbolic tangent of <i>x</i>.
267 */
268
269 static VALUE
270 math_atanh(VALUE obj, VALUE x)
271 {
272 double d;
273
274 Need_Float(x);
275 errno = 0;
276 d = atanh(RFLOAT_VALUE(x));
277 domain_check(d, "atanh");
278 return DOUBLE2NUM(d);
279 }
280
281 /*
282 * call-seq:
283 * Math.exp(x) => float
284 *
285 * Returns e**x.
286 */
287
288 static VALUE
289 math_exp(VALUE obj, VALUE x)
290 {
291 Need_Float(x);
292 return DOUBLE2NUM(exp(RFLOAT_VALUE(x)));
293 }
294
295 #if defined __CYGWIN__
296 # include <cygwin/version.h>
297 # if CYGWIN_VERSION_DLL_MAJOR < 1005
298 # define nan(x) nan()
299 # endif
300 # define log(x) ((x) < 0.0 ? nan("") : log(x))
301 # define log10(x) ((x) < 0.0 ? nan("") : log10(x))
302 #endif
303
304 /*
305 * call-seq:
306 * Math.log(numeric) => float
307 * Math.log(num,base) => float
308 *
309 * Returns the natural logarithm of <i>numeric</i>.
310 * If additional second argument is given, it will be the base
311 * of logarithm.
312 */
313
314 static VALUE
315 math_log(int argc, VALUE *argv)
316 {
317 VALUE x, base;
318 double d;
319
320 rb_scan_args(argc, argv, "11", &x, &base);
321 Need_Float(x);
322 errno = 0;
323 d = log(RFLOAT_VALUE(x));
324 if (!NIL_P(base)) {
325 Need_Float(base);
326 d /= log(RFLOAT_VALUE(base));
327 }
328 domain_check(d, "log");
329 return DOUBLE2NUM(d);
330 }
331
332 #ifndef log2
333 #ifndef HAVE_LOG2
334 double
335 log2(double x)
336 {
337 return log10(x)/log10(2.0);
338 }
339 #else
340 extern double log2(double);
341 #endif
342 #endif
343
344 /*
345 * call-seq:
346 * Math.log2(numeric) => float
347 *
348 * Returns the base 2 logarithm of <i>numeric</i>.
349 */
350
351 static VALUE
352 math_log2(VALUE obj, VALUE x)
353 {
354 double d;
355
356 Need_Float(x);
357 errno = 0;
358 d = log2(RFLOAT_VALUE(x));
359 if (errno) {
360 rb_sys_fail("log2");
361 }
362 return DOUBLE2NUM(d);
363 }
364
365 /*
366 * call-seq:
367 * Math.log10(numeric) => float
368 *
369 * Returns the base 10 logarithm of <i>numeric</i>.
370 */
371
372 static VALUE
373 math_log10(VALUE obj, VALUE x)
374 {
375 double d;
376
377 Need_Float(x);
378 errno = 0;
379 d = log10(RFLOAT_VALUE(x));
380 domain_check(d, "log10");
381 return DOUBLE2NUM(d);
382 }
383
384 /*
385 * call-seq:
386 * Math.sqrt(numeric) => float
387 *
388 * Returns the non-negative square root of <i>numeric</i>.
389 */
390
391 static VALUE
392 math_sqrt(VALUE obj, VALUE x)
393 {
394 double d;
395
396 Need_Float(x);
397 errno = 0;
398 d = sqrt(RFLOAT_VALUE(x));
399 domain_check(d, "sqrt");
400 return DOUBLE2NUM(d);
401 }
402
403 /*
404 * call-seq:
405 * Math.cbrt(numeric) => float
406 *
407 * Returns the cube root of <i>numeric</i>.
408 */
409
410 static VALUE
411 math_cbrt(VALUE obj, VALUE x)
412 {
413 Need_Float(x);
414 return DOUBLE2NUM(cbrt(RFLOAT_VALUE(x)));
415 }
416
417 /*
418 * call-seq:
419 * Math.frexp(numeric) => [ fraction, exponent ]
420 *
421 * Returns a two-element array containing the normalized fraction (a
422 * <code>Float</code>) and exponent (a <code>Fixnum</code>) of
423 * <i>numeric</i>.
424 *
425 * fraction, exponent = Math.frexp(1234) #=> [0.6025390625, 11]
426 * fraction * 2**exponent #=> 1234.0
427 */
428
429 static VALUE
430 math_frexp(VALUE obj, VALUE x)
431 {
432 double d;
433 int exp;
434
435 Need_Float(x);
436
437 d = frexp(RFLOAT_VALUE(x), &exp);
438 return rb_assoc_new(DOUBLE2NUM(d), INT2NUM(exp));
439 }
440
441 /*
442 * call-seq:
443 * Math.ldexp(flt, int) -> float
444 *
445 * Returns the value of <i>flt</i>*(2**<i>int</i>).
446 *
447 * fraction, exponent = Math.frexp(1234)
448 * Math.ldexp(fraction, exponent) #=> 1234.0
449 */
450
451 static VALUE
452 math_ldexp(VALUE obj, VALUE x, VALUE n)
453 {
454 Need_Float(x);
455 return DOUBLE2NUM(ldexp(RFLOAT_VALUE(x), NUM2INT(n)));
456 }
457
458 /*
459 * call-seq:
460 * Math.hypot(x, y) => float
461 *
462 * Returns sqrt(x**2 + y**2), the hypotenuse of a right-angled triangle
463 * with sides <i>x</i> and <i>y</i>.
464 *
465 * Math.hypot(3, 4) #=> 5.0
466 */
467
468 static VALUE
469 math_hypot(VALUE obj, VALUE x, VALUE y)
470 {
471 Need_Float2(x, y);
472 return DOUBLE2NUM(hypot(RFLOAT_VALUE(x), RFLOAT_VALUE(y)));
473 }
474
475 /*
476 * call-seq:
477 * Math.erf(x) => float
478 *
479 * Calculates the error function of x.
480 */
481
482 static VALUE
483 math_erf(VALUE obj, VALUE x)
484 {
485 Need_Float(x);
486 return DOUBLE2NUM(erf(RFLOAT_VALUE(x)));
487 }
488
489 /*
490 * call-seq:
491 * Math.erfc(x) => float
492 *
493 * Calculates the complementary error function of x.
494 */
495
496 static VALUE
497 math_erfc(VALUE obj, VALUE x)
498 {
499 Need_Float(x);
500 return DOUBLE2NUM(erfc(RFLOAT_VALUE(x)));
501 }
502
503 /*
504 * call-seq:
505 * Math.gamma(x) => float
506 *
507 * Calculates the gamma function of x.
508 *
509 * Note that gamma(n) is same as fact(n-1) for integer n >= 0.
510 * However gamma(n) returns float and possibly has error in calculation.
511 *
512 * def fact(n) (1..n).inject(1) {|r,i| r*i } end
513 * 0.upto(25) {|i| p [i, Math.gamma(i+1), fact(i)] }
514 * =>
515 * [0, 1.0, 1]
516 * [1, 1.0, 1]
517 * [2, 2.0, 2]
518 * [3, 6.0, 6]
519 * [4, 24.0, 24]
520 * [5, 120.0, 120]
521 * [6, 720.0, 720]
522 * [7, 5040.0, 5040]
523 * [8, 40320.0, 40320]
524 * [9, 362880.0, 362880]
525 * [10, 3628800.0, 3628800]
526 * [11, 39916800.0, 39916800]
527 * [12, 479001599.999999, 479001600]
528 * [13, 6227020800.00001, 6227020800]
529 * [14, 87178291199.9998, 87178291200]
530 * [15, 1307674368000.0, 1307674368000]
531 * [16, 20922789888000.0, 20922789888000]
532 * [17, 3.55687428096001e+14, 355687428096000]
533 * [18, 6.40237370572799e+15, 6402373705728000]
534 * [19, 1.21645100408832e+17, 121645100408832000]
535 * [20, 2.43290200817664e+18, 2432902008176640000]
536 * [21, 5.10909421717094e+19, 51090942171709440000]
537 * [22, 1.12400072777761e+21, 1124000727777607680000]
538 * [23, 2.58520167388851e+22, 25852016738884976640000]
539 * [24, 6.20448401733239e+23, 620448401733239439360000]
540 * [25, 1.5511210043331e+25, 15511210043330985984000000]
541 *
542 */
543
544 static VALUE
545 math_gamma(VALUE obj, VALUE x)
546 {
547 double d;
548 Need_Float(x);
549 errno = 0;
550 d = tgamma(RFLOAT_VALUE(x));
551 domain_check(d, "gamma");
552 return DOUBLE2NUM(d);
553 }
554
555 /*
556 * call-seq:
557 * Math.lgamma(x) => [float, -1 or 1]
558 *
559 * Calculates the logarithmic gamma of x and
560 * the sign of gamma of x.
561 *
562 * Math.lgamma(x) is same as
563 * [Math.log(Math.gamma(x).abs), Math.gamma(x) < 0 ? -1 : 1]
564 * but avoid overflow by Math.gamma(x) for large x.
565 */
566
567 static VALUE
568 math_lgamma(VALUE obj, VALUE x)
569 {
570 double d;
571 int sign;
572 VALUE v;
573 Need_Float(x);
574 errno = 0;
575 d = lgamma_r(RFLOAT_VALUE(x), &sign);
576 domain_check(d, "lgamma");
577 v = DOUBLE2NUM(d);
578 return rb_assoc_new(v, INT2FIX(sign));
579 }
580
581 /*
582 * The <code>Math</code> module contains module functions for basic
583 * trigonometric and transcendental functions. See class
584 * <code>Float</code> for a list of constants that
585 * define Ruby's floating point accuracy.
586 */
587
588
589 void
590 Init_Math(void)
591 {
592 rb_mMath = rb_define_module("Math");
593
594 #ifdef M_PI
595 rb_define_const(rb_mMath, "PI", DOUBLE2NUM(M_PI));
596 #else
597 rb_define_const(rb_mMath, "PI", DOUBLE2NUM(atan(1.0)*4.0));
598 #endif
599
600 #ifdef M_E
601 rb_define_const(rb_mMath, "E", DOUBLE2NUM(M_E));
602 #else
603 rb_define_const(rb_mMath, "E", DOUBLE2NUM(exp(1.0)));
604 #endif
605
606 rb_define_module_function(rb_mMath, "atan2", math_atan2, 2);
607 rb_define_module_function(rb_mMath, "cos", math_cos, 1);
608 rb_define_module_function(rb_mMath, "sin", math_sin, 1);
609 rb_define_module_function(rb_mMath, "tan", math_tan, 1);
610
611 rb_define_module_function(rb_mMath, "acos", math_acos, 1);
612 rb_define_module_function(rb_mMath, "asin", math_asin, 1);
613 rb_define_module_function(rb_mMath, "atan", math_atan, 1);
614
615 rb_define_module_function(rb_mMath, "cosh", math_cosh, 1);
616 rb_define_module_function(rb_mMath, "sinh", math_sinh, 1);
617 rb_define_module_function(rb_mMath, "tanh", math_tanh, 1);
618
619 rb_define_module_function(rb_mMath, "acosh", math_acosh, 1);
620 rb_define_module_function(rb_mMath, "asinh", math_asinh, 1);
621 rb_define_module_function(rb_mMath, "atanh", math_atanh, 1);
622
623 rb_define_module_function(rb_mMath, "exp", math_exp, 1);
624 rb_define_module_function(rb_mMath, "log", math_log, -1);
625 rb_define_module_function(rb_mMath, "log2", math_log2, 1);
626 rb_define_module_function(rb_mMath, "log10", math_log10, 1);
627 rb_define_module_function(rb_mMath, "sqrt", math_sqrt, 1);
628 rb_define_module_function(rb_mMath, "cbrt", math_cbrt, 1);
629
630 rb_define_module_function(rb_mMath, "frexp", math_frexp, 1);
631 rb_define_module_function(rb_mMath, "ldexp", math_ldexp, 2);
632
633 rb_define_module_function(rb_mMath, "hypot", math_hypot, 2);
634
635 rb_define_module_function(rb_mMath, "erf", math_erf, 1);
636 rb_define_module_function(rb_mMath, "erfc", math_erfc, 1);
637
638 rb_define_module_function(rb_mMath, "gamma", math_gamma, 1);
639 rb_define_module_function(rb_mMath, "lgamma", math_lgamma, 1);
640 }
Something went wrong with that request. Please try again.