path: root/src/algebra/herm.as.pamphlet

\documentclass{article}
\usepackage{axiom}
\begin{document}
\title{\$SPAD/src/algebra herm.as}
\author{Michael Richardson}
\maketitle
\begin{abstract}
\end{abstract}
\eject
\tableofcontents
\eject
\begin{verbatim}
-- N.B. ndftip.as inlines this, must be recompiled if this is.

-- To test:
-- sed -ne '1,/^#if NeverAssertThis/d;/#endif/d;p' < herm.as > herm.input    
-- axiom
-- )set nag host <some machine running nagd>
-- )r herm.input
\end{verbatim}
\section{PackedHermitianSequence}
<<PackedHermitianSequence>>=
#include "axiom.as"

INT ==> Integer ;
NNI ==> NonNegativeInteger ;
PHS ==> PackedHermitianSequence ;
 
+++ Author: M.G. Richardson
+++ Date Created: 1995 Nov. 24
+++ Date Last Updated:
+++ Basic Functions:
+++ Related Constructors: Vector
+++ Also See:
+++ AMS Classifications:
+++ Keywords:
+++ References:
+++ Description:
+++ This type represents packed Hermitian sequences - that is, complex
+++ sequences s, whose tails ("rest s", in Axiom terms) are conjugate to
+++ themselves reversed - by real sequences, in the "standard" manner;
+++ in this, the real parts of the elements of the first "half" of the
+++ tail are stored there and the imaginary parts are stored in reverse
+++ order in the second "half" of the tail.
+++ (If the tail has an odd number of elements, its middle element is
+++ real and is stored unchanged.  The "halves" mentioned above then
+++ refer to the elements before and after the middle, respectively.)

PackedHermitianSequence(R : CommutativeRing) : LinearAggregate(R) with{

  pHS : List R -> % ;
++ pHS(l) converts the list l to a packedHermitianSequence.

  expand : % -> Vector Complex R ;
++ expand(h) converts the packedHermitianSequence h to a Hermitian
++ sequence (a complex vector).

  packHS : Vector Complex R -> % ;
++ packHS(v) checks that the complex vector v represents a Hermitian
++ sequences and, if so, converts it to a packedHermitianSequence;
++ otherwise an error message is printed.

  conjHerm : % -> % ;
++ conjHerm(h) returns the packedHermitianSequence which represents the
++ Hermitian sequence conjugate to that represented by h.

 coerce : % -> OutputForm -- shouldn't need this, should be inherited
			  -- from Vector.

} == Vector(R) add {
  Rep ==> Vector R;
  import from Rep;
  import from INT ;
  import from R ;
  import from Vector R ;
  import from Complex R ;
  import from Vector Complex R ;
  import from ErrorFunctions ;
  import from String ;
  import from List String ;

  local (..)(a:INT,b:INT):Generator INT == {
    generate {
      t := a ;
      while (t <= b) repeat {
        yield t ;
        t := t + 1 ;
        }
      }
    }

  pHS(l : List R) : % == (vector l) pretend PHS R ;

  expand(h : %) : Vector Complex R == {
 
    local len : NNI ;
    local nvals, npairs, n1, n2 : INT ;
    local fullh : Vector Complex R ;

    {
    len := # h ;
    nvals := len pretend INT ; -- pretend since :: fails
    npairs := (nvals - 1) quo 2 ;
    fullh := new(len, 0) ;
    (nvals = 0) => () ;
    fullh.1 := complex(h.1,0) ;
    (nvals = 1) => () ;
    fullh.(npairs+2) := complex(h.(npairs+2),0) ; -- need if even length
                                                  -- (not worth testing)
    for j in 1 .. npairs repeat {
      n1 := j + 1 ;
      n2 := nvals - j + 1 ;
      fullh.n1 := complex(h.n1, h.n2) ;
      fullh.n2 := complex(h.n1, -(h.n2)) ;
      }

    }
      
    fullh

  }

  packHS(v : Vector Complex R) : % == {

    local len : NNI ;
    local nonhs : String ==  "The argument of packHS is not Hermitian" ;
    local nvals, testprs, n1, n2 : INT ;
    local hpacked : Vector R ;
    local v1, v2 : Complex R ;
    local r1, i1, r2, i2 : R ;

    {
    len := # v ;
    nvals := len pretend INT ; -- pretend since :: fails
    testprs := nvals quo 2 ;
    hpacked := new(len, 0) ;
    (nvals = 0) => () ;
    if imag(v.1) ~= 0 
    then error [nonhs, " - the first element must be real."]
    else {
      hpacked.1 := real(v.1) ;
      (nvals = 1) => () ;
      for j in 1 .. testprs repeat {
        n1 := j + 1 ;
        n2 := nvals - j + 1 ;
        v1 := v.n1 ;
        v2 := v.n2 ;
        r1 := real v1 ;
        i1 := imag v1 ;
        r2 := real v2 ;
        i2 := imag v2 ;
        if r1 ~= r2 or i1 ~= -i2
        then if n1 = n2
	     then error [nonhs,
			 " - element ",
			 string(n1),
			 " must be real to be self-conjugate."]
             else error [nonhs,
                         " - elements ",
			 string(n1),
			 " and ",
			 string(n2),
			 " are not conjugate."]
        else {
	  hpacked.n2 := i1 ; -- This order means that when the tail of v
	  hpacked.n1 := r1 ; -- has odd length, the (real part) of its
			     -- middle element ends up in that position.
	  }
	}
      }
    }

    hpacked pretend %

  }

  local set!(x: %, i: INT, v: R): () == {
	(rep x).i := v;
  }
  conjHerm(h : %) : %  == {

    local len : NNI ;
    local nvals, npairs : INT ;
    local ch : % ;

    ch := copy h ;
    len := # h ;
    (len < 3) => ch ; -- these Hermitian sequences are self-conjugate.
    nvals := len pretend INT ; -- pretend since :: fails
    npairs := (nvals - 1) quo 2 ;
    for j in (nvals - npairs + 1) .. nvals repeat ch.j := - h.j ;
    ch

  }
   
  import from List OutputForm ;
  
  coerce(h : %) : OutputForm ==
    bracket commaSeparate [
       qelt(h, k) :: OutputForm for k in minIndex h .. maxIndex h]

}

#if NeverAssertThis

)lib herm

h0 := pHS([] :: List INT)

--       []

h1 := pHS [1]

--       [1]

h2 := pHS [1,2]

--       [1,2]

h3 := pHS [1,2,3]

--       [1,2,3]

h4 := pHS [1,2,3,4]

--       [1,2,3,4]

h5 := pHS [1,2,3,4,5]

--       [1,2,3,4,5]


f0 := expand h0

--       []

f1 := expand h1

--       [1]

f2 := expand h2

--       [1,2]

f3 := expand h3

--       [1,2 + 3%i,2 - 3%i]

f4 := expand h4

--       [1,2 + 4%i,3,2 - 4%i]

f5 := expand h5

--       [1,2 + 5%i,3 + 4%i,3 - 4%i,2 - 5%i]
      
packHS f0

--       []

packHS f1

--       [1]

packHS f2

--       [1,2]

packHS f3

--       [1,2,3]

packHS f4

--       [1,2,3,4]

packHS f5

--       [1,2,3,4,5]

packHS vector[%i,3,3,3]

-- Error signalled from user code:
--    The argument of packHS is not Hermitian - the first element must
--    be real.

packHS vector [1, 3, 5, 7]

-- Error signalled from user code:
--    The argument of packHS is not Hermitian - elements 2 and 4 are 
--    not conjugate.

packHS [1, 3, %i, 3]

-- Error signalled from user code:
--    The argument of packHS is not Hermitian - element 3 must be real
--    to be self-conjugate.

conjHerm h0

--       []

conjHerm h1

--       [1]

conjHerm h2

--       [1,2]

conjHerm h3

--       [1,2,- 3]

conjHerm h4

--       [1,2,3,- 4]

conjHerm h5

--       [1,2,3,- 4,- 5]

output "End of tests"

#endif
@
\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>>

<<PackedHermitianSequence>>
@
\eject
\begin{thebibliography}{99}
\bibitem{1} nothing
\end{thebibliography}
\end{document}