Построение неявных уравнений в 3d
Я хотел бы построить неявное уравнение F (x, y, z) = 0 в 3D. Возможно ли это в Matplotlib?
7 ответов
Вы можете обмануть 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)
Вы можете упростить визуализацию, добавив подсказки глубины с помощью креативной цветовой карты:
Вот как выглядит сюжет ОП:
def hyp_part1(x,y,z): return -(x**2) - (y**2) + (z**2) - 1 plot_implicit(hyp_part1, bbox=(-100.,100.))
Бонус: вы можете использовать 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.))
Мотивация
Поздний ответ, мне просто нужно было сделать то же самое, и я нашел другой способ сделать это в какой-то степени. Так что я разделяю это с другой точки зрения.
Этот пост не отвечает: (1) Как построить любую неявную функцию F(x,y,z)=0 ? Но отвечает ли: (2) Как построить параметрические поверхности (не все неявные функции, но некоторые из них), используя сетку с matplotlib ?
Метод Пола имеет то преимущество, что он не параметрический, поэтому мы можем построить практически все, что захотим, используя контурный метод на каждом топоре, он полностью обращается к (1). Но matplotlib не может легко построить меш из этого метода, поэтому мы не можем напрямую получить поверхность из него, вместо этого мы получаем плоские кривые во всех направлениях. Это то, что мотивировало мой ответ, я хотел обратиться (2).
Рендеринг сетки
Если мы можем параметризовать (это может быть сложно или невозможно), используя не более 2 параметров, поверхность, которую мы хотим построить, то мы можем построить ее с помощью matplotlib.plot_trisurf .
То есть из неявного уравнения F(x,y,z)=0 , если мы можем получить параметрическую систему S= , то мы можем легко построить ее с помощью 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))
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))
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))
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))
Ограничение
В большинстве случаев для координации требуется Triangulation построение сетки метода plot_trisurf , и этот объект принимает только два параметра, поэтому мы ограничены двумерными параметрическими поверхностями. Вряд ли мы могли бы представить Путаницу Гурса этим методом.
Насколько я знаю, это невозможно. Вы должны решить это уравнение самостоятельно. Использование 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)
Наконец, я сделал это (я обновил свой 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)
Matplotlib ожидает ряд очков; это сделает построение графика, если вы сможете выяснить, как визуализировать ваше уравнение.
Ссылаясь на Можно ли строить неявные уравнения с помощью Matplotlib? В ответе Майка Грэма предлагается использовать scipy.optimize для численного исследования неявной функции.
На есть интересная галерея http://xrt.wikidot.com/gallery:implicit. множество неявных функций с трассировкой лучей — если ваше уравнение соответствует одной из них, это может дать вам лучшее представление о том, на что вы смотрите.
В противном случае, если вы хотите поделиться фактическим уравнением, возможно, кто-то может предложить более простой подход.
MathGL (библиотека черчения GPL) может легко ее отобразить. Просто создайте сетку данных со значениями функций f [i, j, k] и используйте функцию Surf3 (), чтобы создать изоповерхность со значением f [i, j, k] = 0. Посмотрите этот пример.
Можно ли построить неявные уравнения с помощью Matplotlib?
Я хотел бы построить неявные уравнения (вида f(x, y)=g (x, y), например. X^y=y^x) в Matplotlib. Возможно ли это?
7 ответов
Я не верю, что для этого есть очень хорошая поддержка, но вы можете попробовать что-то вроде
import matplotlib.pyplot from numpy import arange from numpy import meshgrid delta = 0.025 xrange = arange(-5.0, 20.0, delta) yrange = arange(-5.0, 20.0, delta) X, Y = meshgrid(xrange,yrange) # F is one side of the equation, G is the other F = Y**X G = X**Y matplotlib.pyplot.contour(X, Y, (F - G), [0]) matplotlib.pyplot.show()
посмотреть API docs на contour : если четвертый аргумент является последовательностью, он указывает, какие контурные линии следует построить. Но сюжет будет так же хорош, как разрешение ваших диапазонов, и есть определенные функции, которые он может никогда не получить правильно, часто в точках пересечения.
поскольку вы пометили этот вопрос с симпатией, я приведу такой пример.
from sympy import var, plot_implicit var('x y') plot_implicit(x*y**3 - y*x**3)
matplotlib не строит уравнения; он строит серии точек. Вы можете использовать такой инструмент, как scipy.optimize численно вычислить точки y из значений x (или наоборот) неявных уравнений численно или любое количество других инструментов по мере необходимости.
например, вот пример, где я строю неявное уравнение x ** 2 + x * y + y ** 2 = 10 в определенной области.
from functools import partial import numpy import scipy.optimize import matplotlib.pyplot as pp def z(x, y): return x ** 2 + x * y + y ** 2 - 10 x_window = 0, 5 y_window = 0, 5 xs = [] ys = [] for x in numpy.linspace(*x_window, num=200): try: # A more efficient technique would use the last-found-y-value as a # starting point y = scipy.optimize.brentq(partial(z, x), *y_window) except ValueError: # Should we not be able to find a solution in this window. pass else: xs.append(x) ys.append(y) pp.plot(xs, ys) pp.xlim(*x_window) pp.ylim(*y_window) pp.show()
существует неявное уравнение (и неравенство) плоттера в sympy. Он создается как часть GSoC и создает графики как экземпляры фигур matplotlib.
начиная с версии sympy 0.7.2-это:
>>> from sympy.plotting import plot_implicit >>> p = plot_implicit(x < sin(x)) # also creates a window with the plot >>> the_matplotlib_axes_instance = p._backend._ax
Если вы хотите использовать что-то другое, кроме matplotlib (но все же python), есть sage:
большое спасибо Стив, Майк, Алекс. Я согласился с решением Стива (см. код ниже). Моя единственная оставшаяся проблема заключается в том, что контурный график появляется за моими линиями сетки, в отличие от обычного графика, который я могу заставить спереди с помощью zorder. Еще один халп очень признателен.
import matplotlib.pyplot as plt from matplotlib.ticker import MultipleLocator, FormatStrFormatter import numpy as np fig = plt.figure(1) ax = fig.add_subplot(111) # set up axis ax.spines['left'].set_position('zero') ax.spines['right'].set_color('none') ax.spines['bottom'].set_position('zero') ax.spines['top'].set_color('none') ax.xaxis.set_ticks_position('bottom') ax.yaxis.set_ticks_position('left') # setup x and y ranges and precision x = np.arange(-0.5,5.5,0.01) y = np.arange(-0.5,5.5,0.01) # draw a curve line, = ax.plot(x, x**2,zorder=100) # draw a contour X,Y=np.meshgrid(x,y) F=X**Y G=Y**X ax.contour(X,Y,(F-G),[0],zorder=100) #set bounds ax.set_xbound(-1,7) ax.set_ybound(-1,7) #produce gridlines of different colors/widths ax.xaxis.set_minor_locator(MultipleLocator(0.2)) ax.yaxis.set_minor_locator(MultipleLocator(0.2)) ax.xaxis.grid(True,'minor',linestyle='-') ax.yaxis.grid(True,'minor',linestyle='-') minor_grid_lines = [tick.gridline for tick in ax.xaxis.get_minor_ticks()] for idx,loc in enumerate(ax.xaxis.get_minorticklocs()): if loc % 2.0 == 0: minor_grid_lines[idx].set_color('0.3') minor_grid_lines[idx].set_linewidth(2) elif loc % 1.0 == 0: minor_grid_lines[idx].set_c('0.5') minor_grid_lines[idx].set_linewidth(1) else: minor_grid_lines[idx].set_c('0.7') minor_grid_lines[idx].set_linewidth(1) minor_grid_lines = [tick.gridline for tick in ax.yaxis.get_minor_ticks()] for idx,loc in enumerate(ax.yaxis.get_minorticklocs()): if loc % 2.0 == 0: minor_grid_lines[idx].set_color('0.3') minor_grid_lines[idx].set_linewidth(2) elif loc % 1.0 == 0: minor_grid_lines[idx].set_c('0.5') minor_grid_lines[idx].set_linewidth(1) else: minor_grid_lines[idx].set_c('0.7') minor_grid_lines[idx].set_linewidth(1) plt.show()