使用 Numpy 和 Pyplot 进行条件绘图
Conditional Plotting with Numpy and Pyplot
我正在尝试绘制一个有条件定义的函数。具体来说:
U(x) = (2**delta)/((D-d)**delta)*(D/2 - (x-x0))**delta ,对于 abs(x-x0) 小于 D/2 否则为 0。
但我的问题是我想将 x、x0 作为 numpy 数组,因为这就是我在其余实际代码中使用它们的方式。
我设置了以下示例:
import numpy as np
import matplotlib.pyplot as plt
AD = 0.2
D = 0.4
delta = 8
def Parabolic(x, delta, D, AD):
x0 = np.round(x)
tempx = np.abs(x-x0)
tempD = D/2*np.ones(len(x))
if tempx<tempD:
return ((2**delta)/(D-AD)**delta)*(D/2 - (x-x0))**delta
else:
return 0
figure = plt.figure(figsize=(10,8), dpi=72)
xmin = -1.0
xmax = 1.0
X = np.linspace(xmin,xmax,1000)
plt.plot(X, Parabolic(X, delta=8, D=0.4, AD=0.2))
显然这个例子不起作用,因为行 tempx<tempD
引发了列表的真值不明确的错误。
我搜索了 numpy 的文档并找到了函数 np.less(tempx, tempD)。但是,如果我用 np.less(tempx, tempD)
替换 tempx < tempD
它仍然不起作用,因为我再次要求整个列表的真值。我知道问题不在于numpy,而是我无法理解如何使用numpy提供的逻辑函数。
很抱歉,如果这在另一个 post 中以某种方式回答,我在这个论坛中搜索但找不到除了 curve()
方法之外的其他方法。但是我想保留我的 numpy.array 格式以用于我的实际代码。我敢打赌,答案一定很简单,我就是想不出来。
试试这个使用 numpy 逻辑数组的方法:
import numpy as np
import matplotlib.pyplot as plt
AD = 0.2
D = 0.4
delta = 8
def Parabolic(x, delta, D, AD):
rtn_arr = np.zeros(len(x))
x0 = np.round(x)
tempx = np.abs(x-x0)
tempD = D/2*np.ones(len(x))
lgc_arr = tempx<tempD
x_cut = x[lgc_arr]
x0_cut = x0[lgc_arr]
rtn_arr[lgc_arr] = ((2**delta)/(D-AD)**delta)*(D/2 - (x_cut-x0_cut))**delta
return rtn_arr
figure = plt.figure(figsize=(10,8), dpi=72)
xmin = -1.0
xmax = 1.0
X = np.linspace(xmin,xmax,1000)
plt.plot(X, Parabolic(X, delta=8, D=0.4, AD=0.2))
Parabolic 必须是一个 ufunc,所以你不能把 python 测试放在你的代码中。
一个简单的解决方法是:
def Parabolic(x, delta, D, AD):
x0 = np.round(x)
tempx = np.abs(x-x0)
tempD = D/2*np.ones(len(x))
u=(((2**delta)/(D-AD)**delta)*(D/2 - (x-x0))**delta)
u[tempx>=tempD]=0
return u
或者为了避免不必要的计算:
def Parabolic2(x, delta, D, AD):
x0 = np.round(x)
tempx = np.abs(x-x0)
tempD = D/2*np.ones(len(x))
u= zeros_like(x)
valid=tempx<tempD
u[valid]=(((2**delta)/(D-AD)**delta)*(D/2 - (x-x0)[valid])**delta)
return u
第二个稍微快一点:
In [141]: %timeit Parabolic(x,8,.4,.2)
1000 loops, best of 3: 310 µs per loop
In [142]: %timeit Parabolic2(x,8,.4,.2)
1000 loops, best of 3: 218 µs per loop
我正在尝试绘制一个有条件定义的函数。具体来说: U(x) = (2**delta)/((D-d)**delta)*(D/2 - (x-x0))**delta ,对于 abs(x-x0) 小于 D/2 否则为 0。
但我的问题是我想将 x、x0 作为 numpy 数组,因为这就是我在其余实际代码中使用它们的方式。
我设置了以下示例:
import numpy as np
import matplotlib.pyplot as plt
AD = 0.2
D = 0.4
delta = 8
def Parabolic(x, delta, D, AD):
x0 = np.round(x)
tempx = np.abs(x-x0)
tempD = D/2*np.ones(len(x))
if tempx<tempD:
return ((2**delta)/(D-AD)**delta)*(D/2 - (x-x0))**delta
else:
return 0
figure = plt.figure(figsize=(10,8), dpi=72)
xmin = -1.0
xmax = 1.0
X = np.linspace(xmin,xmax,1000)
plt.plot(X, Parabolic(X, delta=8, D=0.4, AD=0.2))
显然这个例子不起作用,因为行 tempx<tempD
引发了列表的真值不明确的错误。
我搜索了 numpy 的文档并找到了函数 np.less(tempx, tempD)。但是,如果我用 np.less(tempx, tempD)
替换 tempx < tempD
它仍然不起作用,因为我再次要求整个列表的真值。我知道问题不在于numpy,而是我无法理解如何使用numpy提供的逻辑函数。
很抱歉,如果这在另一个 post 中以某种方式回答,我在这个论坛中搜索但找不到除了 curve()
方法之外的其他方法。但是我想保留我的 numpy.array 格式以用于我的实际代码。我敢打赌,答案一定很简单,我就是想不出来。
试试这个使用 numpy 逻辑数组的方法:
import numpy as np
import matplotlib.pyplot as plt
AD = 0.2
D = 0.4
delta = 8
def Parabolic(x, delta, D, AD):
rtn_arr = np.zeros(len(x))
x0 = np.round(x)
tempx = np.abs(x-x0)
tempD = D/2*np.ones(len(x))
lgc_arr = tempx<tempD
x_cut = x[lgc_arr]
x0_cut = x0[lgc_arr]
rtn_arr[lgc_arr] = ((2**delta)/(D-AD)**delta)*(D/2 - (x_cut-x0_cut))**delta
return rtn_arr
figure = plt.figure(figsize=(10,8), dpi=72)
xmin = -1.0
xmax = 1.0
X = np.linspace(xmin,xmax,1000)
plt.plot(X, Parabolic(X, delta=8, D=0.4, AD=0.2))
Parabolic 必须是一个 ufunc,所以你不能把 python 测试放在你的代码中。
一个简单的解决方法是:
def Parabolic(x, delta, D, AD):
x0 = np.round(x)
tempx = np.abs(x-x0)
tempD = D/2*np.ones(len(x))
u=(((2**delta)/(D-AD)**delta)*(D/2 - (x-x0))**delta)
u[tempx>=tempD]=0
return u
或者为了避免不必要的计算:
def Parabolic2(x, delta, D, AD):
x0 = np.round(x)
tempx = np.abs(x-x0)
tempD = D/2*np.ones(len(x))
u= zeros_like(x)
valid=tempx<tempD
u[valid]=(((2**delta)/(D-AD)**delta)*(D/2 - (x-x0)[valid])**delta)
return u
第二个稍微快一点:
In [141]: %timeit Parabolic(x,8,.4,.2)
1000 loops, best of 3: 310 µs per loop
In [142]: %timeit Parabolic2(x,8,.4,.2)
1000 loops, best of 3: 218 µs per loop