--* From RMC%WATSON.vnet.ibm.com@yktvmv.watson.ibm.com  Wed May 25 19:55:33 1994
--* Received: from yktvmv-ob.watson.ibm.com by asharp.watson.ibm.com (AIX 3.2/UCB 5.64/930311)
--*           id AA16736; Wed, 25 May 1994 19:55:33 -0400
--* Received: from watson.vnet.ibm.com by yktvmv.watson.ibm.com (IBM VM SMTP V2R3)
--*    with BSMTP id 2959; Wed, 25 May 94 19:55:34 EDT
--* Received: from YKTVMV by watson.vnet.ibm.com with "VAGENT.V1.0"
--*           id <A.RMC.NOTE.VAGENT2.2019.May.25.19:55:33.-0400>
--*           for asbugs@watson; Wed, 25 May 94 19:55:33 -0400
--* Received: from YKTVMV by watson.vnet.ibm.com with "VAGENT.V1.0"
--*           id 2015; Wed, 25 May 1994 19:55:33 EDT
--* Received: from matteo.watson.ibm.com by yktvmv.watson.ibm.com (IBM VM SMTP V2R3)
--*    with TCP; Wed, 25 May 94 19:55:32 EDT
--* Received: by matteo.watson.ibm.com (AIX 3.2/UCB 5.64/920123)
--*           id AA28239; Wed, 25 May 1994 19:53:13 -0400
--* Date: Wed, 25 May 1994 19:53:13 -0400
--* From: rmc@matteo.watson.ibm.com
--* Message-Id: <9405252353.AA28239@matteo.watson.ibm.com>
--* To: asbugs@watson.ibm.com
--* Subject: [2] Bug: Bad case 3 in terrorAssignOrSetBang [AXRR.as][0.35.3]

--@ Fixed  by:  PI   Tue Aug 9 15:47:40 EDT 1994 
--@ Tested by:  none 
--@ Summary:    I fixed this bug times ago. 


--
-- First new type: Rounded Rationals.  RMC 25/5/1994
--
-- Basic idea: use continued fractions to keep rationals short.
--

-- Lessons.
-- 1) Rep is a type (capitalizations indicate types).
--    RoundedRatio is also a type, by the convention.
--    rep, on the other hand, is a procedure to cast something
--    into its representation.
--
-- 2) Records are easy.
--
-- 3) Exports:  some are automatic, for example if you are building an Ordered Ring or Field.
--    In that case, things like %^Integer are necessary to make.  (Though you may make them
--    give an error message, as Ratio does).
--
-- 4)

#include "axiom.as"


-- The parameter MAXDENOM ensures that denominators
-- are not too big.

RoundedRatio( I: IntegerNumberSystem,
              MAXDENOM: I):
         Join(OrderedRing, Field) with {
	        *:       (I, %) -> %;
        	/:       (I, I) -> %;
--		^:	 (%,I)  -> %;              -- This is a bad idea if I = SingleInteger.
--		^:	 (%,SingleInteger) -> %;   -- This is a bad idea if Integers are around.
	        numer:   % -> I;
	        denom:   % -> I;
	        coerce:  I -> %;

	     export from I;

  }

== add {
        Rep ==> Record(numer: I, denom: I);

        import from Rep;

        default a, b: %;


        ratio(n: I, d: I): % == per [n, d];

	gcd(a: %, b: %): % == 1;
	(a: %) quo (b: %): % == a/b;
	(a: %) rem (b: %): % == 0;
	divide(a: %, b: %): (%, %) == (a/b, 0);

	-- Public Part --
        numer(a: %): I == rep(a).numer;
        denom(a: %): I == rep(a).denom;

--        apply(p: OutPort, z: %): OutPort ==
--                p(numer z)("/")(denom z);

        coerce(z:%):OutputForm == {
           rn:Fraction I := numer(rn)/denom(rn);
           rn::OutputForm}

        (a: %) =  (b: %): Boolean == numer a = numer b and denom a = denom b;
        (a: %) ~= (b: %): Boolean == ~(a = b);
        (a: %) <  (b: %): Boolean == denom b * numer a  <  denom a * numer b;
        (a: %) <= (b: %): Boolean == denom b * numer a  <= denom a * numer b;
        (a: %) >  (b: %): Boolean == denom b * numer a  >  denom a * numer b;
        (a: %) >= (b: %): Boolean == denom b * numer a  >= denom a * numer b;
	max(a: %, b: %): %        == if a > b then a else b;
	min(a: %, b: %): %        == if a < b then a else b;

	sign     (a: %): %       == ratio(sign numer a, 1);
	abs      (a: %): %       == if negative? a then -a else a;
        zero?    (a: %): Boolean == zero?     numer a;
	negative?(a: %): Boolean == negative? numer a;
	positive?(a: %): Boolean == positive? numer a;
	
        coerce(n: I): % == ratio(n, 1);

        inv(a: %): % == round( [denom a, numer a] );

        0: % == ratio(0,1);
        1: % == ratio(1,1);

	+(a: %): % == a;
        -(a: %): % == ratio(-numer a, denom a);

        (a: %) + (b: %): % == {
		z := [numer a * denom b + numer b * denom a,
			denom a * denom b ];
		round z;
	}
        (a: %) - (b: %): % == a + (-b);

        (a: %) * (b: %):% == {
		round [ numer a * numer b, denom a * denom b]
	}
        (n: I) * (b: %): % == round( [n * numer b, denom b] );
        (n: I) / (d: I): % == round( [n, d] ) ;
        (a: %) / (b: %): % == a * inv b;
        (a: %) \ (b: %): % == inv a * b;
#if 0
        (alpha: %) ^ (n: I): % == {
		local a,c:%;
		local m:I;
		a := alpha;
		m := abs n;
		c := 1;
		while m > 0 repeat {
			{ odd? m => {
				c := c * a;
				m := m - 1;
			            }
			a := a * a;
			m := m quo 2;
			}
		}
		c
	}
        (alpha: %) ^ (n: SingleInteger): % == {
		local a,c:%;
		local m:SingleInteger;
		a := alpha;
		m := abs n;
		c := 1;
		while m > 0 repeat {
			{ odd? m => {
				c := c * a;
				m := m - 1;
			            }
			a := a * a;
			m := m quo 2;   -- bit is supposed to be more efficient.
			}
		}
		c
	}
#endif

        (alpha: %) ^ (n: Integer): % == {
		local a,c:%;
		local m:Integer;
		a := alpha;
		m := abs n;
		c := 1;
		while m > 0 repeat {
			{ odd? m => {
				c := c * a;
				m := m - 1;
			            }
			a := a * a;
			m := m quo 2;
			}
		}
		c
	}


	round( x: Rep) : % == {
	  local a0,a1,q,r,p0,p1,p2,q0,q1,q2 : I;
	  a0 := x.numer;
	  a1 := x.denom;

	  (q,r) := divide(a0,a1);

 	  a0 := a1; a1 := r;
	
	  p0 := 1; p1 := q;
	  q0 := 0; q1 := 1;
	
	
	  while r ~= 0 repeat {

            (q,r) := divide(a0,a1);

	    a0 := a1; a1 := r;
	
	    p2 := q*p1 + p0;
	    q2 := q*q1 + q0;
	
	    q2 > MAXDENOM => break;

	    p0 := p1; p1 := p2;
	    q0 := q1; q1 := q2;
	    }
	
	per [p1, q1]
	}


}