9
$\begingroup$

Building on from my previous question, how can I plot something more complex than a sphere (and not as easy to split into a set of parametric equations)?

Say, for example

$$2x^2+y^2-5z^2+z-7x=16$$

Graphed on Desmos (use the z-slider to change the z value) for those who cannot as easily visualize the 3D surface.

$\endgroup$
  • $\begingroup$ You could potentially use Python script to build a mesh by layers using successive values of Z. $\endgroup$ – Rich Sedman Sep 16 '17 at 9:41
  • 2
    $\begingroup$ Welcome to the site :) It looks like you have three separate accounts. If you would care to register one of them fully, we can merge them together and you'll be able to edit your question, accept answers, and do all those other fun things. Thanks $\endgroup$ – gandalf3 Sep 16 '17 at 22:24
8
$\begingroup$

As far as I know, there isn't an add-on or method that can graph more complex mathematical functions in the way you want. As an alternative, you can put your equation into nodes and then render it with a Volume Scatter shader.

Begin with the default cube. Add a new material, and add a Texture Coordinate node (Shift + A > Input > Texture Coordinate) and Separate XYZ node (Shift + A > Converter > Separate XYZ):

enter image description here

Then, create your equation using Math nodes (Shift + A > Converter > Math). Here, I created your example equation:

enter image description here

Now, add a Color Mix node (Shift + A > Color > MixRGB) and set it to Multiply. Plug the texture coordinates in to one slot, and a Value node in to the other. This will scale our window:

enter image description here

Finally, multiply the output by 10 to increase the value, then plug it into the Density of a Volume Scatter and a Volume Absorption shader added together:

enter image description here

This will nicely display the function when rendered, though the shading isn't ideal:

enter image description here
click for larger image

To optimize rendering, change the following settings:

enter image description here

Although it doesn't give you a workable mesh, this can definitely give you a renderable graph.

We can also improve the look of it by coloring it:

enter image description here
click for larger image


|improve this answer|||||
$\endgroup$
  • $\begingroup$ Need 50 rep to comment... @ScottMilner I was hoping there was a way to preview in viewport directly (not viewport shading: rendered, but viewport shading: solid). The less CPU power, the better I am actually looking only to visualize the shape (so quality is not important. Quantity is) $\endgroup$ – math Sep 16 '17 at 5:25
7
$\begingroup$

Discrete approach with Python

A possible approch that avoids solving the equation would consist in interating through all possible points in the 3D space, check whether or not in each one of them the equation is verified, and if it is create a vertex there.

As in the 3D space there are infinite points, we must discretize it and limits the boundaries of the investigation domain. For the sake of semplicity, I choose to consider, for example, only integer numbers from -100 to +100 for each axis. Here's the code for such iteration:

for x in range(-100,100):
    for y in range(-100,100):
        for z in range(-100,100):

Then, for each vector we'll calculate the left member of the equation and check if it's equal to the right member of the equation previusly defined. As we discretized the domain, it's best to give to the result a bit of "tollerance" in verifying the equation (that's why I'm using the "+/-toll" trick upon the right member verification).

        function=2*x**2+y**2-5*z**2+z-7*x #function here

        if((function<(goal+toll))and(function>(goal-toll))):
          bm.verts.new([x,y,z])

Here's a possible code with all the lines needed to create an object in the scene:

import bpy
import bmesh

mesh = bpy.data.meshes.new("mesh")
obj = bpy.data.objects.new("Function", mesh)  

scene = bpy.context.scene
scene.objects.link(obj)
scene.objects.active = obj
obj.select = True 

mesh = bpy.context.object.data
bm = bmesh.new()

goal=16 #right member of the equation
toll=5 #threshold

for x in range(-100,100):
    for y in range(-100,100):
        for z in range(-100,100):
            function=2*x**2+y**2-5*z**2+z-7*x #function here
            if((function<(goal+toll))and(function>(goal-toll))):
                bm.verts.new([x,y,z])

bm.to_mesh(mesh)  
bm.free()

Examples of usage:

    function=2*x**2+y**2-5*z**2+z-7*x

enter image description here Note: as you can see, it takes a while to iterate through all the 8 000 000 positions...


This script only provides some of the vertices due to the approximation, but, with a proper resolution, you can obtain a good representation of the surface.

The meshing of the vertices into a surface needs it's own approach. This "brute force plotting" I'm proposing does not go well with the "on the fly" definitions of faces. The set of point is not structured in any way. You'll need something to mesh the pointcloud.

Here's another version of the code where you can better control the domain size and resolution using floats to inspect the region closer to the origin in high definition:

import bpy
import bmesh

max=7
min=-7
prec=10 #precision of the calc (1=1 units, 10 = 0.1 units, 100 = 0.01 units...)
tol=0.5 #tollerance

def f(x, y, z): return 2*x**2+y**2-5*z**2+z-7*x -16 #the function = 0

mesh = bpy.data.meshes.new("mesh")
obj = bpy.data.objects.new("XYZ Function", mesh)  
scene = bpy.context.scene
scene.objects.link(obj)
scene.objects.active = obj
obj.select = True 
mesh = bpy.context.object.data
bm = bmesh.new()

for x in [float(j)/prec for j in range(min*prec, max*prec+1,1)]:
    for y in [float(j)/prec for j in range(min*prec, max*prec+1,1)]:
        for z in [float(j)/prec for j in range(min*prec, max*prec+1,1)]:
            if(abs(f(x,y,z))<tol):
                bm.verts.new([x,y,z])

bm.to_mesh(mesh)  
bm.free()

Here's a possible result:

enter image description here

You can for example export points to Meshlab and take advantage of the built in tools. Here's an example using the previusly obtained set of points.

enter image description here

|improve this answer|||||
$\endgroup$
6
$\begingroup$

Quadric Surfaces

The equation in the form in the question is a quadric surface. By completing the square we can translate the axis and reduce to form

(x-7/4)² - 49/16 + 1/2y² - 5/2(z - 1/10)² + 5/200 - 8 = 0

or

X² + 0.5Y² - 2.5Z² - - 11.0375 = 0 

enter image description here

A sample of the surface created by iterating through z (an xy plane) then finding the 0, 1 or 2 roots for each y using the quadratic formula to solve for x. Notice the missing faces where x=0. Can run thru solving for y and have missing faces around y=0. Overlaying the two gives the shape, but not a pretty mesh.

I've put in both the original equation f(x, y, z) and the translated form f2(x, y, z).

import bpy
import bmesh
from math import sqrt

def f(x, y, z):
    return 2 * x * x + y * y - 5 * z * z + z -7 * x -16
# f transformed to have origin (0, 0, 0)
def f2(x, y, z):
    x * x  + 0.5 * y * y - 2.5 * z * z - 11.0375

def quad_root(a, b, c):
    axis = -b / (2 * a)
    check = b * b - 4 * a * c
    if abs(check) < 0.0001:
        return True, axis, axis
    if check < 0:
        return False, None, None
    dist = sqrt(check) / (2 * a)
    return True, axis - dist, axis + dist

bm = bmesh.new()

dom = [m for m in range(-10, 11)]
use_f2 = False # use f2()
planes = []
# coefficents of x
if use_f2:
    a, b = 1, 0 # f2()
else:
    a, b = 2, -7 # f()


for z in dom:
    # for each xy plane
    neg_line = []
    pos_line = []
    for y in dom:
        # for each line in plane
        if use_f2:
            # constant f2()
            c = 0.5 * y * y - 2.5 * z * z -11.0375
        else:
            # constant f()
            c = y * y - 5 * z * z + z -16
        # find roots
        has_solution, neg_root, pos_root = quad_root(a, b, c)
        if has_solution:
            neg_root = bm.verts.new((neg_root, y, z))
            pos_root = bm.verts.new((pos_root, y, z))
        pos_line.append(pos_root)
        neg_line.append(neg_root)
    planes.append([neg_line, pos_line])
# skin it
if len(planes) > 1:
    for i in range(len(planes) - 1):
        n0, p0 = planes[i]
        n1, p1 = planes[i + 1]       
        # pos 
        for j in range(len(n0) - 1):
            verts = [n0[j], n1[j], n1[j+1], n0[j + 1]]
            if None not in verts:
                f = bm.faces.new(verts)
            verts = [p0[j], p1[j], p1[j+1], p0[j + 1]]
            if None not in verts:
                f = bm.faces.new(verts)        
# make an object        
scene = bpy.context.scene
mesh = bpy.data.meshes.new("Quadric")
bm.to_mesh(mesh)
ob = bpy.data.objects.new("Quadric", mesh)
scene.objects.link(ob)
scene.objects.active = ob
ob.scale *= 0.01

To solve for y for f(x, y, z) edit

        a, b = 1, 0
        c = 2 * x * x - 5 * z * z + z -7 * x - 16
        has_solution, neg_root, pos_root = quad_root(a, b, c)
        if has_solution:
            #print(neg_root, pos_root)
            neg_root = bm.verts.new((x, neg_root, z))
            pos_root = bm.verts.new((x, pos_root, z))

Result after running in both x and y and joining

enter image description here

|improve this answer|||||
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.