path: root/src/algebra/radix.spad.pamphlet

\documentclass{article}
\usepackage{open-axiom}
\begin{document}
\title{\$SPAD/src/algebra radix.spad}
\author{Stephen M. Watt, Clifton J. Williamson}
\maketitle
\begin{abstract}
\end{abstract}
\eject
\tableofcontents
\eject

\section{domain RADIX RadixExpansion}

<<domain RADIX RadixExpansion>>=
import Integer
import Fraction
import List
import Stream
)abbrev domain RADIX RadixExpansion
++ Author: Stephen M. Watt
++ Date Created: October 1986
++ Date Last Updated: May 15, 1991
++ Basic Operations: wholeRadix, fractRadix, wholeRagits, fractRagits
++ Related Domains: BinaryExpansion, DecimalExpansion, HexadecimalExpansion,
++    RadixUtilities
++ Also See:
++ AMS Classifications:
++ Keywords: radix, base, repeating decimal
++ Examples:
++ References:
++ Description:
++   This domain allows rational numbers to be presented as repeating
++   decimal expansions or more generally as repeating expansions in any base.

RadixExpansion(bb): Exports == Implementation where
  bb   :  Integer
  I   ==> Integer
  NNI ==> NonNegativeInteger
  OUT ==> OutputForm
  RN  ==> Fraction Integer
  ST  ==> Stream Integer
  QuoRem ==> Record(quotient: Integer, remainder: Integer)

  Exports == Join(QuotientFieldCategory(Integer),_
                    CoercibleTo Fraction Integer) with
    fractionPart: % -> Fraction Integer
      ++ fractionPart(rx) returns the fractional part of a radix expansion.
    wholeRagits: % -> List Integer
      ++ wholeRagits(rx) returns the ragits of the integer part
      ++ of a radix expansion.
    fractRagits: % -> Stream Integer
      ++ fractRagits(rx) returns the ragits of the fractional part
      ++ of a radix expansion.
    prefixRagits: % -> List Integer
      ++ prefixRagits(rx) returns the non-cyclic part of the ragits
      ++ of the fractional part of a radix expansion.
      ++ For example, if \spad{x = 3/28 = 0.10 714285 714285 ...},
      ++ then \spad{prefixRagits(x)=[1,0]}.
    cycleRagits: % -> List Integer
      ++ cycleRagits(rx) returns the cyclic part of the ragits of the
      ++ fractional part of a radix expansion.
      ++ For example, if \spad{x = 3/28 = 0.10 714285 714285 ...},
      ++ then \spad{cycleRagits(x) = [7,1,4,2,8,5]}.
    wholeRadix: List Integer -> %
      ++ wholeRadix(l) creates an integral radix expansion from a list
      ++ of ragits.
      ++ For example, \spad{wholeRadix([1,3,4])} will return \spad{134}.
    fractRadix: (List Integer, List Integer) -> %
      ++ fractRadix(pre,cyc) creates a fractional radix expansion
      ++ from a list of prefix ragits and a list of cyclic ragits.
      ++ For example, \spad{fractRadix([1],[6])} will return \spad{0.16666666...}.

  Implementation ==> add
    -- The efficiency of arithmetic operations is poor.
    -- Could use a lazy eval where either rational rep
    -- or list of ragit rep (the current) or both are kept
    -- as demanded.

    bb < 2 => error "Radix base must be at least 2"
    Rep := Record(sgn: Integer,      int: List Integer,
                  pfx: List Integer, cyc: List Integer)

    q:     RN
    qr:    QuoRem
    a,b:   %
    n:     I

    radixInt:    (I, I)    -> List I
    radixFrac:   (I, I, I) -> Record(pfx: List I, cyc: List I)
    checkRagits: List I    -> Boolean

    -- Arithmetic operations
    characteristic == 0
    differentiate a == 0

    0     == [1, nil(),  nil(), nil()]
    1     == [1, [1], nil(), nil()]
    - a   == (a = 0 => 0; [-a.sgn, a.int, a.pfx, a.cyc])
    a + b == (a::RN + b::RN)::%
    a - b == (a::RN - b::RN)@RN::%
    n * a == (n     * a::RN)::%
    a * b == (a::RN * b::RN)::%
    a / b == (a::RN / b::RN)::%
    (i:I) / (j:I) == (i/j)@RN :: %
    a < b == a::RN < b::RN
    a = b == a.sgn = b.sgn and a.int = b.int and
             a.pfx = b.pfx and a.cyc = b.cyc
    numer a == numer(a::RN)
    denom a == denom(a::RN)

    -- Algebraic coercions
    coerce(a):RN == (wholePart a) :: RN + fractionPart a
    coerce(n):%  == n :: RN :: %
    coerce(q):%  ==
      s := 1; if negative? q then (s := -1; q := -q)
      qr      := divide(numer q,denom q)
      whole   := radixInt (qr.quotient,bb)
      fractn  := radixFrac(qr.remainder,denom q,bb)
      cycle   := (fractn.cyc = [0] => nil(); fractn.cyc)
      [s,whole,fractn.pfx,cycle]

    retractIfCan(a):Union(RN,"failed") == a::RN
    retractIfCan(a):Union(I,"failed") ==
      empty?(a.pfx) and empty?(a.cyc) => wholePart a
      "failed"

    -- Exported constructor/destructors
    ceiling a == ceiling(a::RN)
    floor a == floor(a::RN)

    wholePart a ==
      n0 := 0
      for r in a.int repeat n0 := bb*n0 + r
      a.sgn*n0
    fractionPart(a: %): Fraction Integer ==
      n0 := 0
      for r in a.pfx repeat n0 := bb*n0 + r
      null a.cyc =>
          a.sgn*n0/bb**((#a.pfx)::NNI)
      n1 := n0
      for r in a.cyc repeat n1 := bb*n1 + r
      n := n1 - n0
      d := (bb**((#a.cyc)::NNI) - 1) * bb**((#a.pfx)::NNI)
      a.sgn*n/d

    wholeRagits  a == a.int
    fractRagits  a == concat(construct(a.pfx)@ST,repeating a.cyc)
    prefixRagits a == a.pfx
    cycleRagits  a == a.cyc

    wholeRadix li ==
      checkRagits li
      [1, li, nil(), nil()]
    fractRadix(lpfx, lcyc) ==
      checkRagits lpfx; checkRagits lcyc
      [1, nil(), lpfx, lcyc]

    -- Output

    ALPHAS : String := "ABCDEFGHIJKLMNOPQRSTUVWXYZ"

    intToExpr(i:I): OUT ==
      -- computes a digit for bases between 11 and 36
      i < 10 => i :: OUT
      elt(ALPHAS,(i-10) + minIndex(ALPHAS)) :: OUT

    exprgroup(le: List OUT): OUT ==
      empty? le      => error "exprgroup needs non-null list"
      empty? rest le => first le
      abs bb <= 36 => hconcat le
      blankSeparate le

    intgroup(li: List I): OUT ==
      empty? li      => error "intgroup needs non-null list"
      empty? rest li => 
        abs bb <= 36 => intToExpr first(li)
        first(li)::OUT
      abs bb <= 10 => hconcat [i :: OUT for i in li]
      abs bb <= 36 => hconcat [intToExpr(i) for i in li]
      blankSeparate [i :: OUT for i in li]

    overBar(li: List I): OUT == overbar intgroup li

    coerce(a): OUT ==
      le : List OUT := nil()
      if not null a.cyc then le := concat(overBar  a.cyc,le)
      if not null a.pfx then le := concat(intgroup a.pfx,le)
      if not null le    then le := concat("." :: OUT,le)
      if not null a.int then le := concat(intgroup a.int,le)
      else le := concat(0 :: OUT,le)
      rex := exprgroup le
      if negative? a.sgn then -rex else rex

    -- Construction utilities
    checkRagits li ==
      for i in li repeat if negative? i or i >= bb then
        error "Each ragit (digit) must be between 0 and base-1"
      true

    radixInt(n,bas) ==
      rits: List I := nil()
      while abs n ~= 0 repeat
        qr   := divide(n,bas)
        n    := qr.quotient
        rits := concat(qr.remainder,rits)
      rits

    radixFrac(num,den,bas) ==
      -- Rits is the sequence of quotient/remainder pairs
      -- in calculating the radix expansion of the rational number.
      -- We wish to find p and c such that
      --    rits.i are distinct    for 0<=i<=p+c-1
      --    rits.i = rits.(i+p)    for i>p
      -- I.e. p is the length of the non-periodic prefix and c is
      -- the length of the cycle.

      -- Compute p and c using Floyd's algorithm.
      -- 1. Find smallest n s.t. rits.n = rits.(2*n)
      qr    := divide(bas * num, den)
      i : I := 0
      qr1i  := qr2i := qr
      rits: List QuoRem := [qr]
      until qr1i = qr2i repeat
        qr1i := divide(bas * qr1i.remainder,den)
        qrt  := divide(bas * qr2i.remainder,den)
        qr2i := divide(bas * qrt.remainder,den)
        rits := concat(qr2i, concat(qrt, rits))
        i    := i + 1
      rits := reverse! rits
      n    := i
      -- 2. Find p = first i such that rits.i = rits.(i+n)
      ritsi := rits
      ritsn := rits; for i: local in 1..n repeat ritsn := rest ritsn
      i := 0
      while first(ritsi) ~= first(ritsn) repeat
        ritsi := rest ritsi
        ritsn := rest ritsn
        i     := i + 1
      p := i
      -- 3. Find c = first i such that rits.p = rits.(p+i)
      ritsn := rits; for i: local in 1..n repeat ritsn := rest ritsn
      rn    := first ritsn
      cfound:= false
      c : I := 0
      for i: local in 1..p while not cfound repeat
        ritsn := rest ritsn
        if rn = first(ritsn) then
          c := i
          cfound := true
      if not cfound then c := n
      -- 4. Now produce the lists of ragits.
      ritspfx: List I := nil()
      ritscyc: List I := nil()
      for i: local in 1..p repeat
        ritspfx := concat(first(rits).quotient, ritspfx)
        rits    := rest rits
      for i: local in 1..c repeat
        ritscyc := concat(first(rits).quotient, ritscyc)
        rits    := rest rits
      [reverse! ritspfx, reverse! ritscyc]

@
\section{domain BINARY BinaryExpansion}
<<domain BINARY BinaryExpansion>>=
)abbrev domain BINARY BinaryExpansion
++ Author: Clifton J. Williamson
++ Date Created: April 26, 1990
++ Date Last Updated: May 15, 1991
++ Basic Operations:
++ Related Domains: RadixExpansion
++ Also See:
++ AMS Classifications:
++ Keywords: radix, base, binary
++ Examples:
++ References:
++ Description:
++   This domain allows rational numbers to be presented as repeating
++   binary expansions.

BinaryExpansion(): Exports == Implementation where
  Exports == Join(QuotientFieldCategory(Integer),_
               CoercibleTo Fraction Integer,CoercibleTo RadixExpansion(2)) with
    fractionPart: % -> Fraction Integer
      ++ fractionPart(b) returns the fractional part of a binary expansion.
    binary: Fraction Integer -> %
      ++ binary(r) converts a rational number to a binary expansion.

  Implementation ==> RadixExpansion(2) add
    binary r == r :: %
    coerce(x:%): RadixExpansion(2) == rep x

@
\section{domain DECIMAL DecimalExpansion}
<<domain DECIMAL DecimalExpansion>>=
)abbrev domain DECIMAL DecimalExpansion
++ Author: Stephen M. Watt
++ Date Created: October, 1986
++ Date Last Updated: May 15, 1991
++ Basic Operations:
++ Related Domains: RadixExpansion
++ Also See:
++ AMS Classifications:
++ Keywords: radix, base, repeating decimal
++ Examples:
++ References:
++ Description:
++   This domain allows rational numbers to be presented as repeating
++   decimal expansions.
DecimalExpansion(): Exports == Implementation where
  Exports == Join(QuotientFieldCategory(Integer),_
              CoercibleTo Fraction Integer,CoercibleTo RadixExpansion 10) with
    fractionPart: % -> Fraction Integer
      ++ fractionPart(d) returns the fractional part of a decimal expansion.
    decimal: Fraction Integer -> %
      ++ decimal(r) converts a rational number to a decimal expansion.

  Implementation ==> RadixExpansion(10) add
    decimal r == r :: %
    coerce(x:%): RadixExpansion(10) == rep x

@
\section{domain HEXADEC HexadecimalExpansion}
<<domain HEXADEC HexadecimalExpansion>>=
)abbrev domain HEXADEC HexadecimalExpansion
++ Author: Clifton J. Williamson
++ Date Created: April 26, 1990
++ Date Last Updated: May 15, 1991
++ Basic Operations:
++ Related Domains: RadixExpansion
++ Also See:
++ AMS Classifications:
++ Keywords: radix, base, hexadecimal
++ Examples:
++ References:
++ Description:
++   This domain allows rational numbers to be presented as repeating
++   hexadecimal expansions.

HexadecimalExpansion(): Exports == Implementation where
  Exports == Join(QuotientFieldCategory(Integer),_
               CoercibleTo Fraction Integer,_
               CoercibleTo RadixExpansion 16)  with
    fractionPart: % -> Fraction Integer
      ++ fractionPart(h) returns the fractional part of a hexadecimal expansion.
    hex: Fraction Integer -> %
      ++ hex(r) converts a rational number to a hexadecimal expansion.

  Implementation ==> RadixExpansion(16) add
    hex r == r :: %
    coerce(x:%): RadixExpansion(16) == rep x

@
\section{package RADUTIL RadixUtilities}
<<package RADUTIL RadixUtilities>>=
)abbrev package RADUTIL RadixUtilities
++ Author: Stephen M. Watt
++ Date Created: October 1986
++ Date Last Updated: May 15, 1991
++ Basic Operations:
++ Related Domains: RadixExpansion
++ Also See:
++ AMS Classifications:
++ Keywords: radix, base, repeading decimal
++ Examples:
++ References:
++ Description:
++   This package provides tools for creating radix expansions.
RadixUtilities: Exports == Implementation where
  Exports ==> with
    radix: (Fraction Integer,Integer) -> Any
      ++ radix(x,b) converts x to a radix expansion in base b.
  Implementation ==> add
    radix(q, b) ==
      coerce(q :: RadixExpansion(b))$AnyFunctions1(RadixExpansion b)

@
\section{License}
<<license>>=
--Copyright (c) 1991-2002, The Numerical ALgorithms Group Ltd.
--All rights reserved.
--
--Redistribution and use in source and binary forms, with or without
--modification, are permitted provided that the following conditions are
--met:
--
--    - Redistributions of source code must retain the above copyright
--      notice, this list of conditions and the following disclaimer.
--
--    - Redistributions in binary form must reproduce the above copyright
--      notice, this list of conditions and the following disclaimer in
--      the documentation and/or other materials provided with the
--      distribution.
--
--    - Neither the name of The Numerical ALgorithms Group Ltd. nor the
--      names of its contributors may be used to endorse or promote products
--      derived from this software without specific prior written permission.
--
--THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS
--IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
--TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A
--PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER
--OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
--EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
--PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
--PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
--LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
--NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
--SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
@
<<*>>=
<<license>>

<<domain RADIX RadixExpansion>>
<<domain BINARY BinaryExpansion>>
<<domain DECIMAL DecimalExpansion>>
<<domain HEXADEC HexadecimalExpansion>>
<<package RADUTIL RadixUtilities>>
@
\eject
\begin{thebibliography}{99}
\bibitem{1} nothing
\end{thebibliography}
\end{document}