Задача построения трехмерного графика интеграла от f (x, y, t) dt - PullRequest
Новая
       31

Задача построения трехмерного графика интеграла от f (x, y, t) dt

0 голосов
/ 13 января 2020

Я определил функцию I (a, b) = интеграл f (a, b, t) dt и хочу построить ее, чтобы увидеть, как она зависит от переменных a и b. Сначала я написал программу с графиком y = I (k, x), и она работала просто отлично, но я хотел посмотреть, как она зависит от обеих переменных, поэтому я попытался написать программу, которая отображает ее в 3D.

Программа работала для более простых функций, таких как trigonometri c и полиномов, но когда я пытаюсь построить график I (x, y), она просто выдает мне ошибку «Значение истинности массива с более чем одним элементом неоднозначно. Используйте .any () или a.all () "

Это код, я изначально написал свою собственную программу для аппроксимации интеграла, но затем использовал scipy 

from mpl_toolkits.mplot3d import Axes3D
import matplotlib.pyplot as plt
from matplotlib import cm
from matplotlib.ticker import LinearLocator, FormatStrFormatter
import numpy as np
import scipy.integrate as integrate

def integral(x,y):
    return integrate.quad(lambda t: np.sqrt((x**2 + y**2 - 2*x*y*np.cos(np.pi*t*(np.sqrt(1/x**3) - np.sqrt(1/y**3))))/(x**3*y**3)), 0,  np.sqrt(x**3*y**3))

X = np.arange(0.1,5,0.1)
Y = np.arange(0.1,5,0.1)
X,Y = np.meshgrid(X, Y)
Z = integral(X,Y)


fig = plt.figure()
ax = plt.axes(projection="3d")
ax.plot_wireframe(X, Y, Z, color='green')
ax.set_xlabel('x')
ax.set_ylabel('y')
ax.set_zlabel('z')


ax = plt.axes(projection='3d')
ax.plot_surface(X, Y, Z, rstride=1, cstride=1,
                cmap='winter', edgecolor='none')
ax.set_title('copper');

plt.show()
'''

1 Ответ

0 голосов
/ 13 января 2020

scipy.integrate.quad возвращает кортеж. Вы хотите только первое значение этого. Также вам нужно векторизовать функцию. 

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

def integral(x,y):
    return integrate.quad(lambda t: np.sqrt((x**2 + y**2 - 2*x*y*np.cos(np.pi*t*(np.sqrt(1/x**3) - np.sqrt(1/y**3))))/(x**3*y**3)), 0,  np.sqrt(x**3*y**3))[0]

X = np.arange(0.1,5,0.1)
Y = np.arange(0.1,5,0.1)
X,Y = np.meshgrid(X, Y)
Z = np.vectorize(integral)(X,Y)

fig = plt.figure()
ax = plt.axes(projection="3d")
ax.plot_wireframe(X, Y, Z, color='green')
ax.set_xlabel('x')
ax.set_ylabel('y')
ax.set_zlabel('z')


ax = plt.axes(projection='3d')
ax.plot_surface(X, Y, Z, rstride=1, cstride=1,
                cmap='winter', norm=plt.Normalize(np.nanmin(Z), np.nanmax(Z)), edgecolor='none')

plt.show()

enter image description here

...