Построение неявных уравнений в 3d - PullRequest
32 голосов
/ 13 января 2011

Я бы хотел построить неявное уравнение F (x, y, z) = 0 в 3D. Возможно ли это в Matplotlib?

Ответы [ 7 ]

48 голосов
/ 14 января 2011

Вы можете обмануть matplotlib при построении неявных уравнений в 3D.Просто сделайте одноуровневую контурную диаграмму уравнения для каждого значения z в желаемых пределах.Вы можете повторить процесс вдоль осей Y и Z, чтобы получить более твердую форму.

from mpl_toolkits.mplot3d import axes3d
import matplotlib.pyplot as plt
import numpy as np

def plot_implicit(fn, bbox=(-2.5,2.5)):
    ''' create a plot of an implicit function
    fn  ...implicit function (plot where fn==0)
    bbox ..the x,y,and z limits of plotted interval'''
    xmin, xmax, ymin, ymax, zmin, zmax = bbox*3
    fig = plt.figure()
    ax = fig.add_subplot(111, projection='3d')
    A = np.linspace(xmin, xmax, 100) # resolution of the contour
    B = np.linspace(xmin, xmax, 15) # number of slices
    A1,A2 = np.meshgrid(A,A) # grid on which the contour is plotted

    for z in B: # plot contours in the XY plane
        X,Y = A1,A2
        Z = fn(X,Y,z)
        cset = ax.contour(X, Y, Z+z, [z], zdir='z')
        # [z] defines the only level to plot for this contour for this value of z

    for y in B: # plot contours in the XZ plane
        X,Z = A1,A2
        Y = fn(X,y,Z)
        cset = ax.contour(X, Y+y, Z, [y], zdir='y')

    for x in B: # plot contours in the YZ plane
        Y,Z = A1,A2
        X = fn(x,Y,Z)
        cset = ax.contour(X+x, Y, Z, [x], zdir='x')

    # must set plot limits because the contour will likely extend
    # way beyond the displayed level.  Otherwise matplotlib extends the plot limits
    # to encompass all values in the contour.
    ax.set_zlim3d(zmin,zmax)
    ax.set_xlim3d(xmin,xmax)
    ax.set_ylim3d(ymin,ymax)

    plt.show()

Вот график Гусатского Клубка:

def goursat_tangle(x,y,z):
    a,b,c = 0.0,-5.0,11.8
    return x**4+y**4+z**4+a*(x**2+y**2+z**2)**2+b*(x**2+y**2+z**2)+c

plot_implicit(goursat_tangle)

alt text

Вы можете упростить визуализацию, добавив подсказки глубины с креативным цветовым отображением:

alt text

Вот как выглядит график OP:

def hyp_part1(x,y,z):
    return -(x**2) - (y**2) + (z**2) - 1

plot_implicit(hyp_part1, bbox=(-100.,100.))

alt text

Бонус: Вы можете использовать python для функционального объединения этих неявных функций:

def sphere(x,y,z):
    return x**2 + y**2 + z**2 - 2.0**2

def translate(fn,x,y,z):
    return lambda a,b,c: fn(x-a,y-b,z-c)

def union(*fns):
    return lambda x,y,z: np.min(
        [fn(x,y,z) for fn in fns], 0)

def intersect(*fns):
    return lambda x,y,z: np.max(
        [fn(x,y,z) for fn in fns], 0)

def subtract(fn1, fn2):
    return intersect(fn1, lambda *args:-fn2(*args))

plot_implicit(union(sphere,translate(sphere, 1.,1.,1.)), (-2.,3.))

alt text

2 голосов
/ 13 января 2011

Насколько я знаю, это невозможно. Вы должны решить это уравнение самостоятельно. Использование scipy.optimize - хорошая идея. В простейшем случае вы знаете диапазон поверхности, которую хотите построить, и просто создаете регулярную сетку по x и y и пытаетесь решить уравнение F (xi, yi, z) = 0 для z, давая точка г. Ниже приведен очень грязный код, который может помочь вам

from scipy import *
from scipy import optimize

xrange = (0,1)
yrange = (0,1)
density = 100
startz = 1

def F(x,y,z):
    return x**2+y**2+z**2-10

x = linspace(xrange[0],xrange[1],density)
y = linspace(yrange[0],yrange[1],density)

points = []
for xi in x:
    for yi in y:
        g = lambda z:F(xi,yi,z)
        res = optimize.fsolve(g, startz, full_output=1)
        if res[2] == 1:
            zi = res[0]
            points.append([xi,yi,zi])

points = array(points)
2 голосов
/ 13 января 2011

Матплотлиб ожидает ряд очков;он будет выполнять построение графика, если вы сможете выяснить, как отобразить ваше уравнение.

Ссылаясь на Можно ли построить неявные уравнения, используя Matplotlib? В ответе Майка Грэма предлагается численно использовать scipy.optimizeисследовать неявную функцию.

В http://xrt.wikidot.com/gallery:implicit есть интересная галерея, показывающая различные неявные функции с трассировкой лучей - если ваше уравнение соответствует одной из них, это может дать вам лучшее представление о том, что вы ищетеat.

В противном случае, если вы захотите поделиться фактическим уравнением, возможно, кто-то может предложить более простой подход.

1 голос
/ 26 декабря 2018

Мотивация

Поздний ответ, мне просто нужно было сделать то же самое, и я нашел другой способ сделать это в некоторой степени. Так что я разделяю эту другую точку зрения.

Этот пост не отвечает: (1) Как построить любую неявную функцию F(x,y,z)=0? Но отвечает ли: (2) Как построить параметрические поверхности (не все неявные функции, но некоторые из них), используя сетку с matplotlib?

@ Метод Пола имеет то преимущество, что он непараметрический, поэтому мы можем построить практически все, что захотим, используя контурный метод на каждом топоре, он полностью учитывает (1). Но matplotlib не может легко построить меш из этого метода, поэтому мы не можем напрямую получить поверхность из него, вместо этого мы получаем плоские кривые во всех направлениях. Это то, что мотивировало мой ответ, я хотел обратиться (2).

Рендеринг сетки

Если мы можем параметризовать (это может быть сложно или невозможно), используя не более 2 параметров, поверхность, которую мы хотим построить, то мы можем построить ее с помощью метода matplotlib.plot_trisurf.

То есть из неявного уравнения F(x,y,z)=0, если нам удастся получить параметрическую систему S={x=f(u,v), y=g(u,v), z=h(u,v)}, то мы можем легко построить ее с помощью matplotlib, не прибегая к contour.

Затем рендеринг такой трехмерной поверхности сводится к:

# Render:
ax = plt.axes(projection='3d')
ax.plot_trisurf(x, y, z, triangles=tri.triangles, cmap='jet', antialiased=True) 

Где (x, y, z) - векторы (не meshgrid, см. ravel), функционально вычисленные из параметров (u, v), а параметр triangles - это триангуляция, полученная из (u,v) параметров для построения сетки.

Импорт

Требуемый импорт:

import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits import mplot3d
from matplotlib.tri import Triangulation

Некоторые поверхности

Позволяет параметризовать некоторые поверхности ...

сфера
# Parameters:
theta = np.linspace(0, 2*np.pi, 20)
phi = np.linspace(0, np.pi, 20)
theta, phi = np.meshgrid(theta, phi)
rho = 1

# Parametrization:
x = np.ravel(rho*np.cos(theta)*np.sin(phi))
y = np.ravel(rho*np.sin(theta)*np.sin(phi))
z = np.ravel(rho*np.cos(phi))

# Triangulation:
tri = Triangulation(np.ravel(theta), np.ravel(phi))

enter image description here

шишка
theta = np.linspace(0, 2*np.pi, 20)
rho = np.linspace(-2, 2, 20)
theta, rho = np.meshgrid(theta, rho)

x = np.ravel(rho*np.cos(theta))
y = np.ravel(rho*np.sin(theta))
z = np.ravel(rho)

tri = Triangulation(np.ravel(theta), np.ravel(rho))

enter image description here

торус
a, c = 1, 4
u = np.linspace(0, 2*np.pi, 20)
v = u.copy()
u, v = np.meshgrid(u, v)

x = np.ravel((c + a*np.cos(v))*np.cos(u))
y = np.ravel((c + a*np.cos(v))*np.sin(u))
z = np.ravel(a*np.sin(v))

tri = Triangulation(np.ravel(u), np.ravel(v))

enter image description here

Лента Мёбиуса
u = np.linspace(0, 2*np.pi, 20)
v = np.linspace(-1, 1, 20)
u, v = np.meshgrid(u, v)

x = np.ravel((2 + (v/2)*np.cos(u/2))*np.cos(u))
y = np.ravel((2 + (v/2)*np.cos(u/2))*np.sin(u))
z = np.ravel(v/2*np.sin(u/2))

tri = Triangulation(np.ravel(u), np.ravel(v))

enter image description here

Ограничение

В большинстве случаев Triangulation требуется для координации построения сетки метода plot_trisurf, и этот объект принимает только два параметра, поэтому мы ограничены двумерными параметрическими поверхностями. Маловероятно, что мы могли бы представить Путаницу Гурса этим методом.

1 голос
/ 13 января 2011

Вы смотрели на mplot3d на matplotlib?

0 голосов
/ 08 июня 2011

MathGL (библиотека черчения GPL) может легко построить его. Просто создайте сетку данных со значениями функций f [i, j, k] и используйте функцию Surf3 (), чтобы создать изоповерхность со значением f [i, j, k] = 0. Смотрите этот образец .

0 голосов
/ 21 января 2011

Наконец, я сделал это (я обновил свой matplotlib до 1.0.1). Вот код:

import matplotlib.pyplot as plt
import numpy as np
from mpl_toolkits.mplot3d import Axes3D

def hyp_part1(x,y,z):
    return -(x**2) - (y**2) + (z**2) - 1

fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')

x_range = np.arange(-100,100,10) 
y_range = np.arange(-100,100,10)
X,Y = np.meshgrid(x_range,y_range)
A = np.linspace(-100, 100, 15)

A1,A2 = np.meshgrid(A,A)    

for z in A: 
    X,Y = A1, A2
    Z = hyp_part1(X,Y,z)
    ax.contour(X, Y, Z+z, [z], zdir='z')

for y in A: 
    X,Z= A1, A2
    Y = hyp_part1(X,y,Z)
    ax.contour(X, Y+y, Z, [y], zdir='y')

for x in A:
    Y,Z = A1, A2 
    X = hyp_part1(x,Y,Z)
    ax.contour(X+x, Y, Z, [x], zdir='x')

ax.set_zlim3d(-100,100)
ax.set_xlim3d(-100,100)
ax.set_ylim3d(-100,100)

Вот результат: alt text

Спасибо, Пол!

Добро пожаловать на сайт PullRequest, где вы можете задавать вопросы и получать ответы от других членов сообщества.
...