# For CalculusWithJulia, is there Riemann' Midpoint method?

**URL:** <https://discourse.julialang.org/t/for-calculuswithjulia-is-there-riemann-midpoint-method/91001>\
**Category:** General Usage\
**Tags:** package\
**Created:** [November 29, 2022, 3:08pm UTC](https://discourse.julialang.org/t/for-calculuswithjulia-is-there-riemann-midpoint-method/91001 "2022-11-29T15:08:02Z")\
**Posts on this page:** 11\
**Page:** 1

<div class="post-metadata">

**Author:** ![Freya\_the\_Goddess](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/freya_the_goddess/32/36835_2.png) [@Freya\_the\_Goddess](https://discourse.julialang.org/u/Freya_the_Goddess)\
**Post date:** [November 29, 2022, 3:08pm UTC](https://discourse.julialang.org/t/for-calculuswithjulia-is-there-riemann-midpoint-method/91001/1 "2022-11-29T15:08:03Z")

</div>

I read this:

[https://docs.juliahub.com/CalculusWithJulia/AZHbv/0.0.5/integrals/area.html](https://docs.juliahub.com/CalculusWithJulia/AZHbv/0.0.5/integrals/area.html)

and I want to calculate Riemann sum for several functions, comparing the

1. Left
2. Right
3. Midpoint

methods. Is the midpoint method already available? If not then how to add it to the current Riemann Sum function?

---

<div class="post-metadata">

**Author:** ![j\_verzani](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/j_verzani/32/8551_2.png) [@j\_verzani](https://discourse.julialang.org/u/j_verzani)\
**Post date:** [November 29, 2022, 6:50pm UTC](https://discourse.julialang.org/t/for-calculuswithjulia-is-there-riemann-midpoint-method/91001/2 "2022-11-29T18:50:53Z")

</div>

No, that function is for illustrative purposes only and doesn’t allow modification of the provided methods. You could certainly suggest some modification to the package, but for now, you might follow this simple enough pattern:

```julia
function Riemann(f, a, b, n; meth=(a,b) -> (a+b)/2)
    xs = range(a, b, n+1)
    pairs = zip(xs[begin:end-1], xs[begin+1:end]) # (x₀,x₁), …, (xₙ₋₁,xₙ)
    sum(f(meth(xs...)) for xs in pairs)
end

```

Pass in `meth=(a,b) -> a` for left sums, `meth=(a,b)->b` for right…

BTW, that link is out of date. Try [https://calculuswithjulia.github.io](https://calculuswithjulia.github.io) for a more current version.

---

<div class="post-metadata">

**Author:** ![Freya\_the\_Goddess](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/freya_the_goddess/32/36835_2.png) [@Freya\_the\_Goddess](https://discourse.julialang.org/u/Freya_the_Goddess)\
**Post date:** [December 17, 2022, 8:03am UTC](https://discourse.julialang.org/t/for-calculuswithjulia-is-there-riemann-midpoint-method/91001/3 "2022-12-17T08:03:50Z")

</div>

Hi @j_verzani ,

I was using this function code:

```julia
function riemann(f, a, b, n; method="left")
  xs = a:(b-a)/n:b
  deltas = diff(xs) # forms x2-x1, x3-x2, ..., xn-xn-1
  if method == "left"
    cs = xs[1:end-1]
  elseif method == "right"
    cs = xs[2:end]
  elseif method == "midpoint"
    cs = [(xs[i] + xs[i+1])/2 for i in 1:length(deltas)]
  end
  sum(f(cs[i]) * deltas[i] for i in 1:length(deltas))
end

```

and it can calculate the midpoint method correctly.

I think it is another version of your suggested function:

```julia
function Riemann(f, a, b, n; meth=(a,b) -> (a+b)/2)
    xs = range(a, b, n+1)
    pairs = zip(xs[begin:end-1], xs[begin+1:end]) # (x₀,x₁), …, (xₙ₋₁,xₙ)
    sum(f(meth(xs...)) for xs in pairs)
end

```

How to check which function works better? with btime?

---

<div class="post-metadata">

**Author:** ![Freya\_the\_Goddess](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/freya_the_goddess/32/36835_2.png) [@Freya\_the\_Goddess](https://discourse.julialang.org/u/Freya_the_Goddess)\
**Post date:** [December 17, 2022, 8:20am UTC](https://discourse.julialang.org/t/for-calculuswithjulia-is-there-riemann-midpoint-method/91001/4 "2022-12-17T08:20:07Z")

</div>

If all the materials in ([Calculus with Julia](https://jverzani.github.io/CalculusWithJuliaNotes.jl)) printed to become a book I will look forward to it, it is really helping me to study better than in classes.

---

<div class="post-metadata">

**Author:** ![Freya\_the\_Goddess](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/freya_the_goddess/32/36835_2.png) [@Freya\_the\_Goddess](https://discourse.julialang.org/u/Freya_the_Goddess)\
**Post date:** [December 18, 2022, 10:20am UTC](https://discourse.julialang.org/t/for-calculuswithjulia-is-there-riemann-midpoint-method/91001/5 "2022-12-18T10:20:22Z")

</div>

I have been looking at this function from the old blog of yours with midpoint method:

```julia
function riemann(f, a, b, n; method="left")
  xs = a:(b-a)/n:b
  deltas = diff(xs) # forms x2-x1, x3-x2, ..., xn-xn-1
  if method == "left"
    cs = xs[1:end-1]
  elseif method == "right"
    cs = xs[2:end]
  elseif method == "midpoint"
    cs = [(xs[i] + xs[i+1])/2 for i in 1:length(deltas)]
  end
  sum(f(cs[i]) * deltas[i] for i in 1:length(deltas))
end

```

I wonder how to convert the `midpoint method` into the new function in the updated blog, so all methods will be available:

```julia
function riemann(f::Function, a::Real, b::Real, n::Int; method="right")
    if method == "right"
        meth = f -> (lr -> begin l,r = lr; f(r) * (r-l) end)
    elseif method == "left"
        meth = f -> (lr -> begin l,r = lr; f(l) * (r-l) end)
    elseif method == "trapezoid"
        meth = f -> (lr -> begin l,r = lr; (1/2) * (f(l) + f(r)) * (r-l) end)
    elseif method == "simpsons"
        meth = f -> (lr -> begin l,r = lr; (1/6) * (f(l) + 4*(f((l+r)/2)) + f(r)) * (r-l) end)
    end

    xs = range(a, b, n+1)
    pairs = zip(xs[begin:end-1], xs[begin+1:end]) # (x₀,x₁), …, (xₙ₋₁,xₙ)
    sum(meth(f), pairs)

end

```

---

<div class="post-metadata">

**Author:** ![Freya\_the\_Goddess](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/freya_the_goddess/32/36835_2.png) [@Freya\_the\_Goddess](https://discourse.julialang.org/u/Freya_the_Goddess)\
**Post date:** [December 18, 2022, 10:24am UTC](https://discourse.julialang.org/t/for-calculuswithjulia-is-there-riemann-midpoint-method/91001/6 "2022-12-18T10:24:07Z")

</div>

My inexperience resulting in this to make all methods available:

```julia
function riemann(f::Function, a::Real, b::Real, n::Int; method="right")
    if method == "right"
        meth = f -> (lr -> begin l,r = lr; f(r) * (r-l) end)
	xs = range(a, b, n+1)
        pairs = zip(xs[begin:end-1], xs[begin+1:end]) # (x₀,x₁), …, (xₙ₋₁,xₙ)
        sum(meth(f), pairs)
    elseif method == "left"
        meth = f -> (lr -> begin l,r = lr; f(l) * (r-l) end)
	xs = range(a, b, n+1)
        pairs = zip(xs[begin:end-1], xs[begin+1:end]) # (x₀,x₁), …, (xₙ₋₁,xₙ)
        sum(meth(f), pairs)
    elseif method == "trapezoid"
        meth = f -> (lr -> begin l,r = lr; (1/2) * (f(l) + f(r)) * (r-l) end)
	xs = range(a, b, n+1)
        pairs = zip(xs[begin:end-1], xs[begin+1:end]) # (x₀,x₁), …, (xₙ₋₁,xₙ)
        sum(meth(f), pairs)
    elseif method == "simpsons"
        meth = f -> (lr -> begin l,r = lr; (1/6) * (f(l) + 4*(f((l+r)/2)) + f(r)) * (r-l) end)
	xs = range(a, b, n+1)
        pairs = zip(xs[begin:end-1], xs[begin+1:end]) # (x₀,x₁), …, (xₙ₋₁,xₙ)
        sum(meth(f), pairs)
    elseif method == "midpoint"
	xs = a:(b-a)/n:b
        deltas = diff(xs) # forms x2-x1, x3-x2, ..., xn-xn-1
        cs = [(xs[i] + xs[i+1])/2 for i in 1:length(deltas)]
        sum(f(cs[i]) * deltas[i] for i in 1:length(deltas))
    end
end

```

but it works!!!

![ch1-calculus-5thedefiniteintegral-59](https://global.discourse-cdn.com/julialang/original/3X/b/4/b485177b7240da0291106bce21a2ba8e469cf54c.png)

---

<div class="post-metadata">

**Author:** ![j\_verzani](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/j_verzani/32/8551_2.png) [@j\_verzani](https://discourse.julialang.org/u/j_verzani)\
**Post date:** [December 18, 2022, 11:04pm UTC](https://discourse.julialang.org/t/for-calculuswithjulia-is-there-riemann-midpoint-method/91001/7 "2022-12-18T23:04:37Z")

</div>

Maybe this is a bit more clear with how you can illustrate differences with various methods accompanying Riemann sums if performance isn’t the point:

```julia
function riemann(M, f, a, b, n)
    xs = range(a,b,n)
    riemann(M, f, xs[1:end-1], xs[2:end])
end
riemann(M, f, as, bs) = sum(M(f,a,b) * (b-a) for (a,b) in zip(as, bs))

Ms = (left = (f,a,b) -> f(a),
      right= (f,a,b) -> f(b),
      mid = (f,a,b) -> f(a/2 + b/2),
      m̃ = (f,a,b) -> first(findmin(f, range(a,b,10))),
      M̃ = (f,a,b) -> first(findmax(f, range(a,b,10))),
      trapezoid = (f,a,b) -> (f(a) + f(b))/2,
      simpsons = (f,a,b) -> (c = a/2 + b/2;(1/6) * (f(a) + 4*f(c) + f(b)))
      )

ns = (5, 50, 500)

using DataFrames
d = DataFrame(n = Int[], M = Symbol[], value = Float64[])
for (m,M) in pairs(Ms)
    for n in ns
        v = riemann(M, sin, 0, pi, n)
        push!(d, (n=n, M=m, value=v))
    end
end
  
unstack(d, :n, :value)                                                  

```

Yielding with this f, a,b:

```julia
7×4 DataFrame
 Row │ M 5 50 500      
     │ Symbol Float64? Float64? Float64? 
─────┼─────────────────────────────────────────
   1 │ left 1.89612 1.99931 1.99999
   2 │ right 1.89612 1.99931 1.99999
   3 │ mid 2.05234 2.00034 2.0
   4 │ m̃ 1.11072 1.93523 1.9937
   5 │ M̃ 2.68152 2.06343 2.00629
   6 │ trapezoid 1.89612 1.99931 1.99999
   7 │ simpsons 2.00027 2.0 2.0

```

---

<div class="post-metadata">

**Author:** ![Freya\_the\_Goddess](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/freya_the_goddess/32/36835_2.png) [@Freya\_the\_Goddess](https://discourse.julialang.org/u/Freya_the_Goddess)\
**Post date:** [December 19, 2022, 4:17am UTC](https://discourse.julialang.org/t/for-calculuswithjulia-is-there-riemann-midpoint-method/91001/8 "2022-12-19T04:17:56Z")

</div>

I try to change ns = (5, 9, 17)

the ns actually represents the number of partitions + 1

when I see the solution manual for function 1/x it matches with the solution:

 ![Capture d’écran_2022-12-19_11-16-06](https://global.discourse-cdn.com/julialang/original/3X/b/3/b33cd8e4d4c079367534b3218af5f089b8757da4.png)

---

<div class="post-metadata">

**Author:** ![Freya\_the\_Goddess](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/freya_the_goddess/32/36835_2.png) [@Freya\_the\_Goddess](https://discourse.julialang.org/u/Freya_the_Goddess)\
**Post date:** [December 19, 2022, 11:22am UTC](https://discourse.julialang.org/t/for-calculuswithjulia-is-there-riemann-midpoint-method/91001/9 "2022-12-19T11:22:58Z")

</div>

I have tried to add a new method in the Ms, basically I want to follow this formula:

T - \frac{(f'(b) - f'(a))h^{2}}{12}

with T is the result from the Trapezoidal, and h = \frac{b-a}{n}, n is the number of partition, and the interval is [a,b]

it seems the code still not produces the correct calculation, I shall try again:

```julia
...
modifiedtrapezoid = (f,a,b) -> (f(a) + f(b))/2 - (((f'(b)-f'(a))*h^2)/12),
...

```

Edited, this code works:

```julia
using CalculusWithJulia

# For partition of 8
ns = (9)

# Define the function, the partitions, and the interval 
f(x) = x^4
a = 1
b = 3
n = 8
h = (b-a)/n
g = (((f'(b)-f'(a))*h^2)/12)
function riemann(M, f, a, b, n)
    xs = range(a,b,n)
    riemann(M, f, xs[1:end-1], xs[2:end])
end
riemann(M, f, as, bs) = sum(M(f,a,b) * (b-a) for (a,b) in zip(as, bs))

Ms = (left = (f,a,b) -> f(a),
      right= (f,a,b) -> f(b),
      midpoint = (f,a,b) -> f(a/2 + b/2),
      #m̃ = (f,a,b) -> first(findmin(f, range(a,b,10))),
      #M̃ = (f,a,b) -> first(findmax(f, range(a,b,10))),
      trapezoid = (f,a,b) -> (f(a) + f(b))/2,
      modifiedtrapezoid = (f,a,b) -> (f(a) + f(b))/2 - g/2,
      simpsons = (f,a,b) -> (c = a/2 + b/2;(1/6) * (f(a) + 4*f(c) + f(b)))
      )

using DataFrames
d = DataFrame(n = Int[], M = Symbol[], value = Float64[])
for (m,M) in pairs(Ms)
    for n in ns
        v = riemann(M, f, a, b, n)
        push!(d, (n=n-1, M=m, value=v))
    end
end
  
unstack(d, :n, :value)

```

---

<div class="post-metadata">

**Author:** ![j\_verzani](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/j_verzani/32/8551_2.png) [@j\_verzani](https://discourse.julialang.org/u/j_verzani)\
**Post date:** [December 19, 2022, 1:16pm UTC](https://discourse.julialang.org/t/for-calculuswithjulia-is-there-riemann-midpoint-method/91001/10 "2022-12-19T13:16:50Z")

</div>

> [@Freya\_the\_Goddess](#):
>
> ```julia
> modifiedtrapezoid = (f,a,b) -> (f(a) + f(b))/2 - (((f'(b)-f'(a))*h^2)/12),
> 
> ```

The formula used to add the terms us: `sum(M(f,a,b) * (b-a) for (a,b) in zip(as, bs))` has a factor `h` in the `(b-a)` term, so you only have one `h` in the formula, not `h^2`. Try

```julia
ct = (f,a,b)-> (c = a/2+b/2; (f(a)+f(b))/2 + 1/12 * (f'(b)-f'(a))*(b-a))

```

---

<div class="post-metadata">

**Author:** ![Freya\_the\_Goddess](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/freya_the_goddess/32/36835_2.png) [@Freya\_the\_Goddess](https://discourse.julialang.org/u/Freya_the_Goddess)\
**Post date:** [December 19, 2022, 1:47pm UTC](https://discourse.julialang.org/t/for-calculuswithjulia-is-there-riemann-midpoint-method/91001/11 "2022-12-19T13:47:34Z")

</div>

Alright, thanks !