R + plotly: solid of revolution

Question

I have a function r(x) that I want to rotate around the x axis to get a solid of revolution that I want to add to an existing plot_ly plot using add_surface (colored by x).

Here is an example:

library(dplyr)
library(plotly)

# radius depends on x
r <- function(x) x^2

# interval of interest
int <- c(1, 3)

# number of points along the x-axis
nx <- 20

# number of points along the rotation
ntheta <- 36

# set x points and get corresponding radii
coords <- data_frame(x = seq(int[1], int[2], length.out = nx), r = r(x))

# for each x: rotate r to get y and z coordinates
# edit: ensure 0 and pi are both amongst the angles used
coords %<>%
  rowwise() %>%
  do(data_frame(x = .$x, r = .$r,
                theta = seq(0, pi, length.out = ntheta / 2 + 1) %>%
                c(pi + .[-c(1, length(.))]))) %>%

  ungroup %>%
  mutate(y = r * cos(theta), z = r * sin(theta))

# plot points to make sure the coordinates define the desired shape
coords %>%
  plot_ly(x = ~x, y = ~y, z = ~z, color = ~x) %>%
  add_markers()

How can I generate the shape indicated by the points above as a plotly surface (ideally open on both ends)?

---

edit (1):

Here is my best attempt so far:

# get all x & y values used (sort to connect halves on the side)
xs <-
  unique(coords$x) %>%
  sort
ys <-
  unique(coords$y) %>%
  sort

# for each possible x/y pair: get z^2 value
coords <-
  expand.grid(x = xs, y = ys) %>%
  as_data_frame %>%
  mutate(r = r(x), z2 = r^2 - y^2)

# format z coordinates above x/y plane as matrix where columns
# represent x and rows y
zs <- matrix(sqrt(coords$z2), ncol = length(xs), byrow = TRUE)

# format x coordiantes as matrix as above (for color gradient)
gradient <-
  rep(xs, length(ys)) %>%
  matrix(ncol = length(xs), byrow = TRUE)
  
# plot upper half of shape as surface
p <- plot_ly(x = xs, y = ys, z = zs, surfacecolor = gradient,
             type = "surface", colorbar = list(title = 'x'))

# plot lower have of shape as second surface
p %>%
  add_surface(z = -zs, showscale = FALSE)

While this gives the desired shape,

1. It has 'razor teeth' close to the x/y plane.
2. The halves parts don't touch. (resolved by including 0 and pi in the theta vectors)
3. I didn't figure out how to color it by x instead of z (though I didn't look much into this so far). (resolved by gradient matrix)

---

edit (2):

Here is an attempt using a single surface:

# close circle in y-direction
ys <- c(ys, rev(ys), ys[1])

# get corresponding z-values
zs <- rbind(zs, -zs[nrow(zs):1, ], zs[1, ])

# as above, but for color gradient
gradient <-
  rbind(gradient, gradient[nrow(gradient):1, ], gradient[1, ])

# plot single surface
plot_ly(x = xs, y = ys, z = zs, surfacecolor = gradient,
        type = "surface", colorbar = list(title = 'x'))

Surprisingly, while this should connect the two halves orthogonal to the x / y plane to create the full shape, it still suffers from the same 'razor teeth' effect as the above solution:

---

edit (3):

It turns out the missing parts result from z-values being NaN when close to 0:

# color points 'outside' the solid purple
gradient[is.nan(zs)] <- -1

# show those previously hidden points
zs[is.nan(zs)] <- 0

# plot exactly as before
plot_ly(x = xs, y = ys, z = zs, surfacecolor = gradient,
        type = "surface", colorbar = list(title = 'x'))

This could be caused by numerical instability of the substraction when r^2 and y get too close, resulting in negative input for sqrt where the actual input is still non-negative.

This seams unrelated to numerical issues as even when considering +-4 'close' to zero, the 'razor teeth' effect can not be avoided completely:

# re-calculate z-values rounding to zero if 'close'
eps <- 4
zs <- with(coords, ifelse(abs(z2) < eps, 0, sqrt(z2))) %>%
      matrix(ncol = length(xs), byrow = TRUE) %>%
      rbind(-.[nrow(.):1, ], .[1, ])

# plot exactly as before
plot_ly(x = xs, y = ys, z = zs, surfacecolor = gradient,
        type = "surface", colorbar = list(title = 'x'))

Answer 1

One solution would be to flip your axes so that you are rotating around the z axis rather than the x axis. I don't know if this is feasible, given the existing chart that you are adding this figure to, but it does easily solve the 'teeth' problem.

xs <- seq(-9,9,length.out = 20)
ys <- seq(-9,9,length.out = 20)

coords <-
  expand.grid(x = xs, y = ys) %>%
  mutate(z2 = (x^2 + y^2)^(1/4))

zs <- matrix(coords$z2, ncol = length(xs), byrow = TRUE)

plot_ly(x = xs, y = ys, z = zs, surfacecolor = zs,
             type = "surface", colorbar = list(title = 'x')) %>% 
  layout(scene = list(zaxis = list(range = c(1,3))))