Re: [Axiom-developer] Math Types inclusion

[William Sit]

Dear Bertfried:

Date: Tue, 6 Apr 2004 00:12:49 +0200 (CEST)
From: Bertfried Fauser <address@hidden>

> Most helpful, and I think mathematical modest, would be the development of
> a category in AXIOM which could handle polynomial algebra, where the
> "variables can be as complex as possible data structures (to be able to
> iterate, at least one should have
> 
> A polynomial algebra in several variables, where the variables are signed
> (commutative, anticommutative (two types) and neutral) every "variable"
> should be replaceable by an element from (another instance of) the
> polynomial algebra itself.

Axiom already has these categories, domains, and packages. Check out
MonoidRing(R,M) category constructor (R a commutative ring, M a monoid), which
can be used to construct for example polynomial rings in which the variables
need not be commuting (but the coefficient ring is). On the other hand,
Polynomial(R) can construct polynomial rings where R is not commutative, but the
variables commute. MonoidRing can probably be generalized to take a
non-commutative ring without too much trouble. I think what you mean ny
"replaceable" is a tower construction (combining or iterating constructors, such
as Polynomial (R) where R is a ring and "replacing" can change such a polynomial
into some domain of the category Algebra(R)) and this can also be done on the
fly from the categorical constructors. 
> 
> I think I can do this, but I will take my time (more likely years than
> month). To be able to get started, a simple test domain for polynomial
> algebra would be a great learning filed. I have to face the fact, that I
> am a bad programmer, so my problem really is more like "what is the syntax
> of this and that" and questions like "how to program this and that awkward
> index set" etc. 

Too bad that the open source Axiom has not revived the hyperdoc, which would
help you with this. For the moment you will have to use 
)show       [e.g. )show MonoidRing], or
)display op [e.g. )display op leadingMonomial] to find the syntax.

> And -- efficiency is a major topic since all actual
> computations are rather pretty longwinded, even beyond present day
> compuert power, see the quantum simulator at the FIRST Fraunhofer
> institute in Berlin (soory no url, but easily found via goole)

Axiom is better for testing conceptual algorithms and is not known for its
efficiency. Efficiency can only be gained with better algorithms and fine tuning
in implementation. It is not uncommon that an axiom programmer uses recursion
(carelessly) when a non-recursive version is far more efficient and not that
much harder to code. The justification is that recursion is easier to be proved
correct, and compiler technology and computer hardware will catch up to save the
day.

> Regarding the discussion about Haskell, I am not sure how "functional" one
> can be. My personal experience is, that to abstract mathematics often
> fails the computer algebra texst. I found rather a couple of "theorems"
> which didn't hold true on being tested via Clifford/Bigebra. AXIOM is the
> only tool, which really seems to allow to program nearly functional by
> specifying only types and not actual "elements", but...

Beware that sometimes Axiom is correct in subtle ways even if the result
surprises. It is very strongly typed, and in the interpreter, the type and hence
the operation chosen (guessed?) by the interpreter need not be the one the user
intends. Such "wrong" theorems should be reported as possible bugs though.

> Look at commutativity, given a multiplication m which is commutative
> m(a,b)=m(b,a). A CAS needs _further_ information on term ordering to make
> sense from this. Eg assume a>b then a "simplify" would not affect m(a,b),
> but would change the arguments of m(b,a) into m(a,b). Such issues come up
> in Groebner basis methods, and they really become a point in
> noncommutative algebra. Actually such things as the Euclidean algorithm
> breaks down in such domains and computability is very weak. There math is
> an issue, not only the implementation.

> Hence, (see Singular), even if you have a commutative ring, (monoid, with
> action of an Abelian group [beware, there are rings having no monoid
> structure, ie are not build over a "vector space"]) you need to specify a
> term order to be able to compute. There lays the main difference between
> categorial mathematics and computer algebra needs.

Axiom has also included all sorts of ordering required for computation. Indeed,
different term-orderings in a polynomial ring are implemented for different
polynomial domains (DMP, HDMP, GDMP). Axiom distinguishes an ordered ring from a
ring so there are categories for OrderedRing, OrderedMonoid, etc.

> Hence while eg.
>          monoid < group,
> you need representations of these groups (monoids), which increases the
> computational complexity by magnitudes. To have an SU(2) variable would not
> be entirely satisfying. An SU(2) variable would not even allow you to say
> what spin it has, since the spin value is given to a representation, an
> SU(2) element would however _act_ on such representations. [A group
> representation could be projective, over complex numbers, over a finite
> field etc...]

You can construct as many (parametrized or not) representations of SU(2) as you
like. You can first construct the domain SU(2) as a group of matrices (data
representation) over whatever ground field (that is, the ground field R would be
a parameter of the domain constructor) together with generators and relations.
The group representations can then be constructed using the built-in domain
constructor MonoidRing(R,SU(2)) -- which mathematically speaking, is the algebra
of all maps from SU(2) to R with finite support. You can substitute R with R^n
or other direct product of rings if the representation is more involved (so for
the spin property, you can use as one factor FiniteField(2,1) -- 0 or 1). In
MonoidRing(R, SU(2)), the action need only be defined for the generators and
extended to all of SU(2). Now, may be this is not exactly what you want (and I
am not a physicist, nor have I tested this!), but I think this would be a good
start.

> (AXIOM can already deal rudimentarily with symmetric functions, and hence
> with group representations, see that chapter in the book, alas, the
> current AXIOM does not allow to load this package as described there, so I
> was not able to test what can be done actually)

Axiom has two packages on representation of finite groups:
RepresentationPackage1(R),
RepresentationPackage2(R), the first is for permutation groups, and the second
for modular representtion of finite groups and algebra (using the meat-Axe
algorithm). While these may not fit your application, the source (rep1.spad and
rep2.spad) may be useful as reference.

> Since the present discussion very fast could come into an idle discussion,
> the main point recently has to be to bring AXIOM into life as complete as
> possible, and to provide documentation which allows silly people as I am
> to create a category and a domain, and to understand what actual types are
> already there (I wont try to hack in bad code for things in AXIOM done
> already much neater)

You can take a look at the source code catdef.spad. I am also attaching two old
files monoid.spad and monoid.input as simple examples of constructing your own
monoid categories and domains (with some brief comments comparing with C++). It
also illustrates on-the-fly construction of new domains (direct product of two
monoid domains -- this would be similar to repeated tensor products of SU(2)
constructed on the fly). As long as you put these in your own private directory,
they will not interfere with the system libraries.

William
-- 
William Sit
Department of Mathematics..............Email: address@hidden
City College of New York..........................Tel: 212-650-5179
Convent Ave at West 138th Street..................Fax: 212-862-0004
New York, NY 10031............Axiom, A Scientific Computation Sytem
USA..........................http://www.nongnu.org/axiom/index.html

)load INTPAIR MCAT MCAT- MONOID1 MONOID2 DPM
x:IntegerPair:=makePair(2,3)
x1:Monoid1:=x
x1:MonoidOne:=x
x2:MonoidTwo:=x
x1 * x1
x2 * x2
power(x1,4)
power(x2,4)
M12:=DirectProductOfMonoids(MonoidOne, MonoidTwo)
x12 := makeProduct(x1,x2)
x12 * x12
power(x12,4)
-- you can even form new direct product of the direct product!
-- with just a single line of code
M1212:=DirectProductOfMonoids(M12,M12)
x1212:=makeProduct(x12,x12)
x1212 * x1212
power(x1212,4)

)abb category MCAT MonoidCategory
)abb domain   MONOID1 MonoidOne
)abb domain   INTPAIR IntegerPair
)abb domain   MONOID2 MonoidTwo
)abb domain   DPM DirectProductOfMonoids

--% IntegerPair
IntegerPair(): Header == Code where
   Header ==> SetCategory with -- derived Class
     makePair: (Integer, Integer) -> $ -- C++ constructor, public
     FirstCoordinate:$->Integer -- C++ member function, public
     SecondCoordinate:$->Integer
   Code ==> add
     Rep := Record(x:Integer, y:Integer)  -- encapsulated data, private
     p,q:$
     (p = q):Boolean == (p.x = q.x) and (p.y = q.y)
     coerce(p):OutputForm == bracket([p.x::OutputForm, p.y::OutputForm])
     makePair(a, b) == [a,b]$Rep  
     FirstCoordinate(p) == p.x
     SecondCoordinate(p) == p.y

--% MonoidCategory
-- example of abstract base class
MonoidCategory():Category == SetCategory with
-- specify class member functions (similar to header files)
    "*": ($,$) -> $  -- pure virtual C++ functions
    1: constant -> $
    power: ($, Integer) -> $  -- virtual C++ function
  add
    power(x, n) ==   -- with default implementation
      n < 0 => error("Only positive powers in monoids")
      n = 0 => 1
      n = 1 => x
      x * power(x, n-1)

--% MonoidOne
MonoidOne(): Header == Code where
  Header ==> MonoidCategory with  -- derived from abstract base class
    coerce:IntegerPair -> $ -- new member function, a constructor
  Code ==> IntegerPair add -- derived from real base class
    import IntegerPair
    Rep:=IntegerPair -- layered inheritance
    coerce x == x
    (p:$ * q:$):$ ==
       x := FirstCoordinate(p) * FirstCoordinate(q)
       y := SecondCoordinate(p) * SecondCoordinate(q)
       makePair(x,y)
    1 == makePair(1$Integer, 1$Integer)

--% MonoidTwo
MonoidTwo(): Header == Code where
  Header ==> MonoidCategory with
    coerce:IntegerPair -> $
  Code ==> IntegerPair add
    import IntegerPair
    Rep:=IntegerPair
    coerce x == x  -- equivalent to C++ constructor
    (p:$ * q:$):$ ==
       x := FirstCoordinate(p) * SecondCoordinate(q)
       y := FirstCoordinate(p) * SecondCoordinate(q)
       makePair(x,y)
    1 == makePair(1$Integer, 1$Integer)

--% DirectProductOfMonoids
DirectProductOfMonoids(A:MonoidCategory,B:MonoidCategory):Header == Code where
    Header ==> MonoidCategory with -- derived from abstract base class
      makeProduct:(A, B) -> $
      FirstProjection:$ -> A
      SecondProjection:$-> B
    Code ==> add
      Rep:= Record(first:A, second:B)  -- layered inheritance
      makeProduct(a,b) == [a,b]$Rep
      x, y:$
      x * y == makeProduct(x.first *  y.first, x.second * y.second)
      FirstProjection(x) == x.first
      SecondProjection(x) == x.second
      coerce(x:$):OutputForm ==  -- ostream equivalent
        bracket([x.first::OutputForm, x.second::OutputForm])
      1 == makeProduct(1$A, 1$B)

-- note that we can now form many direct products of monoids
-- and do so "on the fly", all without any notion of pointers.
-- To do so in C++ would require so many levels of pointers
-- only the diehard would love it (every abstraction level
-- generally require at least one level of pointers).