0

実際、私はこのコードを持っています:

from mpl_toolkits.mplot3d import Axes3D
from matplotlib import cm
import matplotlib.pyplot as plt
from numpy import *

def f(x,y):
    r=x**2 + y**2
    return r

n=4.
b=1.
a=-b
h=(2*b)/n
print h
hx=h ##This line##
fig = plt.figure()
ax = Axes3D(fig)
X = arange(a, b+hx, hx)
Y = arange(a, b+h, h)
n = len(X)
m = len(Y)
Z = zeros([n,m])

for i in arange(n):
    for j in arange(m):
        Z[i,j] = f(X[i],Y[j])
X, Y = meshgrid(X, Y)
ax.plot_surface(Y, X, Z, rstride=1, cstride=1, cmap=cm.jet)
ax.set_xlabel("X Axis")
ax.set_ylabel("Y Axis")
ax.set_zlabel("Z Axis")
plt.show()

これは正常に実行され、探しているグラフが表示されます。しかし、##This line## を hx=h/2 に変更すると。実行すると、グラフは地獄に落ち、恐ろしく、理解することは不可能です。Y 軸よりも X の方が近いグリッドが必要です。どうすればこれを行うことができますか??

もちろん、これは私が偏微分方程式を解いている例であり、数値的な安定性を得るには、一方の軸がもう一方の軸よりも近いグリッドを持つ必要があります。

4

1 に答える 1

1

次元を反転しました

Z = zeros([m,n])

for i in arange(n):
    for j in arange(m):
        Z[j,i] = f(X[i],Y[j])
X, Y = meshgrid(X, Y)

の任意の比率で機能nしmます。

あなたが持っている機能を使えば、numpyのブロードキャストを使用して、このセクション全体を次のように書くことができます。

X, Y = meshgrid(X, Y)
Z = f(X,Y)

これは読みやすく、高速です。

このコード ブロック全体を次のように書き直します。

from mpl_toolkits.mplot3d import Axes3D
from matplotlib import cm
import matplotlib.pyplot as plt
from numpy import *

def f(x,y):
    r=x**2 + y**2
    return r


n = 5
m = 10
b = 1.
a = -b

fig = plt.figure()
ax = Axes3D(fig)
X = linspace(a,b,n)
Y = linspace(a,b,m)

X, Y = meshgrid(X, Y)
Z = f(X,Y)

ax.plot_surface(Y, X, Z, rstride=1, cstride=1, cmap=cm.jet)
ax.set_xlabel("X Axis")
ax.set_ylabel("Y Axis")
ax.set_zlabel("Z Axis")
plt.show()
于 2012-12-04T04:49:58.197 に答える