我有几个变量x，y，z

# Variables
x = np.linspace(-1.0, 1.0, num=5, endpoint=True)
y = np.linspace(-1.0, 1.0, num=5, endpoint=True)
z = np.linspace(-1.0, 1.0, num=5, endpoint=True)

# And the corresponding meshgrid
x_v, y_v, z_v = np.meshgrid(x, y, z)

字符串
假设我想计算一个二次多项式的值，其系数为C，在网格的每个点上。为此，我有一个函数来返回一组参数的多项式项（注：这里的函数是通用的，使其适用于任何数量的参数）：

def pol_terms(*args):
    # -- Calculate terms of the polynomial --
    # Linear part (b0, b11*x1, b12*x2, ..., b1k*xk, where k=len(args))
    entries = [1] + [pt for pt in args]
    
    # Combinations part (b12*x12, b13*x13, ..., bk-1k*xk-1k)
    n_args = len(args)
    for i in range(n_args):
        for j in range(i+1, n_args):
            entries.append(args[i]*args[j])
    
    # Quadratic part
    entries += [pt**2 for pt in args]
    
    return np.array([entries])

型
因此，对于一组变量的值，我可以得到我的多项式的值如下：

X = pol_terms(1,2,3)
X.dot(C)

型
我不知道如何为meshgrid x_v，y_v，z_v做同样的事情。在pol_terms上应用np.vectorize将不起作用，因为它返回一个数组。我也不想手动遍历所有维度，因为我希望解决方案是通用的（因此它适用于2变量和5变量）。
提前感谢！

我同意你不能使用np.vectorize，但不是因为你给出的原因。
np.vectorize可以很好地与返回数组的函数一起使用。只是，你一定要说出来。

ptv = np.vectorize(pol_terms, signature='(),(),()->(m)')
ptv([1,10],[2,20],[3,30])

字符串
退货

array([[  1,   1,   2,   3,   2,   3,   6,   1,   4,   9],
       [  1,  10,  20,  30, 200, 300, 600, 100, 400, 900]])

型
不能使用vectorize的原因是，首先，它并不比自己循环好多少（正如文档中所述，vectorize不是为了性能）。第二，你可以在我的例子中看到：它需要恒定数量的参数。至少如果您要指定签名
但是您的代码大部分都已经矢量化了。你所做的一切就是加法和乘法。
如果使用数组（我指的是非0维数组，因为numpy中的scalar只是0维数组），我能看到的唯一问题是1项。该项意味着所有其他项的维数相同（因此维数为0）。
但那是很容易解决的。
替换行

entries = [1] + [pt for pt in args]

型
带线

entries = [np.ones_like(args[0])] + [pt for pt in args]

型
（我假设至少有一个arg。如果没有，我让你添加特殊的“无参数”情况所需的测试）
然后，

pol_terms(np.array([1,10]),np.array([2,20]),np.array([3,30]))

型
退货

array([[[  1,   1],
        [  1,  10],
        [  2,  20],
        [  3,  30],
        [  2, 200],
        [  3, 300],
        [  6, 600],
        [  1, 100],
        [  4, 400],
        [  9, 900]]])

型
(Note，即额外的括号是由于您返回的np.array([entries])而不是np.array(entries)。我不知道这是不是心甘情愿的。但这意味着，无论条目的形状如何，返回值的形状都是1×形状。所以在你的返回中总是有一个大小为1的轴0）
我不会复制pol_terms(x_v,y_v,z_v)的结果。但是（对于那个奇怪的第一轴）他们似乎是我所期望的。
比如说

pt=pol_terms(x_v,y_v,z_v)
c123=pt[0,:,1,2,3]

型
是位置（x[1],x[2],x[3]）处的项。0坐标是因为我提到的那个额外的轴。
所以，除非我不明白你想做什么，否则几乎没有什么可做的。您的函数已经能够处理向量。但是对于硬编码的标量1。

一行字

def pol_terms(*args):
    entries = np.array([np.ones_like(args[0])] + [pt for pt in args])
    return (entries[None,...]*entries[:,None,...])[np.triu_indices(len(entries))]

型
但它只是更紧凑，而不是更快。您的方法已经矢量化了真正重要的东西（网格的轴）。这样做所增加的只是最外层轴的矢量化，即args本身的一个轴（您的2个嵌套for循环，它只是在参数上迭代，而不是在这些参数的轴上迭代）
所以，为了一个毫无意义的矢量化，我付出了代价：我计算了x*y和y*x，然后放弃了几乎一半的计算（只保留x*y，而忽略y*x，因为它是相同的）
因此，从不对称的Angular 来看，这个版本应该比你的版本慢两倍左右（建议修正）。这只是为了好玩的准一行。

numpy 网格上的多项式函数

1条答案

一行字

相关问题

热门标签

最新问答