diff options
Diffstat (limited to 'src/algebra/numquad.spad.pamphlet')
| -rw-r--r-- | src/algebra/numquad.spad.pamphlet | 600 |
1 files changed, 600 insertions, 0 deletions
diff --git a/src/algebra/numquad.spad.pamphlet b/src/algebra/numquad.spad.pamphlet new file mode 100644 index 00000000..79987429 --- /dev/null +++ b/src/algebra/numquad.spad.pamphlet @@ -0,0 +1,600 @@ +\documentclass{article} +\usepackage{axiom} +\begin{document} +\title{\$SPAD/src/algebra numquad.spad} +\author{Yurij Baransky} +\maketitle +\begin{abstract} +\end{abstract} +\eject +\tableofcontents +\eject +\section{package NUMQUAD NumericalQuadrature} +<<package NUMQUAD NumericalQuadrature>>= +)abbrev package NUMQUAD NumericalQuadrature +++ Author: Yurij A. Baransky +++ Date Created: October 90 +++ Date Last Updated: October 90 +++ Basic Operations: +++ Related Constructors: +++ Also See: +++ AMS Classifications: +++ Keywords: +++ References: +++ Description: +++ This suite of routines performs numerical quadrature using +++ algorithms derived from the basic trapezoidal rule. Because +++ the error term of this rule contains only even powers of the +++ step size (for open and closed versions), fast convergence +++ can be obtained if the integrand is sufficiently smooth. +++ +++ Each routine returns a Record of type TrapAns, which contains\indent{3} +++ \newline value (\spadtype{Float}):\tab{20} estimate of the integral +++ \newline error (\spadtype{Float}):\tab{20} estimate of the error in the computation +++ \newline totalpts (\spadtype{Integer}):\tab{20} total number of function evaluations +++ \newline success (\spadtype{Boolean}):\tab{20} if the integral was computed within the user specified error criterion +++ \indent{0}\indent{0} +++ To produce this estimate, each routine generates an internal +++ sequence of sub-estimates, denoted by {\em S(i)}, depending on the +++ routine, to which the various convergence criteria are applied. +++ The user must supply a relative accuracy, \spad{eps_r}, and an absolute +++ accuracy, \spad{eps_a}. Convergence is obtained when either +++ \center{\spad{ABS(S(i) - S(i-1)) < eps_r * ABS(S(i-1))}} +++ \center{or \spad{ABS(S(i) - S(i-1)) < eps_a}} +++ are true statements. +++ +++ The routines come in three families and three flavors: +++ \newline\tab{3} closed:\tab{20}romberg,\tab{30}simpson,\tab{42}trapezoidal +++ \newline\tab{3} open: \tab{20}rombergo,\tab{30}simpsono,\tab{42}trapezoidalo +++ \newline\tab{3} adaptive closed:\tab{20}aromberg,\tab{30}asimpson,\tab{42}atrapezoidal +++ \par +++ The {\em S(i)} for the trapezoidal family is the value of the +++ integral using an equally spaced absicca trapezoidal rule for +++ that level of refinement. +++ \par +++ The {\em S(i)} for the simpson family is the value of the integral +++ using an equally spaced absicca simpson rule for that level of +++ refinement. +++ \par +++ The {\em S(i)} for the romberg family is the estimate of the integral +++ using an equally spaced absicca romberg method. For +++ the \spad{i}-th level, this is an appropriate combination of all the +++ previous trapezodial estimates so that the error term starts +++ with the \spad{2*(i+1)} power only. +++ \par +++ The three families come in a closed version, where the formulas +++ include the endpoints, an open version where the formulas do not +++ include the endpoints and an adaptive version, where the user +++ is required to input the number of subintervals over which the +++ appropriate closed family integrator will apply with the usual +++ convergence parmeters for each subinterval. This is useful +++ where a large number of points are needed only in a small fraction +++ of the entire domain. +++ \par +++ Each routine takes as arguments: +++ \newline f\tab{10} integrand +++ \newline a\tab{10} starting point +++ \newline b\tab{10} ending point +++ \newline \spad{eps_r}\tab{10} relative error +++ \newline \spad{eps_a}\tab{10} absolute error +++ \newline \spad{nmin} \tab{10} refinement level when to start checking for convergence (> 1) +++ \newline \spad{nmax} \tab{10} maximum level of refinement +++ \par +++ The adaptive routines take as an additional parameter +++ \newline \spad{nint}\tab{10} the number of independent intervals to apply a closed +++ family integrator of the same name. +++ \par Notes: +++ \newline Closed family level i uses \spad{1 + 2**i} points. +++ \newline Open family level i uses \spad{1 + 3**i} points. +NumericalQuadrature(): Exports == Implementation where + L ==> List + V ==> Vector + I ==> Integer + B ==> Boolean + E ==> OutputForm + F ==> Float + PI ==> PositiveInteger + OFORM ==> OutputForm + TrapAns ==> Record(value:F, error:F, totalpts:I, success:B ) + + Exports ==> with + aromberg : (F -> F,F,F,F,F,I,I,I) -> TrapAns + ++ aromberg(fn,a,b,epsrel,epsabs,nmin,nmax,nint) + ++ uses the adaptive romberg method to numerically integrate function + ++ \spad{fn} over the closed interval from \spad{a} to \spad{b}, + ++ with relative accuracy \spad{epsrel} and absolute accuracy + ++ \spad{epsabs}, with the refinement levels for convergence checking + ++ vary from \spad{nmin} to \spad{nmax}, and where \spad{nint} + ++ is the number of independent intervals to apply the integrator. + ++ The value returned is a record containing the value of the integral, + ++ the estimate of the error in the computation, the total number of + ++ function evaluations, and either a boolean value which is true if + ++ the integral was computed within the user specified error criterion. + ++ See \spadtype{NumericalQuadrature} for details. + asimpson : (F -> F,F,F,F,F,I,I,I) -> TrapAns + ++ asimpson(fn,a,b,epsrel,epsabs,nmin,nmax,nint) uses the + ++ adaptive simpson method to numerically integrate function \spad{fn} + ++ over the closed interval from \spad{a} to \spad{b}, with relative + ++ accuracy \spad{epsrel} and absolute accuracy \spad{epsabs}, with the + ++ refinement levels for convergence checking vary from \spad{nmin} + ++ to \spad{nmax}, and where \spad{nint} is the number of independent + ++ intervals to apply the integrator. The value returned is a record + ++ containing the value of the integral, the estimate of the error in + ++ the computation, the total number of function evaluations, and + ++ either a boolean value which is true if the integral was computed + ++ within the user specified error criterion. + ++ See \spadtype{NumericalQuadrature} for details. + atrapezoidal : (F -> F,F,F,F,F,I,I,I) -> TrapAns + ++ atrapezoidal(fn,a,b,epsrel,epsabs,nmin,nmax,nint) uses the + ++ adaptive trapezoidal method to numerically integrate function + ++ \spad{fn} over the closed interval from \spad{a} to \spad{b}, with + ++ relative accuracy \spad{epsrel} and absolute accuracy \spad{epsabs}, + ++ with the refinement levels for convergence checking vary from + ++ \spad{nmin} to \spad{nmax}, and where \spad{nint} is the number + ++ of independent intervals to apply the integrator. The value returned + ++ is a record containing the value of the integral, the estimate of + ++ the error in the computation, the total number of function + ++ evaluations, and either a boolean value which is true if + ++ the integral was computed within the user specified error criterion. + ++ See \spadtype{NumericalQuadrature} for details. + romberg : (F -> F,F,F,F,F,I,I) -> TrapAns + ++ romberg(fn,a,b,epsrel,epsabs,nmin,nmax) uses the romberg + ++ method to numerically integrate function \spadvar{fn} over the closed + ++ interval \spad{a} to \spad{b}, with relative accuracy \spad{epsrel} + ++ and absolute accuracy \spad{epsabs}, with the refinement levels + ++ for convergence checking vary from \spad{nmin} to \spad{nmax}. + ++ The value returned is a record containing the value + ++ of the integral, the estimate of the error in the computation, the + ++ total number of function evaluations, and either a boolean value + ++ which is true if the integral was computed within the user specified + ++ error criterion. See \spadtype{NumericalQuadrature} for details. + simpson : (F -> F,F,F,F,F,I,I) -> TrapAns + ++ simpson(fn,a,b,epsrel,epsabs,nmin,nmax) uses the simpson + ++ method to numerically integrate function \spad{fn} over the closed + ++ interval \spad{a} to \spad{b}, with + ++ relative accuracy \spad{epsrel} and absolute accuracy \spad{epsabs}, + ++ with the refinement levels for convergence checking vary from + ++ \spad{nmin} to \spad{nmax}. The value returned + ++ is a record containing the value of the integral, the estimate of + ++ the error in the computation, the total number of function + ++ evaluations, and either a boolean value which is true if + ++ the integral was computed within the user specified error criterion. + ++ See \spadtype{NumericalQuadrature} for details. + trapezoidal : (F -> F,F,F,F,F,I,I) -> TrapAns + ++ trapezoidal(fn,a,b,epsrel,epsabs,nmin,nmax) uses the + ++ trapezoidal method to numerically integrate function \spadvar{fn} over + ++ the closed interval \spad{a} to \spad{b}, with relative accuracy + ++ \spad{epsrel} and absolute accuracy \spad{epsabs}, with the + ++ refinement levels for convergence checking vary + ++ from \spad{nmin} to \spad{nmax}. The value + ++ returned is a record containing the value of the integral, the + ++ estimate of the error in the computation, the total number of + ++ function evaluations, and either a boolean value which is true + ++ if the integral was computed within the user specified error criterion. + ++ See \spadtype{NumericalQuadrature} for details. + rombergo : (F -> F,F,F,F,F,I,I) -> TrapAns + ++ rombergo(fn,a,b,epsrel,epsabs,nmin,nmax) uses the romberg + ++ method to numerically integrate function \spad{fn} over + ++ the open interval from \spad{a} to \spad{b}, with + ++ relative accuracy \spad{epsrel} and absolute accuracy \spad{epsabs}, + ++ with the refinement levels for convergence checking vary from + ++ \spad{nmin} to \spad{nmax}. The value returned + ++ is a record containing the value of the integral, the estimate of + ++ the error in the computation, the total number of function + ++ evaluations, and either a boolean value which is true if + ++ the integral was computed within the user specified error criterion. + ++ See \spadtype{NumericalQuadrature} for details. + simpsono : (F -> F,F,F,F,F,I,I) -> TrapAns + ++ simpsono(fn,a,b,epsrel,epsabs,nmin,nmax) uses the + ++ simpson method to numerically integrate function \spad{fn} over + ++ the open interval from \spad{a} to \spad{b}, with + ++ relative accuracy \spad{epsrel} and absolute accuracy \spad{epsabs}, + ++ with the refinement levels for convergence checking vary from + ++ \spad{nmin} to \spad{nmax}. The value returned + ++ is a record containing the value of the integral, the estimate of + ++ the error in the computation, the total number of function + ++ evaluations, and either a boolean value which is true if + ++ the integral was computed within the user specified error criterion. + ++ See \spadtype{NumericalQuadrature} for details. + trapezoidalo : (F -> F,F,F,F,F,I,I) -> TrapAns + ++ trapezoidalo(fn,a,b,epsrel,epsabs,nmin,nmax) uses the + ++ trapezoidal method to numerically integrate function \spad{fn} + ++ over the open interval from \spad{a} to \spad{b}, with + ++ relative accuracy \spad{epsrel} and absolute accuracy \spad{epsabs}, + ++ with the refinement levels for convergence checking vary from + ++ \spad{nmin} to \spad{nmax}. The value returned + ++ is a record containing the value of the integral, the estimate of + ++ the error in the computation, the total number of function + ++ evaluations, and either a boolean value which is true if + ++ the integral was computed within the user specified error criterion. + ++ See \spadtype{NumericalQuadrature} for details. + + Implementation ==> add + trapclosed : (F -> F,F,F,F,I) -> F + trapopen : (F -> F,F,F,F,I) -> F + import OutputPackage + +--------------------------------------------------- + + aromberg(func,a,b,epsrel,epsabs,nmin,nmax,nint) == + ans : TrapAns + sum : F := 0.0 + err : F := 0.0 + pts : I := 1 + done : B := true + hh : F := (b-a) / nint + x1 : F := a + x2 : F := a + hh + io : L OFORM := [x1::E,x2::E] + i : I + for i in 1..nint repeat + ans := romberg(func,x1,x2,epsrel,epsabs,nmin,nmax) + if (not ans.success) then + io.1 := x1::E + io.2 := x2::E + print blankSeparate cons("accuracy not reached in interval"::E,io) + sum := sum + ans.value + err := err + abs(ans.error) + pts := pts + ans.totalpts-1 + done := (done and ans.success) + x1 := x2 + x2 := x2 + hh + return( [sum , err , pts , done] ) + +--------------------------------------------------- + + asimpson(func,a,b,epsrel,epsabs,nmin,nmax,nint) == + ans : TrapAns + sum : F := 0.0 + err : F := 0.0 + pts : I := 1 + done : B := true + hh : F := (b-a) / nint + x1 : F := a + x2 : F := a + hh + io : L OFORM := [x1::E,x2::E] + i : I + for i in 1..nint repeat + ans := simpson(func,x1,x2,epsrel,epsabs,nmin,nmax) + if (not ans.success) then + io.1 := x1::E + io.2 := x2::E + print blankSeparate cons("accuracy not reached in interval"::E,io) + sum := sum + ans.value + err := err + abs(ans.error) + pts := pts + ans.totalpts-1 + done := (done and ans.success) + x1 := x2 + x2 := x2 + hh + return( [sum , err , pts , done] ) + +--------------------------------------------------- + + atrapezoidal(func,a,b,epsrel,epsabs,nmin,nmax,nint) == + ans : TrapAns + sum : F := 0.0 + err : F := 0.0 + pts : I := 1 + i : I + done : B := true + hh : F := (b-a) / nint + x1 : F := a + x2 : F := a + hh + io : L OFORM := [x1::E,x2::E] + for i in 1..nint repeat + ans := trapezoidal(func,x1,x2,epsrel,epsabs,nmin,nmax) + if (not ans.success) then + io.1 := x1::E + io.2 := x2::E + print blankSeparate cons("accuracy not reached in interval"::E,io) + sum := sum + ans.value + err := err + abs(ans.error) + pts := pts + ans.totalpts-1 + done := (done and ans.success) + x1 := x2 + x2 := x2 + hh + return( [sum , err , pts , done] ) + +--------------------------------------------------- + + romberg(func,a,b,epsrel,epsabs,nmin,nmax) == + length : F := (b-a) + delta : F := length + newsum : F := 0.5 * length * (func(a)+func(b)) + newest : F := 0.0 + oldsum : F := 0.0 + oldest : F := 0.0 + change : F := 0.0 + qx1 : F := newsum + table : V F := new((nmax+1)::PI,0.0) + n : I := 1 + pts : I := 1 + four : I + j : I + i : I + if (nmin < 2) then + output("romberg: nmin to small (nmin > 1) nmin = ",nmin::E) + return([0.0,0.0,0,false]) + if (nmax < nmin) then + output("romberg: nmax < nmin : nmax = ",nmax::E) + output(" nmin = ",nmin::E) + return([0.0,0.0,0,false]) + if (a = b) then + output("romberg: integration limits are equal = ",a::E) + return([0.0,0.0,1,true]) + if (epsrel < 0.0) then + output("romberg: eps_r < 0.0 eps_r = ",epsrel::E) + return([0.0,0.0,0,false]) + if (epsabs < 0.0) then + output("romberg: eps_a < 0.0 eps_a = ",epsabs::E) + return([0.0,0.0,0,false]) + for n in 1..nmax repeat + oldsum := newsum + newsum := trapclosed(func,a,delta,oldsum,pts) + newest := (4.0 * newsum - oldsum) / 3.0 + four := 4 + table(n) := newest + for j in 2..n repeat + i := n+1-j + four := four * 4 + table(i) := table(i+1) + (table(i+1)-table(i)) / (four-1) + if n > nmin then + change := abs(table(1) - qx1) + if change < abs(epsrel*qx1) then + return( [table(1) , change , 2*pts+1 , true] ) + if change < epsabs then + return( [table(1) , change , 2*pts+1 , true] ) + oldsum := newsum + oldest := newest + delta := 0.5*delta + pts := 2*pts + qx1 := table(1) + return( [table(1) , 1.25*change , pts+1 ,false] ) + +--------------------------------------------------- + + simpson(func,a,b,epsrel,epsabs,nmin,nmax) == + length : F := (b-a) + delta : F := length + newsum : F := 0.5*(b-a)*(func(a)+func(b)) + newest : F := 0.0 + oldsum : F := 0.0 + oldest : F := 0.0 + change : F := 0.0 + n : I := 1 + pts : I := 1 + if (nmin < 2) then + output("simpson: nmin to small (nmin > 1) nmin = ",nmin::E) + return([0.0,0.0,0,false]) + if (nmax < nmin) then + output("simpson: nmax < nmin : nmax = ",nmax::E) + output(" nmin = ",nmin::E) + return([0.0,0.0,0,false]) + if (a = b) then + output("simpson: integration limits are equal = ",a::E) + return([0.0,0.0,1,true]) + if (epsrel < 0.0) then + output("simpson: eps_r < 0.0 : eps_r = ",epsrel::E) + return([0.0,0.0,0,false]) + if (epsabs < 0.0) then + output("simpson: eps_a < 0.0 : eps_a = ",epsabs::E) + return([0.0,0.0,0,false]) + for n in 1..nmax repeat + oldsum := newsum + newsum := trapclosed(func,a,delta,oldsum,pts) + newest := (4.0 * newsum - oldsum) / 3.0 + if n > nmin then + change := abs(newest-oldest) + if change < abs(epsrel*oldest) then + return( [newest , 1.25*change , 2*pts+1 , true] ) + if change < epsabs then + return( [newest , 1.25*change , 2*pts+1 , true] ) + oldsum := newsum + oldest := newest + delta := 0.5*delta + pts := 2*pts + return( [newest , 1.25*change , pts+1 ,false] ) + +--------------------------------------------------- + + trapezoidal(func,a,b,epsrel,epsabs,nmin,nmax) == + length : F := (b-a) + delta : F := length + newsum : F := 0.5*(b-a)*(func(a)+func(b)) + change : F := 0.0 + oldsum : F + n : I := 1 + pts : I := 1 + if (nmin < 2) then + output("trapezoidal: nmin to small (nmin > 1) nmin = ",nmin::E) + return([0.0,0.0,0,false]) + if (nmax < nmin) then + output("trapezoidal: nmax < nmin : nmax = ",nmax::E) + output(" nmin = ",nmin::E) + return([0.0,0.0,0,false]) + if (a = b) then + output("trapezoidal: integration limits are equal = ",a::E) + return([0.0,0.0,1,true]) + if (epsrel < 0.0) then + output("trapezoidal: eps_r < 0.0 : eps_r = ",epsrel::E) + return([0.0,0.0,0,false]) + if (epsabs < 0.0) then + output("trapezoidal: eps_a < 0.0 : eps_a = ",epsabs::E) + return([0.0,0.0,0,false]) + for n in 1..nmax repeat + oldsum := newsum + newsum := trapclosed(func,a,delta,oldsum,pts) + if n > nmin then + change := abs(newsum-oldsum) + if change < abs(epsrel*oldsum) then + return( [newsum , 1.25*change , 2*pts+1 , true] ) + if change < epsabs then + return( [newsum , 1.25*change , 2*pts+1 , true] ) + delta := 0.5*delta + pts := 2*pts + return( [newsum , 1.25*change , pts+1 ,false] ) + +--------------------------------------------------- + + rombergo(func,a,b,epsrel,epsabs,nmin,nmax) == + length : F := (b-a) + delta : F := length / 3.0 + newsum : F := length * func( 0.5*(a+b) ) + newest : F := 0.0 + oldsum : F := 0.0 + oldest : F := 0.0 + change : F := 0.0 + qx1 : F := newsum + table : V F := new((nmax+1)::PI,0.0) + four : I + j : I + i : I + n : I := 1 + pts : I := 1 + for n in 1..nmax repeat + oldsum := newsum + newsum := trapopen(func,a,delta,oldsum,pts) + newest := (9.0 * newsum - oldsum) / 8.0 + table(n) := newest + nine := 9 + output(newest::E) + for j in 2..n repeat + i := n+1-j + nine := nine * 9 + table(i) := table(i+1) + (table(i+1)-table(i)) / (nine-1) + if n > nmin then + change := abs(table(1) - qx1) + if change < abs(epsrel*qx1) then + return( [table(1) , 1.5*change , 3*pts , true] ) + if change < epsabs then + return( [table(1) , 1.5*change , 3*pts , true] ) + output(table::E) + oldsum := newsum + oldest := newest + delta := delta / 3.0 + pts := 3*pts + qx1 := table(1) + return( [table(1) , 1.5*change , pts ,false] ) + +--------------------------------------------------- + + simpsono(func,a,b,epsrel,epsabs,nmin,nmax) == + length : F := (b-a) + delta : F := length / 3.0 + newsum : F := length * func( 0.5*(a+b) ) + newest : F := 0.0 + oldsum : F := 0.0 + oldest : F := 0.0 + change : F := 0.0 + n : I := 1 + pts : I := 1 + for n in 1..nmax repeat + oldsum := newsum + newsum := trapopen(func,a,delta,oldsum,pts) + newest := (9.0 * newsum - oldsum) / 8.0 + output(newest::E) + if n > nmin then + change := abs(newest - oldest) + if change < abs(epsrel*oldest) then + return( [newest , 1.5*change , 3*pts , true] ) + if change < epsabs then + return( [newest , 1.5*change , 3*pts , true] ) + oldsum := newsum + oldest := newest + delta := delta / 3.0 + pts := 3*pts + return( [newest , 1.5*change , pts ,false] ) + +--------------------------------------------------- + + trapezoidalo(func,a,b,epsrel,epsabs,nmin,nmax) == + length : F := (b-a) + delta : F := length/3.0 + newsum : F := length*func( 0.5*(a+b) ) + change : F := 0.0 + pts : I := 1 + oldsum : F + n : I + for n in 1..nmax repeat + oldsum := newsum + newsum := trapopen(func,a,delta,oldsum,pts) + output(newsum::E) + if n > nmin then + change := abs(newsum-oldsum) + if change < abs(epsrel*oldsum) then + return([newsum , 1.5*change , 3*pts , true] ) + if change < epsabs then + return([newsum , 1.5*change , 3*pts , true] ) + delta := delta / 3.0 + pts := 3*pts + return([newsum , 1.5*change , pts ,false] ) + +--------------------------------------------------- + + trapclosed(func,start,h,oldsum,numpoints) == + x : F := start + 0.5*h + sum : F := 0.0 + i : I + for i in 1..numpoints repeat + sum := sum + func(x) + x := x + h + return( 0.5*(oldsum + sum*h) ) + +--------------------------------------------------- + + trapopen(func,start,del,oldsum,numpoints) == + ddel : F := 2.0*del + x : F := start + 0.5*del + sum : F := 0.0 + i : I + for i in 1..numpoints repeat + sum := sum + func(x) + x := x + ddel + sum := sum + func(x) + x := x + del + return( (oldsum/3.0 + sum*del) ) + +@ +\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>> + +<<package NUMQUAD NumericalQuadrature>> +@ +\eject +\begin{thebibliography}{99} +\bibitem{1} nothing +\end{thebibliography} +\end{document} |