StopReadingThisUsername - 1 year ago 101

R Question

`integral <- function(x) {2.393794315*((3320)*(x/2.24581)^5+(-1613880/11)*(x/2.24581)^4+(171163181/66)*(x/2.24581)^3+(-7563546913/330)*(x/2.24581)^2+(835541173981/8250)*(x/2.24581)+(-2953570085669/16500))*(((483793.161846485)*x^8+(-76823340.9717028)*x^7+(5337025908.822)*x^6+(-211866341077.587)*x^5+(5256530719898.47)*x^4+(-83466263852549.1)*x^3+(828318375700455)*x^2+(-4697211251008830)*x+(11653475160809900))^0.5)}`

integrate(integral, lower=19.538547, upper=20.3245805)

This gave me

Error in integrate(integral, lower = 19.538547, upper = 20.3245805) :

non-finite function value

I am not sure what to do. It'd be great if somebody could input the integral into a software like Maple to see if it works - or, giving an answer to how to get around this error would be great, too :P

Thanks in advance!

EDIT:

The function that I am trying to revolve around the x-axis is:

y=2.393794315*((3320)

When plotted with limits {19.538547 < x < 20.3245805}:

As in the link provided by @StephenWade (in the comments), to find the surface area, I needed (f'(x))^2 + 1. Using Wolfram, Alpha, I computed this to be:

(f'(x))^2 + 1 = (483793.161846485)*x^8+(-76823340.9717028)*x^7+(5337025908.822)*x^6+(-211866341077.587)*x^5+(5256530719898.47)*x^4+(-83466263852549.1)*x^3+(828318375700455)*x^2+(-4697211251008830)*x+(11653475160809900)

Plotting this, without limits, gives:

However, computing this on R, I got an error:

My input into R was the following:

2.393794315*((3320)

I just need a reasonable estimate of the surface area...Any ideas?

Answer Source

I would approach it using the `polynom`

package, which can calculate the derivative for us, and then apply the formula from Wolfram.

```
revol_coef <- 2.393794315 * c(-2953570085669 / 16500,
835541173981 / (8250 * 2.24581),
-7563546913/ (330 * 2.24581^2),
171163181 / (66 * 2.24581^3),
-1613880 / (11 * 2.24581^4),
3320 / 2.24581^5)
y <- polynomial(revol_coef)
y_d <- deriv(y)
f <- function(x) {
2 * pi * predict(y, x) * sqrt(1 + predict(y_d, x)^2)
}
integrate(f, lower = 19.538547, upper = 20.3245805)
```

The output I get is

```
-45.71118 with absolute error < 0.0017
```

The answer is negative, due to the surface being specified in the lower-right quadrant of the x-y plane (as per your plot). To fix this, take the negative of the number supplied to get a surface area.

Translating formulas from Wolfram Alpha, manipulating them, and then putting them into R can easily involve making a mistake.

In this case, I would recommend doing as much of the work within R, as the tools/functions to do the necessary calculations and simplifications are readily available.