增长一个numpy数组

Growing a numpy array

我有一个 3x256 规则的规则集。每个规则都映射到一个 3x3 的值网格,而这些值本身就是规则。

示例规则:

0 -> [[0,0,0],[0,1,0],[0,0,0]]  
1 -> [[1,1,1],[0,0,0],[1,1,1]]

种子:

[[0]]

1 次迭代后:

[[0,0,0],
[0,1,0],
[0,0,0]]

2次迭代后:

[[0,0,0,0,0,0,0,0,0],
[0,1,0,0,1,0,0,1,0],
[0,0,0,0,0,0,0,0,0],
[0,0,0,1,1,1,0,0,0],
[0,1,0,0,0,0,0,1,0],
[0,0,0,1,1,1,0,0,0],
[0,0,0,0,0,0,0,0,0],
[0,1,0,0,1,0,0,1,0],
[0,0,0,0,0,0,0,0,0]]

现在我有了一个可用的实现,但是,它是我的脚本中最慢的函数。我想知道是否有更 pythonic 和更有效的方法来重写这个函数。

def decode(rules,fractal_iterations,seed):
    final_seed_matrix = np.zeros((3,3**fractal_iterations,3**fractal_iterations))

    for i in range(dimensions):
        seed_matrix = np.array([[seed]])
        for j in range(fractal_iterations):
            size_y = seed_matrix.shape[0]
            size_x = seed_matrix.shape[1]
            new_matrix = np.zeros((size_y*rule_size_sqrt,size_x*rule_size_sqrt))

            for y in range(size_y):
                for x in range(size_x):
                    seed_value = seed_matrix[y,x]
                    new_matrix[y*rule_size_sqrt : y*rule_size_sqrt+rule_size_sqrt, x*rule_size_sqrt : x*rule_size_sqrt+rule_size_sqrt] = rules[int(seed_value),i]

            seed_matrix = new_matrix
        final_seed_matrix[i] = seed_matrix

    return np.moveaxis(final_seed_matrix,0,-1)

这是一个优化版本,它使用高级索引 select 并在一个索引步骤中将所有规则拼凑在一起。这会在它替换的像素位置创建一个具有适当规则的 4D 数组。将其展平为 2D 就是交换中轴和重塑的问题。它似乎给出了与您相同的结果,但速度明显更快(到目前为止仅针对整数规则进行了测试):

results equal: True
OP       : 24.883304461836815 ms
optimized: 1.093490980565548 ms

代码:

import numpy as np

dimensions = 3
rule_size_sqrt = 3

def decode(rules,fractal_iterations,seed):
    final_seed_matrix = np.zeros((3,3**fractal_iterations,3**fractal_iterations))

    for i in range(dimensions):
        seed_matrix = np.array([[seed]])
        for j in range(fractal_iterations):
            size_y = seed_matrix.shape[0]
            size_x = seed_matrix.shape[1]
            new_matrix = np.zeros((size_y*rule_size_sqrt,size_x*rule_size_sqrt))

            for y in range(size_y):
                for x in range(size_x):
                    seed_value = seed_matrix[y,x]
                    new_matrix[y*rule_size_sqrt : y*rule_size_sqrt+rule_size_sqrt, x*rule_size_sqrt : x*rule_size_sqrt+rule_size_sqrt] = rules[int(seed_value),i]

            seed_matrix = new_matrix
        final_seed_matrix[i] = seed_matrix

    return np.moveaxis(final_seed_matrix,0,-1)

def decode_fast(rules, fractal_iterations, seed):
    rules_int = rules.astype(int)
    seed = np.array([[seed]], dtype=int)
    res = np.empty((3**fractal_iterations, 3**fractal_iterations, dimensions),
                   dtype=rules.dtype)
    for i in range(dimensions):
        grow = seed
        for j in range(1, fractal_iterations):
            grow = rules_int[grow, i].swapaxes(1, 2).reshape(3**j, -1)
        grow = rules[grow, i].swapaxes(1, 2).reshape(3**fractal_iterations, -1)
        res[..., i] = grow
    return res

rules = np.random.randint(0, 4, (4, dimensions, 3, 3))
seed = 1
fractal_iterations = 5
print('results equal:', np.all(decode(rules, fractal_iterations, seed) == decode_fast(rules, fractal_iterations, seed)))

from timeit import repeat
print('OP       :', min(repeat('decode(rules, fractal_iterations, seed)', globals=globals(), number=50))*20, 'ms')
print('optimized:', min(repeat('decode_fast(rules, fractal_iterations, seed)', globals=globals(), number=50))*20, 'ms')