\documentclass{article} \usepackage{axiom} \begin{document} \title{\$SPAD/src/algebra d01Package.spad} \author{Brian Dupee} \maketitle \begin{abstract} \end{abstract} \eject \tableofcontents \eject \section{package INTPACK AnnaNumericalIntegrationPackage} <>= )abbrev package INTPACK AnnaNumericalIntegrationPackage ++ Author: Brian Dupee ++ Date Created: August 1994 ++ Date Last Updated: December 1997 ++ Basic Operations: integrate, measure ++ Related Constructors: Result, RoutinesTable ++ Also See: ++ AMS Classifications: ++ Keywords: ++ References: ++ Description: ++ \axiomType{AnnaNumericalIntegrationPackage} is a \axiom{package} ++ of functions for the \axiom{category} \axiomType{NumericalIntegrationCategory} ++ with \axiom{measure}, and \axiom{integrate}. EDF ==> Expression DoubleFloat DF ==> DoubleFloat EF ==> Expression Float F ==> Float INT ==> Integer SOCDF ==> Segment OrderedCompletion DoubleFloat OCDF ==> OrderedCompletion DoubleFloat SBOCF ==> SegmentBinding OrderedCompletion Float LSOCF ==> List Segment OrderedCompletion Float SOCF ==> Segment OrderedCompletion Float OCF ==> OrderedCompletion Float LS ==> List Symbol S ==> Symbol LST ==> List String ST ==> String RT ==> RoutinesTable NIA ==> Record(var:S, fn:EDF, range:SOCDF, abserr:DF, relerr:DF) MDNIA ==> Record(fn:EDF,range:List SOCDF,abserr:DF,relerr:DF) IFL ==> List(Record(ifail:Integer,instruction:String)) Entry ==> Record(chapter:String, type:String, domainName: String, defaultMin:F, measure:F, failList:IFL, explList:List String) Measure ==> Record(measure:F, name:ST, explanations:LST, extra:Result) AnnaNumericalIntegrationPackage(): with integrate: (EF,SOCF,F,F,RT) -> Result ++ integrate(exp, a..b, epsrel, routines) is a top level ANNA function ++ to integrate an expression, {\tt exp}, over a given range {\tt a} ++ to {\tt b} to the required absolute and relative accuracy using ++ the routines available in the RoutinesTable provided. ++ ++ It iterates over the \axiom{domains} of ++ \axiomType{NumericalIntegrationCategory} ++ to get the name and other ++ relevant information of the the (domain of the) numerical ++ routine likely to be the most appropriate, ++ i.e. have the best \axiom{measure}. ++ ++ It then performs the integration of the given expression ++ on that \axiom{domain}. integrate: NumericalIntegrationProblem -> Result ++ integrate(IntegrationProblem) is a top level ANNA function ++ to integrate an expression over a given range or ranges ++ to the required absolute and relative accuracy. ++ ++ It iterates over the \axiom{domains} of ++ \axiomType{NumericalIntegrationCategory} to get the name and other ++ relevant information of the the (domain of the) numerical ++ routine likely to be the most appropriate, ++ i.e. have the best \axiom{measure}. ++ ++ It then performs the integration of the given expression ++ on that \axiom{domain}. integrate: (EF,SOCF,F,F) -> Result ++ integrate(exp, a..b, epsabs, epsrel) is a top level ANNA function ++ to integrate an expression, {\tt exp}, over a given range {\tt a} ++ to {\tt b} to the required absolute and relative accuracy. ++ ++ It iterates over the \axiom{domains} of ++ \axiomType{NumericalIntegrationCategory} to get the name and other ++ relevant information of the the (domain of the) numerical ++ routine likely to be the most appropriate, ++ i.e. have the best \axiom{measure}. ++ ++ It then performs the integration of the given expression ++ on that \axiom{domain}. integrate: (EF,SOCF,F) -> Result ++ integrate(exp, a..b, epsrel) is a top level ANNA ++ function to integrate an expression, {\tt exp}, over a given ++ range {\tt a} to {\tt b} to the required relative accuracy. ++ ++ It iterates over the \axiom{domains} of ++ \axiomType{NumericalIntegrationCategory} to get the name and other ++ relevant information of the the (domain of the) numerical ++ routine likely to be the most appropriate, ++ i.e. have the best \axiom{measure}. ++ ++ It then performs the integration of the given expression ++ on that \axiom{domain}. ++ ++ If epsrel = 0, a default absolute accuracy is used. integrate: (EF,SOCF) -> Result ++ integrate(exp, a..b) is a top ++ level ANNA function to integrate an expression, {\tt exp}, ++ over a given range {\tt a} to {\tt b}. ++ ++ It iterates over the \axiom{domains} of ++ \axiomType{NumericalIntegrationCategory} to get the name and other ++ relevant information of the the (domain of the) numerical ++ routine likely to be the most appropriate, ++ i.e. have the best \axiom{measure}. ++ ++ It then performs the integration of the given expression ++ on that \axiom{domain}. ++ ++ Default values for the absolute and relative error are used. integrate:(EF,LSOCF) -> Result ++ integrate(exp, [a..b,c..d,...]) is a top ++ level ANNA function to integrate a multivariate expression, {\tt exp}, ++ over a given set of ranges. ++ ++ It iterates over the \axiom{domains} of ++ \axiomType{NumericalIntegrationCategory} to get the name and other ++ relevant information of the the (domain of the) numerical ++ routine likely to be the most appropriate, ++ i.e. have the best \axiom{measure}. ++ ++ It then performs the integration of the given expression ++ on that \axiom{domain}. ++ ++ Default values for the absolute and relative error are used. integrate:(EF,LSOCF,F) -> Result ++ integrate(exp, [a..b,c..d,...], epsrel) is a top ++ level ANNA function to integrate a multivariate expression, {\tt exp}, ++ over a given set of ranges to the required relative ++ accuracy. ++ ++ It iterates over the \axiom{domains} of ++ \axiomType{NumericalIntegrationCategory} to get the name and other ++ relevant information of the the (domain of the) numerical ++ routine likely to be the most appropriate, ++ i.e. have the best \axiom{measure}. ++ ++ It then performs the integration of the given expression ++ on that \axiom{domain}. ++ ++ If epsrel = 0, a default absolute accuracy is used. integrate:(EF,LSOCF,F,F) -> Result ++ integrate(exp, [a..b,c..d,...], epsabs, epsrel) is a top ++ level ANNA function to integrate a multivariate expression, {\tt exp}, ++ over a given set of ranges to the required absolute and relative ++ accuracy. ++ ++ It iterates over the \axiom{domains} of ++ \axiomType{NumericalIntegrationCategory} to get the name and other ++ relevant information of the the (domain of the) numerical ++ routine likely to be the most appropriate, ++ i.e. have the best \axiom{measure}. ++ ++ It then performs the integration of the given expression ++ on that \axiom{domain}. integrate:(EF,LSOCF,F,F,RT) -> Result ++ integrate(exp, [a..b,c..d,...], epsabs, epsrel, routines) is a top ++ level ANNA function to integrate a multivariate expression, {\tt exp}, ++ over a given set of ranges to the required absolute and relative ++ accuracy, using the routines available in the RoutinesTable provided. ++ ++ It iterates over the \axiom{domains} of ++ \axiomType{NumericalIntegrationCategory} to get the name and other ++ relevant information of the the (domain of the) numerical ++ routine likely to be the most appropriate, ++ i.e. have the best \axiom{measure}. ++ ++ It then performs the integration of the given expression ++ on that \axiom{domain}. measure:NumericalIntegrationProblem -> Measure ++ measure(prob) is a top level ANNA function for identifying the most ++ appropriate numerical routine for solving the numerical integration ++ problem defined by \axiom{prob}. ++ ++ It calls each \axiom{domain} of \axiom{category} ++ \axiomType{NumericalIntegrationCategory} in turn to calculate all measures ++ and returns the best ++ i.e. the name of the most appropriate domain and any other relevant ++ information. measure:(NumericalIntegrationProblem,RT) -> Measure ++ measure(prob,R) is a top level ANNA function for identifying the most ++ appropriate numerical routine from those in the routines table ++ provided for solving the numerical integration ++ problem defined by \axiom{prob}. ++ ++ It calls each \axiom{domain} listed in \axiom{R} of \axiom{category} ++ \axiomType{NumericalIntegrationCategory} in turn to calculate all measures ++ and returns the best ++ i.e. the name of the most appropriate domain and any other relevant ++ information. integrate:(EF,SBOCF,ST) -> Union(Result,"failed") ++ integrate(exp, x = a..b, "numerical") is a top level ANNA function to ++ integrate an expression, {\tt exp}, over a given range, {\tt a} ++ to {\tt b}. ++ ++ It iterates over the \axiom{domains} of ++ \axiomType{NumericalIntegrationCategory} to get the name and other ++ relevant information of the the (domain of the) numerical ++ routine likely to be the most appropriate, ++ i.e. have the best \axiom{measure}. ++ ++ It then performs the integration of the given expression ++ on that \axiom{domain}.\newline ++ ++ Default values for the absolute and relative error are used. ++ ++ It is an error of the last argument is not {\tt "numerical"}. integrate:(EF,SBOCF,S) -> Union(Result,"failed") ++ integrate(exp, x = a..b, numerical) is a top level ANNA function to ++ integrate an expression, {\tt exp}, over a given range, {\tt a} ++ to {\tt b}. ++ ++ It iterates over the \axiom{domains} of ++ \axiomType{NumericalIntegrationCategory} to get the name and other ++ relevant information of the the (domain of the) numerical ++ routine likely to be the most appropriate, ++ i.e. have the best \axiom{measure}. ++ ++ It then performs the integration of the given expression ++ on that \axiom{domain}.\newline ++ ++ Default values for the absolute and relative error are used. ++ ++ It is an error if the last argument is not {\tt numerical}. == add zeroMeasure: Measure -> Result scriptedVariables?: MDNIA -> Boolean preAnalysis:(Union(nia:NIA,mdnia:MDNIA),RT) -> RT measureSpecific:(ST,RT,Union(nia:NIA,mdnia:MDNIA)) -> Record(measure:F,explanations:LST,extra:Result) changeName:(Result,ST) -> Result recoverAfterFail:(Union(nia:NIA,mdnia:MDNIA),RT,Measure,INT,Result) -> Record(a:Result,b:Measure) better?:(Result,Result) -> Boolean integrateConstant:(EF,SOCF) -> Result integrateConstantList: (EF,LSOCF) -> Result integrateArgs:(NumericalIntegrationProblem,RT) -> Result integrateSpecific:(Union(nia:NIA,mdnia:MDNIA),ST,Result) -> Result import ExpertSystemToolsPackage integrateConstantList(exp:EF,ras:LSOCF):Result == c:OCF := ((retract(exp)@F)$EF)::OCF b := [hi(j)-lo(j) for j in ras] c := c*reduce((#1)*(#2),b) a := coerce(c)$AnyFunctions1(OCF) text := coerce("Constant Function")$AnyFunctions1(ST) construct([[result@S,a],[method@S,text]])$Result integrateConstant(exp:EF,ra:SOCF):Result == c := (retract(exp)@F)$EF r:OCF := (c::OCF)*(hi(ra)-lo(ra)) a := coerce(r)$AnyFunctions1(OCF) text := coerce("Constant Function")$AnyFunctions1(ST) construct([[result@S,a],[method@S,text]])$Result zeroMeasure(m:Measure):Result == a := coerce(0$DF)$AnyFunctions1(DF) text := coerce("Constant Function")$AnyFunctions1(String) r := construct([[result@Symbol,a],[method@Symbol,text]])$Result concat(measure2Result m,r)$ExpertSystemToolsPackage scriptedVariables?(mdnia:MDNIA):Boolean == vars:List Symbol := variables(mdnia.fn)$EDF var1 := first(vars)$(List Symbol) not scripted?(var1) => false name1 := name(var1)$Symbol for i in 2..# vars repeat not ((scripted?(vars.i)$Symbol) and (name1 = name(vars.i)$Symbol)) => return false true preAnalysis(args:Union(nia:NIA,mdnia:MDNIA),t:RT):RT == import RT r:RT := selectIntegrationRoutines t args case nia => arg:NIA := args.nia rangeIsFinite(arg)$d01AgentsPackage case finite => selectFiniteRoutines r selectNonFiniteRoutines r selectMultiDimensionalRoutines r changeName(ans:Result,name:ST):Result == sy:S := coerce(name "Answer")$S anyAns:Any := coerce(ans)$AnyFunctions1(Result) construct([[sy,anyAns]])$Result measureSpecific(name:ST,R:RT,args:Union(nia:NIA,mdnia:MDNIA)): Record(measure:F,explanations:ST,extra:Result) == args case nia => arg:NIA := args.nia name = "d01ajfAnnaType" => measure(R,arg)$d01ajfAnnaType name = "d01akfAnnaType" => measure(R,arg)$d01akfAnnaType name = "d01alfAnnaType" => measure(R,arg)$d01alfAnnaType name = "d01amfAnnaType" => measure(R,arg)$d01amfAnnaType name = "d01anfAnnaType" => measure(R,arg)$d01anfAnnaType name = "d01apfAnnaType" => measure(R,arg)$d01apfAnnaType name = "d01aqfAnnaType" => measure(R,arg)$d01aqfAnnaType name = "d01asfAnnaType" => measure(R,arg)$d01asfAnnaType name = "d01TransformFunctionType" => measure(R,arg)$d01TransformFunctionType error("measureSpecific","invalid type name: " name)$ErrorFunctions args case mdnia => arg2:MDNIA := args.mdnia name = "d01gbfAnnaType" => measure(R,arg2)$d01gbfAnnaType name = "d01fcfAnnaType" => measure(R,arg2)$d01fcfAnnaType error("measureSpecific","invalid type name: " name)$ErrorFunctions error("measureSpecific","invalid type name")$ErrorFunctions measure(a:NumericalIntegrationProblem,R:RT):Measure == args:Union(nia:NIA,mdnia:MDNIA) := retract(a)$NumericalIntegrationProblem sofar := 0$F best := "none" :: ST routs := copy R routs := preAnalysis(args,routs) empty?(routs)$RT => error("measure", "no routines found")$ErrorFunctions rout := inspect(routs)$RT e := retract(rout.entry)$AnyFunctions1(Entry) meth:LST := ["Trying " e.type " integration routines"] ext := empty()$Result for i in 1..# routs repeat rout := extract!(routs)$RT e := retract(rout.entry)$AnyFunctions1(Entry) n := e.domainName if e.defaultMin > sofar then m := measureSpecific(n,R,args) if m.measure > sofar then sofar := m.measure best := n ext := concat(m.extra,ext)$ExpertSystemToolsPackage str:LST := [string(rout.key)$S "measure: " outputMeasure(m.measure) " - " m.explanations] else str:LST := [string(rout.key)$S " is no better than other routines"] meth := append(meth,str)$LST [sofar,best,meth,ext] measure(a:NumericalIntegrationProblem):Measure == measure(a,routines()$RT) integrateSpecific(args:Union(nia:NIA,mdnia:MDNIA),n:ST,ex:Result):Result == args case nia => arg:NIA := args.nia n = "d01ajfAnnaType" => numericalIntegration(arg,ex)$d01ajfAnnaType n = "d01TransformFunctionType" => numericalIntegration(arg,ex)$d01TransformFunctionType n = "d01amfAnnaType" => numericalIntegration(arg,ex)$d01amfAnnaType n = "d01apfAnnaType" => numericalIntegration(arg,ex)$d01apfAnnaType n = "d01aqfAnnaType" => numericalIntegration(arg,ex)$d01aqfAnnaType n = "d01alfAnnaType" => numericalIntegration(arg,ex)$d01alfAnnaType n = "d01akfAnnaType" => numericalIntegration(arg,ex)$d01akfAnnaType n = "d01anfAnnaType" => numericalIntegration(arg,ex)$d01anfAnnaType n = "d01asfAnnaType" => numericalIntegration(arg,ex)$d01asfAnnaType error("integrateSpecific","invalid type name: " n)$ErrorFunctions args case mdnia => arg2:MDNIA := args.mdnia n = "d01gbfAnnaType" => numericalIntegration(arg2,ex)$d01gbfAnnaType n = "d01fcfAnnaType" => numericalIntegration(arg2,ex)$d01fcfAnnaType error("integrateSpecific","invalid type name: " n)$ErrorFunctions error("integrateSpecific","invalid type name: " n)$ErrorFunctions better?(r:Result,s:Result):Boolean == a1 := search("abserr"::S,r)$Result a1 case "failed" => false abserr1 := retract(a1)$AnyFunctions1(DF) negative?(abserr1) => false a2 := search("abserr"::S,s)$Result a2 case "failed" => true abserr2 := retract(a2)$AnyFunctions1(DF) negative?(abserr2) => true (abserr1 < abserr2) -- true if r.abserr better than s.abserr recoverAfterFail(n:Union(nia:NIA,mdnia:MDNIA),routs:RT,m:Measure,iint:INT, r:Result):Record(a:Result,b:Measure) == bestName := m.name while positive?(iint) repeat routineName := m.name s := recoverAfterFail(routs,routineName(1..6),iint)$RoutinesTable s case "failed" => iint := 0 if s = "changeEps" then nn := n.nia zero?(nn.abserr) => nn.abserr := 1.0e-8 :: DF m := measure(n::NumericalIntegrationProblem,routs) zero?(m.measure) => iint := 0 r := integrateSpecific(n,m.name,m.extra) iint := 0 rn := routineName(1..6) buttVal := getButtonValue(rn,"functionEvaluations")$AttributeButtons if (s = "incrFunEvals") and (buttVal < 0.8) then increase(rn,"functionEvaluations")$AttributeButtons if s = "increase tolerance" then (n.nia).relerr := (n.nia).relerr*(10.0::DF) if s = "decrease tolerance" then (n.nia).relerr := (n.nia).relerr/(10.0::DF) fl := coerce(s)$AnyFunctions1(ST) flrec:Record(key:S,entry:Any):=[failure@S,fl] m2 := measure(n::NumericalIntegrationProblem,routs) zero?(m2.measure) => iint := 0 r2:Result := integrateSpecific(n,m2.name,m2.extra) better?(r,r2) => m.name := m2.name insert!(flrec,r)$Result bestName := m2.name m := m2 insert!(flrec,r2)$Result r := concat(r2,changeName(r,routineName))$ExpertSystemToolsPackage iany := search(ifail@S,r2)$Result iany case "failed" => iint := 0 iint := retract(iany)$AnyFunctions1(INT) m.name := bestName [r,m] integrateArgs(prob:NumericalIntegrationProblem,t:RT):Result == args:Union(nia:NIA,mdnia:MDNIA) := retract(prob)$NumericalIntegrationProblem routs := copy(t)$RT if args case mdnia then arg := args.mdnia v := (# variables(arg.fn)) not scriptedVariables?(arg) => error("MultiDimensionalNumericalIntegrationPackage", "invalid variable names")$ErrorFunctions (v ~= # arg.range)@Boolean => error("MultiDimensionalNumericalIntegrationPackage", "number of variables do not match number of ranges")$ErrorFunctions m := measure(prob,routs) zero?(m.measure) => zeroMeasure m r := integrateSpecific(args,m.name,m.extra) iany := search(ifail@S,r)$Result iint := 0$INT if (iany case Any) then iint := retract(iany)$AnyFunctions1(INT) if positive?(iint) then tu:Record(a:Result,b:Measure) := recoverAfterFail(args,routs,m,iint,r) r := tu.a m := tu.b r := concat(measure2Result m,r)$ExpertSystemToolsPackage n := m.name nn:ST := (# n > 14) => "d01transform" n(1..6) expl := getExplanations(routs,nn)$RoutinesTable expla := coerce(expl)$AnyFunctions1(LST) explaa:Record(key:Symbol,entry:Any) := ["explanations"::Symbol,expla] r := concat(construct([explaa]),r) args case nia => att := showAttributes(args.nia)$IntegrationFunctionsTable att case "failed" => r concat(att2Result att,r)$ExpertSystemToolsPackage r integrate(args:NumericalIntegrationProblem):Result == integrateArgs(args,routines()$RT) integrate(exp:EF,ra:SOCF,epsabs:F,epsrel:F,r:RT):Result == Var:LS := variables(exp)$EF empty?(Var)$LS => integrateConstant(exp,ra) args:NIA := [first(Var)$LS,ef2edf exp,socf2socdf ra,f2df epsabs,f2df epsrel] integrateArgs(args::NumericalIntegrationProblem,r) integrate(exp:EF,ra:SOCF,epsabs:F,epsrel:F):Result == integrate(exp,ra,epsabs,epsrel,routines()$RT) integrate(exp:EF,ra:SOCF,err:F):Result == positive?(err)$F => integrate(exp,ra,0$F,err) integrate(exp,ra,1.0E-5,err) integrate(exp:EF,ra:SOCF):Result == integrate(exp,ra,0$F,1.0E-5) integrate(exp:EF,sb:SBOCF, st:ST) == st = "numerical" => integrate(exp,segment sb) "failed" integrate(exp:EF,sb:SBOCF, s:S) == s = (numerical::Symbol) => integrate(exp,segment sb) "failed" integrate(exp:EF,ra:LSOCF,epsabs:F,epsrel:F,r:RT):Result == vars := variables(exp)$EF empty?(vars)$LS => integrateConstantList(exp,ra) args:MDNIA := [ef2edf exp,convert ra,f2df epsabs,f2df epsrel] integrateArgs(args::NumericalIntegrationProblem,r) integrate(exp:EF,ra:LSOCF,epsabs:F,epsrel:F):Result == integrate(exp,ra,epsabs,epsrel,routines()$RT) integrate(exp:EF,ra:LSOCF,epsrel:F):Result == zero? epsrel => integrate(exp,ra,1.0e-6,epsrel) integrate(exp,ra,0$F,epsrel) integrate(exp:EF,ra:LSOCF):Result == integrate(exp,ra,1.0e-4) @ \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}