Python 一维不匹配两个数组上的广义ufunc

首先要知道的是，您的函数应该在每次调用中多次计算结果.计算次数存储在dimensions[0]中.您指定了gufunc签名(i),(j)->().当调用该函数时，dimensions[1]和dimensions[2]将保持i和j的GUFUC签名(即，核心操作的1-D输入的长度).

例如，在此代码片段中，

x = np.array([[0, 1, 2, 3, 4], [0, 0, 0, 1, 1.5]])  # shape (2, 5)
y = np.array([[-1, -3, 7], [2, 2, 2]])              # shape (2, 3)
result = mean_diff(x, y)                            # shape (2,)

您的函数将被称为once，将dimensions[0]设置为2.dimensions[1]和dimensions[2]将分别为5和3.您的外部循环将计算([0, 1, 2, 3, 4], [-1, -3, 7])对的结果，然后计算([0, 0, 0, 1, 1.5], [2, 2, 2])对的结果.

args[0]和args[1]将指向输入数组的第一个元素，args[2]将指向输出数组的第一个元素.您的代码执行in1 = args[0]、in2 = args[1]和out = args[2]，因此它必须计算in1和in2指向的数据的结果，并将结果存储在*out中.然后，它必须递增in1、in2和out以指向下一组数据并再次进行计算.因此，您将实现一个外部循环，该循环 for each 输入对执行计算，并递增指针以移动到下一组输入和输出.增量in1、in2和out的量分别存储在steps[0]、steps[1]和steps[2]中.(例如，steps[0]是从x的第一行到第二行在内存中跳转的字节数.)

然后，在这个外部循环中，对当前从x到y(由in1和in2指向)的数组对进行计算.通常，您不能假设数组的元素是连续的，因此在计算平均值时，您必须按适当的步长递增指向数据的指针；在本例中，这些步长是steps[3]和steps[4].

这里有一种方法可以实现这一点.我已经注释了我从dimensions和steps中提取的变量，以帮助澄清我上面的简洁描述.

static void mean_diff(
    char **args,
    const npy_intp *dimensions,
    const npy_intp *steps,
    void *extra)
{
    char *in1 = args[0];
    char *in2 = args[1];
    char *out = args[2];

    npy_intp nloops = dimensions[0];  // Number of outer loops
    npy_intp len1 = dimensions[1];    // Core dimension i
    npy_intp len2 = dimensions[2];    // Core dimension j

    npy_intp step1 = steps[0];        // Outer loop step size for the first input
    npy_intp step2 = steps[1];        // Outer loop step size for the second input
    npy_intp step_out = steps[2];     // Outer loop step size for the output
    npy_intp innerstep1 = steps[3];   // Step size of elements within the first input
    npy_intp innerstep2 = steps[4];   // Step size of elements within the second input

    for (npy_intp i = 0; i < nloops;
         i++, in1 += step1, in2 += step2, out += step_out) {

        // The core calculation.  in1 and in2 point to the 1-d input arrays,
        // and out points to the output array.

        double s1 = 0.0;
        double s2 = 0.0;
        npy_intp n1 = 0;
        npy_intp n2 = 0;

        for (npy_intp j = 0; j < len1; ++j) {
            double val = *(double *)(in1 + j*innerstep1);
            if (!isnan(val)) {
                s1 += val;
                n1++;
            }
        }

        for (npy_intp j = 0; j < len2; ++j) {
            double val = *(double *)(in2 + j*innerstep2);
            if (!isnan(val)) {
                s2 += val;
                n2++;
            }
        }

        double mean1 = (n1 > 0) ? s1 / n1 : 0.0;
        double mean2 = (n2 > 0) ? s2 / n2 : 0.0;

        *((double *)out) = mean1 - mean2;
    }
}

顺便说一句，如果像这样实现的gufunc实际上比用mean方法计算平均值并减go 结果要少slower，也不要感到惊讶.NumPy正在不断改进其性能，对于最近的SIMD实现，uuncf和guuncf的朴素实现可能比几个优化的NumPy方法调用的组合要慢.