Векторизованная функция, возвращающая список векторов в Python

Question

Я пытаюсь написать функцию в Python, которая принимает список параметров a в качестве ввода и возвращает список векторов u, которые зависят от этих параметров. Затем я хотел бы иметь возможность действовать по u со стандартными векторными операциями (такими как скалярное произведение). В приведенном ниже (вдохновленный от второго ответа this post) Например, функция возвращает комбинацию двух исходных векторов v и w:

import numpy as np 
v = np.array([1,0,0]) 
w = np.array([0,1,0]) 
a = np.linspace(0,np.pi,4) 
def u(a): 
    u = np.cos(a)*v + np.sin(a)*w 
    return u 
uvec = np.vectorize(u, otypes=[np.ndarray])  
ufin = np.array(uvec(a).tolist()) 
print "ufin =", ufin 
print "ufin.v =", np.dot(ufin,v) 

Это возвращает:

ufin = [[ 1.00000000e+00 0.00000000e+00 0.00000000e+00] 
[ 5.00000000e-01 8.66025404e-01 0.00000000e+00] 
[ -5.00000000e-01 8.66025404e-01 0.00000000e+00] 
[ -1.00000000e+00 1.22464680e-16 0.00000000e+00]] 
ufin.v = [ 1. 0.5 -0.5 -1. ] 

, который я хотел получить с ufin ведет себя как вектор. Не могли бы вы рассказать мне, есть ли другие более простые способы реализовать это? Мне нужно написать код, в котором должны быть определены многочисленные векторы и векторные операции, и хотелось бы, чтобы он был как можно более компактным.

Спасибо заранее!

---

Edit:

Я нашел другой (по-видимому, более компактное) решение, основанное на последнем ответе на this post. Идея состоит в том, чтобы изменить список параметров в массив столбцов, так что вывод функции автоматически (без необходимости векторизации) возвращает список векторов как двухмерный массив. Это делается таким образом:

import numpy as np 
from numpy.core.umath_tests import inner1d 

v = np.array([1,0,0]) 
w = np.array([0,1,0]) 
a = np.linspace(0,np.pi,4).reshape((4,1)) 
b = np.linspace(0,np.pi/2,4).reshape((4,1)) 
def u(a): 
    u = np.cos(a)*v + np.sin(a)*w 
    return u 
print "u(a) =",u(a) 
print "u(b) =",u(b) 
print "u(a).v =",np.dot(u(a),v) 
print "u(a)^v =",np.cross(u(a),v) 
# print "u(a).u(b) =",np.dot(u(a),u(b)) # does not work 
print "u(a).u(b) =",inner1d(u(a),u(b)) # works 
print "u(a)^u(b) =",np.cross(u(a),u(b)) 

Это возвращает:

u(a) = [[ 1.00000000e+00 0.00000000e+00 0.00000000e+00] 
[ 5.00000000e-01 8.66025404e-01 0.00000000e+00] 
[ -5.00000000e-01 8.66025404e-01 0.00000000e+00] 
[ -1.00000000e+00 1.22464680e-16 0.00000000e+00]] 
u(b) = [[ 1.00000000e+00 0.00000000e+00 0.00000000e+00] 
[ 8.66025404e-01 5.00000000e-01 0.00000000e+00] 
[ 5.00000000e-01 8.66025404e-01 0.00000000e+00] 
[ 6.12323400e-17 1.00000000e+00 0.00000000e+00]] 
u(a).v = [ 1. 0.5 -0.5 -1. ] 
u(a)^v = [[ 0.00000000e+00 0.00000000e+00 0.00000000e+00] 
[ 0.00000000e+00 0.00000000e+00 -8.66025404e-01] 
[ 0.00000000e+00 0.00000000e+00 -8.66025404e-01] 
[ 0.00000000e+00 0.00000000e+00 -1.22464680e-16]] 
u(a).u(b) = [ 1.00000000e+00 8.66025404e-01 5.00000000e-01 6.12323400e-17] 
u(a)^u(b) = [[ 0.   0.   0.  ] 
[ 0.   0.  -0.5  ] 
[ 0.   0.  -0.8660254] 
[ 0.   0.  -1.  ]] 

, который является правильным поведением как для операций, связанных с начальными (u и v) и выходные списки векторов (u(a) и u(b)) и для операций с двумя выходными списками векторов. Единственное предостережение (для операций между выходным списком векторов) состоит в том, что вместо стандартного np.dot следует использовать функцию inner1d, потому что последняя интерпретируется как матричный продукт, который не может быть выполнен, поскольку две матрицы имеют несогласованные размеры.

Answer 1

Если вы не возражаете против получения 2d массива вместо массива объектов (бывший является более удобной, чем последняя для большинства применений, во всяком случае), то внешнего продукт полезен для этой задачи :

>>> np.outer(np.cos(a), v) + np.outer(np.sin(a), w) 
array([[ 1.00000000e+00, 0.00000000e+00, 0.00000000e+00], 
     [ 5.00000000e-01, 8.66025404e-01, 0.00000000e+00], 
     [ -5.00000000e-01, 8.66025404e-01, 0.00000000e+00], 
     [ -1.00000000e+00, 1.22464680e-16, 0.00000000e+00]]) 
>>> _ @ v 
array([ 1. , 0.5, -0.5, -1. ]) 

оператор @ инфикс доступен только на питона 3.5 или более поздней версии. В более старых версиях вы можете использовать np.dot. Смысл либо в этом контексте - matrix multiplication.

Answer 2

Итак, вы сообщаете vectorze, что функция u возвращает массив, или фактически object, так как нет dtype=array.

In [475]: uvec=np.vectorize(u, otypes=[object]) 
In [476]: uvec(a) 
Out[476]: 
array([array([ 1., 0., 0.]), array([ 0.5  , 0.8660254, 0.  ]), 
     array([-0.5  , 0.8660254, 0.  ]), 
     array([ -1.00000000e+00, 1.22464680e-16, 0.00000000e+00])], dtype=object) 

Так что, если вы хотите, чтобы «разгадать», что и создать 2d массив, который вы должны преобразовать его в список массивов, а затем использовать np.array или (np.stack), чтобы присоединиться к ним, а также новое измерение:

In [477]: uvec(a).tolist() 
Out[477]: 
[array([ 1., 0., 0.]), 
array([ 0.5  , 0.8660254, 0.  ]), 
array([-0.5  , 0.8660254, 0.  ]), 
array([ -1.00000000e+00, 1.22464680e-16, 0.00000000e+00])] 

Но для итерации на одном массиве 1d, таком как a, vectorize является чрезмерным. Простой список понимание работает так же хорошо

In [483]: np.array([u(i) for i in a]) 
Out[483]: 
array([[ 1.00000000e+00, 0.00000000e+00, 0.00000000e+00], 
     [ 5.00000000e-01, 8.66025404e-01, 0.00000000e+00], 
     [ -5.00000000e-01, 8.66025404e-01, 0.00000000e+00], 
     [ -1.00000000e+00, 1.22464680e-16, 0.00000000e+00]] 

vectorize является более полезным, когда существует несколько входных массивов, и они могут иметь несколько измерений. Он заботится о трансляции и повторении нескольких уровней. Вы могли бы также использовать np.frompyfunc(u,1,1), который всегда возвращает массив объектов, и немного быстрее.

Но вы также можете вычислить это без зацикливания:

def ufunc(a,v,w): 
    a = a[:,None] 
    return np.cos(a)*v + np.sin(a)*w 

In [498]: ufunc(a,v,w) 
Out[498]: 
array([[ 1.00000000e+00, 0.00000000e+00, 0.00000000e+00], 
     [ 5.00000000e-01, 8.66025404e-01, 0.00000000e+00], 
     [ -5.00000000e-01, 8.66025404e-01, 0.00000000e+00], 
     [ -1.00000000e+00, 1.22464680e-16, 0.00000000e+00]]) 

Это превращает a в (4,1) массив, который транслирует против (3,)v произвести (4,3) массив.

---

Другой подход заключается в определении функции, которая принимает скаляры (все вокруг)

import math 
def u1(a,v,w): 
    return math.cos(a)*v+math.sin(a)*w 
u1v = np.vectorize(u1, otypes=[float]) 

In [505]: u1v(a[:,None], v, w) 
Out[505]: 
array([[ 1.00000000e+00, 0.00000000e+00, 0.00000000e+00], 
     [ 5.00000000e-01, 8.66025404e-01, 0.00000000e+00], 
     [ -5.00000000e-01, 8.66025404e-01, 0.00000000e+00], 
     [ -1.00000000e+00, 1.22464680e-16, 0.00000000e+00]]) 

Это имеет преимущество полного вещания способности vectorize. Я, хотя это будет медленнее, чем ваш код, потому что он должен зацикливаться в обоих измерениях, но он на самом деле проверяет быстрее (хотя и не так быстро, как мой ufunc). Пытаясь провести половину векторизации функции u, вы не сможете купить вам какую-либо скорость.

---

Вы могли бы также определили свою первоначальную функцию, чтобы взять v и w в качестве аргументов, и использовать параметр vectorizeexcluded передать их как есть. Я не думаю, что у него есть преимущество в производительности.

In [521]: def u2(a,v,w): 
    ...:  return np.cos(a)*v + np.sin(a)*w 
    ...: 
In [522]: u2vec=np.vectorize(u2, otypes=[object],excluded=[1,2]) 
In [523]: u2vec(a,v,w) 
Out[523]: 
array([array([ 1., 0., 0.]), array([ 0.5  , 0.8660254, 0.  ]), 
     array([-0.5  , 0.8660254, 0.  ]), 
     array([ -1.00000000e+00, 1.22464680e-16, 0.00000000e+00])], dtype=object)