\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 -- )r herm.input \end{verbatim} \section{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} <>= --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. @ <<*>>= <> <> @ \eject \begin{thebibliography}{99} \bibitem{1} nothing \end{thebibliography} \end{document}