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

\documentclass{article}
\usepackage{axiom}
\begin{document}
\title{\$SPAD/src/algebra opalg.spad}
\author{Manuel Bronstein}
\maketitle
\begin{abstract}
\end{abstract}
\eject
\tableofcontents
\eject
\section{domain MODOP ModuleOperator}
<<domain MODOP ModuleOperator>>=
)abbrev domain MODOP ModuleOperator
++ Author: Manuel Bronstein
++ Date Created: 15 May 1990
++ Date Last Updated: 17 June 1993
++ Description:
++ Algebra of ADDITIVE operators on a module.
ModuleOperator(R: Ring, M:LeftModule(R)): Exports == Implementation where
  O    ==> OutputForm
  OP   ==> BasicOperator
  FG   ==> FreeGroup OP
  RM   ==> Record(coef:R, monom:FG)
  TERM ==> List RM
  FAB  ==> FreeAbelianGroup TERM
  OPADJ   ==> "%opAdjoint"
  OPEVAL  ==> "%opEval"
  INVEVAL ==> "%invEval"

  Exports ==> Join(Ring, RetractableTo R, RetractableTo OP,
                   Eltable(M, M)) with
    if R has CharacteristicZero then CharacteristicZero
    if R has CharacteristicNonZero then CharacteristicNonZero
    if R has CommutativeRing then
      Algebra(R)
      adjoint: $ -> $
        ++ adjoint(op) returns the adjoint of the operator \spad{op}.
      adjoint: ($, $) -> $
        ++ adjoint(op1, op2) sets the adjoint of op1 to be op2.
        ++ op1 must be a basic operator
      conjug  : R -> R
        ++ conjug(x)should be local but conditional
    evaluate: ($, M -> M) -> $
      ++ evaluate(f, u +-> g u) attaches the map g to f.
      ++ f must be a basic operator
      ++ g MUST be additive, i.e. \spad{g(a + b) = g(a) + g(b)} for
      ++ any \spad{a}, \spad{b} in M.
      ++ This implies that \spad{g(n a) = n g(a)} for
      ++ any \spad{a} in M and integer \spad{n > 0}.
    evaluateInverse: ($, M -> M) -> $
	++ evaluateInverse(x,f) \undocumented
    "**": (OP, Integer) -> $
	++ op**n \undocumented
    "**": ($, Integer) -> $
	++ op**n \undocumented
    opeval  : (OP, M) -> M
      ++ opeval should be local but conditional
    makeop   : (R, FG) -> $
      ++ makeop should be local but conditional

  Implementation ==> FAB add
    import NoneFunctions1($)
    import BasicOperatorFunctions1(M)

    Rep := FAB

    inv      : TERM -> $
    termeval : (TERM, M) -> M
    rmeval   : (RM, M) -> M
    monomeval: (FG, M) -> M
    opInvEval: (OP, M) -> M
    mkop     : (R, FG) -> $
    termprod0: (Integer, TERM, TERM) -> $
    termprod : (Integer, TERM, TERM) -> TERM
    termcopy : TERM -> TERM
    trm2O    : (Integer, TERM) -> O
    term2O   : TERM -> O
    rm2O     : (R, FG) -> O
    nocopy   : OP -> $

    1                   == makeop(1, 1)
    coerce(n:Integer):$ == n::R::$
    coerce(r:R):$       == (zero? r => 0; makeop(r, 1))
    coerce(op:OP):$     == nocopy copy op
    nocopy(op:OP):$     == makeop(1, op::FG)
    elt(x:$, r:M)       == +/[t.exp * termeval(t.gen, r) for t in terms x]
    rmeval(t, r)        == t.coef * monomeval(t.monom, r)
    termcopy t          == [[rm.coef, rm.monom] for rm in t]
    characteristic == characteristic$R
    mkop(r, fg)         == [[r, fg]$RM]$TERM :: $
    evaluate(f, g)   == nocopy setProperty(retract(f)@OP,OPEVAL,g pretend None)

    if R has OrderedSet then
      makeop(r, fg) == (r >= 0 => mkop(r, fg); - mkop(-r, fg))
    else makeop(r, fg) == mkop(r, fg)

    inv(t:TERM):$ ==
      empty? t => 1
      c := first(t).coef
      m := first(t).monom
      inv(rest t) * makeop(1, inv m) * (recip(c)::R::$)

    x:$ ** i:Integer ==
      i = 0 => 1
      i > 0 => expt(x,i pretend PositiveInteger)$RepeatedSquaring($)
      (inv(retract(x)@TERM)) ** (-i)

    evaluateInverse(f, g) ==
      nocopy setProperty(retract(f)@OP, INVEVAL, g pretend None)

    coerce(x:$):O ==
      zero? x => (0$R)::O
      reduce(_+, [trm2O(t.exp, t.gen) for t in terms x])$List(O)

    trm2O(c, t) ==
      one? c => term2O t
      c = -1 => - term2O t
      c::O * term2O t

    term2O t ==
      reduce(_*, [rm2O(rm.coef, rm.monom) for rm in t])$List(O)

    rm2O(c, m) ==
      one? c => m::O
      one? m => c::O
      c::O * m::O

    x:$ * y:$ ==
      +/[ +/[termprod0(t.exp * s.exp, t.gen, s.gen) for s in terms y]
          for t in terms x]

    termprod0(n, x, y) ==
      n >= 0 => termprod(n, x, y)::$
      - (termprod(-n, x, y)::$)

    termprod(n, x, y) ==
      lc := first(xx := termcopy x)
      lc.coef := n * lc.coef
      rm := last xx
      one?(first(y).coef) =>
        rm.monom := rm.monom * first(y).monom
        concat_!(xx, termcopy rest y)
      one?(rm.monom) =>
        rm.coef := rm.coef * first(y).coef
        rm.monom := first(y).monom
        concat_!(xx, termcopy rest y)
      concat_!(xx, termcopy y)

    if M has ExpressionSpace then
      opeval(op, r) ==
        (func := property(op, OPEVAL)) case "failed" => kernel(op, r)
        ((func::None) pretend (M -> M)) r

    else
      opeval(op, r) ==
        (func := property(op, OPEVAL)) case "failed" =>
          error "eval: operator has no evaluation function"
        ((func::None) pretend (M -> M)) r

    opInvEval(op, r) ==
      (func := property(op, INVEVAL)) case "failed" =>
         error "eval: operator has no inverse evaluation function"
      ((func::None) pretend (M -> M)) r

    termeval(t, r)  ==
      for rm in reverse t repeat r := rmeval(rm, r)
      r

    monomeval(m, r) ==
      for rec in reverse_! factors m repeat
        e := rec.exp
        g := rec.gen
        e > 0 =>
          for i in 1..e repeat r := opeval(g, r)
        e < 0 =>
          for i in 1..(-e) repeat r := opInvEval(g, r)
      r

    recip x ==
      (r := retractIfCan(x)@Union(R, "failed")) case "failed" => "failed"
      (r1 := recip(r::R)) case "failed" => "failed"
      r1::R::$

    retractIfCan(x:$):Union(R, "failed") ==
      (r:= retractIfCan(x)@Union(TERM,"failed")) case "failed" => "failed"
      empty?(t := r::TERM) => 0$R
      empty? rest t =>
        rm := first t
        one?(rm.monom) => rm.coef
        "failed"
      "failed"

    retractIfCan(x:$):Union(OP, "failed") ==
      (r:= retractIfCan(x)@Union(TERM,"failed")) case "failed" => "failed"
      empty?(t := r::TERM) => "failed"
      empty? rest t =>
        rm := first t
        one?(rm.coef) => retractIfCan(rm.monom)
        "failed"
      "failed"

    if R has CommutativeRing then
      termadj  : TERM -> $
      rmadj    : RM -> $
      monomadj : FG -> $
      opadj    : OP -> $

      r:R * x:$        == r::$ * x
      x:$ * r:R        == x * (r::$)
      adjoint x        == +/[t.exp * termadj(t.gen) for t in terms x]
      rmadj t          == conjug(t.coef) * monomadj(t.monom)
      adjoint(op, adj) == nocopy setProperty(retract(op)@OP, OPADJ, adj::None)

      termadj t ==
        ans:$ := 1
        for rm in t repeat ans := rmadj(rm) * ans
        ans

      monomadj m ==
        ans:$ := 1
        for rec in factors m repeat ans := (opadj(rec.gen) ** rec.exp) * ans
        ans

      opadj op ==
        (adj := property(op, OPADJ)) case "failed" =>
           error "adjoint: operator does not have a defined adjoint"
        (adj::None) pretend $

      if R has conjugate:R -> R then conjug r == conjugate r else conjug r == r

@
\section{domain OP Operator}
<<domain OP Operator>>=
)abbrev domain OP Operator
++ Author: Manuel Bronstein
++ Date Created: 15 May 1990
++ Date Last Updated: 12 February 1993
++ Description:
++ Algebra of ADDITIVE operators over a ring.
Operator(R: Ring) == ModuleOperator(R,R)

@
\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 MODOP ModuleOperator>>
<<domain OP Operator>>
@
\eject
\begin{thebibliography}{99}
\bibitem{1} nothing
\end{thebibliography}
\end{document}