\documentclass{article} \usepackage{open-axiom} \begin{document} \title{\$SPAD/src/algebra taylor.spad} \author{Clifton J. Williamson} \maketitle \begin{abstract} \end{abstract} \eject \tableofcontents \eject \section{domain ITAYLOR InnerTaylorSeries} <>= )abbrev domain ITAYLOR InnerTaylorSeries ++ Author: Clifton J. Williamson ++ Date Created: 21 December 1989 ++ Date Last Updated: October 26, 2009 ++ Basic Operations: ++ Related Domains: UnivariateTaylorSeries(Coef,var,cen) ++ Also See: ++ AMS Classifications: ++ Keywords: stream, dense Taylor series ++ Examples: ++ References: ++ Description: Internal package for dense Taylor series. ++ This is an internal Taylor series type in which Taylor series ++ are represented by a \spadtype{Stream} of \spadtype{Ring} elements. ++ For univariate series, the \spad{Stream} elements are the Taylor ++ coefficients. For multivariate series, the \spad{n}th Stream element ++ is a form of degree n in the power series variables. InnerTaylorSeries(Coef): Exports == Implementation where Coef : Ring I ==> Integer NNI ==> NonNegativeInteger ST ==> Stream Coef STT ==> StreamTaylorSeriesOperations Coef Exports == Join(Ring, BiModule(Coef,Coef)) with coefficients: % -> Stream Coef ++\spad{coefficients(x)} returns a stream of ring elements. ++ When x is a univariate series, this is a stream of Taylor ++ coefficients. When x is a multivariate series, the ++ \spad{n}th element of the stream is a form of ++ degree n in the power series variables. series: Stream Coef -> % ++\spad{series(s)} creates a power series from a stream of ++ ring elements. ++ For univariate series types, the stream s should be a stream ++ of Taylor coefficients. For multivariate series types, the ++ stream s should be a stream of forms the \spad{n}th element ++ of which is a ++ form of degree n in the power series variables. pole?: % -> Boolean ++\spad{pole?(x)} tests if the series x has a pole. ++ Note: this is false when x is a Taylor series. order: % -> NNI ++\spad{order(x)} returns the order of a power series x, ++ i.e. the degree of the first non-zero term of the series. order: (%,NNI) -> NNI ++\spad{order(x,n)} returns the minimum of n and the order of x. * : (%,Integer)->% ++\spad{x*i} returns the product of integer i and the series x. if Coef has IntegralDomain then IntegralDomain --++ An IntegralDomain provides 'exquo' Implementation == add Rep == ST -- In what follows, we will be calling operations on Streams -- which are NOT defined in the package Stream. Thus, it is -- necessary to explicitly pass back and forth between Rep and %. series st == per st 0 == per coerce(0)$STT 1 == per coerce(1)$STT x = y == -- tests if two power series are equal -- difference must be a finite stream of zeroes of length <= n + 1, -- where n = $streamCount$Lisp st : ST := rep(x - y) n : I := _$streamCount$Lisp for i in 0..n repeat empty? st => return true frst st ~= 0 => return false st := rst st empty? st coefficients x == rep x x + y == per(rep x +$STT rep y) x - y == per(rep x -$STT rep y) (x:%) * (y:%) == per(rep x *$STT rep y) - x == per(-$STT rep x) (i:I) * (x:%) == per((i::Coef) *$STT rep x) (x:%) * (i:I) == per(rep x *$STT (i::Coef)) (c:Coef) * (x:%) == per(c *$STT rep x) (x:%) * (c:Coef) == per(rep x *$STT c) recip x == (rec := recip$STT rep x) case "failed" => "failed" per(rec :: ST) if Coef has IntegralDomain then x exquo y == (quot := rep x exquo$STT rep y) case "failed" => "failed" per(quot :: ST) x:% ** n:NNI == n = 0 => 1 expt(x,n :: PositiveInteger)$RepeatedSquaring(%) characteristic == characteristic$Coef pole? x == false iOrder: (ST,NNI,NNI) -> NNI iOrder(st,n,n0) == (n = n0) or (empty? st) => n0 zero? frst st => iOrder(rst st,n + 1,n0) n order(x,n) == iOrder(rep x,0,n) iOrder2: (ST,NNI) -> NNI iOrder2(st,n) == empty? st => error "order: series has infinite order" zero? frst st => iOrder2(rst st,n + 1) n order x == iOrder2(rep x,0) zero? x == empty? rep x @ \section{domain UTS UnivariateTaylorSeries} <>= )abbrev domain UTS UnivariateTaylorSeries ++ Author: Clifton J. Williamson ++ Date Created: 21 December 1989 ++ Date Last Updated: June 18, 2010 ++ Basic Operations: ++ Related Domains: UnivariateLaurentSeries(Coef,var,cen), UnivariatePuiseuxSeries(Coef,var,cen) ++ Also See: ++ AMS Classifications: ++ Keywords: dense, Taylor series ++ Examples: ++ References: ++ Description: Dense Taylor series in one variable ++ \spadtype{UnivariateTaylorSeries} is a domain representing Taylor ++ series in ++ one variable with coefficients in an arbitrary ring. The parameters ++ of the type specify the coefficient ring, the power series variable, ++ and the center of the power series expansion. For example, ++ \spadtype{UnivariateTaylorSeries}(Integer,x,3) represents ++ Taylor series in ++ \spad{(x - 3)} with \spadtype{Integer} coefficients. UnivariateTaylorSeries(Coef,var,cen): Exports == Implementation where Coef : Ring var : Symbol cen : Coef I ==> Integer NNI ==> NonNegativeInteger P ==> Polynomial Coef RN ==> Fraction Integer ST ==> Stream STT ==> StreamTaylorSeriesOperations Coef TERM ==> Record(k:NNI,c:Coef) UP ==> UnivariatePolynomial(var,Coef) Exports ==> Join(UnivariateTaylorSeriesCategory(Coef),_ PartialDifferentialDomain(%,Variable var)) with coerce: UP -> % ++\spad{coerce(p)} converts a univariate polynomial p in the variable ++\spad{var} to a univariate Taylor series in \spad{var}. univariatePolynomial: (%,NNI) -> UP ++\spad{univariatePolynomial(f,k)} returns a univariate polynomial ++ consisting of the sum of all terms of f of degree \spad{<= k}. coerce: Variable(var) -> % ++\spad{coerce(var)} converts the series variable \spad{var} into a ++ Taylor series. lagrange: % -> % ++\spad{lagrange(g(x))} produces the Taylor series for \spad{f(x)} ++ where \spad{f(x)} is implicitly defined as \spad{f(x) = x*g(f(x))}. lambert: % -> % ++\spad{lambert(f(x))} returns \spad{f(x) + f(x^2) + f(x^3) + ...}. ++ This function is used for computing infinite products. ++ \spad{f(x)} should have zero constant coefficient. ++ If \spad{f(x)} is a Taylor series with constant term 1, then ++ \spad{product(n = 1..infinity,f(x^n)) = exp(log(lambert(f(x))))}. oddlambert: % -> % ++\spad{oddlambert(f(x))} returns \spad{f(x) + f(x^3) + f(x^5) + ...}. ++ \spad{f(x)} should have a zero constant coefficient. ++ This function is used for computing infinite products. ++ If \spad{f(x)} is a Taylor series with constant term 1, then ++ \spad{product(n=1..infinity,f(x^(2*n-1)))=exp(log(oddlambert(f(x))))}. evenlambert: % -> % ++\spad{evenlambert(f(x))} returns \spad{f(x^2) + f(x^4) + f(x^6) + ...}. ++ \spad{f(x)} should have a zero constant coefficient. ++ This function is used for computing infinite products. ++ If \spad{f(x)} is a Taylor series with constant term 1, then ++ \spad{product(n=1..infinity,f(x^(2*n))) = exp(log(evenlambert(f(x))))}. generalLambert: (%,I,I) -> % ++\spad{generalLambert(f(x),a,d)} returns \spad{f(x^a) + f(x^(a + d)) + ++ f(x^(a + 2 d)) + ... }. \spad{f(x)} should have zero constant ++ coefficient and \spad{a} and d should be positive. revert: % -> % ++ \spad{revert(f(x))} returns a Taylor series \spad{g(x)} such that ++ \spad{f(g(x)) = g(f(x)) = x}. Series \spad{f(x)} should have constant ++ coefficient 0 and invertible 1st order coefficient. multisect: (I,I,%) -> % ++\spad{multisect(a,b,f(x))} selects the coefficients of ++ \spad{x^((a+b)*n+a)}, and changes this monomial to \spad{x^n}. invmultisect: (I,I,%) -> % ++\spad{invmultisect(a,b,f(x))} substitutes \spad{x^((a+b)*n)} ++ for \spad{x^n} and multiples by \spad{x^b}. if Coef has Algebra Fraction Integer then integrate: (%,Variable(var)) -> % ++ \spad{integrate(f(x),x)} returns an anti-derivative of the power ++ series \spad{f(x)} with constant coefficient 0. ++ We may integrate a series when we can divide coefficients ++ by integers. Implementation ==> InnerTaylorSeries(Coef) add Rep := Stream Coef --% creation and destruction of series stream: % -> Stream Coef stream x == x pretend Stream(Coef) coerce(v:Variable(var)) == zero? cen => monomial(1,1) monomial(1,1) + monomial(cen,0) coerce(n:I) == n :: Coef :: % coerce(r:Coef) == coerce(r)$STT monomial(c,n) == monom(c,n)$STT getExpon: TERM -> NNI getExpon term == term.k getCoef: TERM -> Coef getCoef term == term.c rec: (NNI,Coef) -> TERM rec(expon,coef) == [expon,coef] recs: (ST Coef,NNI) -> ST TERM recs(st,n) == delay$ST(TERM) empty? st => empty() zero? (coef := frst st) => recs(rst st,n + 1) concat(rec(n,coef),recs(rst st,n + 1)) terms x == recs(stream x,0) recsToCoefs: (ST TERM,NNI) -> ST Coef recsToCoefs(st,n) == delay empty? st => empty() term := frst st; expon := getExpon term n = expon => concat(getCoef term,recsToCoefs(rst st,n + 1)) concat(0,recsToCoefs(st,n + 1)) series(st: ST TERM) == recsToCoefs(st,0) stToPoly: (ST Coef,P,NNI,NNI) -> P stToPoly(st,term,n,n0) == (n > n0) or (empty? st) => 0 frst(st) * term ** n + stToPoly(rst st,term,n + 1,n0) polynomial(x,n) == stToPoly(stream x,(var :: P) - (cen :: P),0,n) polynomial(x,n1,n2) == if n1 > n2 then (n1,n2) := (n2,n1) stToPoly(rest(stream x,n1),(var :: P) - (cen :: P),n1,n2) stToUPoly: (ST Coef,UP,NNI,NNI) -> UP stToUPoly(st,term,n,n0) == (n > n0) or (empty? st) => 0 frst(st) * term ** n + stToUPoly(rst st,term,n + 1,n0) univariatePolynomial(x,n) == stToUPoly(stream x,monomial(1,1)$UP - monomial(cen,0)$UP,0,n) coerce(p:UP) == zero? p => 0 if not zero? cen then p := p(monomial(1,1)$UP + monomial(cen,0)$UP) st : ST Coef := empty() oldDeg : NNI := degree(p) + 1 while not zero? p repeat deg := degree p delta := (oldDeg - deg - 1) :: NNI for i in 1..delta repeat st := concat(0$Coef,st) st := concat(leadingCoefficient p,st) oldDeg := deg; p := reductum p for i in 1..oldDeg repeat st := concat(0$Coef,st) st if Coef has coerce: Symbol -> Coef then if Coef has "**": (Coef,NNI) -> Coef then stToCoef: (ST Coef,Coef,NNI,NNI) -> Coef stToCoef(st,term,n,n0) == (n > n0) or (empty? st) => 0 frst(st) * term ** n + stToCoef(rst st,term,n + 1,n0) approximate(x,n) == stToCoef(stream x,(var :: Coef) - cen,0,n) --% values variable x == var center s == cen coefficient(x,n) == -- Cannot use elt! Should return 0 if stream doesn't have it. u := stream x while not empty? u and positive? n repeat u := rst u n := (n - 1) :: NNI empty? u or n ~= 0 => 0 frst u elt(x:%,n:NNI) == coefficient(x,n) --% functions map(f,x) == map(f,x)$Rep eval(x:%,r:Coef) == eval(stream x,r-cen)$STT differentiate x == deriv(stream x)$STT differentiate(x:%,v:Variable(var)) == differentiate x if Coef has PartialDifferentialRing(Symbol) then differentiate(x:%,s:Symbol) == (s = variable(x)) => differentiate x map(differentiate(#1,s),x) - differentiate(center x,s)*differentiate(x) multiplyCoefficients(f,x) == gderiv(f,stream x)$STT lagrange x == lagrange(stream x)$STT lambert x == lambert(stream x)$STT oddlambert x == oddlambert(stream x)$STT evenlambert x == evenlambert(stream x)$STT generalLambert(x:%,a:I,d:I) == generalLambert(stream x,a,d)$STT extend(x,n) == extend(x,n+1)$Rep complete x == complete(x)$Rep truncate(x,n) == first(stream x,n + 1)$Rep truncate(x,n1,n2) == if n2 < n1 then (n1,n2) := (n2,n1) m := (n2 - n1) :: NNI st := first(rest(stream x,n1)$Rep,m + 1)$Rep for i in 1..n1 repeat st := concat(0$Coef,st) st elt(x:%,y:%) == compose(stream x,stream y)$STT revert x == revert(stream x)$STT multisect(a,b,x) == multisect(a,b,stream x)$STT invmultisect(a,b,x) == invmultisect(a,b,stream x)$STT multiplyExponents(x,n) == invmultisect(n,0,x) quoByVar x == (empty? x => 0; rst x) if Coef has IntegralDomain then unit? x == unit? coefficient(x,0) if Coef has Field then if Coef is RN then (x:%) ** (s:Coef) == powern(s,stream x)$STT else (x:%) ** (s:Coef) == power(s,stream x)$STT if Coef has Algebra Fraction Integer then coerce(r:RN) == r :: Coef :: % integrate x == integrate(0,stream x)$STT integrate(x:%,v:Variable(var)) == integrate x if Coef has integrate: (Coef,Symbol) -> Coef and _ Coef has variables: Coef -> List Symbol then integrate(x:%,s:Symbol) == (s = variable(x)) => integrate x not entry?(s,variables center x) => map(integrate(#1,s),x) error "integrate: center is a function of variable of integration" if Coef has TranscendentalFunctionCategory and _ Coef has PrimitiveFunctionCategory and _ Coef has AlgebraicallyClosedFunctionSpace Integer then integrateWithOneAnswer: (Coef,Symbol) -> Coef integrateWithOneAnswer(f,s) == res := integrate(f,s)$FunctionSpaceIntegration(I,Coef) res case Coef => res :: Coef first(res :: List Coef) integrate(x:%,s:Symbol) == (s = variable(x)) => integrate x not entry?(s,variables center x) => map(integrateWithOneAnswer(#1,s),x) error "integrate: center is a function of variable of integration" --% OutputForms -- We use the default coerce: % -> OutputForm in UTSCAT& @ \section{package UTS2 UnivariateTaylorSeriesFunctions2} <>= )abbrev package UTS2 UnivariateTaylorSeriesFunctions2 ++ Author: Clifton J. Williamson ++ Date Created: 9 February 1990 ++ Date Last Updated: 9 February 1990 ++ Basic Operations: ++ Related Domains: UnivariateTaylorSeries(Coef1,var,cen) ++ Also See: ++ AMS Classifications: ++ Keywords: Taylor series, map ++ Examples: ++ References: ++ Description: Mapping package for univariate Taylor series. ++ This package allows one to apply a function to the coefficients of ++ a univariate Taylor series. UnivariateTaylorSeriesFunctions2(Coef1,Coef2,UTS1,UTS2):_ Exports == Implementation where Coef1 : Ring Coef2 : Ring UTS1 : UnivariateTaylorSeriesCategory Coef1 UTS2 : UnivariateTaylorSeriesCategory Coef2 ST2 ==> StreamFunctions2(Coef1,Coef2) Exports ==> with map: (Coef1 -> Coef2,UTS1) -> UTS2 ++\spad{map(f,g(x))} applies the map f to the coefficients of ++ the Taylor series \spad{g(x)}. Implementation ==> add map(f,uts) == series map(f,coefficients uts)$ST2 @ \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}