(* :Name: NumericalMath`ListIntegrate` *)

(* :Title:  Approximate Integration of Functions That Are Known Only at a Few
		Distinct Points *)

(* :Author: Jerry B. Keiper *)

(* :Summary:
Given a function in the form of a sorted list of ordered pairs
(or simply as a list of ordinates, assuming equally spaced
abscissas), ListIntegrate gives an approximation to the definite
integral of that function from the first to the last abscissa.
*)

(*:Context: NumericalMath`ListIntegrate` *)

(*:Package Version: Mathematica 3.0 *)

(* :Copyright: Copyright 1990-1999, Wolfram Research, Inc.
*)

(* :History:
	Original version by Jerry B. Keiper, March 1989.
	Revised by Jerry B. Keiper, November 1990.
	Revised by Robert J. Knapp, January, 1995.
*)

(* :Keywords: integration, numerical integration, interpolation *)

(* :Source: any elementary numerical analysis textbook. *)

(* :Warnings:
	ListIntegrate should not be used if the data contain
	significant error.  In such cases the error should be smoothed
	using Fit, and the fitted function should be integrated
	using Integrate or NIntegrate.
*)

(*  :Comments
	This package is really superceded by modifications made to 
	Interpolation and InterpolatingFunction objects.  It is included
	for compatibility.
*)
BeginPackage["NumericalMath`ListIntegrate`"]

ListIntegrate::usage =
"ListIntegrate[{y0, y1, ..., yn}, h, k] uses an InterpolatingFunction
object to give an approximation to the integral of a function with
values equal to  y0,...,yn at points equally spaced a distance h
apart.  If k is odd it is increased by 1.  The default value of k is 4.
ListIntegrate[{{x0,y0}, {x1,y1}, ..., {xn,yn}}, k] can be used for
variable stepsize data.  If the data are known to contain errors, you
may be better off performing Integrate on the result of Fit
applied to the data."

Begin["NumericalMath`ListIntegrate`Private`"]

ListIntegrate[cl_?VectorQ, h_, k_Integer:4] :=
	Module[{ans},
		ans /; (ans = ListIntegrate0[h, cl, k]) =!= $Failed
	] /; k > 0;

ListIntegrate[cl_?MatrixQ, k_Integer:4] :=
	Module[{ans, dim = Length[Dimensions[cl]]},
		ans /; ((dim == 2) &&
			(Length[cl[[1]]] == 2) &&
			((ans = ListIntegrate1[cl, k]) =!= $Failed))
	] /; k > 0;


ListIntegrate0[h_, cl_, k_] := 
Module[{ifun,order = If[EvenQ[k],k-1,k],x},
	ifun = Interpolation[cl,InterpolationOrder->order];
	If[Head[ifun] === Interpolation,Return[$Failed]];
	h*Integrate[ifun[x],{x,1,Length[cl]}]]

ListIntegrate1[cl_, k_] := 
Module[{a = cl[[1,1]],b = (Last[cl])[[1]],order = If[EvenQ[k],k-1,k],ifun,x},
	ifun = Interpolation[cl,InterpolationOrder->order];
	If[Head[ifun] === Interpolation,Return[$Failed]];
	Integrate[ifun[x],{x,a,b}]]

End[] (* NumericalMath`ListIntegrate`Private` *)

Protect[ListIntegrate];

EndPackage[] (* NumericalMath`ListIntegrate` *)

