1/*
2  rational.c: Coded by Tadayoshi Funaba 2008-2012
3
4  This implementation is based on Keiju Ishitsuka's Rational library
5  which is written in ruby.
6*/
7
8#include "ruby.h"
9#include "internal.h"
10#include <math.h>
11#include <float.h>
12
13#ifdef HAVE_IEEEFP_H
14#include <ieeefp.h>
15#endif
16
17#define NDEBUG
18#include <assert.h>
19
20#define ZERO INT2FIX(0)
21#define ONE INT2FIX(1)
22#define TWO INT2FIX(2)
23
24VALUE rb_cRational;
25
26static ID id_abs, id_cmp, id_convert, id_eqeq_p, id_expt, id_fdiv,
27    id_floor, id_idiv, id_inspect, id_integer_p, id_negate, id_to_f,
28    id_to_i, id_to_s, id_truncate, id_i_num, id_i_den;
29
30#define f_boolcast(x) ((x) ? Qtrue : Qfalse)
31
32#define binop(n,op) \
33inline static VALUE \
34f_##n(VALUE x, VALUE y)\
35{\
36  return rb_funcall(x, (op), 1, y);\
37}
38
39#define fun1(n) \
40inline static VALUE \
41f_##n(VALUE x)\
42{\
43    return rb_funcall(x, id_##n, 0);\
44}
45
46#define fun2(n) \
47inline static VALUE \
48f_##n(VALUE x, VALUE y)\
49{\
50    return rb_funcall(x, id_##n, 1, y);\
51}
52
53inline static VALUE
54f_add(VALUE x, VALUE y)
55{
56    if (FIXNUM_P(y) && FIX2LONG(y) == 0)
57	return x;
58    else if (FIXNUM_P(x) && FIX2LONG(x) == 0)
59	return y;
60    return rb_funcall(x, '+', 1, y);
61}
62
63inline static VALUE
64f_cmp(VALUE x, VALUE y)
65{
66    if (FIXNUM_P(x) && FIXNUM_P(y)) {
67	long c = FIX2LONG(x) - FIX2LONG(y);
68	if (c > 0)
69	    c = 1;
70	else if (c < 0)
71	    c = -1;
72	return INT2FIX(c);
73    }
74    return rb_funcall(x, id_cmp, 1, y);
75}
76
77inline static VALUE
78f_div(VALUE x, VALUE y)
79{
80    if (FIXNUM_P(y) && FIX2LONG(y) == 1)
81	return x;
82    return rb_funcall(x, '/', 1, y);
83}
84
85inline static VALUE
86f_gt_p(VALUE x, VALUE y)
87{
88    if (FIXNUM_P(x) && FIXNUM_P(y))
89	return f_boolcast(FIX2LONG(x) > FIX2LONG(y));
90    return rb_funcall(x, '>', 1, y);
91}
92
93inline static VALUE
94f_lt_p(VALUE x, VALUE y)
95{
96    if (FIXNUM_P(x) && FIXNUM_P(y))
97	return f_boolcast(FIX2LONG(x) < FIX2LONG(y));
98    return rb_funcall(x, '<', 1, y);
99}
100
101binop(mod, '%')
102
103inline static VALUE
104f_mul(VALUE x, VALUE y)
105{
106    if (FIXNUM_P(y)) {
107	long iy = FIX2LONG(y);
108	if (iy == 0) {
109	    if (FIXNUM_P(x) || RB_TYPE_P(x, T_BIGNUM))
110		return ZERO;
111	}
112	else if (iy == 1)
113	    return x;
114    }
115    else if (FIXNUM_P(x)) {
116	long ix = FIX2LONG(x);
117	if (ix == 0) {
118	    if (FIXNUM_P(y) || RB_TYPE_P(y, T_BIGNUM))
119		return ZERO;
120	}
121	else if (ix == 1)
122	    return y;
123    }
124    return rb_funcall(x, '*', 1, y);
125}
126
127inline static VALUE
128f_sub(VALUE x, VALUE y)
129{
130    if (FIXNUM_P(y) && FIX2LONG(y) == 0)
131	return x;
132    return rb_funcall(x, '-', 1, y);
133}
134
135fun1(abs)
136fun1(floor)
137fun1(inspect)
138fun1(integer_p)
139fun1(negate)
140
141inline static VALUE
142f_to_i(VALUE x)
143{
144    if (RB_TYPE_P(x, T_STRING))
145	return rb_str_to_inum(x, 10, 0);
146    return rb_funcall(x, id_to_i, 0);
147}
148inline static VALUE
149f_to_f(VALUE x)
150{
151    if (RB_TYPE_P(x, T_STRING))
152	return DBL2NUM(rb_str_to_dbl(x, 0));
153    return rb_funcall(x, id_to_f, 0);
154}
155
156fun1(to_s)
157fun1(truncate)
158
159inline static VALUE
160f_eqeq_p(VALUE x, VALUE y)
161{
162    if (FIXNUM_P(x) && FIXNUM_P(y))
163	return f_boolcast(FIX2LONG(x) == FIX2LONG(y));
164    return rb_funcall(x, id_eqeq_p, 1, y);
165}
166
167fun2(expt)
168fun2(fdiv)
169fun2(idiv)
170
171#define f_expt10(x) f_expt(INT2FIX(10), x)
172
173inline static VALUE
174f_negative_p(VALUE x)
175{
176    if (FIXNUM_P(x))
177	return f_boolcast(FIX2LONG(x) < 0);
178    return rb_funcall(x, '<', 1, ZERO);
179}
180
181#define f_positive_p(x) (!f_negative_p(x))
182
183inline static VALUE
184f_zero_p(VALUE x)
185{
186    switch (TYPE(x)) {
187      case T_FIXNUM:
188	return f_boolcast(FIX2LONG(x) == 0);
189      case T_BIGNUM:
190	return Qfalse;
191      case T_RATIONAL:
192      {
193	  VALUE num = RRATIONAL(x)->num;
194
195	  return f_boolcast(FIXNUM_P(num) && FIX2LONG(num) == 0);
196      }
197    }
198    return rb_funcall(x, id_eqeq_p, 1, ZERO);
199}
200
201#define f_nonzero_p(x) (!f_zero_p(x))
202
203inline static VALUE
204f_one_p(VALUE x)
205{
206    switch (TYPE(x)) {
207      case T_FIXNUM:
208	return f_boolcast(FIX2LONG(x) == 1);
209      case T_BIGNUM:
210	return Qfalse;
211      case T_RATIONAL:
212      {
213	  VALUE num = RRATIONAL(x)->num;
214	  VALUE den = RRATIONAL(x)->den;
215
216	  return f_boolcast(FIXNUM_P(num) && FIX2LONG(num) == 1 &&
217			    FIXNUM_P(den) && FIX2LONG(den) == 1);
218      }
219    }
220    return rb_funcall(x, id_eqeq_p, 1, ONE);
221}
222
223inline static VALUE
224f_minus_one_p(VALUE x)
225{
226    switch (TYPE(x)) {
227      case T_FIXNUM:
228	return f_boolcast(FIX2LONG(x) == -1);
229      case T_BIGNUM:
230	return Qfalse;
231      case T_RATIONAL:
232      {
233	  VALUE num = RRATIONAL(x)->num;
234	  VALUE den = RRATIONAL(x)->den;
235
236	  return f_boolcast(FIXNUM_P(num) && FIX2LONG(num) == -1 &&
237			    FIXNUM_P(den) && FIX2LONG(den) == 1);
238      }
239    }
240    return rb_funcall(x, id_eqeq_p, 1, INT2FIX(-1));
241}
242
243inline static VALUE
244f_kind_of_p(VALUE x, VALUE c)
245{
246    return rb_obj_is_kind_of(x, c);
247}
248
249inline static VALUE
250k_numeric_p(VALUE x)
251{
252    return f_kind_of_p(x, rb_cNumeric);
253}
254
255inline static VALUE
256k_integer_p(VALUE x)
257{
258    return f_kind_of_p(x, rb_cInteger);
259}
260
261inline static VALUE
262k_float_p(VALUE x)
263{
264    return f_kind_of_p(x, rb_cFloat);
265}
266
267inline static VALUE
268k_rational_p(VALUE x)
269{
270    return f_kind_of_p(x, rb_cRational);
271}
272
273#define k_exact_p(x) (!k_float_p(x))
274#define k_inexact_p(x) k_float_p(x)
275
276#define k_exact_zero_p(x) (k_exact_p(x) && f_zero_p(x))
277#define k_exact_one_p(x) (k_exact_p(x) && f_one_p(x))
278
279#ifndef NDEBUG
280#define f_gcd f_gcd_orig
281#endif
282
283inline static long
284i_gcd(long x, long y)
285{
286    if (x < 0)
287	x = -x;
288    if (y < 0)
289	y = -y;
290
291    if (x == 0)
292	return y;
293    if (y == 0)
294	return x;
295
296    while (x > 0) {
297	long t = x;
298	x = y % x;
299	y = t;
300    }
301    return y;
302}
303
304inline static VALUE
305f_gcd(VALUE x, VALUE y)
306{
307    VALUE z;
308
309    if (FIXNUM_P(x) && FIXNUM_P(y))
310	return LONG2NUM(i_gcd(FIX2LONG(x), FIX2LONG(y)));
311
312    if (f_negative_p(x))
313	x = f_negate(x);
314    if (f_negative_p(y))
315	y = f_negate(y);
316
317    if (f_zero_p(x))
318	return y;
319    if (f_zero_p(y))
320	return x;
321
322    for (;;) {
323	if (FIXNUM_P(x)) {
324	    if (FIX2LONG(x) == 0)
325		return y;
326	    if (FIXNUM_P(y))
327		return LONG2NUM(i_gcd(FIX2LONG(x), FIX2LONG(y)));
328	}
329	z = x;
330	x = f_mod(y, x);
331	y = z;
332    }
333    /* NOTREACHED */
334}
335
336#ifndef NDEBUG
337#undef f_gcd
338
339inline static VALUE
340f_gcd(VALUE x, VALUE y)
341{
342    VALUE r = f_gcd_orig(x, y);
343    if (f_nonzero_p(r)) {
344	assert(f_zero_p(f_mod(x, r)));
345	assert(f_zero_p(f_mod(y, r)));
346    }
347    return r;
348}
349#endif
350
351inline static VALUE
352f_lcm(VALUE x, VALUE y)
353{
354    if (f_zero_p(x) || f_zero_p(y))
355	return ZERO;
356    return f_abs(f_mul(f_div(x, f_gcd(x, y)), y));
357}
358
359#define get_dat1(x) \
360    struct RRational *dat;\
361    dat = ((struct RRational *)(x))
362
363#define get_dat2(x,y) \
364    struct RRational *adat, *bdat;\
365    adat = ((struct RRational *)(x));\
366    bdat = ((struct RRational *)(y))
367
368inline static VALUE
369nurat_s_new_internal(VALUE klass, VALUE num, VALUE den)
370{
371    NEWOBJ_OF(obj, struct RRational, klass, T_RATIONAL);
372
373    obj->num = num;
374    obj->den = den;
375
376    return (VALUE)obj;
377}
378
379static VALUE
380nurat_s_alloc(VALUE klass)
381{
382    return nurat_s_new_internal(klass, ZERO, ONE);
383}
384
385#define rb_raise_zerodiv() rb_raise(rb_eZeroDivError, "divided by 0")
386
387#if 0
388static VALUE
389nurat_s_new_bang(int argc, VALUE *argv, VALUE klass)
390{
391    VALUE num, den;
392
393    switch (rb_scan_args(argc, argv, "11", &num, &den)) {
394      case 1:
395	if (!k_integer_p(num))
396	    num = f_to_i(num);
397	den = ONE;
398	break;
399      default:
400	if (!k_integer_p(num))
401	    num = f_to_i(num);
402	if (!k_integer_p(den))
403	    den = f_to_i(den);
404
405	switch (FIX2INT(f_cmp(den, ZERO))) {
406	  case -1:
407	    num = f_negate(num);
408	    den = f_negate(den);
409	    break;
410	  case 0:
411	    rb_raise_zerodiv();
412	    break;
413	}
414	break;
415    }
416
417    return nurat_s_new_internal(klass, num, den);
418}
419#endif
420
421inline static VALUE
422f_rational_new_bang1(VALUE klass, VALUE x)
423{
424    return nurat_s_new_internal(klass, x, ONE);
425}
426
427inline static VALUE
428f_rational_new_bang2(VALUE klass, VALUE x, VALUE y)
429{
430    assert(f_positive_p(y));
431    assert(f_nonzero_p(y));
432    return nurat_s_new_internal(klass, x, y);
433}
434
435#ifdef CANONICALIZATION_FOR_MATHN
436#define CANON
437#endif
438
439#ifdef CANON
440static int canonicalization = 0;
441
442RUBY_FUNC_EXPORTED void
443nurat_canonicalization(int f)
444{
445    canonicalization = f;
446}
447#endif
448
449inline static void
450nurat_int_check(VALUE num)
451{
452    switch (TYPE(num)) {
453      case T_FIXNUM:
454      case T_BIGNUM:
455	break;
456      default:
457	if (!k_numeric_p(num) || !f_integer_p(num))
458	    rb_raise(rb_eTypeError, "not an integer");
459    }
460}
461
462inline static VALUE
463nurat_int_value(VALUE num)
464{
465    nurat_int_check(num);
466    if (!k_integer_p(num))
467	num = f_to_i(num);
468    return num;
469}
470
471inline static VALUE
472nurat_s_canonicalize_internal(VALUE klass, VALUE num, VALUE den)
473{
474    VALUE gcd;
475
476    switch (FIX2INT(f_cmp(den, ZERO))) {
477      case -1:
478	num = f_negate(num);
479	den = f_negate(den);
480	break;
481      case 0:
482	rb_raise_zerodiv();
483	break;
484    }
485
486    gcd = f_gcd(num, den);
487    num = f_idiv(num, gcd);
488    den = f_idiv(den, gcd);
489
490#ifdef CANON
491    if (f_one_p(den) && canonicalization)
492	return num;
493#endif
494    return nurat_s_new_internal(klass, num, den);
495}
496
497inline static VALUE
498nurat_s_canonicalize_internal_no_reduce(VALUE klass, VALUE num, VALUE den)
499{
500    switch (FIX2INT(f_cmp(den, ZERO))) {
501      case -1:
502	num = f_negate(num);
503	den = f_negate(den);
504	break;
505      case 0:
506	rb_raise_zerodiv();
507	break;
508    }
509
510#ifdef CANON
511    if (f_one_p(den) && canonicalization)
512	return num;
513#endif
514    return nurat_s_new_internal(klass, num, den);
515}
516
517static VALUE
518nurat_s_new(int argc, VALUE *argv, VALUE klass)
519{
520    VALUE num, den;
521
522    switch (rb_scan_args(argc, argv, "11", &num, &den)) {
523      case 1:
524	num = nurat_int_value(num);
525	den = ONE;
526	break;
527      default:
528	num = nurat_int_value(num);
529	den = nurat_int_value(den);
530	break;
531    }
532
533    return nurat_s_canonicalize_internal(klass, num, den);
534}
535
536inline static VALUE
537f_rational_new1(VALUE klass, VALUE x)
538{
539    assert(!k_rational_p(x));
540    return nurat_s_canonicalize_internal(klass, x, ONE);
541}
542
543inline static VALUE
544f_rational_new2(VALUE klass, VALUE x, VALUE y)
545{
546    assert(!k_rational_p(x));
547    assert(!k_rational_p(y));
548    return nurat_s_canonicalize_internal(klass, x, y);
549}
550
551inline static VALUE
552f_rational_new_no_reduce1(VALUE klass, VALUE x)
553{
554    assert(!k_rational_p(x));
555    return nurat_s_canonicalize_internal_no_reduce(klass, x, ONE);
556}
557
558inline static VALUE
559f_rational_new_no_reduce2(VALUE klass, VALUE x, VALUE y)
560{
561    assert(!k_rational_p(x));
562    assert(!k_rational_p(y));
563    return nurat_s_canonicalize_internal_no_reduce(klass, x, y);
564}
565
566/*
567 * call-seq:
568 *    Rational(x[, y])  ->  numeric
569 *
570 * Returns x/y;
571 *
572 *    Rational(1, 2)   #=> (1/2)
573 *    Rational('1/2')  #=> (1/2)
574 *
575 * Syntax of string form:
576 *
577 *   string form = extra spaces , rational , extra spaces ;
578 *   rational = [ sign ] , unsigned rational ;
579 *   unsigned rational = numerator | numerator , "/" , denominator ;
580 *   numerator = integer part | fractional part | integer part , fractional part ;
581 *   denominator = digits ;
582 *   integer part = digits ;
583 *   fractional part = "." , digits , [ ( "e" | "E" ) , [ sign ] , digits ] ;
584 *   sign = "-" | "+" ;
585 *   digits = digit , { digit | "_" , digit } ;
586 *   digit = "0" | "1" | "2" | "3" | "4" | "5" | "6" | "7" | "8" | "9" ;
587 *   extra spaces = ? \s* ? ;
588 *
589 * See String#to_r.
590 */
591static VALUE
592nurat_f_rational(int argc, VALUE *argv, VALUE klass)
593{
594    return rb_funcall2(rb_cRational, id_convert, argc, argv);
595}
596
597/*
598 * call-seq:
599 *    rat.numerator  ->  integer
600 *
601 * Returns the numerator.
602 *
603 *    Rational(7).numerator        #=> 7
604 *    Rational(7, 1).numerator     #=> 7
605 *    Rational(9, -4).numerator    #=> -9
606 *    Rational(-2, -10).numerator  #=> 1
607 */
608static VALUE
609nurat_numerator(VALUE self)
610{
611    get_dat1(self);
612    return dat->num;
613}
614
615/*
616 * call-seq:
617 *    rat.denominator  ->  integer
618 *
619 * Returns the denominator (always positive).
620 *
621 *    Rational(7).denominator             #=> 1
622 *    Rational(7, 1).denominator          #=> 1
623 *    Rational(9, -4).denominator         #=> 4
624 *    Rational(-2, -10).denominator       #=> 5
625 *    rat.numerator.gcd(rat.denominator)  #=> 1
626 */
627static VALUE
628nurat_denominator(VALUE self)
629{
630    get_dat1(self);
631    return dat->den;
632}
633
634#ifndef NDEBUG
635#define f_imul f_imul_orig
636#endif
637
638inline static VALUE
639f_imul(long a, long b)
640{
641    VALUE r;
642
643    if (a == 0 || b == 0)
644	return ZERO;
645    else if (a == 1)
646	return LONG2NUM(b);
647    else if (b == 1)
648	return LONG2NUM(a);
649
650    if (MUL_OVERFLOW_LONG_P(a, b))
651	r = rb_big_mul(rb_int2big(a), rb_int2big(b));
652    else
653        r = LONG2NUM(a * b);
654    return r;
655}
656
657#ifndef NDEBUG
658#undef f_imul
659
660inline static VALUE
661f_imul(long x, long y)
662{
663    VALUE r = f_imul_orig(x, y);
664    assert(f_eqeq_p(r, f_mul(LONG2NUM(x), LONG2NUM(y))));
665    return r;
666}
667#endif
668
669inline static VALUE
670f_addsub(VALUE self, VALUE anum, VALUE aden, VALUE bnum, VALUE bden, int k)
671{
672    VALUE num, den;
673
674    if (FIXNUM_P(anum) && FIXNUM_P(aden) &&
675	FIXNUM_P(bnum) && FIXNUM_P(bden)) {
676	long an = FIX2LONG(anum);
677	long ad = FIX2LONG(aden);
678	long bn = FIX2LONG(bnum);
679	long bd = FIX2LONG(bden);
680	long ig = i_gcd(ad, bd);
681
682	VALUE g = LONG2NUM(ig);
683	VALUE a = f_imul(an, bd / ig);
684	VALUE b = f_imul(bn, ad / ig);
685	VALUE c;
686
687	if (k == '+')
688	    c = f_add(a, b);
689	else
690	    c = f_sub(a, b);
691
692	b = f_idiv(aden, g);
693	g = f_gcd(c, g);
694	num = f_idiv(c, g);
695	a = f_idiv(bden, g);
696	den = f_mul(a, b);
697    }
698    else {
699	VALUE g = f_gcd(aden, bden);
700	VALUE a = f_mul(anum, f_idiv(bden, g));
701	VALUE b = f_mul(bnum, f_idiv(aden, g));
702	VALUE c;
703
704	if (k == '+')
705	    c = f_add(a, b);
706	else
707	    c = f_sub(a, b);
708
709	b = f_idiv(aden, g);
710	g = f_gcd(c, g);
711	num = f_idiv(c, g);
712	a = f_idiv(bden, g);
713	den = f_mul(a, b);
714    }
715    return f_rational_new_no_reduce2(CLASS_OF(self), num, den);
716}
717
718/*
719 * call-seq:
720 *    rat + numeric  ->  numeric
721 *
722 * Performs addition.
723 *
724 *    Rational(2, 3)  + Rational(2, 3)   #=> (4/3)
725 *    Rational(900)   + Rational(1)      #=> (900/1)
726 *    Rational(-2, 9) + Rational(-9, 2)  #=> (-85/18)
727 *    Rational(9, 8)  + 4                #=> (41/8)
728 *    Rational(20, 9) + 9.8              #=> 12.022222222222222
729 */
730static VALUE
731nurat_add(VALUE self, VALUE other)
732{
733    switch (TYPE(other)) {
734      case T_FIXNUM:
735      case T_BIGNUM:
736	{
737	    get_dat1(self);
738
739	    return f_addsub(self,
740			    dat->num, dat->den,
741			    other, ONE, '+');
742	}
743      case T_FLOAT:
744	return f_add(f_to_f(self), other);
745      case T_RATIONAL:
746	{
747	    get_dat2(self, other);
748
749	    return f_addsub(self,
750			    adat->num, adat->den,
751			    bdat->num, bdat->den, '+');
752	}
753      default:
754	return rb_num_coerce_bin(self, other, '+');
755    }
756}
757
758/*
759 * call-seq:
760 *    rat - numeric  ->  numeric
761 *
762 * Performs subtraction.
763 *
764 *    Rational(2, 3)  - Rational(2, 3)   #=> (0/1)
765 *    Rational(900)   - Rational(1)      #=> (899/1)
766 *    Rational(-2, 9) - Rational(-9, 2)  #=> (77/18)
767 *    Rational(9, 8)  - 4                #=> (23/8)
768 *    Rational(20, 9) - 9.8              #=> -7.577777777777778
769 */
770static VALUE
771nurat_sub(VALUE self, VALUE other)
772{
773    switch (TYPE(other)) {
774      case T_FIXNUM:
775      case T_BIGNUM:
776	{
777	    get_dat1(self);
778
779	    return f_addsub(self,
780			    dat->num, dat->den,
781			    other, ONE, '-');
782	}
783      case T_FLOAT:
784	return f_sub(f_to_f(self), other);
785      case T_RATIONAL:
786	{
787	    get_dat2(self, other);
788
789	    return f_addsub(self,
790			    adat->num, adat->den,
791			    bdat->num, bdat->den, '-');
792	}
793      default:
794	return rb_num_coerce_bin(self, other, '-');
795    }
796}
797
798inline static VALUE
799f_muldiv(VALUE self, VALUE anum, VALUE aden, VALUE bnum, VALUE bden, int k)
800{
801    VALUE num, den;
802
803    if (k == '/') {
804	VALUE t;
805
806	if (f_negative_p(bnum)) {
807	    anum = f_negate(anum);
808	    bnum = f_negate(bnum);
809	}
810	t = bnum;
811	bnum = bden;
812	bden = t;
813    }
814
815    if (FIXNUM_P(anum) && FIXNUM_P(aden) &&
816	FIXNUM_P(bnum) && FIXNUM_P(bden)) {
817	long an = FIX2LONG(anum);
818	long ad = FIX2LONG(aden);
819	long bn = FIX2LONG(bnum);
820	long bd = FIX2LONG(bden);
821	long g1 = i_gcd(an, bd);
822	long g2 = i_gcd(ad, bn);
823
824	num = f_imul(an / g1, bn / g2);
825	den = f_imul(ad / g2, bd / g1);
826    }
827    else {
828	VALUE g1 = f_gcd(anum, bden);
829	VALUE g2 = f_gcd(aden, bnum);
830
831	num = f_mul(f_idiv(anum, g1), f_idiv(bnum, g2));
832	den = f_mul(f_idiv(aden, g2), f_idiv(bden, g1));
833    }
834    return f_rational_new_no_reduce2(CLASS_OF(self), num, den);
835}
836
837/*
838 * call-seq:
839 *    rat * numeric  ->  numeric
840 *
841 * Performs multiplication.
842 *
843 *    Rational(2, 3)  * Rational(2, 3)   #=> (4/9)
844 *    Rational(900)   * Rational(1)      #=> (900/1)
845 *    Rational(-2, 9) * Rational(-9, 2)  #=> (1/1)
846 *    Rational(9, 8)  * 4                #=> (9/2)
847 *    Rational(20, 9) * 9.8              #=> 21.77777777777778
848 */
849static VALUE
850nurat_mul(VALUE self, VALUE other)
851{
852    switch (TYPE(other)) {
853      case T_FIXNUM:
854      case T_BIGNUM:
855	{
856	    get_dat1(self);
857
858	    return f_muldiv(self,
859			    dat->num, dat->den,
860			    other, ONE, '*');
861	}
862      case T_FLOAT:
863	return f_mul(f_to_f(self), other);
864      case T_RATIONAL:
865	{
866	    get_dat2(self, other);
867
868	    return f_muldiv(self,
869			    adat->num, adat->den,
870			    bdat->num, bdat->den, '*');
871	}
872      default:
873	return rb_num_coerce_bin(self, other, '*');
874    }
875}
876
877/*
878 * call-seq:
879 *    rat / numeric     ->  numeric
880 *    rat.quo(numeric)  ->  numeric
881 *
882 * Performs division.
883 *
884 *    Rational(2, 3)  / Rational(2, 3)   #=> (1/1)
885 *    Rational(900)   / Rational(1)      #=> (900/1)
886 *    Rational(-2, 9) / Rational(-9, 2)  #=> (4/81)
887 *    Rational(9, 8)  / 4                #=> (9/32)
888 *    Rational(20, 9) / 9.8              #=> 0.22675736961451246
889 */
890static VALUE
891nurat_div(VALUE self, VALUE other)
892{
893    switch (TYPE(other)) {
894      case T_FIXNUM:
895      case T_BIGNUM:
896	if (f_zero_p(other))
897	    rb_raise_zerodiv();
898	{
899	    get_dat1(self);
900
901	    return f_muldiv(self,
902			    dat->num, dat->den,
903			    other, ONE, '/');
904	}
905      case T_FLOAT:
906	{
907	    double x = RFLOAT_VALUE(other), den;
908	    get_dat1(self);
909
910	    if (isnan(x)) return DBL2NUM(NAN);
911	    if (isinf(x)) return INT2FIX(0);
912	    if (x != 0.0 && modf(x, &den) == 0.0) {
913		return rb_rational_raw2(dat->num, f_mul(rb_dbl2big(den), dat->den));
914	    }
915	}
916	return rb_funcall(f_to_f(self), '/', 1, other);
917      case T_RATIONAL:
918	if (f_zero_p(other))
919	    rb_raise_zerodiv();
920	{
921	    get_dat2(self, other);
922
923	    if (f_one_p(self))
924		return f_rational_new_no_reduce2(CLASS_OF(self),
925						 bdat->den, bdat->num);
926
927	    return f_muldiv(self,
928			    adat->num, adat->den,
929			    bdat->num, bdat->den, '/');
930	}
931      default:
932	return rb_num_coerce_bin(self, other, '/');
933    }
934}
935
936/*
937 * call-seq:
938 *    rat.fdiv(numeric)  ->  float
939 *
940 * Performs division and returns the value as a float.
941 *
942 *    Rational(2, 3).fdiv(1)       #=> 0.6666666666666666
943 *    Rational(2, 3).fdiv(0.5)     #=> 1.3333333333333333
944 *    Rational(2).fdiv(3)          #=> 0.6666666666666666
945 */
946static VALUE
947nurat_fdiv(VALUE self, VALUE other)
948{
949    if (f_zero_p(other))
950	return f_div(self, f_to_f(other));
951    return f_to_f(f_div(self, other));
952}
953
954inline static VALUE
955f_odd_p(VALUE integer)
956{
957    if (rb_funcall(integer, '%', 1, INT2FIX(2)) != INT2FIX(0)) {
958	return Qtrue;
959    }
960    return Qfalse;
961
962}
963
964/*
965 * call-seq:
966 *    rat ** numeric  ->  numeric
967 *
968 * Performs exponentiation.
969 *
970 *    Rational(2)    ** Rational(3)    #=> (8/1)
971 *    Rational(10)   ** -2             #=> (1/100)
972 *    Rational(10)   ** -2.0           #=> 0.01
973 *    Rational(-4)   ** Rational(1,2)  #=> (1.2246063538223773e-16+2.0i)
974 *    Rational(1, 2) ** 0              #=> (1/1)
975 *    Rational(1, 2) ** 0.0            #=> 1.0
976 */
977static VALUE
978nurat_expt(VALUE self, VALUE other)
979{
980    if (k_numeric_p(other) && k_exact_zero_p(other))
981	return f_rational_new_bang1(CLASS_OF(self), ONE);
982
983    if (k_rational_p(other)) {
984	get_dat1(other);
985
986	if (f_one_p(dat->den))
987	    other = dat->num; /* c14n */
988    }
989
990    /* Deal with special cases of 0**n and 1**n */
991    if (k_numeric_p(other) && k_exact_p(other)) {
992	get_dat1(self);
993	if (f_one_p(dat->den))
994	    if (f_one_p(dat->num))
995		return f_rational_new_bang1(CLASS_OF(self), ONE);
996	    else if (f_minus_one_p(dat->num) && k_integer_p(other))
997		return f_rational_new_bang1(CLASS_OF(self), INT2FIX(f_odd_p(other) ? -1 : 1));
998	    else if (f_zero_p(dat->num))
999		if (FIX2INT(f_cmp(other, ZERO)) == -1)
1000		    rb_raise_zerodiv();
1001		else
1002		    return f_rational_new_bang1(CLASS_OF(self), ZERO);
1003    }
1004
1005    /* General case */
1006    switch (TYPE(other)) {
1007      case T_FIXNUM:
1008	{
1009	    VALUE num, den;
1010
1011	    get_dat1(self);
1012
1013	    switch (FIX2INT(f_cmp(other, ZERO))) {
1014	      case 1:
1015		num = f_expt(dat->num, other);
1016		den = f_expt(dat->den, other);
1017		break;
1018	      case -1:
1019		num = f_expt(dat->den, f_negate(other));
1020		den = f_expt(dat->num, f_negate(other));
1021		break;
1022	      default:
1023		num = ONE;
1024		den = ONE;
1025		break;
1026	    }
1027	    return f_rational_new2(CLASS_OF(self), num, den);
1028	}
1029      case T_BIGNUM:
1030	rb_warn("in a**b, b may be too big");
1031	/* fall through */
1032      case T_FLOAT:
1033      case T_RATIONAL:
1034	return f_expt(f_to_f(self), other);
1035      default:
1036	return rb_num_coerce_bin(self, other, id_expt);
1037    }
1038}
1039
1040/*
1041 * call-seq:
1042 *    rational <=> numeric  ->  -1, 0, +1 or nil
1043 *
1044 * Performs comparison and returns -1, 0, or +1.
1045 *
1046 * +nil+ is returned if the two values are incomparable.
1047 *
1048 *    Rational(2, 3)  <=> Rational(2, 3)  #=> 0
1049 *    Rational(5)     <=> 5               #=> 0
1050 *    Rational(2,3)   <=> Rational(1,3)   #=> 1
1051 *    Rational(1,3)   <=> 1               #=> -1
1052 *    Rational(1,3)   <=> 0.3             #=> 1
1053 */
1054static VALUE
1055nurat_cmp(VALUE self, VALUE other)
1056{
1057    switch (TYPE(other)) {
1058      case T_FIXNUM:
1059      case T_BIGNUM:
1060	{
1061	    get_dat1(self);
1062
1063	    if (FIXNUM_P(dat->den) && FIX2LONG(dat->den) == 1)
1064		return f_cmp(dat->num, other); /* c14n */
1065	    return f_cmp(self, f_rational_new_bang1(CLASS_OF(self), other));
1066	}
1067      case T_FLOAT:
1068	return f_cmp(f_to_f(self), other);
1069      case T_RATIONAL:
1070	{
1071	    VALUE num1, num2;
1072
1073	    get_dat2(self, other);
1074
1075	    if (FIXNUM_P(adat->num) && FIXNUM_P(adat->den) &&
1076		FIXNUM_P(bdat->num) && FIXNUM_P(bdat->den)) {
1077		num1 = f_imul(FIX2LONG(adat->num), FIX2LONG(bdat->den));
1078		num2 = f_imul(FIX2LONG(bdat->num), FIX2LONG(adat->den));
1079	    }
1080	    else {
1081		num1 = f_mul(adat->num, bdat->den);
1082		num2 = f_mul(bdat->num, adat->den);
1083	    }
1084	    return f_cmp(f_sub(num1, num2), ZERO);
1085	}
1086      default:
1087	return rb_num_coerce_cmp(self, other, id_cmp);
1088    }
1089}
1090
1091/*
1092 * call-seq:
1093 *    rat == object  ->  true or false
1094 *
1095 * Returns true if rat equals object numerically.
1096 *
1097 *    Rational(2, 3)  == Rational(2, 3)   #=> true
1098 *    Rational(5)     == 5                #=> true
1099 *    Rational(0)     == 0.0              #=> true
1100 *    Rational('1/3') == 0.33             #=> false
1101 *    Rational('1/2') == '1/2'            #=> false
1102 */
1103static VALUE
1104nurat_eqeq_p(VALUE self, VALUE other)
1105{
1106    switch (TYPE(other)) {
1107      case T_FIXNUM:
1108      case T_BIGNUM:
1109	{
1110	    get_dat1(self);
1111
1112	    if (f_zero_p(dat->num) && f_zero_p(other))
1113		return Qtrue;
1114
1115	    if (!FIXNUM_P(dat->den))
1116		return Qfalse;
1117	    if (FIX2LONG(dat->den) != 1)
1118		return Qfalse;
1119	    if (f_eqeq_p(dat->num, other))
1120		return Qtrue;
1121	    return Qfalse;
1122	}
1123      case T_FLOAT:
1124	return f_eqeq_p(f_to_f(self), other);
1125      case T_RATIONAL:
1126	{
1127	    get_dat2(self, other);
1128
1129	    if (f_zero_p(adat->num) && f_zero_p(bdat->num))
1130		return Qtrue;
1131
1132	    return f_boolcast(f_eqeq_p(adat->num, bdat->num) &&
1133			      f_eqeq_p(adat->den, bdat->den));
1134	}
1135      default:
1136	return f_eqeq_p(other, self);
1137    }
1138}
1139
1140/* :nodoc: */
1141static VALUE
1142nurat_coerce(VALUE self, VALUE other)
1143{
1144    switch (TYPE(other)) {
1145      case T_FIXNUM:
1146      case T_BIGNUM:
1147	return rb_assoc_new(f_rational_new_bang1(CLASS_OF(self), other), self);
1148      case T_FLOAT:
1149	return rb_assoc_new(other, f_to_f(self));
1150      case T_RATIONAL:
1151	return rb_assoc_new(other, self);
1152      case T_COMPLEX:
1153	if (k_exact_zero_p(RCOMPLEX(other)->imag))
1154	    return rb_assoc_new(f_rational_new_bang1
1155				(CLASS_OF(self), RCOMPLEX(other)->real), self);
1156	else
1157	    return rb_assoc_new(other, rb_Complex(self, INT2FIX(0)));
1158    }
1159
1160    rb_raise(rb_eTypeError, "%s can't be coerced into %s",
1161	     rb_obj_classname(other), rb_obj_classname(self));
1162    return Qnil;
1163}
1164
1165#if 0
1166/* :nodoc: */
1167static VALUE
1168nurat_idiv(VALUE self, VALUE other)
1169{
1170    return f_idiv(self, other);
1171}
1172
1173/* :nodoc: */
1174static VALUE
1175nurat_quot(VALUE self, VALUE other)
1176{
1177    return f_truncate(f_div(self, other));
1178}
1179
1180/* :nodoc: */
1181static VALUE
1182nurat_quotrem(VALUE self, VALUE other)
1183{
1184    VALUE val = f_truncate(f_div(self, other));
1185    return rb_assoc_new(val, f_sub(self, f_mul(other, val)));
1186}
1187#endif
1188
1189#if 0
1190/* :nodoc: */
1191static VALUE
1192nurat_true(VALUE self)
1193{
1194    return Qtrue;
1195}
1196#endif
1197
1198static VALUE
1199nurat_floor(VALUE self)
1200{
1201    get_dat1(self);
1202    return f_idiv(dat->num, dat->den);
1203}
1204
1205static VALUE
1206nurat_ceil(VALUE self)
1207{
1208    get_dat1(self);
1209    return f_negate(f_idiv(f_negate(dat->num), dat->den));
1210}
1211
1212/*
1213 * call-seq:
1214 *    rat.to_i  ->  integer
1215 *
1216 * Returns the truncated value as an integer.
1217 *
1218 * Equivalent to
1219 *    rat.truncate.
1220 *
1221 *    Rational(2, 3).to_i   #=> 0
1222 *    Rational(3).to_i      #=> 3
1223 *    Rational(300.6).to_i  #=> 300
1224 *    Rational(98,71).to_i  #=> 1
1225 *    Rational(-30,2).to_i  #=> -15
1226 */
1227static VALUE
1228nurat_truncate(VALUE self)
1229{
1230    get_dat1(self);
1231    if (f_negative_p(dat->num))
1232	return f_negate(f_idiv(f_negate(dat->num), dat->den));
1233    return f_idiv(dat->num, dat->den);
1234}
1235
1236static VALUE
1237nurat_round(VALUE self)
1238{
1239    VALUE num, den, neg;
1240
1241    get_dat1(self);
1242
1243    num = dat->num;
1244    den = dat->den;
1245    neg = f_negative_p(num);
1246
1247    if (neg)
1248	num = f_negate(num);
1249
1250    num = f_add(f_mul(num, TWO), den);
1251    den = f_mul(den, TWO);
1252    num = f_idiv(num, den);
1253
1254    if (neg)
1255	num = f_negate(num);
1256
1257    return num;
1258}
1259
1260static VALUE
1261f_round_common(int argc, VALUE *argv, VALUE self, VALUE (*func)(VALUE))
1262{
1263    VALUE n, b, s;
1264
1265    if (argc == 0)
1266	return (*func)(self);
1267
1268    rb_scan_args(argc, argv, "01", &n);
1269
1270    if (!k_integer_p(n))
1271	rb_raise(rb_eTypeError, "not an integer");
1272
1273    b = f_expt10(n);
1274    s = f_mul(self, b);
1275
1276    if (k_float_p(s)) {
1277	if (f_lt_p(n, ZERO))
1278	    return ZERO;
1279	return self;
1280    }
1281
1282    if (!k_rational_p(s)) {
1283	s = f_rational_new_bang1(CLASS_OF(self), s);
1284    }
1285
1286    s = (*func)(s);
1287
1288    s = f_div(f_rational_new_bang1(CLASS_OF(self), s), b);
1289
1290    if (f_lt_p(n, ONE))
1291	s = f_to_i(s);
1292
1293    return s;
1294}
1295
1296/*
1297 * call-seq:
1298 *    rat.floor               ->  integer
1299 *    rat.floor(precision=0)  ->  rational
1300 *
1301 * Returns the truncated value (toward negative infinity).
1302 *
1303 *    Rational(3).floor      #=> 3
1304 *    Rational(2, 3).floor   #=> 0
1305 *    Rational(-3, 2).floor  #=> -1
1306 *
1307 *           decimal      -  1  2  3 . 4  5  6
1308 *                          ^  ^  ^  ^   ^  ^
1309 *          precision      -3 -2 -1  0  +1 +2
1310 *
1311 *    '%f' % Rational('-123.456').floor(+1)  #=> "-123.500000"
1312 *    '%f' % Rational('-123.456').floor(-1)  #=> "-130.000000"
1313 */
1314static VALUE
1315nurat_floor_n(int argc, VALUE *argv, VALUE self)
1316{
1317    return f_round_common(argc, argv, self, nurat_floor);
1318}
1319
1320/*
1321 * call-seq:
1322 *    rat.ceil               ->  integer
1323 *    rat.ceil(precision=0)  ->  rational
1324 *
1325 * Returns the truncated value (toward positive infinity).
1326 *
1327 *    Rational(3).ceil      #=> 3
1328 *    Rational(2, 3).ceil   #=> 1
1329 *    Rational(-3, 2).ceil  #=> -1
1330 *
1331 *           decimal      -  1  2  3 . 4  5  6
1332 *                          ^  ^  ^  ^   ^  ^
1333 *          precision      -3 -2 -1  0  +1 +2
1334 *
1335 *    '%f' % Rational('-123.456').ceil(+1)  #=> "-123.400000"
1336 *    '%f' % Rational('-123.456').ceil(-1)  #=> "-120.000000"
1337 */
1338static VALUE
1339nurat_ceil_n(int argc, VALUE *argv, VALUE self)
1340{
1341    return f_round_common(argc, argv, self, nurat_ceil);
1342}
1343
1344/*
1345 * call-seq:
1346 *    rat.truncate               ->  integer
1347 *    rat.truncate(precision=0)  ->  rational
1348 *
1349 * Returns the truncated value (toward zero).
1350 *
1351 *    Rational(3).truncate      #=> 3
1352 *    Rational(2, 3).truncate   #=> 0
1353 *    Rational(-3, 2).truncate  #=> -1
1354 *
1355 *           decimal      -  1  2  3 . 4  5  6
1356 *                          ^  ^  ^  ^   ^  ^
1357 *          precision      -3 -2 -1  0  +1 +2
1358 *
1359 *    '%f' % Rational('-123.456').truncate(+1)  #=>  "-123.400000"
1360 *    '%f' % Rational('-123.456').truncate(-1)  #=>  "-120.000000"
1361 */
1362static VALUE
1363nurat_truncate_n(int argc, VALUE *argv, VALUE self)
1364{
1365    return f_round_common(argc, argv, self, nurat_truncate);
1366}
1367
1368/*
1369 * call-seq:
1370 *    rat.round               ->  integer
1371 *    rat.round(precision=0)  ->  rational
1372 *
1373 * Returns the truncated value (toward the nearest integer;
1374 * 0.5 => 1; -0.5 => -1).
1375 *
1376 *    Rational(3).round      #=> 3
1377 *    Rational(2, 3).round   #=> 1
1378 *    Rational(-3, 2).round  #=> -2
1379 *
1380 *           decimal      -  1  2  3 . 4  5  6
1381 *                          ^  ^  ^  ^   ^  ^
1382 *          precision      -3 -2 -1  0  +1 +2
1383 *
1384 *    '%f' % Rational('-123.456').round(+1)  #=> "-123.500000"
1385 *    '%f' % Rational('-123.456').round(-1)  #=> "-120.000000"
1386 */
1387static VALUE
1388nurat_round_n(int argc, VALUE *argv, VALUE self)
1389{
1390    return f_round_common(argc, argv, self, nurat_round);
1391}
1392
1393/*
1394 * call-seq:
1395 *    rat.to_f  ->  float
1396 *
1397 * Return the value as a float.
1398 *
1399 *    Rational(2).to_f      #=> 2.0
1400 *    Rational(9, 4).to_f   #=> 2.25
1401 *    Rational(-3, 4).to_f  #=> -0.75
1402 *    Rational(20, 3).to_f  #=> 6.666666666666667
1403 */
1404static VALUE
1405nurat_to_f(VALUE self)
1406{
1407    get_dat1(self);
1408    return f_fdiv(dat->num, dat->den);
1409}
1410
1411/*
1412 * call-seq:
1413 *    rat.to_r  ->  self
1414 *
1415 * Returns self.
1416 *
1417 *    Rational(2).to_r      #=> (2/1)
1418 *    Rational(-8, 6).to_r  #=> (-4/3)
1419 */
1420static VALUE
1421nurat_to_r(VALUE self)
1422{
1423    return self;
1424}
1425
1426#define id_ceil rb_intern("ceil")
1427#define f_ceil(x) rb_funcall((x), id_ceil, 0)
1428
1429#define id_quo rb_intern("quo")
1430#define f_quo(x,y) rb_funcall((x), id_quo, 1, (y))
1431
1432#define f_reciprocal(x) f_quo(ONE, (x))
1433
1434/*
1435  The algorithm here is the method described in CLISP.  Bruno Haible has
1436  graciously given permission to use this algorithm.  He says, "You can use
1437  it, if you present the following explanation of the algorithm."
1438
1439  Algorithm (recursively presented):
1440    If x is a rational number, return x.
1441    If x = 0.0, return 0.
1442    If x < 0.0, return (- (rationalize (- x))).
1443    If x > 0.0:
1444      Call (integer-decode-float x). It returns a m,e,s=1 (mantissa,
1445      exponent, sign).
1446      If m = 0 or e >= 0: return x = m*2^e.
1447      Search a rational number between a = (m-1/2)*2^e and b = (m+1/2)*2^e
1448      with smallest possible numerator and denominator.
1449      Note 1: If m is a power of 2, we ought to take a = (m-1/4)*2^e.
1450        But in this case the result will be x itself anyway, regardless of
1451        the choice of a. Therefore we can simply ignore this case.
1452      Note 2: At first, we need to consider the closed interval [a,b].
1453        but since a and b have the denominator 2^(|e|+1) whereas x itself
1454        has a denominator <= 2^|e|, we can restrict the search to the open
1455        interval (a,b).
1456      So, for given a and b (0 < a < b) we are searching a rational number
1457      y with a <= y <= b.
1458      Recursive algorithm fraction_between(a,b):
1459        c := (ceiling a)
1460        if c < b
1461          then return c       ; because a <= c < b, c integer
1462          else
1463            ; a is not integer (otherwise we would have had c = a < b)
1464            k := c-1          ; k = floor(a), k < a < b <= k+1
1465            return y = k + 1/fraction_between(1/(b-k), 1/(a-k))
1466                              ; note 1 <= 1/(b-k) < 1/(a-k)
1467
1468  You can see that we are actually computing a continued fraction expansion.
1469
1470  Algorithm (iterative):
1471    If x is rational, return x.
1472    Call (integer-decode-float x). It returns a m,e,s (mantissa,
1473      exponent, sign).
1474    If m = 0 or e >= 0, return m*2^e*s. (This includes the case x = 0.0.)
1475    Create rational numbers a := (2*m-1)*2^(e-1) and b := (2*m+1)*2^(e-1)
1476    (positive and already in lowest terms because the denominator is a
1477    power of two and the numerator is odd).
1478    Start a continued fraction expansion
1479      p[-1] := 0, p[0] := 1, q[-1] := 1, q[0] := 0, i := 0.
1480    Loop
1481      c := (ceiling a)
1482      if c >= b
1483        then k := c-1, partial_quotient(k), (a,b) := (1/(b-k),1/(a-k)),
1484             goto Loop
1485    finally partial_quotient(c).
1486    Here partial_quotient(c) denotes the iteration
1487      i := i+1, p[i] := c*p[i-1]+p[i-2], q[i] := c*q[i-1]+q[i-2].
1488    At the end, return s * (p[i]/q[i]).
1489    This rational number is already in lowest terms because
1490    p[i]*q[i-1]-p[i-1]*q[i] = (-1)^i.
1491*/
1492
1493static void
1494nurat_rationalize_internal(VALUE a, VALUE b, VALUE *p, VALUE *q)
1495{
1496    VALUE c, k, t, p0, p1, p2, q0, q1, q2;
1497
1498    p0 = ZERO;
1499    p1 = ONE;
1500    q0 = ONE;
1501    q1 = ZERO;
1502
1503    while (1) {
1504	c = f_ceil(a);
1505	if (f_lt_p(c, b))
1506	    break;
1507	k = f_sub(c, ONE);
1508	p2 = f_add(f_mul(k, p1), p0);
1509	q2 = f_add(f_mul(k, q1), q0);
1510	t = f_reciprocal(f_sub(b, k));
1511	b = f_reciprocal(f_sub(a, k));
1512	a = t;
1513	p0 = p1;
1514	q0 = q1;
1515	p1 = p2;
1516	q1 = q2;
1517    }
1518    *p = f_add(f_mul(c, p1), p0);
1519    *q = f_add(f_mul(c, q1), q0);
1520}
1521
1522/*
1523 * call-seq:
1524 *    rat.rationalize       ->  self
1525 *    rat.rationalize(eps)  ->  rational
1526 *
1527 * Returns a simpler approximation of the value if the optional
1528 * argument eps is given (rat-|eps| <= result <= rat+|eps|), self
1529 * otherwise.
1530 *
1531 *    r = Rational(5033165, 16777216)
1532 *    r.rationalize                    #=> (5033165/16777216)
1533 *    r.rationalize(Rational('0.01'))  #=> (3/10)
1534 *    r.rationalize(Rational('0.1'))   #=> (1/3)
1535 */
1536static VALUE
1537nurat_rationalize(int argc, VALUE *argv, VALUE self)
1538{
1539    VALUE e, a, b, p, q;
1540
1541    if (argc == 0)
1542	return self;
1543
1544    if (f_negative_p(self))
1545	return f_negate(nurat_rationalize(argc, argv, f_abs(self)));
1546
1547    rb_scan_args(argc, argv, "01", &e);
1548    e = f_abs(e);
1549    a = f_sub(self, e);
1550    b = f_add(self, e);
1551
1552    if (f_eqeq_p(a, b))
1553	return self;
1554
1555    nurat_rationalize_internal(a, b, &p, &q);
1556    return f_rational_new2(CLASS_OF(self), p, q);
1557}
1558
1559/* :nodoc: */
1560static VALUE
1561nurat_hash(VALUE self)
1562{
1563    st_index_t v, h[2];
1564    VALUE n;
1565
1566    get_dat1(self);
1567    n = rb_hash(dat->num);
1568    h[0] = NUM2LONG(n);
1569    n = rb_hash(dat->den);
1570    h[1] = NUM2LONG(n);
1571    v = rb_memhash(h, sizeof(h));
1572    return LONG2FIX(v);
1573}
1574
1575static VALUE
1576f_format(VALUE self, VALUE (*func)(VALUE))
1577{
1578    VALUE s;
1579    get_dat1(self);
1580
1581    s = (*func)(dat->num);
1582    rb_str_cat2(s, "/");
1583    rb_str_concat(s, (*func)(dat->den));
1584
1585    return s;
1586}
1587
1588/*
1589 * call-seq:
1590 *    rat.to_s  ->  string
1591 *
1592 * Returns the value as a string.
1593 *
1594 *    Rational(2).to_s      #=> "2/1"
1595 *    Rational(-8, 6).to_s  #=> "-4/3"
1596 *    Rational('1/2').to_s  #=> "1/2"
1597 */
1598static VALUE
1599nurat_to_s(VALUE self)
1600{
1601    return f_format(self, f_to_s);
1602}
1603
1604/*
1605 * call-seq:
1606 *    rat.inspect  ->  string
1607 *
1608 * Returns the value as a string for inspection.
1609 *
1610 *    Rational(2).inspect      #=> "(2/1)"
1611 *    Rational(-8, 6).inspect  #=> "(-4/3)"
1612 *    Rational('1/2').inspect  #=> "(1/2)"
1613 */
1614static VALUE
1615nurat_inspect(VALUE self)
1616{
1617    VALUE s;
1618
1619    s = rb_usascii_str_new2("(");
1620    rb_str_concat(s, f_format(self, f_inspect));
1621    rb_str_cat2(s, ")");
1622
1623    return s;
1624}
1625
1626/* :nodoc: */
1627static VALUE
1628nurat_dumper(VALUE self)
1629{
1630    return self;
1631}
1632
1633/* :nodoc: */
1634static VALUE
1635nurat_loader(VALUE self, VALUE a)
1636{
1637    get_dat1(self);
1638
1639    dat->num = rb_ivar_get(a, id_i_num);
1640    dat->den = rb_ivar_get(a, id_i_den);
1641
1642    return self;
1643}
1644
1645/* :nodoc: */
1646static VALUE
1647nurat_marshal_dump(VALUE self)
1648{
1649    VALUE a;
1650    get_dat1(self);
1651
1652    a = rb_assoc_new(dat->num, dat->den);
1653    rb_copy_generic_ivar(a, self);
1654    return a;
1655}
1656
1657/* :nodoc: */
1658static VALUE
1659nurat_marshal_load(VALUE self, VALUE a)
1660{
1661    rb_check_frozen(self);
1662    rb_check_trusted(self);
1663
1664    Check_Type(a, T_ARRAY);
1665    if (RARRAY_LEN(a) != 2)
1666	rb_raise(rb_eArgError, "marshaled rational must have an array whose length is 2 but %ld", RARRAY_LEN(a));
1667    if (f_zero_p(RARRAY_PTR(a)[1]))
1668	rb_raise_zerodiv();
1669
1670    rb_ivar_set(self, id_i_num, RARRAY_PTR(a)[0]);
1671    rb_ivar_set(self, id_i_den, RARRAY_PTR(a)[1]);
1672
1673    return self;
1674}
1675
1676/* --- */
1677
1678VALUE
1679rb_rational_reciprocal(VALUE x)
1680{
1681    get_dat1(x);
1682    return f_rational_new_no_reduce2(CLASS_OF(x), dat->den, dat->num);
1683}
1684
1685/*
1686 * call-seq:
1687 *    int.gcd(int2)  ->  integer
1688 *
1689 * Returns the greatest common divisor (always positive).  0.gcd(x)
1690 * and x.gcd(0) return abs(x).
1691 *
1692 *    2.gcd(2)                    #=> 2
1693 *    3.gcd(-7)                   #=> 1
1694 *    ((1<<31)-1).gcd((1<<61)-1)  #=> 1
1695 */
1696VALUE
1697rb_gcd(VALUE self, VALUE other)
1698{
1699    other = nurat_int_value(other);
1700    return f_gcd(self, other);
1701}
1702
1703/*
1704 * call-seq:
1705 *    int.lcm(int2)  ->  integer
1706 *
1707 * Returns the least common multiple (always positive).  0.lcm(x) and
1708 * x.lcm(0) return zero.
1709 *
1710 *    2.lcm(2)                    #=> 2
1711 *    3.lcm(-7)                   #=> 21
1712 *    ((1<<31)-1).lcm((1<<61)-1)  #=> 4951760154835678088235319297
1713 */
1714VALUE
1715rb_lcm(VALUE self, VALUE other)
1716{
1717    other = nurat_int_value(other);
1718    return f_lcm(self, other);
1719}
1720
1721/*
1722 * call-seq:
1723 *    int.gcdlcm(int2)  ->  array
1724 *
1725 * Returns an array; [int.gcd(int2), int.lcm(int2)].
1726 *
1727 *    2.gcdlcm(2)                    #=> [2, 2]
1728 *    3.gcdlcm(-7)                   #=> [1, 21]
1729 *    ((1<<31)-1).gcdlcm((1<<61)-1)  #=> [1, 4951760154835678088235319297]
1730 */
1731VALUE
1732rb_gcdlcm(VALUE self, VALUE other)
1733{
1734    other = nurat_int_value(other);
1735    return rb_assoc_new(f_gcd(self, other), f_lcm(self, other));
1736}
1737
1738VALUE
1739rb_rational_raw(VALUE x, VALUE y)
1740{
1741    return nurat_s_new_internal(rb_cRational, x, y);
1742}
1743
1744VALUE
1745rb_rational_new(VALUE x, VALUE y)
1746{
1747    return nurat_s_canonicalize_internal(rb_cRational, x, y);
1748}
1749
1750static VALUE nurat_s_convert(int argc, VALUE *argv, VALUE klass);
1751
1752VALUE
1753rb_Rational(VALUE x, VALUE y)
1754{
1755    VALUE a[2];
1756    a[0] = x;
1757    a[1] = y;
1758    return nurat_s_convert(2, a, rb_cRational);
1759}
1760
1761#define id_numerator rb_intern("numerator")
1762#define f_numerator(x) rb_funcall((x), id_numerator, 0)
1763
1764#define id_denominator rb_intern("denominator")
1765#define f_denominator(x) rb_funcall((x), id_denominator, 0)
1766
1767#define id_to_r rb_intern("to_r")
1768#define f_to_r(x) rb_funcall((x), id_to_r, 0)
1769
1770/*
1771 * call-seq:
1772 *    num.numerator  ->  integer
1773 *
1774 * Returns the numerator.
1775 */
1776static VALUE
1777numeric_numerator(VALUE self)
1778{
1779    return f_numerator(f_to_r(self));
1780}
1781
1782/*
1783 * call-seq:
1784 *    num.denominator  ->  integer
1785 *
1786 * Returns the denominator (always positive).
1787 */
1788static VALUE
1789numeric_denominator(VALUE self)
1790{
1791    return f_denominator(f_to_r(self));
1792}
1793
1794/*
1795 * call-seq:
1796 *    int.numerator  ->  self
1797 *
1798 * Returns self.
1799 */
1800static VALUE
1801integer_numerator(VALUE self)
1802{
1803    return self;
1804}
1805
1806/*
1807 * call-seq:
1808 *    int.denominator  ->  1
1809 *
1810 * Returns 1.
1811 */
1812static VALUE
1813integer_denominator(VALUE self)
1814{
1815    return INT2FIX(1);
1816}
1817
1818/*
1819 * call-seq:
1820 *    flo.numerator  ->  integer
1821 *
1822 * Returns the numerator.  The result is machine dependent.
1823 *
1824 *    n = 0.3.numerator    #=> 5404319552844595
1825 *    d = 0.3.denominator  #=> 18014398509481984
1826 *    n.fdiv(d)            #=> 0.3
1827 */
1828static VALUE
1829float_numerator(VALUE self)
1830{
1831    double d = RFLOAT_VALUE(self);
1832    if (isinf(d) || isnan(d))
1833	return self;
1834    return rb_call_super(0, 0);
1835}
1836
1837/*
1838 * call-seq:
1839 *    flo.denominator  ->  integer
1840 *
1841 * Returns the denominator (always positive).  The result is machine
1842 * dependent.
1843 *
1844 * See numerator.
1845 */
1846static VALUE
1847float_denominator(VALUE self)
1848{
1849    double d = RFLOAT_VALUE(self);
1850    if (isinf(d) || isnan(d))
1851	return INT2FIX(1);
1852    return rb_call_super(0, 0);
1853}
1854
1855/*
1856 * call-seq:
1857 *    nil.to_r  ->  (0/1)
1858 *
1859 * Returns zero as a rational.
1860 */
1861static VALUE
1862nilclass_to_r(VALUE self)
1863{
1864    return rb_rational_new1(INT2FIX(0));
1865}
1866
1867/*
1868 * call-seq:
1869 *    nil.rationalize([eps])  ->  (0/1)
1870 *
1871 * Returns zero as a rational.  The optional argument eps is always
1872 * ignored.
1873 */
1874static VALUE
1875nilclass_rationalize(int argc, VALUE *argv, VALUE self)
1876{
1877    rb_scan_args(argc, argv, "01", NULL);
1878    return nilclass_to_r(self);
1879}
1880
1881/*
1882 * call-seq:
1883 *    int.to_r  ->  rational
1884 *
1885 * Returns the value as a rational.
1886 *
1887 *    1.to_r        #=> (1/1)
1888 *    (1<<64).to_r  #=> (18446744073709551616/1)
1889 */
1890static VALUE
1891integer_to_r(VALUE self)
1892{
1893    return rb_rational_new1(self);
1894}
1895
1896/*
1897 * call-seq:
1898 *    int.rationalize([eps])  ->  rational
1899 *
1900 * Returns the value as a rational.  The optional argument eps is
1901 * always ignored.
1902 */
1903static VALUE
1904integer_rationalize(int argc, VALUE *argv, VALUE self)
1905{
1906    rb_scan_args(argc, argv, "01", NULL);
1907    return integer_to_r(self);
1908}
1909
1910static void
1911float_decode_internal(VALUE self, VALUE *rf, VALUE *rn)
1912{
1913    double f;
1914    int n;
1915
1916    f = frexp(RFLOAT_VALUE(self), &n);
1917    f = ldexp(f, DBL_MANT_DIG);
1918    n -= DBL_MANT_DIG;
1919    *rf = rb_dbl2big(f);
1920    *rn = INT2FIX(n);
1921}
1922
1923#if 0
1924static VALUE
1925float_decode(VALUE self)
1926{
1927    VALUE f, n;
1928
1929    float_decode_internal(self, &f, &n);
1930    return rb_assoc_new(f, n);
1931}
1932#endif
1933
1934#define id_lshift rb_intern("<<")
1935#define f_lshift(x,n) rb_funcall((x), id_lshift, 1, (n))
1936
1937/*
1938 * call-seq:
1939 *    flt.to_r  ->  rational
1940 *
1941 * Returns the value as a rational.
1942 *
1943 * NOTE: 0.3.to_r isn't the same as '0.3'.to_r.  The latter is
1944 * equivalent to '3/10'.to_r, but the former isn't so.
1945 *
1946 *    2.0.to_r    #=> (2/1)
1947 *    2.5.to_r    #=> (5/2)
1948 *    -0.75.to_r  #=> (-3/4)
1949 *    0.0.to_r    #=> (0/1)
1950 *
1951 * See rationalize.
1952 */
1953static VALUE
1954float_to_r(VALUE self)
1955{
1956    VALUE f, n;
1957
1958    float_decode_internal(self, &f, &n);
1959#if FLT_RADIX == 2
1960    {
1961	long ln = FIX2LONG(n);
1962
1963	if (ln == 0)
1964	    return f_to_r(f);
1965	if (ln > 0)
1966	    return f_to_r(f_lshift(f, n));
1967	ln = -ln;
1968	return rb_rational_new2(f, f_lshift(ONE, INT2FIX(ln)));
1969    }
1970#else
1971    return f_to_r(f_mul(f, f_expt(INT2FIX(FLT_RADIX), n)));
1972#endif
1973}
1974
1975/*
1976 * call-seq:
1977 *    flt.rationalize([eps])  ->  rational
1978 *
1979 * Returns a simpler approximation of the value (flt-|eps| <= result
1980 * <= flt+|eps|).  if the optional eps is not given, it will be chosen
1981 * automatically.
1982 *
1983 *    0.3.rationalize          #=> (3/10)
1984 *    1.333.rationalize        #=> (1333/1000)
1985 *    1.333.rationalize(0.01)  #=> (4/3)
1986 *
1987 * See to_r.
1988 */
1989static VALUE
1990float_rationalize(int argc, VALUE *argv, VALUE self)
1991{
1992    VALUE e, a, b, p, q;
1993
1994    if (f_negative_p(self))
1995	return f_negate(float_rationalize(argc, argv, f_abs(self)));
1996
1997    rb_scan_args(argc, argv, "01", &e);
1998
1999    if (argc != 0) {
2000	e = f_abs(e);
2001	a = f_sub(self, e);
2002	b = f_add(self, e);
2003    }
2004    else {
2005	VALUE f, n;
2006
2007	float_decode_internal(self, &f, &n);
2008	if (f_zero_p(f) || f_positive_p(n))
2009	    return rb_rational_new1(f_lshift(f, n));
2010
2011#if FLT_RADIX == 2
2012	{
2013	    VALUE two_times_f, den;
2014
2015	    two_times_f = f_mul(TWO, f);
2016	    den = f_lshift(ONE, f_sub(ONE, n));
2017
2018	    a = rb_rational_new2(f_sub(two_times_f, ONE), den);
2019	    b = rb_rational_new2(f_add(two_times_f, ONE), den);
2020	}
2021#else
2022	{
2023	    VALUE radix_times_f, den;
2024
2025	    radix_times_f = f_mul(INT2FIX(FLT_RADIX), f);
2026	    den = f_expt(INT2FIX(FLT_RADIX), f_sub(ONE, n));
2027
2028	    a = rb_rational_new2(f_sub(radix_times_f, INT2FIX(FLT_RADIX - 1)), den);
2029	    b = rb_rational_new2(f_add(radix_times_f, INT2FIX(FLT_RADIX - 1)), den);
2030	}
2031#endif
2032    }
2033
2034    if (f_eqeq_p(a, b))
2035	return f_to_r(self);
2036
2037    nurat_rationalize_internal(a, b, &p, &q);
2038    return rb_rational_new2(p, q);
2039}
2040
2041#include <ctype.h>
2042
2043inline static int
2044issign(int c)
2045{
2046    return (c == '-' || c == '+');
2047}
2048
2049static int
2050read_sign(const char **s)
2051{
2052    int sign = '?';
2053
2054    if (issign(**s)) {
2055	sign = **s;
2056	(*s)++;
2057    }
2058    return sign;
2059}
2060
2061inline static int
2062isdecimal(int c)
2063{
2064    return isdigit((unsigned char)c);
2065}
2066
2067static int
2068read_digits(const char **s, int strict,
2069	    VALUE *num, int *count)
2070{
2071    char *b, *bb;
2072    int us = 1, ret = 1;
2073
2074    if (!isdecimal(**s)) {
2075	*num = ZERO;
2076	return 0;
2077    }
2078
2079    bb = b = ALLOCA_N(char, strlen(*s) + 1);
2080
2081    while (isdecimal(**s) || **s == '_') {
2082	if (**s == '_') {
2083	    if (strict) {
2084		if (us) {
2085		    ret = 0;
2086		    goto conv;
2087		}
2088	    }
2089	    us = 1;
2090	}
2091	else {
2092	    if (count)
2093		(*count)++;
2094	    *b++ = **s;
2095	    us = 0;
2096	}
2097	(*s)++;
2098    }
2099    if (us)
2100	do {
2101	    (*s)--;
2102	} while (**s == '_');
2103  conv:
2104    *b = '\0';
2105    *num = rb_cstr_to_inum(bb, 10, 0);
2106    return ret;
2107}
2108
2109inline static int
2110islettere(int c)
2111{
2112    return (c == 'e' || c == 'E');
2113}
2114
2115static int
2116read_num(const char **s, int numsign, int strict,
2117	 VALUE *num)
2118{
2119    VALUE ip, fp, exp;
2120
2121    *num = rb_rational_new2(ZERO, ONE);
2122    exp = Qnil;
2123
2124    if (**s != '.') {
2125	if (!read_digits(s, strict, &ip, NULL))
2126	    return 0;
2127	*num = rb_rational_new2(ip, ONE);
2128    }
2129
2130    if (**s == '.') {
2131	int count = 0;
2132
2133	(*s)++;
2134	if (!read_digits(s, strict, &fp, &count))
2135	    return 0;
2136	{
2137	    VALUE l = f_expt10(INT2NUM(count));
2138	    *num = f_mul(*num, l);
2139	    *num = f_add(*num, fp);
2140	    *num = f_div(*num, l);
2141	}
2142    }
2143
2144    if (islettere(**s)) {
2145	int expsign;
2146
2147	(*s)++;
2148	expsign = read_sign(s);
2149	if (!read_digits(s, strict, &exp, NULL))
2150	    return 0;
2151	if (expsign == '-')
2152	    exp = f_negate(exp);
2153    }
2154
2155    if (numsign == '-')
2156	*num = f_negate(*num);
2157    if (!NIL_P(exp)) {
2158	VALUE l = f_expt10(exp);
2159	*num = f_mul(*num, l);
2160    }
2161    return 1;
2162}
2163
2164inline static int
2165read_den(const char **s, int strict,
2166	 VALUE *num)
2167{
2168    if (!read_digits(s, strict, num, NULL))
2169	return 0;
2170    return 1;
2171}
2172
2173static int
2174read_rat_nos(const char **s, int sign, int strict,
2175	     VALUE *num)
2176{
2177    VALUE den;
2178
2179    if (!read_num(s, sign, strict, num))
2180	return 0;
2181    if (**s == '/') {
2182	(*s)++;
2183	if (!read_den(s, strict, &den))
2184	    return 0;
2185	if (!(FIXNUM_P(den) && FIX2LONG(den) == 1))
2186	    *num = f_div(*num, den);
2187    }
2188    return 1;
2189}
2190
2191static int
2192read_rat(const char **s, int strict,
2193	 VALUE *num)
2194{
2195    int sign;
2196
2197    sign = read_sign(s);
2198    if (!read_rat_nos(s, sign, strict, num))
2199	return 0;
2200    return 1;
2201}
2202
2203inline static void
2204skip_ws(const char **s)
2205{
2206    while (isspace((unsigned char)**s))
2207	(*s)++;
2208}
2209
2210static int
2211parse_rat(const char *s, int strict,
2212	  VALUE *num)
2213{
2214    skip_ws(&s);
2215    if (!read_rat(&s, strict, num))
2216	return 0;
2217    skip_ws(&s);
2218
2219    if (strict)
2220	if (*s != '\0')
2221	    return 0;
2222    return 1;
2223}
2224
2225static VALUE
2226string_to_r_strict(VALUE self)
2227{
2228    char *s;
2229    VALUE num;
2230
2231    rb_must_asciicompat(self);
2232
2233    s = RSTRING_PTR(self);
2234
2235    if (!s || memchr(s, '\0', RSTRING_LEN(self)))
2236	rb_raise(rb_eArgError, "string contains null byte");
2237
2238    if (s && s[RSTRING_LEN(self)]) {
2239	rb_str_modify(self);
2240	s = RSTRING_PTR(self);
2241	s[RSTRING_LEN(self)] = '\0';
2242    }
2243
2244    if (!s)
2245	s = (char *)"";
2246
2247    if (!parse_rat(s, 1, &num)) {
2248	VALUE ins = f_inspect(self);
2249	rb_raise(rb_eArgError, "invalid value for convert(): %s",
2250		 StringValuePtr(ins));
2251    }
2252
2253    if (RB_TYPE_P(num, T_FLOAT))
2254	rb_raise(rb_eFloatDomainError, "Infinity");
2255    return num;
2256}
2257
2258/*
2259 * call-seq:
2260 *    str.to_r  ->  rational
2261 *
2262 * Returns a rational which denotes the string form.  The parser
2263 * ignores leading whitespaces and trailing garbage.  Any digit
2264 * sequences can be separated by an underscore.  Returns zero for null
2265 * or garbage string.
2266 *
2267 * NOTE: '0.3'.to_r isn't the same as 0.3.to_r.  The former is
2268 * equivalent to '3/10'.to_r, but the latter isn't so.
2269 *
2270 *    '  2  '.to_r       #=> (2/1)
2271 *    '300/2'.to_r       #=> (150/1)
2272 *    '-9.2'.to_r        #=> (-46/5)
2273 *    '-9.2e2'.to_r      #=> (-920/1)
2274 *    '1_234_567'.to_r   #=> (1234567/1)
2275 *    '21 june 09'.to_r  #=> (21/1)
2276 *    '21/06/09'.to_r    #=> (7/2)
2277 *    'bwv 1079'.to_r    #=> (0/1)
2278 *
2279 * See Kernel.Rational.
2280 */
2281static VALUE
2282string_to_r(VALUE self)
2283{
2284    char *s;
2285    VALUE num;
2286
2287    rb_must_asciicompat(self);
2288
2289    s = RSTRING_PTR(self);
2290
2291    if (s && s[RSTRING_LEN(self)]) {
2292	rb_str_modify(self);
2293	s = RSTRING_PTR(self);
2294	s[RSTRING_LEN(self)] = '\0';
2295    }
2296
2297    if (!s)
2298	s = (char *)"";
2299
2300    (void)parse_rat(s, 0, &num);
2301
2302    if (RB_TYPE_P(num, T_FLOAT))
2303	rb_raise(rb_eFloatDomainError, "Infinity");
2304    return num;
2305}
2306
2307VALUE
2308rb_cstr_to_rat(const char *s, int strict) /* for complex's internal */
2309{
2310    VALUE num;
2311
2312    (void)parse_rat(s, strict, &num);
2313
2314    if (RB_TYPE_P(num, T_FLOAT))
2315	rb_raise(rb_eFloatDomainError, "Infinity");
2316    return num;
2317}
2318
2319static VALUE
2320nurat_s_convert(int argc, VALUE *argv, VALUE klass)
2321{
2322    VALUE a1, a2, backref;
2323
2324    rb_scan_args(argc, argv, "11", &a1, &a2);
2325
2326    if (NIL_P(a1) || (argc == 2 && NIL_P(a2)))
2327	rb_raise(rb_eTypeError, "can't convert nil into Rational");
2328
2329    switch (TYPE(a1)) {
2330      case T_COMPLEX:
2331	if (k_exact_zero_p(RCOMPLEX(a1)->imag))
2332	    a1 = RCOMPLEX(a1)->real;
2333    }
2334
2335    switch (TYPE(a2)) {
2336      case T_COMPLEX:
2337	if (k_exact_zero_p(RCOMPLEX(a2)->imag))
2338	    a2 = RCOMPLEX(a2)->real;
2339    }
2340
2341    backref = rb_backref_get();
2342    rb_match_busy(backref);
2343
2344    switch (TYPE(a1)) {
2345      case T_FIXNUM:
2346      case T_BIGNUM:
2347	break;
2348      case T_FLOAT:
2349	a1 = f_to_r(a1);
2350	break;
2351      case T_STRING:
2352	a1 = string_to_r_strict(a1);
2353	break;
2354    }
2355
2356    switch (TYPE(a2)) {
2357      case T_FIXNUM:
2358      case T_BIGNUM:
2359	break;
2360      case T_FLOAT:
2361	a2 = f_to_r(a2);
2362	break;
2363      case T_STRING:
2364	a2 = string_to_r_strict(a2);
2365	break;
2366    }
2367
2368    rb_backref_set(backref);
2369
2370    switch (TYPE(a1)) {
2371      case T_RATIONAL:
2372	if (argc == 1 || (k_exact_one_p(a2)))
2373	    return a1;
2374    }
2375
2376    if (argc == 1) {
2377	if (!(k_numeric_p(a1) && k_integer_p(a1)))
2378	    return rb_convert_type(a1, T_RATIONAL, "Rational", "to_r");
2379    }
2380    else {
2381	if ((k_numeric_p(a1) && k_numeric_p(a2)) &&
2382	    (!f_integer_p(a1) || !f_integer_p(a2)))
2383	    return f_div(a1, a2);
2384    }
2385
2386    {
2387	VALUE argv2[2];
2388	argv2[0] = a1;
2389	argv2[1] = a2;
2390	return nurat_s_new(argc, argv2, klass);
2391    }
2392}
2393
2394/*
2395 * A rational number can be represented as a paired integer number;
2396 * a/b (b>0).  Where a is numerator and b is denominator.  Integer a
2397 * equals rational a/1 mathematically.
2398 *
2399 * In ruby, you can create rational object with Rational, to_r or
2400 * rationalize method.  The return values will be irreducible.
2401 *
2402 *    Rational(1)      #=> (1/1)
2403 *    Rational(2, 3)   #=> (2/3)
2404 *    Rational(4, -6)  #=> (-2/3)
2405 *    3.to_r           #=> (3/1)
2406 *
2407 * You can also create rational object from floating-point numbers or
2408 * strings.
2409 *
2410 *    Rational(0.3)    #=> (5404319552844595/18014398509481984)
2411 *    Rational('0.3')  #=> (3/10)
2412 *    Rational('2/3')  #=> (2/3)
2413 *
2414 *    0.3.to_r         #=> (5404319552844595/18014398509481984)
2415 *    '0.3'.to_r       #=> (3/10)
2416 *    '2/3'.to_r       #=> (2/3)
2417 *    0.3.rationalize  #=> (3/10)
2418 *
2419 * A rational object is an exact number, which helps you to write
2420 * program without any rounding errors.
2421 *
2422 *    10.times.inject(0){|t,| t + 0.1}              #=> 0.9999999999999999
2423 *    10.times.inject(0){|t,| t + Rational('0.1')}  #=> (1/1)
2424 *
2425 * However, when an expression has inexact factor (numerical value or
2426 * operation), will produce an inexact result.
2427 *
2428 *    Rational(10) / 3   #=> (10/3)
2429 *    Rational(10) / 3.0 #=> 3.3333333333333335
2430 *
2431 *    Rational(-8) ** Rational(1, 3)
2432 *                       #=> (1.0000000000000002+1.7320508075688772i)
2433 */
2434void
2435Init_Rational(void)
2436{
2437    VALUE compat;
2438#undef rb_intern
2439#define rb_intern(str) rb_intern_const(str)
2440
2441    assert(fprintf(stderr, "assert() is now active\n"));
2442
2443    id_abs = rb_intern("abs");
2444    id_cmp = rb_intern("<=>");
2445    id_convert = rb_intern("convert");
2446    id_eqeq_p = rb_intern("==");
2447    id_expt = rb_intern("**");
2448    id_fdiv = rb_intern("fdiv");
2449    id_floor = rb_intern("floor");
2450    id_idiv = rb_intern("div");
2451    id_inspect = rb_intern("inspect");
2452    id_integer_p = rb_intern("integer?");
2453    id_negate = rb_intern("-@");
2454    id_to_f = rb_intern("to_f");
2455    id_to_i = rb_intern("to_i");
2456    id_to_s = rb_intern("to_s");
2457    id_truncate = rb_intern("truncate");
2458    id_i_num = rb_intern("@numerator");
2459    id_i_den = rb_intern("@denominator");
2460
2461    rb_cRational = rb_define_class("Rational", rb_cNumeric);
2462
2463    rb_define_alloc_func(rb_cRational, nurat_s_alloc);
2464    rb_undef_method(CLASS_OF(rb_cRational), "allocate");
2465
2466#if 0
2467    rb_define_private_method(CLASS_OF(rb_cRational), "new!", nurat_s_new_bang, -1);
2468    rb_define_private_method(CLASS_OF(rb_cRational), "new", nurat_s_new, -1);
2469#else
2470    rb_undef_method(CLASS_OF(rb_cRational), "new");
2471#endif
2472
2473    rb_define_global_function("Rational", nurat_f_rational, -1);
2474
2475    rb_define_method(rb_cRational, "numerator", nurat_numerator, 0);
2476    rb_define_method(rb_cRational, "denominator", nurat_denominator, 0);
2477
2478    rb_define_method(rb_cRational, "+", nurat_add, 1);
2479    rb_define_method(rb_cRational, "-", nurat_sub, 1);
2480    rb_define_method(rb_cRational, "*", nurat_mul, 1);
2481    rb_define_method(rb_cRational, "/", nurat_div, 1);
2482    rb_define_method(rb_cRational, "quo", nurat_div, 1);
2483    rb_define_method(rb_cRational, "fdiv", nurat_fdiv, 1);
2484    rb_define_method(rb_cRational, "**", nurat_expt, 1);
2485
2486    rb_define_method(rb_cRational, "<=>", nurat_cmp, 1);
2487    rb_define_method(rb_cRational, "==", nurat_eqeq_p, 1);
2488    rb_define_method(rb_cRational, "coerce", nurat_coerce, 1);
2489
2490#if 0 /* NUBY */
2491    rb_define_method(rb_cRational, "//", nurat_idiv, 1);
2492#endif
2493
2494#if 0
2495    rb_define_method(rb_cRational, "quot", nurat_quot, 1);
2496    rb_define_method(rb_cRational, "quotrem", nurat_quotrem, 1);
2497#endif
2498
2499#if 0
2500    rb_define_method(rb_cRational, "rational?", nurat_true, 0);
2501    rb_define_method(rb_cRational, "exact?", nurat_true, 0);
2502#endif
2503
2504    rb_define_method(rb_cRational, "floor", nurat_floor_n, -1);
2505    rb_define_method(rb_cRational, "ceil", nurat_ceil_n, -1);
2506    rb_define_method(rb_cRational, "truncate", nurat_truncate_n, -1);
2507    rb_define_method(rb_cRational, "round", nurat_round_n, -1);
2508
2509    rb_define_method(rb_cRational, "to_i", nurat_truncate, 0);
2510    rb_define_method(rb_cRational, "to_f", nurat_to_f, 0);
2511    rb_define_method(rb_cRational, "to_r", nurat_to_r, 0);
2512    rb_define_method(rb_cRational, "rationalize", nurat_rationalize, -1);
2513
2514    rb_define_method(rb_cRational, "hash", nurat_hash, 0);
2515
2516    rb_define_method(rb_cRational, "to_s", nurat_to_s, 0);
2517    rb_define_method(rb_cRational, "inspect", nurat_inspect, 0);
2518
2519    rb_define_private_method(rb_cRational, "marshal_dump", nurat_marshal_dump, 0);
2520    compat = rb_define_class_under(rb_cRational, "compatible", rb_cObject);
2521    rb_define_private_method(compat, "marshal_load", nurat_marshal_load, 1);
2522    rb_marshal_define_compat(rb_cRational, compat, nurat_dumper, nurat_loader);
2523
2524    /* --- */
2525
2526    rb_define_method(rb_cInteger, "gcd", rb_gcd, 1);
2527    rb_define_method(rb_cInteger, "lcm", rb_lcm, 1);
2528    rb_define_method(rb_cInteger, "gcdlcm", rb_gcdlcm, 1);
2529
2530    rb_define_method(rb_cNumeric, "numerator", numeric_numerator, 0);
2531    rb_define_method(rb_cNumeric, "denominator", numeric_denominator, 0);
2532
2533    rb_define_method(rb_cInteger, "numerator", integer_numerator, 0);
2534    rb_define_method(rb_cInteger, "denominator", integer_denominator, 0);
2535
2536    rb_define_method(rb_cFloat, "numerator", float_numerator, 0);
2537    rb_define_method(rb_cFloat, "denominator", float_denominator, 0);
2538
2539    rb_define_method(rb_cNilClass, "to_r", nilclass_to_r, 0);
2540    rb_define_method(rb_cNilClass, "rationalize", nilclass_rationalize, -1);
2541    rb_define_method(rb_cInteger, "to_r", integer_to_r, 0);
2542    rb_define_method(rb_cInteger, "rationalize", integer_rationalize, -1);
2543    rb_define_method(rb_cFloat, "to_r", float_to_r, 0);
2544    rb_define_method(rb_cFloat, "rationalize", float_rationalize, -1);
2545
2546    rb_define_method(rb_cString, "to_r", string_to_r, 0);
2547
2548    rb_define_private_method(CLASS_OF(rb_cRational), "convert", nurat_s_convert, -1);
2549}
2550
2551/*
2552Local variables:
2553c-file-style: "ruby"
2554End:
2555*/
2556