# QuadGK Computation Question

**URL:** <https://discourse.julialang.org/t/quadgk-computation-question/92983>\
**Category:** General Usage\
**Tags:** question, package, quadgk\
**Created:** [January 15, 2023, 6:17am UTC](https://discourse.julialang.org/t/quadgk-computation-question/92983 "2023-01-15T06:17:41Z")\
**Posts on this page:** 9\
**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:** [January 15, 2023, 6:17am UTC](https://discourse.julialang.org/t/quadgk-computation-question/92983/1 "2023-01-15T06:17:41Z")

</div>

Hi all,

I have tried `QuadGK`, but the result is not the same as the solution manual.

Is there a problem in between?

this is my code:

```julia
using SymPy, QuadGK

x = symbols("x")

f(x) = (x^6 + 2)/(8x^2)
# simplify(diff(f(x),x)
g(x) = sqrt(1 + (x^6 - 2)/(2x^3))

h(x) = f(x)*g(x)
Area(x) = 2pi*h(x)

d = quadgk(x -> Area(x), 1, 3, rtol=1e-10)
println("c = ", d)

```

 ![Capture d’écran_2023-01-15_13-16-16](https://global.discourse-cdn.com/julialang/original/3X/1/2/123dfe8210fabc7d754cf8ce73ec2bc104651545.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:** [January 15, 2023, 6:19am UTC](https://discourse.julialang.org/t/quadgk-computation-question/92983/2 "2023-01-15T06:19:43Z")

</div>

No its okay,

I forgot to power the derivative. Now the problem is solved.

But the result still differ a lot!

**325.4437377507612 from QuadGK**

---

<div class="post-metadata">

**Author:** ![PeterSimon](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/petersimon/32/25193_2.png) [@PeterSimon](https://discourse.julialang.org/u/PeterSimon)\
**Post date:** [January 16, 2023, 3:56am UTC](https://discourse.julialang.org/t/quadgk-computation-question/92983/3 "2023-01-16T03:56:44Z")

</div>

If you use

```julia
g(x) = sqrt(1 + ((x^6 - 1)/(2x^3))^2)

```

the result is

```julia
c = (326.9195614457822, 1.2235672031124523e-8)

```

which compares very favorably to the exact answer

```julia
julia> 8429π/81
326.9195614457823

```

---

<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:** [January 16, 2023, 7:26am UTC](https://discourse.julialang.org/t/quadgk-computation-question/92983/4 "2023-01-16T07:26:45Z")

</div>

This is lovely:

```julia
using SymPy, QuadGK

x = symbols("x")

f(x) = (x^6 + 2)/(8x^2)
# simplify(diff(f(x),x))
# we cannot use diff inside the fd(x) since QuadGK can't comprehend diff
g(x) = sqrt(1 + ((x^6 - 1)/(2x^3))^2)
#g(x) = f(x)*√(1+(fd(x))^2)

h(x) = f(x)*g(x)
Area(x) = 2pi*h(x)

d = quadgk(x -> Area(x), 1, 3, rtol=1e-10)
dq = 2π*quadgk(h, 1, 3)[1]

println("(Area, error) = ", d)

println("Area with QuadGK = ", dq)

```

**(Area, error) = (326.9195614457822, 1.2235672031124523e-8)**  
**Area with QuadGK = 326.91956144578495**

---

<div class="post-metadata">

**Author:** ![DNF](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dnf/32/10191_2.png) [@DNF](https://discourse.julialang.org/u/DNF)\
**Post date:** [January 16, 2023, 8:49am UTC](https://discourse.julialang.org/t/quadgk-computation-question/92983/5 "2023-01-16T08:49:58Z")

</div>

> [@Freya\_the\_Goddess](#):
>
> ```julia
> x = symbols("x")
> 
> f(x) = (x^6 + 2)/(8x^2)
> # simplify(diff(f(x),x))
> # we cannot use diff inside the fd(x) since QuadGK can't comprehend diff
> g(x) = sqrt(1 + ((x^6 - 1)/(2x^3))^2)
> #g(x) = f(x)*√(1+(fd(x))^2)
> 
> ```

I’m not sure if you are aware, but the `x = symbols("x")` doesn’t have any effect here. The `x` inside the definitions of `f(x)`, `g(x)`, etc. are completely unrelated to the symbolic `x` that you create first. Its presence is a bit confusing and it doesn’t seem to belong in a Minimal Working Example.

One thing you could do, is to use automatic differentiation:

```julia
using ForwardDiff: derivative

f(x) = (x^6 + 2)/(8x^2)
A(x) = 2π * f(x) * sqrt(1 + derivative(f, x)^2)

julia> quadgk(A, 1, 3; rtol=1e-10)
(326.9195614457822, 1.2235672031124523e-8)

```

The autodiff version is only _slightly_ slower:

```julia
julia> @btime quadgk(Area, 1, 3, rtol=1e-10)
  897.619 ns (2 allocations: 368 bytes)
(326.9195614457822, 1.2235672031124523e-8)

julia> @btime quadgk(A, 1, 3, rtol=1e-10)
  1.150 μs (2 allocations: 368 bytes)
(326.9195614457822, 1.2235672031124523e-8)

```

(It’s quite incredible how little performance overhead there is, actually.)

BTW, there’s no need to write `x -> Area(x)`, this is just the same thing as `Area`. So in your code you would write: `quadgk(Area, 1, 3; ...)`

---

<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:** [January 16, 2023, 9:05am UTC](https://discourse.julialang.org/t/quadgk-computation-question/92983/6 "2023-01-16T09:05:01Z")

</div>

Well, first I am not aware, I thought using

```julia
x = symbols("x")

```

is needed to calculate the integral.

I never use Automatic Differentiation or try `ForwardDiff`, I will test and try it now.

---

<div class="post-metadata">

**Author:** ![DNF](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dnf/32/10191_2.png) [@DNF](https://discourse.julialang.org/u/DNF)\
**Post date:** [January 16, 2023, 9:07am UTC](https://discourse.julialang.org/t/quadgk-computation-question/92983/7 "2023-01-16T09:07:27Z")

</div>

> [@DNF](#):
>
> The autodiff version is only _slightly_ slower:

In fact, if you re-write your `f(x)` to the equivalent form (taken from your posted equations)

```julia
f(x) = x^4/8 + 1/4x^2 # this is the same
A(x) = 2π * f(x) * sqrt(1 + derivative(f, x)^2)

```

then it’s _even_ faster:

```julia
julia> @btime quadgk(A, 1, 3; rtol=1e-10)
  696.528 ns (2 allocations: 368 bytes)
(326.91956144578234, 1.2235728874543383e-8)

```

The answer is slightly different from before, but still _as_ accurate as the previous version.

---

<div class="post-metadata">

**Author:** ![DNF](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dnf/32/10191_2.png) [@DNF](https://discourse.julialang.org/u/DNF)\
**Post date:** [January 16, 2023, 9:09am UTC](https://discourse.julialang.org/t/quadgk-computation-question/92983/8 "2023-01-16T09:09:22Z")

</div>

> [@Freya\_the\_Goddess](#):
>
> I thought using
> 
> ```julia
> x = symbols("x")
> 
> ```
> 
> is needed to calculate the integral.

No, it does nothing here, and is completely unrelated to quadgk, which does purely numerical integration.

---

<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:** [January 16, 2023, 12:01pm UTC](https://discourse.julialang.org/t/quadgk-computation-question/92983/9 "2023-01-16T12:01:27Z")

</div>

It takes time for me to read and comprehend all your responses, thanks anyway.