\documentclass{article} \usepackage{open-axiom} \begin{document} \title{\$SPAD/src/algebra intaux.spad} \author{Barry Trager, Manuel Bronstein} \maketitle \begin{abstract} \end{abstract} \eject \tableofcontents \eject \section{domain IR IntegrationResult} <>= )abbrev domain IR IntegrationResult ++ The result of a transcendental integration. ++ Author: Barry Trager, Manuel Bronstein ++ Date Created: 1987 ++ Date Last Updated: 12 August 1992 ++ Description: ++ If a function f has an elementary integral g, then g can be written ++ in the form \spad{g = h + c1 log(u1) + c2 log(u2) + ... + cn log(un)} ++ where h, which is in the same field than f, is called the rational ++ part of the integral, and \spad{c1 log(u1) + ... cn log(un)} is called the ++ logarithmic part of the integral. This domain manipulates integrals ++ represented in that form, by keeping both parts separately. The logs ++ are not explicitly computed. ++ Keywords: integration. ++ Examples: )r RATINT INPUT IntegrationResult(F:Field): Exports == Implementation where O ==> OutputForm B ==> Boolean Z ==> Integer Q ==> Fraction Integer SE ==> Symbol UP ==> SparseUnivariatePolynomial F LOG ==> Record(scalar:Q, coeff:UP, logand:UP) NE ==> Record(integrand:F, intvar:F) Exports ==> (Module Q, RetractableTo F) with mkAnswer: (F, List LOG, List NE) -> % ++ mkAnswer(r,l,ne) creates an integration result from ++ a rational part r, a logarithmic part l, and a non-elementary part ne. ratpart : % -> F ++ ratpart(ir) returns the rational part of an integration result logpart : % -> List LOG ++ logpart(ir) returns the logarithmic part of an integration result notelem : % -> List NE ++ notelem(ir) returns the non-elementary part of an integration result elem? : % -> B ++ elem?(ir) tests if an integration result is elementary over F? integral: (F, F) -> % ++ integral(f,x) returns the formal integral of f with respect to x differentiate: (%, F -> F) -> F ++ differentiate(ir,D) differentiates ir with respect to the derivation D. if F has PartialDifferentialRing(SE) then differentiate: (%, Symbol) -> F ++ differentiate(ir,x) differentiates ir with respect to x if F has RetractableTo Symbol then integral: (F, Symbol) -> % ++ integral(f,x) returns the formal integral of f with respect to x Implementation ==> add Rep := Record(ratp: F, logp: List LOG, nelem: List NE) timelog : (Q, LOG) -> LOG timene : (Q, NE) -> NE LOG2O : LOG -> O NE2O : NE -> O Q2F : Q -> F nesimp : List NE -> List NE neselect: (List NE, F) -> F pLogDeriv: (LOG, F -> F) -> F pNeDeriv : (NE, F -> F) -> F alpha:O := new()$Symbol :: O - u == (-1$Z) * u 0 == mkAnswer(0, empty(), empty()) coerce(x:F):% == mkAnswer(x, empty(), empty()) ratpart u == u.ratp logpart u == u.logp notelem u == u.nelem elem? u == empty? notelem u mkAnswer(x, l, n) == [x, l, nesimp n] timelog(r, lg) == [r * lg.scalar, lg.coeff, lg.logand] integral(f:F,x:F) == (zero? f => 0; mkAnswer(0, empty(), [[f, x]])) timene(r, ne) == [Q2F(r) * ne.integrand, ne.intvar] n:Z * u:% == (n::Q) * u Q2F r == numer(r)::F / denom(r)::F neselect(l, x) == +/[ne.integrand for ne in l | ne.intvar = x] if F has RetractableTo Symbol then integral(f:F, x:Symbol):% == integral(f, x::F) LOG2O rec == one? degree rec.coeff => -- deg 1 minimal poly doesn't get sigma lastc := - coefficient(rec.coeff, 0) / coefficient(rec.coeff, 1) lg := (rec.logand) lastc logandp := prefix('log::O, [lg::O]) (cc := Q2F(rec.scalar) * lastc) = 1 => logandp cc = -1 => - logandp cc::O * logandp coeffp:O := (outputForm(rec.coeff, alpha) = 0::Z::O)@O logandp := alpha * prefix('log::O, [outputForm(rec.logand, alpha)]) if not one?(cc := Q2F(rec.scalar)) then logandp := cc::O * logandp sum(logandp, coeffp) nesimp l == [[u,x] for x in removeDuplicates!([ne.intvar for ne in l]$List(F)) | (u := neselect(l, x)) ~= 0] if (F has LiouvillianFunctionCategory) and (F has RetractableTo Symbol) then retractIfCan u == empty? logpart u => ratpart u + +/[integral(ne.integrand, retract(ne.intvar)@Symbol)$F for ne in notelem u] "failed" else retractIfCan u == elem? u and empty? logpart u => ratpart u "failed" r:Q * u:% == r = 0 => 0 mkAnswer(Q2F(r) * ratpart u, map(timelog(r, #1), logpart u), map(timene(r, #1), notelem u)) -- Initial attempt, quick and dirty, no simplification u + v == mkAnswer(ratpart u + ratpart v, concat(logpart u, logpart v), nesimp concat(notelem u, notelem v)) if F has PartialDifferentialRing(Symbol) then differentiate(u:%, x:Symbol):F == differentiate(u, differentiate(#1, x)) differentiate(u:%, derivation:F -> F):F == derivation ratpart u + +/[pLogDeriv(log, derivation) for log in logpart u] + (+/[pNeDeriv(ne, derivation) for ne in notelem u]) pNeDeriv(ne, derivation) == one? derivation(ne.intvar) => ne.integrand zero? derivation(ne.integrand) => 0 error "pNeDeriv: cannot differentiate not elementary part into F" pLogDeriv(log, derivation) == map(derivation, log.coeff) ~= 0 => error "pLogDeriv: can only handle logs with constant coefficients" one?(n := degree(log.coeff)) => c := - (leadingCoefficient reductum log.coeff) / (leadingCoefficient log.coeff) ans := (log.logand) c Q2F(log.scalar) * c * derivation(ans) / ans numlog := map(derivation, log.logand) diflog := extendedEuclidean(log.logand, log.coeff, numlog)::Record(coef1:UP, coef2:UP) algans := diflog.coef1 ans:F := 0 for i in 0..(n-1) repeat algans := algans * monomial(1, 1) rem log.coeff ans := ans + coefficient(algans, i) Q2F(log.scalar) * ans coerce(u:%):O == (r := retractIfCan u) case F => r::F::O l := reverse! [LOG2O f for f in logpart u]$List(O) if ratpart u ~= 0 then l := concat(ratpart(u)::O, l) if not elem? u then l := concat([NE2O f for f in notelem u], l) null l => 0@F::O reduce("+", l) NE2O ne == int((ne.integrand)::O * hconcat ['d::O, (ne.intvar)::O]) @ \section{package IR2 IntegrationResultFunctions2} <>= )abbrev package IR2 IntegrationResultFunctions2 ++ Internally used by the integration packages ++ Author: Manuel Bronstein ++ Date Created: 1987 ++ Date Last Updated: 12 August 1992 ++ Keywords: integration. IntegrationResultFunctions2(E, F): Exports == Implementation where E : Field F : Field SE ==> Symbol Q ==> Fraction Integer IRE ==> IntegrationResult E IRF ==> IntegrationResult F UPE ==> SparseUnivariatePolynomial E UPF ==> SparseUnivariatePolynomial F NEE ==> Record(integrand:E, intvar:E) NEF ==> Record(integrand:F, intvar:F) LGE ==> Record(scalar:Q, coeff:UPE, logand:UPE) LGF ==> Record(scalar:Q, coeff:UPF, logand:UPF) NLE ==> Record(coeff:E, logand:E) NLF ==> Record(coeff:F, logand:F) UFE ==> Union(Record(mainpart:E, limitedlogs:List NLE), "failed") URE ==> Union(Record(ratpart:E, coeff:E), "failed") UE ==> Union(E, "failed") Exports ==> with map: (E -> F, IRE) -> IRF ++ map(f,ire) \undocumented map: (E -> F, URE) -> Union(Record(ratpart:F, coeff:F), "failed") ++ map(f,ure) \undocumented map: (E -> F, UE) -> Union(F, "failed") ++ map(f,ue) \undocumented map: (E -> F, UFE) -> Union(Record(mainpart:F, limitedlogs:List NLF), "failed") ++ map(f,ufe) \undocumented Implementation ==> add import SparseUnivariatePolynomialFunctions2(E, F) NEE2F: (E -> F, NEE) -> NEF LGE2F: (E -> F, LGE) -> LGF NLE2F: (E -> F, NLE) -> NLF NLE2F(func, r) == [func(r.coeff), func(r.logand)] NEE2F(func, n) == [func(n.integrand), func(n.intvar)] map(func:E -> F, u:UE) == (u case "failed" => "failed"; func(u::E)) map(func:E -> F, ir:IRE) == mkAnswer(func ratpart ir, [LGE2F(func, f) for f in logpart ir], [NEE2F(func, g) for g in notelem ir]) map(func:E -> F, u:URE) == u case "failed" => "failed" [func(u.ratpart), func(u.coeff)] map(func:E -> F, u:UFE) == u case "failed" => "failed" [func(u.mainpart), [NLE2F(func, f) for f in u.limitedlogs]] LGE2F(func, lg) == [lg.scalar, map(func, lg.coeff), map(func, lg.logand)] @ \section{License} <>= --Copyright (c) 1991-2002, The Numerical ALgorithms Group Ltd. --All rights reserved. --Copyright (C) 2007-2009, Gabriel Dos Reis. --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. @ <<*>>= <> -- SPAD files for the integration world should be compiled in the -- following order: -- -- INTAUX rderf intrf rdeef intef irexpand integrat <> <> @ \eject \begin{thebibliography}{99} \bibitem{1} nothing \end{thebibliography} \end{document}