How to plot curve over surface with Plots.jl?

I have written some code to generate a gif for the double pendulum problem:

# In another code not shown here, I've:
# - computed solution of double pendulum problem
# - extracted independent variables θ₁, θ₂ from solution

N = length(θ₁)
# plot configuration space, i.e. a torus
Θ₁ = -π:0.01:π
Θ₂ = -π:0.01:π
X_torus = [(1 + cos(Θ₂ᵢ)) * cos(Θ₁ᵢ) for Θ₁ᵢ in Θ₁, Θ₂ᵢ in Θ₂]
Y_torus = [(1 + cos(Θ₂ᵢ)) * sin(Θ₁ᵢ) for Θ₁ᵢ in Θ₁, Θ₂ᵢ in Θ₂]
Z_torus = [sin(Θ₂ᵢ) for Θ₁ᵢ in Θ₁, Θ₂ᵢ in Θ₂]
p₂ = surface(X_torus, Y_torus, Z_torus, colorbar = :none, legend = false, xlims = [-2, 2], ylims = [-2, 2], zlims = [-2, 2])
# create gif
anim = @animate for i = 2:200:N
  # plot first pendulum
  i_ = max(1, i-1000)
  x1 = @. cos(θ₁[i_:i] - π/2)
  y1 = @. sin(θ₁[i_:i] - π/2)
  x₁ = x1[end]
  y₁ = y1[end]
  α = range(0.0, 1.0, length = length(x1))
  p₁ = plot(x1, y1, α = α, color = 1, legend = false, xlims = [-2, 2], ylims = [-2, 2])
  # plot second pendulum
  x2 = @. x1 + cos(θ₂[i_:i] - π/2)
  y2 = @. y1 + sin(θ₂[i_:i] - π/2)
  x₂ = x2[end]
  y₂ = y2[end]
  plot!(p₁, x2, y2, color = 3, α = α, xlabel = L"x", ylabel = L"y")
  plot!(p₁, [0; x₁], [0; y₁], linewidth = 3, color = :black)
  plot!(p₁, [x₁; x₂], [y₁; y₂], linewidth = 3, color = :black)
  scatter!(p₁, (0, 0), markersize = 5, markercolor = :black)
  scatter!(p₁, (x₁, y₁), markersize = 5, markercolor = 1)
  scatter!(p₁, (x₂, y₂), markersize = 5, markercolor = 3)
  # get torus' plot
  p₃ = deepcopy(p₂)
  # plot path of representative particle in configuration space
  X = @. (1 + cos(θ₁[i_:i])) * cos(θ₂[i_:i])
  Y = @. (1 + cos(θ₁[i_:i])) * sin(θ₂[i_:i])
  Z = @. sin(θ₁[i_:i])
  plot!(p₃, X, Y, Z, color = :black)
  scatter!(p₃, (X[end], Y[end], Z[end]), markersize = 5, markercolor = :black)
  plot(p₁, p₃, size = (600, 300))
end
gif(anim, "double-pendulum.gif")

movie

On the left, the actual double pendulum; on the right, the trajectory traced by a representative “particle” over its configuration space (i.e. a torus parametrised by \theta_1 and \theta_2). My problem is that such a trajectory always goes in front of the torus, whereas it should sometimes go behind. What am I doing wrong? Am I missing some option for plot? Is this a Plots or a GR issue?

5 Likes

The plotlyjs() backend doesn’t exhibit this issue.

A simpler MWE to check this, is to plot two antipodal points after the surface command:
scatter!([-2,2], [0, 0], [0, 0], ms=5, color=:black)

MWE
using Plots; gr()  # replace by plotlyjs()
Θ₁ = -π:0.01:π
Θ₂ = -π:0.01:π
X_torus = [(1 + cos(Θ₂ᵢ)) * cos(Θ₁ᵢ) for Θ₁ᵢ in Θ₁, Θ₂ᵢ in Θ₂]
Y_torus = [(1 + cos(Θ₂ᵢ)) * sin(Θ₁ᵢ) for Θ₁ᵢ in Θ₁, Θ₂ᵢ in Θ₂]
Z_torus = [sin(Θ₂ᵢ) for Θ₁ᵢ in Θ₁, Θ₂ᵢ in Θ₂]
surface(X_torus, Y_torus, Z_torus, colorbar = :none, legend = false, xlims = [-2, 2], ylims = [-2, 2], zlims = [-2, 2]);
scatter!([-2,2], [0, 0], [0, 0], ms=5, color=:black)

Thanks! I would really like to get it working with gr(). Should I open an issue over at the GitHub page of GR.jl using your MWE?

Sure.

For reference: https://github.com/jheinen/GR.jl/issues/446