This is a follow up to my question from here: https://math.stackexchange.com/questions/4959810/how-did-they-make-this-graph-on-wikipedia

I am trying to make this graph in R:

I first tried this code:

library(plotly)

x_vals <- seq(-3, 3, length.out = 100)
y_vals <- seq(-3, 3, length.out = 100)

x <- rep(x_vals, each = length(y_vals))
y <- rep(y_vals, times = length(x_vals))

# Define the numerator and denominator
num_real <- x^2 - y^2  # Real part of the numerator
num_imag <- 2 * x * y  # Imaginary part of the numerator

denom_real <- (x^2 - y^2 + 2 * x + 2)  # Real part 
denom_imag <- (2 * x * y + 2 * y)  # Imaginary part

real_part <- (num_real * denom_real + num_imag * denom_imag) / 
    (denom_real^2 + denom_imag^2)

z_matrix <- matrix(real_part, nrow = length(x_vals), ncol = length(y_vals))

plot_ly(x = ~x_vals, y = ~y_vals, z = ~z_matrix) %>%
    add_surface(colorscale = list(
        c(0, 'rgb(255,182,193)'),  # Light pink
        c(1, 'rgb(128,0,128)')     # Purple
    )) %>%
    layout(scene = list(
        xaxis = list(title = "Re(z)"),
        yaxis = list(title = "Im(z)"),
        zaxis = list(title = "Real Part of f(z)")
    ),
    title = "3D Plot of the Real Part of f(z)")`

Based on feedback I got, I then reformatted the code to try a different way:

`library(plotly)

f <- function(x, y) {
    z <- complex(real = x, imaginary = y)
    Re((z^2) / (z^2 + 2*z + 2))
}

# grid
x <- y <- seq(-4, 4, length.out = 200)
grid <- expand.grid(x = x, y = y)

z <- outer(x, y, Vectorize(f))

# Main plot
p <- plot_ly(x = x, y = y, z = z) %>%
    add_surface(
        colorscale = list(
            c(0, "rgb(255,182,193)"),  # Light pink
            c(1, "rgb(128,0,128)")     # Purple
        ),
        lighting = list(ambient = 0.6, diffuse = 0.5, specular = 0.1),
        contours = list(
            z = list(
                show = TRUE,
                usecolormap = TRUE,
                highlightcolor = "#ffffff",
                project = list(z = TRUE)
            )
        )
    )

for (i in seq(-4, 4, length.out = 13)) {
    p <- p %>% add_mesh(
        x = x,
        y = rep(i, length(x)),
        z = outer(x, rep(i, length(x)), Vectorize(f)),
        opacity = 0.2,
        color = I("white"),
        hoverinfo = "none"
    )
    p <- p %>% add_mesh(
        x = rep(i, length(y)),
        y = y,
        z = outer(rep(i, length(y)), y, Vectorize(f)),
        opacity = 0.2,
        color = I("white"),
        hoverinfo = "none"
    )
}

t <- seq(0, 2*pi, length.out = 100)

# Curve 1: 2*exp(i*t)
curve1 <- data.frame(
    x = 2 * cos(t),
    y = 2 * sin(t),
    z = f(2 * cos(t), 2 * sin(t))
)

# Curve 2: exp(i*t)/2 - 1 + i
curve2 <- data.frame(
    x = 0.5 * cos(t) - 1,
    y = 0.5 * sin(t) + 1,
    z = f(0.5 * cos(t) - 1, 0.5 * sin(t) + 1)
)

# Curve 3: exp(i*t)/2 - 1 - i
curve3 <- data.frame(
    x = 0.5 * cos(t) - 1,
    y = 0.5 * sin(t) - 1,
    z = f(0.5 * cos(t) - 1, 0.5 * sin(t) - 1)
)

p <- p %>%
    add_trace(data = curve1, x = ~x, y = ~y, z = ~z, type = 'scatter3d', mode = 'lines', line = list(color = 'black', width = 4), name = 'Curve 1') %>%
    add_trace(data = curve2, x = ~x, y = ~y, z = ~z, type = 'scatter3d', mode = 'lines', line = list(color = 'red', width = 4), name = 'Curve 2') %>%
    add_trace(data = curve3, x = ~x, y = ~y, z = ~z, type = 'scatter3d', mode = 'lines', line = list(color = 'blue', width = 4), name = 'Curve 3')

p <- p %>% layout(
    scene = list(
        aspectmode = "manual",
        aspectratio = list(x = 1, y = 1, z = 0.65),
        xaxis = list(range = c(-4, 4), title = "Re(z)"),
        yaxis = list(range = c(-4, 4), title = "Im(z)"),
        zaxis = list(range = c(-4, 4), title = "Re(f(z))")
    )
)

p`

But in my code the axes are incorrect and the parametrized curves (black, red, and blue loops) are oriented the wrong way and do not lie on the function.

Can someone please help me here?