участок поверхности на нерегулярной сетке в среде питона - PullRequest
0 голосов
/ 24 августа 2018

Можно ли создать поверхностный график и контурный график функции, например, u(x,y) = x^2 + y^2

в следующей области, ограниченной уравнением

r(t) = 1+(cos(4*t))^2, x = r(t)*cos(t),y = r(t)*sin(t), 0 < t < 2*pi

enter image description here

Я хочу вариант поверхности следующего рассеянного графика.enter image description here Я также использовал scipy griddata следующим образом

from matplotlib.pyplot  import *
import numpy as np
from mpl_toolkits.mplot3d import Axes3D
from scipy.interpolate import griddata

data = np.array([[ 2.00000000e+00,  0.00000000e+00],
   [ 1.89614525e+00,  1.51126793e-01],
   [ 1.62613327e+00,  2.60869688e-01],
   [ 1.29639187e+00,  3.15328472e-01],
   [ 1.03183519e+00,  3.39805963e-01],
   [ 9.22330309e-01,  3.87424263e-01],
   [ 9.85968247e-01,  5.09806830e-01],
   [ 1.16496791e+00,  7.25092951e-01],
   [ 1.35481653e+00,  1.00089121e+00],
   [ 1.45215333e+00,  1.26290100e+00],
   [ 1.39897768e+00,  1.42707450e+00],
   [ 1.20351374e+00,  1.44074409e+00],
   [ 9.29854857e-01,  1.31248166e+00],
   [ 6.63453200e-01,  1.11474417e+00],
   [ 4.70396020e-01,  9.55864325e-01],
   [ 3.70174853e-01,  9.32786757e-01],
   [ 3.33831717e-01,  1.08590518e+00],
   [ 3.06157327e-01,  1.37738581e+00],
   [ 2.38732323e-01,  1.70413196e+00],
   [ 1.15900227e-01,  1.94068353e+00],
   [-3.96388724e-02,  1.99329358e+00],
   [-1.83621616e-01,  1.84088591e+00],
   [-2.79344160e-01,  1.54430632e+00],
   [-3.22477316e-01,  1.21965478e+00],
   [-3.47231849e-01,  9.87829012e-01],
   [-4.09520279e-01,  9.23088740e-01],
   [-5.55288737e-01,  1.02352335e+00],
   [-7.90735832e-01,  1.21586141e+00],
   [-1.07111976e+00,  1.39113164e+00],
   [-1.31601256e+00,  1.45381050e+00],
   [-1.44544698e+00,  1.36169834e+00],
   [-1.42013085e+00,  1.13913085e+00],
   [-1.26561597e+00,  8.59376173e-01],
   [-1.06700598e+00,  6.06642746e-01],
   [-9.34540362e-01,  4.37040778e-01],
   [-9.54661396e-01,  3.57053472e-01],
   [-1.14896223e+00,  3.28361061e-01],
   [-1.46070082e+00,  2.94327457e-01],
   [-1.77632959e+00,  2.12929020e-01],
   [-1.97334870e+00,  7.85155418e-02],
   [-1.97334870e+00, -7.85155418e-02],
   [-1.77632959e+00, -2.12929020e-01],
   [-1.46070082e+00, -2.94327457e-01],
   [-1.14896223e+00, -3.28361061e-01],
   [-9.54661396e-01, -3.57053472e-01],
   [-9.34540362e-01, -4.37040778e-01],
   [-1.06700598e+00, -6.06642746e-01],
   [-1.26561597e+00, -8.59376173e-01],
   [-1.42013085e+00, -1.13913085e+00],
   [-1.44544698e+00, -1.36169834e+00],
   [-1.31601256e+00, -1.45381050e+00],
   [-1.07111976e+00, -1.39113164e+00],
   [-7.90735832e-01, -1.21586141e+00],
   [-5.55288737e-01, -1.02352335e+00],
   [-4.09520279e-01, -9.23088740e-01],
   [-3.47231849e-01, -9.87829012e-01],
   [-3.22477316e-01, -1.21965478e+00],
   [-2.79344160e-01, -1.54430632e+00],
   [-1.83621616e-01, -1.84088591e+00],
   [-3.96388724e-02, -1.99329358e+00],
   [ 1.15900227e-01, -1.94068353e+00],
   [ 2.38732323e-01, -1.70413196e+00],
   [ 3.06157327e-01, -1.37738581e+00],
   [ 3.33831717e-01, -1.08590518e+00],
   [ 3.70174853e-01, -9.32786757e-01],
   [ 4.70396020e-01, -9.55864325e-01],
   [ 6.63453200e-01, -1.11474417e+00],
   [ 9.29854857e-01, -1.31248166e+00],
   [ 1.20351374e+00, -1.44074409e+00],
   [ 1.39897768e+00, -1.42707450e+00],
   [ 1.45215333e+00, -1.26290100e+00],
   [ 1.35481653e+00, -1.00089121e+00],
   [ 1.16496791e+00, -7.25092951e-01],
   [ 9.85968247e-01, -5.09806830e-01],
   [ 9.22330309e-01, -3.87424263e-01],
   [ 1.03183519e+00, -3.39805963e-01],
   [ 1.29639187e+00, -3.15328472e-01],
   [ 1.62613327e+00, -2.60869688e-01],
   [ 1.89614525e+00, -1.51126793e-01],
   [ 2.00000000e+00, -4.89858720e-16],
   [ 0.00000000e+00, -1.50000000e+00],
   [-1.00000000e+00, -1.00000000e+00],
   [ 0.00000000e+00, -1.00000000e+00],
   [ 1.00000000e+00, -1.00000000e+00],
   [-5.00000000e-01, -5.00000000e-01],
   [ 0.00000000e+00, -5.00000000e-01],
   [ 5.00000000e-01, -5.00000000e-01],
   [-1.50000000e+00,  0.00000000e+00],
   [-1.00000000e+00,  0.00000000e+00],
   [-5.00000000e-01,  0.00000000e+00],
   [ 0.00000000e+00,  0.00000000e+00],
   [ 5.00000000e-01,  0.00000000e+00],
   [ 1.00000000e+00,  0.00000000e+00],
   [ 1.50000000e+00,  0.00000000e+00],
   [ 2.00000000e+00,  0.00000000e+00],
   [-5.00000000e-01,  5.00000000e-01],
   [ 0.00000000e+00,  5.00000000e-01],
   [ 5.00000000e-01,  5.00000000e-01],
   [-1.00000000e+00,  1.00000000e+00],
   [ 0.00000000e+00,  1.00000000e+00],
   [ 1.00000000e+00,  1.00000000e+00],
   [ 0.00000000e+00,  1.50000000e+00]])   
ua = data[:,0]**2+data[:,1]**2 # u=x^2+y^2


xx,yy = np.meshgrid(np.linspace(-2,2,100),np.linspace(-2,2,100))
Ua = griddata((data[:,0],data[:,1]),ua,(xx,yy),method='cubic') 

fig = figure(1)
plot (data[:,0], data[:,1], '*'); # 
fig = figure(2)
ax = fig.gca(projection='3d')

ax.plot_wireframe(xx,yy,Ua,rstride=1,cstride=1,linewidth=.5) 

, но результат не такой хороший, как показано ниже enter image description here

1 Ответ

0 голосов
/ 25 августа 2018

Просто рассчитайте функцию на прямоугольной сетке!

scipy.interpolate.griddata интерполируется на выпуклой оболочке. Вы можете прочитать это в документах . Это означает, что он также интерполирует между вашими долями. Вот почему вы не получаете правильный домен. Интерполяция также по своей природе менее точна, чем просто вычисление функции по сетке.

plot_wireframe хочет прямоугольную сетку, которую вы уже создали. Все, что вам нужно сделать, это рассчитать значения функции на прямоугольной сетке. Чтобы отображать только значения на границе вашего домена, установите все, что находится за пределами np.NaN (не число).

Вот как это сделать:

import matplotlib.pyplot as plt
import numpy as np

# cartesian coordinates
xx,yy = np.meshgrid(np.linspace(-2,2,100),np.linspace(-2,2,100))

# function value across square domain
Ua = xx**2 + yy**2

# polar coordinates
tt = np.arctan2(yy, xx)
rr = np.sqrt(Ua) # re-using x^2 + y^2 -- only works for this function

# r coordinate of domain boundary
domain_boundary = 1 + (np.cos(4*tt))**2

# function value across actual domain, with rest set to NaN
Ua[rr > domain_boundary] = np.NaN

# plotting
fig = plt.figure(2)
ax = fig.gca(projection='3d')
ax.plot_wireframe(xx,yy,Ua,rstride=1,cstride=1,linewidth=.5)

enter image description here

Это решение не идеально, так как вы можете распознать прямоугольную сетку в результате. Вы можете попробовать работать в полярных координатах, как показано в этом официальном примере с matplotlib .

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