BeginPackage["VisualIntegration`","Graphics`FilledPlot`",
"Graphics`Colors`"]

LeftRSum::usage="LeftRSum[f[x],{x,a,b},n] visualizes and computes the 
left Riemann Sum for the function f[x] on the interval [a,b], 
using n equal-length subintervals.
(The value for n is optional; the default is 20.)"
RightRSum::usage="RightRSum[f[x],{x,a,b},n] visualizes and computes the 
right Riemann Sum for the function f[x] on the interval [a,b], 
using n equal-length subintervals.
(The value for n is optional; the default is 20.)"
DarbouxSum::usage="DarbouxSum[f[x],{x,a,b},n,s] visualizes 
and computes the lower and upper Darboux Sums for the 
function f[x] on the interval [a,b], 
using n equal-length subintervals.
The maximum and the minimum on each subinterval is computed
by sampling f[x] at (s+1) equally-spaced points in each 
subinterval.
(The values for n and s are optional; the defaults are 20,
and 10, respectively.)"
MidRSum::usage="MidRSum[f[x],{x,a,b},n] visualizes and computes the 
midpoint Riemann Sum for the function f[x] on the 
interval [a,b], using n equal-length subintervals.
(The value for n is optional; the default is 20.)"
TrapRule::usage="TrapRule[f[x],{x,a,b},n] visualizes and computes an approximation 
for the integral of the function f[x] on the interval [a,b],
using the Trapezoidal Rule with n equal-length subintervals.  
(The value for n is optional; the default is 20.)"
SimpRule::usage="SimpRule[f[x],{x,a,b},n] visualizes and computes an approximation 
for the integral of the function f[x] on the interval [a,b],
using Simpson's Rule with n equal-length subintervals. 
(The (even) value for n is optional; the default is 10.)"


Begin["`Private`"]

Threshold=50;


LeftRSum[f_,{x_,a_,b_},n_Integer:20,options___Rule]:=
Module[
{h,g,fct},
(h=(b-a)/n;
g[s_]=f/.{x->s};
Compile{x,g[x]};
fct=Plot[g[x],{x,a,b},DisplayFunction->Identity];
Show[Graphics[
Table[{Blue,
Rectangle[{j,0},{j+h,g[j]}]},{j,a,b-h,h}]],fct,
If[n<=Threshold,Graphics[Table[{Black,Line[{{j,0},{j,g[j]}}]},{j,a,b,h}]],{}],
Axes->False,Frame->True,PlotRange->{{a,b},All},options];
Print["The Left Riemann Sum with ",n," intervals is  ",
N[Sum[g[j],{j,a,b-h,h}]*h]])]

RightRSum[f_,{x_,a_,b_},n_Integer:20,options___Rule]:=
Module[
{h,g,fct},
(h=(b-a)/n;
g[s_]=f/.{x->s};
Compile{x,g[x]};
fct=Plot[g[x],{x,a,b},DisplayFunction->Identity];
Show[Graphics[
Table[{Red,
Rectangle[{j,0},{j+h,g[j+h]}]},{j,a,b-h,h}]],fct,
If[n<=Threshold,Graphics[Table[{Black,Line[{{j,0},{j,g[j]}}]},{j,a,b,h}]],{}],
Axes->False,Frame->True,PlotRange->{{a,b},All},options];
Print["The Right Riemann Sum with ",n," intervals is  ",
N[Sum[g[j],{j,a+h,b,h}]*h]])]





DarbouxSum[f_,{x_,a_,b_},n_Integer:20,split_Integer:10,options___Rule]:=
Module[
{h,g,fct,min,max},
(h=(b-a)/n;
g[s_]=f/.{x->s};
Compile{x,g[x]};
fct=Plot[g[x],{x,a,b},DisplayFunction->Identity];
min[j_]=Min[Table[g[j+k*h/split],{k,0,split}]];
max[j_]=Max[Table[g[j+k*h/split],{k,0,split}]];
rect=Flatten[
Table[
If[N[max[j]]<0,
{{Blue,Rectangle[{j,0},{j+h,min[j]}]},
{Red,Rectangle[{j,0},{j+h,max[j]}]}},
{{Red,Rectangle[{j,0},{j+h,max[j]}]},
{Blue,Rectangle[{j,0},{j+h,min[j]}]}}],
{j,a,b-h,h}],1];
Show[Graphics[rect],
fct,
Axes->False,Frame->True,PlotRange->{{a,b},All},options];
Print["The Lower Darboux Sum (in blue) with ",n," intervals is  ",
N[Sum[min[j],{j,a,b-h,h}]*h]];
Print["The Upper Darboux Sum (in red) with ",n," intervals is  ",
N[Sum[max[j],{j,a,b-h,h}]*h]])]

MidRSum[f_,{x_,a_,b_},n_Integer:20,options___Rule]:=
Module[
{h,g,fct},
(h=(b-a)/n;
g[s_]=f/.{x->s};
Compile{x,g[x]};
fct=Plot[g[x],{x,a,b},DisplayFunction->Identity];
Show[Graphics[
Table[{Purple,
Rectangle[{j,0},{j+h,g[j+h/2]}]},{j,a,b-h,h}]],fct,
If[n<=Threshold,Graphics[Table[{Black,Line[{{j,0},{j,g[j]}}]},{j,a,b,h}]],{}],
Axes->False,Frame->True,PlotRange->{{a,b},All},options];
Print["The Midpoint Riemann Sum with ",n, " intervals is  "
,N[Sum[g[j+h/2],{j,a,b-h,h}]*h]])]

TrapRule[f_,{x_,a_,b_},n_Integer:20,options___Rule]:=
Module[
{h,g,fct},
(h=(b-a)/n;
g[s_]=f/.{x->s};
Compile{x,g[x]};
fct=Plot[g[x],{x,a,b},PlotPoints->200,DisplayFunction->Identity];
Show[Graphics[
Table[{Yellow,
Polygon[{{j,0},{j+h,0},{j+h,g[j+h]},{j,g[j]}}]},{j,a,b-h,h}]],fct,
If[n<=Threshold,Graphics[Table[{Black,Line[{{j,0},{j,g[j]}}]},{j,a,b,h}]],{}],
Axes->False,Frame->True,PlotRange->{{a,b},All}];
Print["The Trapezoidal Rule with ",n," intervals"];
Print["yields the result   ",
(2*N[Sum[g[j],{j,a+h,b-h,h}]]+N[g[a]]+N[g[b]])*h/2
])]


SimpRule[f_,{x_,a_,b_},n_Integer:10]:=
Module[
{h,g,fct,m},
(If[EvenQ[n],m=n,(m=n+1;
Print["The number of intervals must be even."];
Print["Number of intervals increased to ",m,"."])];
h=(b-a)/m;
g[s_]=f/.{x->s};
Compile{x,g[x]};
fct=Plot[g[x],{x,a,b},PlotPoints->200,DisplayFunction->Identity];
Show[Table[
FilledPlot[
{0,Evaluate[Fit[{{j,g[j]},{j+h,g[j+h]},{j+2h,g[j+2h]}},{1,x,x^2},x]]},
{x,j,j+2h},DisplayFunction->Identity,
Fills->Green,PlotStyle->{Black,Green}],{j,a,b-2h,2h}],fct,
If[n<=Threshold,Graphics[
Table[{Black,Line[{{j,0},{j,g[j]}}]},{j,a,b,h}]],{}],
Axes->False,Frame->True,
PlotRange->{{a,b},All},
DisplayFunction->$DisplayFunction];
Print["Simpson's Rule with ", m ," intervals yields the result "];
Print[
(N[g[a]]+4N[g[b-h]]+N[g[b]]+N[Sum[4*g[j]+2*g[j+h],
{j,a+h,b-2h,2h}]])*h/3])]





End[ ]

EndPackage[ ]



(* 
Created: March 1996.
Author: Helmut Knaust.
e-Mail: helmut@math.utep.edu            
*)