# Lower and upper Riemann sums

I am trying to create a few images for educational purposes. Therefore I need to be able to illustrate the area under various functions, using lower and upper riemann sums (rectangles)

Reason for not including an MWE is that I do not know which package or tools are best suited for the task, wheter it is tikz, pgf-plots, or something else. Although prefferably not gnuplot, as it is not working for me.

So my goal is to type in the function, specify the start and end value, upper or lower sum, number of rectangles, and produce something like below.

From another question, I was able to discover how to make middle riemann sums. But I absolutely must have lower and upper =(

\documentclass[10pt,a4paper]{minimal}
\usepackage{pgfplots}

\pgfplotsset{
integral segments/.code={\pgfmathsetmacro\integralsegments{#1}},
integral segments=3,
integral/.style args={#1:#2}{
ybar interval,
domain=#1+((#2-#1)/\integralsegments)/2:#2+((#2-#1)/\integralsegments)/2,
samples=\integralsegments+1,
x filter/.code=\pgfmathparse{\pgfmathresult-((#2-#1)/\integralsegments)/2}
}
}

\begin{document}

\pagecolor{black}

\begin{center}
\begin{tikzpicture}[color=white,/pgf/declare function={f=3*e^(-x)*x^3;}]
\begin{axis}[
domain=0:8.1,
samples=100,
axis lines=middle
]
white,
integral segments=4,
integral=0:8
] {f};
\end{axis}
\end{tikzpicture}
\end{center}

\end{document}

-
It might not be very easy to do (unless the function is monotonous, which is a very special case, and including only this case in an introductory material might be quite misleading for students). –  mbork Mar 10 '12 at 14:52
Finding the maximum value of a function on an interval requires (as you know) calculus. And TeX is not a very good CAS. You could cheat by plotting many points in each interval and choosing the largest value. My only solution was not automatic (and not pgfplots-based): find the critical points by hand and plot left/right/interior points as necessary. –  Matthew Leingang Mar 10 '12 at 15:01
As much as I like the question I dislike OP choice of tools to accomplish the task. Can somebody provide a "real" solution in some kind of computer algebra system, Python, Asymptote or PostScript? –  Predrag Punosevac Mar 10 '12 at 15:58
See my answer here tikz-pgf-draw-integral-test-plot/ –  cmhughes Mar 10 '12 at 23:53
I've written a simple package, that works solely by including the pgfplots package, that draws a variety of Riemann sums including upper and lower sums of either fixed width rectangles of any number or of any number of varying width rectangles. It also has routines for combining them on one graph with appropriate order of overlap. See my answer to Placing Coloured Rectangles on a Plot Using Points From the Plot Riemann Sums –  Geoff Pointer Oct 31 '13 at 0:59

A solution with PSTricks. Run it with xelatex or latex->dvips->ps2pdf

\documentclass{article}
\begin{document}

\psset{plotpoints=200,algebraic}
\begin{pspicture}(-0.5,-2.5)(10,3)
\psStep[linecolor=magenta,StepType=upper,
fillstyle=solid,fillcolor=cyan!50](0,9){20}{sqrt(x)*sin(x)}
\psStep[linecolor=blue,fillstyle=solid,fillcolor=blue!30,
opacity=0.4](0,9){20}{sqrt(x)*sin(x)}
\psaxes[labelFontSize=\scriptstyle]{->}(0,0)(0,-2.25)(10,3)
\psplot[linewidth=1.5pt]{0}{10}{sqrt(x)*sin(x)}
\end{pspicture}

\psset{yunit=1.25cm}
\begin{pspicture}(-0.5,-1.75)(10,1.5)
\psaxes[labelFontSize=\scriptstyle]{->}(0,0)(0,-1.5)(10,1.5)
\psStep[StepType=Riemann,fillstyle=solid,
fillcolor=magenta!20](0,10){50}{sqrt(x)*cos(x)*sin(x)}
\psplot[linewidth=1.5pt]{0}{10}{sqrt(x)*cos(x)*sin(x)}
\end{pspicture}
\end{document}


-
Very neat, but it should be pointed out that this example doesn't show the upper and lower Riemann sums, but the left and right Riemann sums. Is there a way to use the maximum/minimum of each segment instead of the left or right value? –  Jake Mar 10 '12 at 22:04
It seems there's a bug in PStricks here. StepType=sup doesn't actually use the supremum over each step, but merely the larger of the left and right bounds of the step. Also, it seems that the terminology is non-standard: The upper sum is commonly understood to be the sum using the supremum over each step, so upper and sup should be equivalent, using the supremum over each step, and what is currently defined as upper should properly be called right (similar for lower, inf and left). –  Jake Mar 11 '12 at 8:52

I don't claim that this is a particularly elegant solution, but this is how I did a small illustration of Riemann sums and upper and lower Darboux sums:

\documentclass[10pt]{article}
\usepackage{mathtools,relsize,tikz}
\usetikzlibrary{calc,intersections,through} % I'm not 100% on which of these are needed
\begin{figure}[h]%\centering
% f(x) = x^3/20 - x^2/5 - x/5 + 3 from -3 to 3 (drawing -3.1 to 2.8).
% f(-3) = 0.45,  f(3) = 1.95
% f has max at approx. -0.4305 and f(-0.4305) ~= 3.0450.
% We're computing the integral from -2.8 to 2.5
% (x0=a,x1,x2,x3,x4=b) = (-2.8,     -1.3,      0.2,      1.3,      2.5)
% (x0*,x1*,x2*,x3*)    = (     -2.1,     -0.8,      0.7,      1.9)
\def\GraphAndAxis{%
\draw[->] (-3.0,0) -- (2.7,0) coordinate (x axis);
\draw[smooth,semithick,fill=none,name path=graph,domain=-2.9:2.6]%
plot[id=poly] (\x,0.05*\x*\x*\x - 0.2*\x*\x - 0.2*\x + 3);} % end GraphAndAxis
\begin{minipage}{0.5\linewidth}
\centering
\begin{tikzpicture}
\def\TagBox#1#2#3#4{% 1=x_i, 2=x*_i, 3=x_(i+1), 4=i
% \useasboundingbox (0,0) rectangle (0,0);
\path[name path=upward at tag] (#2,0.0) -- (#2,3.1);
% In general the meaning of (p) |- (q) is "the intersection of a vertical line through p
% and a horizontal line through q."
\path[name intersections={of=graph and upward at tag, by=fofxstar}];
\draw[fill=gray!40] (#3,0.0 |- fofxstar) rectangle (#1,0.0)
node[below,opacity=1.0] {\smaller\ensuremath{x^{\phantom{*}}_{#4}}};
\draw[thin,dashed] (#2,0.0) node[below] {\smaller\ensuremath{x^*_{#4}}}-- (fofxstar);
} % end TagBox
\GraphAndAxis %draw it so \TagBox can find intersection with it
\TagBox{-2.8}{-2.1}{-1.3}{0}
\TagBox{-1.3}{-0.8}{0.2}{1}
\TagBox{0.2}{0.7}{1.3}{2}
\TagBox{1.3}{1.9}{2.5}{3}
\GraphAndAxis
\begin{scope}
\path (2.5,1.0) -- (2.5,0.0) node [below] {\smaller\ensuremath{x^{\phantom{*}}_{4}}}; %ugly manual hack
\end{scope}
\end{tikzpicture}
\caption{Riemann sum.}\label{Fig::riemann sum}
\end{minipage}
\begin{minipage}{0.5\linewidth}
\centering
\begin{tikzpicture}
\def\UpperAndLowerBoxes#1#2#3#4#5{% #1=x_i, #2=x_(i+1),
% #3=x s.t. f(x) = min_x f(x),
% #4=x s.t. f(x) = max_x f(x),
% #5= index
\path[name path=upward at max x] (#4,0.0) -- (#4,3.1);
\path[name intersections={of=graph and upward at max x, by=maxf}];
\draw[fill=gray!30] (#1,0.0) rectangle (#2,0.0 |- maxf);
\path[name path=upward at min x] (#3,0.0) -- (#3,3.1);
\path[name intersections={of=graph and upward at min x, by=minf}];
\draw[fill=gray!50] (#2,0.0 |- minf) rectangle (#1,0.0)
node[below,opacity=1.0] {\smaller{\ensuremath{x^{\phantom{*}}_{#5}}}};
} % end UpperAndLowerBoxes
\GraphAndAxis
\UpperAndLowerBoxes{-2.8}{-1.3}{-2.8}{-1.3}{0}
\UpperAndLowerBoxes{-1.3}{0.2}{-1.3}{-0.4305}{1}
\UpperAndLowerBoxes{0.2}{1.3}{1.3}{0.2}{2}
\UpperAndLowerBoxes{1.3}{2.5}{2.5}{1.3}{3}
\GraphAndAxis
\begin{scope}
\path (2.5,1.0) -- (2.5,0.0) node [below] {\smaller\ensuremath{x^{\phantom{*}}_{4}}}; %ugly manual hack
\end{scope}
\end{tikzpicture}
\caption{Upper and lower Darboux sums.}\label{Fig::upper and lower sums}
\end{minipage}
\end{figure}


-
Thanks for adding the picture, @egreg –  kahen Mar 10 '12 at 17:23

Now it works as expected. To demonstrate the behavior I plotted two different functions.

How it works:

You have the following keys:

• integral segments

number of the equal parts which the domain of the function is divided into;

• integral samples

number of computation cycles for getting the maximum value in the given interval defined by integral segments

• integral min

draw the lower Riemann sum

• integral max

draw the upper Riemann sum

I tested the example with tikz 2.10 CSV. There is a known bug and I added a fix to the code. For more details have a look at this question. Why am I getting the PGF Math Error: Unknown function getargs'?

Another relevant answer was provided by @Christian Feuersänger: How do I use pgfmathdeclarefunction to create define a new pgf function?

\documentclass[10pt,a4paper]{article}
\usepackage{pgfplots}
\pgfplotsset{
integral segments/.code={\pgfmathsetmacro\integralsegments{#1}},
integral segments=3,
integral samples/.code={\edef\integralsamples{#1}},
integral samples = 10,
integral min/.style args={#1:#2}{
ybar interval,
domain=#1:#2,
samples=\integralsegments+1,
x filter/.code={
\edef\lastx{\pgfmathresult}
\pgfmathresult
},%
y filter/.code={%
\pgfmathparse{(#2/(\integralsegments))/\integralsamples}%
\edef\tempstep{\pgfmathresult}%
\pgfmathparse{f(\lastx)}%
\edef\tempa{\pgfmathresult}%
\edef\tempb{\pgfmathresult}%
\foreach \x in {0,1,...,\integralsamples}%
{%
\pgfmathparse{f(\lastx+\x*\tempstep)}%
\xdef\tempb{\tempb,\pgfmathresult}%
}%
\pgfmathmin{\tempb}{\tempb}
},
},
integral max/.style args={#1:#2}{
ybar interval,
domain=#1:#2,
samples=\integralsegments+1,
x filter/.code={
\edef\lastx{\pgfmathresult}
\pgfmathresult
},%
y filter/.code={%
\pgfmathparse{(#2/(\integralsegments))/\integralsamples}%
\edef\tempstep{\pgfmathresult}%
\pgfmathparse{f(\lastx)}%
\edef\tempa{\pgfmathresult}%
\edef\tempb{\pgfmathresult}%
\foreach \x in {0,1,...,\integralsamples}%
{%
\pgfmathparse{f(\lastx+\x*\tempstep)}%
\xdef\tempb{\tempb,\pgfmathresult}%
}%
\pgfmathmax{\tempb}{\tempb}
},
},
}
\makeatletter
%see http://tex.stackexchange.com/questions/9722/why-am-i-getting-the-pgf-math-error-unknown-function-getargs
\def\pgfmathmax#1#2{%
%   \pgfmathparse{getargs(#1,#2)}%
\pgfmathparse{#1,#2}%
\expandafter\pgfmathmax@\expandafter{\pgfmathresult}%
}
\def\pgfmathmin#1#2{%
\pgfmathparse{#1,#2}%
\expandafter\pgfmathmin@\expandafter{\pgfmathresult}%
}
%see http://tex.stackexchange.com/questions/15435/how-do-i-use-pgfmathdeclarefunction-to-create-define-a-new-pgf-function
\makeatother
\begin{document}
\begin{tikzpicture}
\pgfset{declare function={f(\x)=3*exp(-(\x))*(\x)^3+1;}}
\begin{axis}[
domain=0:8.1,
samples=100,
axis lines=middle
]

red,
integral segments=4,
integral min=0:8
] {f(x)};
\end{axis}
\end{tikzpicture}
\hfill
\begin{tikzpicture}
\pgfset{declare function={f(\x)=3*exp(-(\x))*(\x)^3+1;}}
\begin{axis}[
domain=0:8.1,
samples=100,
axis lines=middle
]

blue,
integral segments=5,
integral max=0:8
] {f(x)};
\end{axis}
\end{tikzpicture}

\begin{tikzpicture}
\pgfset{declare function={f(\x)=sin(2*deg(\x))*exp(0.1*\x)+2;}}
\begin{axis}[
domain=0:8.1,
samples=100,
axis lines=middle
]

red,
integral segments=4,
integral min=0:8
] {f(x)};
\end{axis}
\end{tikzpicture}
\hfill
\begin{tikzpicture}
\pgfset{declare function={f(\x)=sin(2*deg(\x))*exp(0.1*\x)+2;}}
\begin{axis}[
domain=0:8.1,
samples=100,
axis lines=middle
]

blue,
integral segments=5,
integral max=0:8
] {f(x)};
\end{axis}
\end{tikzpicture}

\end{document}


TEST CODE to handle functions with f(x)<0:

\documentclass[10pt,a4paper]{article}
\usepackage{pgfplots}
\pgfplotsset{
integral segments/.code={\pgfmathsetmacro\integralsegments{#1}},
integral segments=3,
integral samples/.code={\edef\integralsamples{#1}},
integral samples = 10,
integral min/.style args={#1:#2}{
ybar interval,
domain=#1:#2,
samples=\integralsegments+1,
x filter/.code={
\edef\lastx{\pgfmathresult}
\pgfmathresult
},%
y filter/.code={%
\pgfmathparse{(#2/(\integralsegments))/\integralsamples}%
\edef\tempstep{\pgfmathresult}%
\pgfmathparse{f(\lastx)}%
\edef\tempa{\pgfmathresult}%
\edef\tempb{\pgfmathresult}%
\foreach \x in {0,1,...,\integralsamples}%
{%
\pgfmathparse{f(\lastx+\x*\tempstep)}%
\xdef\tempb{\tempb,\pgfmathresult}%
}%
\pgfmathmin{\tempb}{\tempb}
\let\savepgfmathresult\pgfmathresult
\pgfmathgreater{\pgfmathresult}{0}
\ifdim\pgfmathresult pt> 0 pt \relax
\let\pgfmathresult\savepgfmathresult
\else
\pgfmathparse{f(\lastx)}%
\edef\tempa{\pgfmathresult}%
\edef\tempb{\pgfmathresult}%
\foreach \x in {0,1,...,\integralsamples}%
{%
\pgfmathparse{f(\lastx+\x*\tempstep)}%
\xdef\tempb{\tempb,\pgfmathresult}%
}%
\pgfmathmax{\tempb}{\tempb}
\fi
},
},
integral max/.style args={#1:#2}{
ybar interval,
domain=#1:#2,
samples=\integralsegments+1,
x filter/.code={
\edef\lastx{\pgfmathresult}
\pgfmathresult
},%
y filter/.code={%
\pgfmathparse{(#2/(\integralsegments))/\integralsamples}%
\edef\tempstep{\pgfmathresult}%
\pgfmathparse{f(\lastx)}%
\edef\tempa{\pgfmathresult}%
\edef\tempb{\pgfmathresult}%
\foreach \x in {0,1,...,\integralsamples}%
{%
\pgfmathparse{f(\lastx+\x*\tempstep)}%
\xdef\tempb{\tempb,\pgfmathresult}%
}%
\pgfmathmax{\tempb}{\tempb,0}
\let\savepgfmathresult\pgfmathresult
%                               \pgfmathparse{ifthenelse(\pgfmathresult>=0,1,0)}
\pgfmathgreater{\pgfmathresult}{0}
\ifdim\pgfmathresult pt> 0 pt \relax
\let\pgfmathresult\savepgfmathresult
\else
\pgfmathparse{f(\lastx)}%
\edef\tempa{\pgfmathresult}%
\edef\tempb{\pgfmathresult}%
\foreach \x in {0,1,...,\integralsamples}%
{%
\pgfmathparse{f(\lastx+\x*\tempstep)}%
\xdef\tempb{\tempb,\pgfmathresult}%
}%
\pgfmathmin{\tempb}{\tempb}
\fi
},
},
}
\makeatletter
%see http://tex.stackexchange.com/questions/9722/why-am-i-getting-the-pgf-math-error-unknown-function-getargs
\def\pgfmathmax#1#2{%
%   \pgfmathparse{getargs(#1,#2)}%
\pgfmathparse{#1,#2}%
\expandafter\pgfmathmax@\expandafter{\pgfmathresult}%
}
\def\pgfmathmin#1#2{%
\pgfmathparse{#1,#2}%
\expandafter\pgfmathmin@\expandafter{\pgfmathresult}%
}
%see http://tex.stackexchange.com/questions/15435/how-do-i-use-pgfmathdeclarefunction-to-create-define-a-new-pgf-function
\makeatother
\begin{document}
\begin{tikzpicture}
\pgfset{declare function={f(\x)=3*exp(-(\x))*(\x)^3+1;}}
\begin{axis}[
domain=0:8.1,
samples=100,
axis lines=middle
]

red,
integral segments=4,
integral min=0:8
] {f(x)};
\end{axis}
\end{tikzpicture}
\hfill
\begin{tikzpicture}
\pgfset{declare function={f(\x)=3*exp(-(\x))*(\x)^3+1;}}
\begin{axis}[
domain=0:8.1,
samples=100,
axis lines=middle
]

blue,
integral segments=5,
integral max=0:8
] {f(x)};
\end{axis}
\end{tikzpicture}

\begin{tikzpicture}
\pgfset{declare function={f(\x)=sin(2*deg(\x))*exp(0.1*\x)+2;}}
\begin{axis}[
domain=0:8.1,
samples=100,
axis lines=middle
]

red,
integral segments=4,
integral min=0:8
] {f(x)};
\end{axis}
\end{tikzpicture}
\hfill
\begin{tikzpicture}
\pgfset{declare function={f(\x)=sin(2*deg(\x))*exp(0.1*\x)+2;}}
\begin{axis}[
domain=0:8.1,
samples=100,
axis lines=middle
]

blue,
integral segments=5,
integral max=0:8,
] {f(x)};
\end{axis}
\end{tikzpicture}

\begin{tikzpicture}
\pgfset{declare function={f(\x)=sqrt(\x)*cos(deg(\x))*sin(deg(\x));}}
\begin{axis}[
domain=0:8.1,
samples=100,
axis lines=middle
]

red,
integral segments=30,
integral min=0:8
] {f(x)};
\end{axis}
\end{tikzpicture}
\hfill
\begin{tikzpicture}
\pgfset{declare function={f(\x)=sqrt(\x)*cos(deg(\x))*sin(deg(\x));}}
\begin{axis}[
domain=0:8.1,
samples=100,
axis lines=middle
]

blue,
integral segments=15,
integral max=0:8
] {1};
\end{axis}
\end{tikzpicture}

\end{document}

-

This seems to work but sort of hackish.

\documentclass{article}
\usepackage{pgfplots}
\begin{document}
\begin{tikzpicture}
\begin{axis}
coordinates
{(0,0.1) (0.1,0.15) (0.2,0.5) (0.3,0.62)
(0.4,0.56) (0.5,0.58) (0.6,0.65) (0.7,0.6)
(0.8,0.58) (0.9,0.55) (1,0.52) }\closedcycle;

coordinates
{(0,0.1) (0.1,0.15) (0.2,0.5) (0.3,0.62)
(0.4,0.56) (0.5,0.58) (0.6,0.65) (0.7,0.6)
(0.8,0.58) (0.9,0.55) (1,0.52) }\closedcycle;

coordinates
{(0,0.1) (0.1,0.15) (0.2,0.5) (0.3,0.62)
(0.4,0.56) (0.5,0.58) (0.6,0.65) (0.7,0.6)
(0.8,0.58) (0.9,0.55) (1,0.52)};
\end{axis}
\end{tikzpicture}
\end{document}


See the manual about const plot mark <right/left/mid>

-

As is mentioned above, TeX is not a good CAS, so I would suggest the following (you'll need 3 tabs open in your browser): First, staying on this page, create another tab, go here and find "Numerical integrals with various rules". Highlight the code from that section and copy it (Control-C). Second, in the third tab go here. Paste the code in over the existing code (Control-V). Press 'Evaluate' on the screen to make sure the code is working properly (i.e. you've copied everything properly [spacing/indentation is important). You should now have a manipulative that does what you've asked for, except there's no PDF output. We'll change that. On the Sage Cell Server page, go to the line of code that says, " show(plot(func,a,b) + rects, xmin = a, xmax = b, ymin = min_y, ymax = max_y)". That 1 line will be replaced by 3 lines.

C= plot(func,a,b) + rects
C.show(xmin = a, xmax = b, ymin = min_y, ymax = max_y)
C.save("MyPicture.pdf")


Again, spacing/indentation is critical. Press 'Evaluate' and you'll get output like what is shown below. I couldn't get everything on the screen, but under the display (below 'Session Notes' there is a link to "MyPicture.pdf.". Right click on it and save it to your computer. Note, it only shows the graph, not the HTML Riemann sum calculation or sliders that you see in the screenshot.

It shouldn't take more than several minutes to run through all the steps. Here is a screenshot which includes the altered code section and image.

-