加速我的Python/Cython代码的提示

Question

我已经编写了Python程序并将其化。Cython获得的加速（30％）不满意。肯定有通过更改代码结构或含糊的方式来优化它的空间。

该程序基本上是采用数字高程模型（DEM）栅格图和具有相同形状的多余水图。对于多余的水图中的每个像素，它搜索四个相邻的像素，并确定像素是否低于周围的邻居，具有相同的水平或高于它们。基于此，它增加了像素的水位，或者将多余的水分配给海拔较低的邻居。该代码一直持续到所有多余的水都分布在地面上。这是代码的Cython版本。

import numpy as np
cimport numpy as np
cdef unravel(np.ndarray[np.int_t, ndim = 1] idx,int shape0, int shape1):
    return idx//shape0, idx%shape1
cdef find_lower_neighbours(int i, int j, np.ndarray[np.double_t, ndim=2] water_level, double friction_head_loss):
    cdef double current_water_el, minlevels, deltav_total, deltav_min
    current_water_el = water_level[i,j]
    cdef np.ndarray[np.double_t, ndim = 1] levels = np.zeros(4, dtype = np.double)
    levels[:] = water_level[i - 1, j], water_level[i, j + 1], water_level[i + 1, j], water_level[i, j - 1]
    minlevels = levels.min()
    if current_water_el - minlevels < 0:
        return 0, minlevels, 0
    elif np.absolute(current_water_el - minlevels) < 0.0001:
        res = np.where(np.absolute(levels - current_water_el) < 0.0001)
        return 1, res[0], 0
    else:
        levels = current_water_el - levels
        low_values_flags = levels < 0
        levels[low_values_flags] = 0
        deltav_total = np.sum(levels)
        deltav_min = levels[levels > 0].min()
        return 2, levels / (deltav_total + deltav_min), deltav_min / (deltav_total + deltav_min)

cpdef np.ndarray[np.double_t, ndim=2] new_algorithm( np.ndarray[np.double_t, ndim=2] DEM, np.ndarray[np.double_t, ndim=2] extra_volume_map, double nodata, double pixel_area, double friction_head_loss , int von_neuman):
    cdef int terminate = 0
    cdef int iteration = 1
    cdef double sum_extra_volume_map

    index_dic_von_neuman = [[-1, 0], [0, 1], [1, 0], [0, -1]]
    cdef int DEMshape0 = DEM.shape[0]
    cdef int DEMshape1 = DEM.shape[1]
    cdef np.ndarray[np.double_t, ndim = 2] water_levels
    water_levels = np.copy(DEM)
    cdef np.ndarray[np.double_t, ndim = 2] temp_water_levels
    temp_water_levels = np.copy(water_levels)
    cdef np.ndarray[np.double_t, ndim = 2] temp_extra_volume_map
    temp_extra_volume_map = np.copy(extra_volume_map)

    cdef np.ndarray[np.int_t, ndim = 2] wetcells
    wetcells = np.zeros((DEMshape0,DEMshape1), dtype= np.int)
    wetcells[extra_volume_map > 0] = 1
    cdef np.ndarray[np.int_t, ndim = 2] temp_wetcells
    temp_wetcells = np.zeros((DEMshape0,DEMshape1), dtype= np.int )
    cdef double min_excess = friction_head_loss * pixel_area
    cdef np.ndarray[np.int_t, ndim = 1] fdx
    cdef int i, j , k, condition, have_any_dry_cells_in_neghbors
    cdef double water_level_difference, n , volume_to_each_neghbour, weight, w0
    if von_neuman == 0:
        index_dic = index_dic_von_neuman

    while terminate != 1:
        fdx = np.flatnonzero(extra_volume_map > min_excess)
        extra_volume_locations = unravel(fdx, DEMshape1,DEMshape1)
        if not extra_volume_locations[0].size:
            terminate = 1
            return water_levels
        for item in zip(*extra_volume_locations):
            i = item[0]
            j = item[1]
            if DEM[i,j] == nodata:
                print "warningggggg", i, j
                temp_extra_volume_map[i,j] = 0.
            else:
                condition, wi, w0 = find_lower_neighbours(i, j, water_levels, friction_head_loss)
                if condition == 0:
                    water_level_difference = wi - water_levels[i,j]
                    temp_water_levels[i,j] = wi
                    temp_extra_volume_map[i,j] -= water_level_difference * pixel_area
                elif condition == 1:
                    n = len(wi)
                    volume_to_each_neghbour = (extra_volume_map[i, j] - friction_head_loss * pixel_area * .02)/ (n * 1. )
                    for itemm in wi:
                        temp_extra_volume_map[i + index_dic[itemm][0], j + index_dic[itemm][1]] += volume_to_each_neghbour
                        temp_wetcells[i + index_dic[itemm][0], j + index_dic[itemm][1]] += 1
                    temp_extra_volume_map[i,j] -= extra_volume_map[i,j]
                    temp_water_levels[i,j] += friction_head_loss * .02
                elif condition == 2:
                    temp_wetcells[i,j] += 1
                    have_any_dry_cells_in_neghbors = 0 # means "no"
                    for k, weight in enumerate(wi):
                        if weight > 0:
                            temp_extra_volume_map[i + index_dic[k][0], j + index_dic[k][1]] += weight * (extra_volume_map[i,j])
                            if wetcells[i + index_dic[k][0], j + index_dic[k][1]] == 0:
                                have_any_dry_cells_in_neghbors = 1
                            temp_wetcells[i + index_dic[k][0], j + index_dic[k][1]] += 1
                    if have_any_dry_cells_in_neghbors == 1:
                        temp_water_levels[i,j] += friction_head_loss
                        temp_extra_volume_map[i, j] -= (1. - w0) * (extra_volume_map[i, j])
                    else:
                        temp_extra_volume_map[i, j] -= (1. - w0) * (extra_volume_map[i, j])

        wetcells = np.copy(temp_wetcells)
        water_levels = np.copy(temp_water_levels)
        extra_volume_map = np.copy(temp_extra_volume_map)
        iteration += 1
        if iteration*1. %500. == 0.:
            sum_extra_volume_map = np.sum(extra_volume_map)
            print "iteration", iteration, "volume left =",  sum_extra_volume_map

    print "finished at iteration =", iteration - 1
    return water_levels

这是分析的结果：

6712150 function calls in 29.005 seconds
Ordered by: standard name
ncalls  tottime  percall  cumtime  percall filename:lineno(function)
      1    0.000    0.000   29.005   29.005 <string>:1(<module>)
2017444    0.497    0.000    4.414    0.000 _methods.py:28(_amin)
289757    0.070    0.000    0.543    0.000 _methods.py:31(_sum)
289757    0.437    0.000    1.098    0.000 fromnumeric.py:1730(sum)
506055    0.203    0.000    0.705    0.000 function_base.py:1453(copy)
168685    0.140    0.000    0.300    0.000 numeric.py:859(flatnonzero)
     1    0.000    0.000    0.000    0.000 test.py:22(print_dem_for_excel)
     1    0.000    0.000   29.005   29.005 test.py:27(test)
289757    0.119    0.000    0.119    0.000 {isinstance}
     1    0.000    0.000    0.000    0.000 {len}
    20    0.000    0.000    0.000    0.000 {map}
    20    0.000    0.000    0.000    0.000 {method 'append' of 'list' objects}
     1    0.000    0.000    0.000    0.000 {method 'disable' of '_lsprof.Profiler' objects}
    20    0.000    0.000    0.000    0.000 {method 'join' of 'str' objects}
168685    0.106    0.000    0.106    0.000 {method 'nonzero' of 'numpy.ndarray' objects}
168685    0.054    0.000    0.054    0.000 {method 'ravel' of 'numpy.ndarray' objects}
2307201    4.390    0.000    4.390    0.000 {method 'reduce' of 'numpy.ufunc' objects}
 506056    0.501    0.000    0.501    0.000 {numpy.core.multiarray.array}
     1    0.000    0.000    0.000    0.000 {numpy.core.multiarray.zeros}
     1    0.000    0.000    0.000    0.000 {time.time}
     1   22.488   22.488   29.005   29.005 {version3Cython.new_algorithm}

Answer 1

您首先需要配置您的程序并精确地测量需要时间的时间（哪些功能在CPU时间最多消耗）。

您可以改善程序中的算法并降低某些功能的时间复杂性。

您也可以手动重写代码的重量（例如CPU消耗）部分，例如如（手写）C中的Python扩展（也许还具有C 中的某些零件，或CO易于从C中调用的任何语言，并且与Python和C的内存模型兼容）。不要指望任何工具（例如Cython）会自动有效地执行此操作。

您的程序足够小。考虑花几个星期（或几个月）在重写其中的某些部分，甚至全部重写和重写。我想静态键入的编译语言（例如OCAML，GO，D，C ）可以提供一些性能的改进。也许您可以将代码完全重新设计为并发应用程序（多线程或opencl，mpi，....）

将现有的Python原型视为开始理解您要解决的问题的一种方式，而是重新设计和重写代码（可能是用不同的语言）。还花一些时间阅读有关该主题的现有文献并与专家讨论。

btw，您的用户可以认为可以将执行时间少于一分钟。那么，您无需重写代码。请记住，您的开发时间也有一些费用！没有银色子弹 ...