Convert Formula to xyz Function (blender python) - Pisot Triaxial Surface
Make a Script.
Seems you have already got there, started answer, thought I'd post as crunching in the equation via the text editor is IMO way simpler

The equations are in parametric form on the Paul Bourke site linked to question.
To make them python expressions only requires adding an asterisk * to represent multiplication.
I find crunching formulas into blender's text boxes a little tedious, so have written as a script.  Enable the Add Mesh Extra Objects addon to have the operator available.
Note: This is the script to run in object mode to create the Pisot Surface (Equations from Paul Bourke link in question) Simply copy script, paste into text editor and press run script.

import bpy
from math import pi

# run in object mode
if context.mode != 'OBJECT':
    bpy.ops.object.mode_set()


bpy.ops.mesh.primitive_xyz_function_surface(
        x_eq=&quot;0.655866 * cos(1.03002 + u) * (2 + cos(v))&quot;,
        y_eq=&quot;0.754878 * cos(1.40772 - u) * (2 + 0.868837 * cos(2.43773 + v))&quot;,
        z_eq=&quot;0.868837 * cos(2.43773 + u) * (2 + 0.495098 * cos(0.377696 - v))&quot;,
        range_u_max=2 * pi, 
        range_u_step=32,
        range_v_max=2 * pi, 
        range_v_step=32,
        edit_mode=False,
        )

Notes. The defaults of the operator on first run can also be seen by bpy.ops.mesh.primitive_add_xyz_function_surface(Tab in blenders python console.  There is no need to add a value to call above that is already a default.
The code below is the (expanded) result of auto completing the operator in the python console. It shows the default settings of the operator that create the sea shell.  It is not the Tisot surface
&gt;&gt;&gt; bpy.ops.mesh.primitive_xyz_function_surface(
primitive_xyz_function_surface()
bpy.ops.mesh.primitive_xyz_function_surface(
    x_eq=&quot;cos(v)*(1+cos(u))*sin(v/8)&quot;, 
    y_eq=&quot;sin(u)*sin(v/8)+cos(v/8)*1.5&quot;, 
    z_eq=&quot;sin(v)*(1+cos(u))*sin(v/8)&quot;, 
    range_u_min=0, 
    range_u_max=6.28319, 
    range_u_step=32, 
    wrap_u=True, 
    range_v_min=0, 
    range_v_max=12.5664, 
    range_v_step=128, 
    wrap_v=False, 
    close_v=False, 
    n_eq=1, 
    a_eq=&quot;0&quot;, 
    b_eq=&quot;0&quot;, 
    c_eq=&quot;0&quot;, 
    f_eq=&quot;0&quot;, 
    g_eq=&quot;0&quot;, 
    h_eq=&quot;0&quot;, 
    show_wire=True, 
    edit_mode=True
)
Add a surface defined defined by 3 functions:
x=F1(u,v), y=F2(u,v) and z=F3(u,v)

Tip to more easily change some of the constants eg 0.868837 could assign them to a helper method.
    x_eq=&quot;A + u&quot;,
    a_eq=&quot;0.868837&quot;