[Tex/LaTex] Plot the cusp catastrophe surface

asymptotediagramsmetapostpgfplotspstricks

I'm trying to draw (plot) the Cusp catastrophe surface into latex using any software logic. I am completely stuck because I'm definitely not a mathematician. The surface has the potential:

V(x) = F(x,u,v) = x^4 + ux^2 + vx

so it should be given by this implicit equation:

4x^3 + 2ux + v = 0

But I don't know how to implement this. Apparently gnuplot can't do this (even if I find something approaching) so I should go with a CAS like Octave or scripting something in Maxima first, then feed gnuplot or export the graphic. Ideally, it would be great if it was as easy as this Mathematica example:

F[x_, u_, v_] :=  x^4 + u x^2 + v*x

Needs["Graphics`ContourPlot3D`"]

ContourPlot3D[Evaluate[D[F[x, u, v], x]],
  {u, -2.5, 2}, {v, -2.2, 3}, {x, -1.4, 1.3}, 
  PlotPoints -> 7, ViewPoint -> {-1.5, 1.5, 1.4}, 
  Axes -> True, ContourStyle -> {EdgeForm[]}, 
  AxesLabel -> TraditionalForm /@ {u, v, x}] // Timing

edit / updates

So I dig I lot, apparently, I can't do it perfectly without the ability of a mathematician (that's for sure) si I decided I should stick with an implicit solving method. I can translate the Mathematica code into Sage:

u, v, x = var('u, v, x')
f(u,v,x) = x^4 + u*x^2 + v*x
g = f.diff(x)
implicit_plot3d(g, (u, -2, 2), (v, -2, 2), (x, -2, 2))

It produces neat results:
cusp surface from Sage

You can also check the result on cloud.sagemath (subscription required).

But unfortunately 3d plot can't be exported in vector format, maybe there's a method for raw (x,y,z) points export, dunno…
I tried an other solution with a the experimental psplotImpIIID feature of ps-tricks, it's ugly and produce a lot of artefacts.

And to end up with pstricks, one could do something like this:

\begin{pspicture}(-6,-4)(6,5)
    \pstThreeDCoor[xMin=-4,xMax=4,yMin=-4,yMax=4,zMin=-4,zMax=4,RotY=90,RotX=25]
    \psplotThreeD[plotstyle=line,hiddenLine,% does not work in my case (xelatex?)
    yPlotpoints=50,xPlotpoints=50,linewidth=.5pt,algebraic](-3,3)(-3,3){4*x^3 - 2*y*x}
\end{pspicture}

pst-plot3d cusp surface

The advantage being that we can easily rotate the axis, the disadvantage (from the perspective of a pstricks noob) being styling, plus I need to clip the result image. There's also a solution that uses a R package which I tried and and it works great and outputs vector graphics. But honestly, I think I loop the web (I even scripted it as a Blender Python) and there is definitely no easy straight solution for this.

Best Answer

\documentclass[pstricks]{standalone}
\usepackage{pst-solides3d}
\begin{document}

\begin{pspicture}(-4,-6)(5,9)
\psset{viewpoint=100 30 40 rtp2xyz,lightsrc=viewpoint, Decran=120}
\psSurface[ngrid=.15 .15,incolor=yellow,hue=0 1,linewidth=0.1\pslinewidth,
  algebraic,axesboxed](-1,-2)(1,2){ 4*x^3 - 2*y*x}
\end{pspicture}

\end{document}

enter image description here

The z max/min values can be handled only on PostScript level (will be changed in the future)

\documentclass[pstricks]{standalone}
\usepackage{pst-solides3d}
\begin{document}

\begin{pspicture}(-4,-6)(5,9)
\psset{viewpoint=100 30 40 rtp2xyz,lightsrc=viewpoint, Decran=120}
\psSurface[ngrid=0.1 0.15,incolor=yellow,hue=0 1,linewidth=0.1\pslinewidth,axesboxed, 
 % algebraic](-2,-2)(2,2){4*x^3 -2*y*x}
  ](-2,-2)(1,3){ 4 x 3 exp mul 2 y mul x mul sub 
    dup 4 gt { pop 4 } if 
    dup -4 lt { pop -4 } if 
    }
\end{pspicture}

\end{document}

However, the values are not ignored, they are only set to the max/min.

enter image description here