Like this?
using Makie, GLMakie
fig = Figure()
range = LinRange(-1, 1, 100) # 100-element LinRange{Float64, Int64}
cube = [ (abs.(x).^10 + abs.(y).^10 + abs.(z).^10) for x = range, y = range, z = range ] # 100×100×100 Array{Float64, 3}
ax = Axis3( fig[1,1], aspect = :data, azimuth = 1.17 * pi, viewmode = :fit, title = "Cube" )
volume!( cube, algorithm = :iso, isorange = 0.05, isovalue = 1 )
fig

Or you can use a parametric surface, a supershape:
using Makie, GLMakie
function r(phi; a, b, m, n1, n2, n3)
return 1 / (abs(cos(m*phi/4)/a)^n2 + (abs(sin(m*phi/4)/b)^n3))^(1/n1)
end
phi = (-pi/2):0.01:(pi/2)
theta = (-pi):0.01:(pi)
x = [
r(theta; a=1, b=1, m=4, n1=10, n2=10, n3=10) * cos(theta) *
r(phi; a=1, b=1, m=4, n1=10, n2=10, n3=10) * cos(phi)
for phi in phi, theta in theta
]
y = [
r(theta; a=1, b=1, m=4, n1=10, n2=10, n3=10) * sin(theta) *
r(phi; a=1, b=1, m=4, n1=10, n2=10, n3=10) * cos(phi)
for phi in phi, theta in theta
]
z = [
r(phi; a=1, b=1, m=4, n1=10, n2=10, n3=10) * sin(phi)
for phi in phi, theta in theta
]
fig, _ = surface(x, y, z)
fig

Better to use a function:
function r(phi; a, b, m, n1, n2, n3)
return 1 / (abs(cos(m*phi/4)/a)^n2 + (abs(sin(m*phi/4)/b)^n3))^(1/n1)
end
function supershape(p1, p2)
phi = (-pi/2):0.01:(pi/2)
theta = (-pi):0.01:(pi)
x = [
r(theta; a=p1.a, b=p1.b, m=p1.m, n1=p1.n1, n2=p1.n2, n3=p1.n3) * cos(theta) *
r(phi; a=p2.a, b=p2.b, m=p2.m, n1=p2.n1, n2=p2.n2, n3=p2.n3) * cos(phi)
for phi in phi, theta in theta
]
y = [
r(theta; a=p1.a, b=p1.b, m=p1.m, n1=p1.n1, n2=p1.n2, n3=p1.n3) * sin(theta) *
r(phi; a=p2.a, b=p2.b, m=p2.m, n1=p2.n1, n2=p2.n2, n3=p2.n3) * cos(phi)
for phi in phi, theta in theta
]
z = [
r(phi; a=p2.a, b=p2.b, m=p2.m, n1=p2.n1, n2=p2.n2, n3=p2.n3) * sin(phi)
for phi in phi, theta in theta
]
return (x = x, y = y, z = z)
end
params1 = (a=1, b=1, m=4, n1=10, n2=10, n3=10)
params2 = params1
x, y, z = supershape(params1, params2)
fig, _ = surface(x, y, z)
EDIT
The rendering is not nice with surface
. It's better to do a mesh:
function parametricMesh(f, umin, umax, vmin, vmax, nu, nv)
u_ = LinRange(umin, umax, nu)
v_ = LinRange(vmin, vmax, nv)
vertices = Array{Float64}(undef, nu*nv, 3)
triangles = Array{Int64}(undef, 2*(nu-1)*(nv-1), 3)
k = 1
for i in 1:nv
v = v_[i]
for j in 1:nu
vertices[k,:] = f(u_[j], v)
k = k+1
end
end
k = 1
for i in 1:(nv-1)
for j in 1:(nu-1)
a = (i-1) * nu + j
b = (i-1) * nu + j + 1
c = i*nu + j + 1
d = i*nu + j
triangles[2*(k-1)+1,:] = [b, a, d]
triangles[2*k,:] = [c, b, d]
k = k+1
end
end
return (vertices = vertices, triangles = triangles)
end
function r(phi; a, b, m, n1, n2, n3)
return 1 / (abs(cos(m*phi/4)/a)^n2 + (abs(sin(m*phi/4)/b)^n3))^(1/n1)
end
p1 = (a=1, b=1, m=4, n1=10, n2=10, n3=10)
p2 = p1
function f(phi, theta)
x = r(theta; a=p1.a, b=p1.b, m=p1.m, n1=p1.n1, n2=p1.n2, n3=p1.n3) * cos(theta) *
r(phi; a=p2.a, b=p2.b, m=p2.m, n1=p2.n1, n2=p2.n2, n3=p2.n3) * cos(phi)
y = r(theta; a=p1.a, b=p1.b, m=p1.m, n1=p1.n1, n2=p1.n2, n3=p1.n3) * sin(theta) *
r(phi; a=p2.a, b=p2.b, m=p2.m, n1=p2.n1, n2=p2.n2, n3=p2.n3) * cos(phi)
z = r(phi; a=p2.a, b=p2.b, m=p2.m, n1=p2.n1, n2=p2.n2, n3=p2.n3) * sin(phi)
return [x, y, z]
end
vertices, triangles = parametricMesh(f, -pi/2, pi/2, -pi, pi, 50, 50)
mesh(vertices, triangles, color = "yellow")
