如何寻找可剔除最多k个连续行区域的最优圆形区域

DevOps 小满

2026-7-1

问题描述

我有一个n×n的整数矩阵，要找到以左上角为圆心的圆形区域，使得区域内元素的和最大。难点在于最多可以剔除k个连续行块（k≤5），需要同时找到最优半径和最优的行块剔除方案，且代码要高效支持大规模矩阵。

示例说明

基础示例（无行剔除）

下面是带圆形覆盖的网格示例：

生成该示例的代码：

import matplotlib.pyplot as plt
from matplotlib.patches import Circle
import numpy as np
plt.yticks(np.arange(0, 10.01, 1))
plt.xticks(np.arange(0, 10.01, 1))
plt.xlim(0,9)
plt.ylim(0,9)
plt.gca().invert_yaxis()
# 设置等比例
plt.gca().set_aspect('equal', adjustable='box')
plt.grid()   
np.random.seed(40)
square = np.empty((10, 10), dtype=np.int_)
for x in np.arange(0, 10, 1):
    for y in np.arange(0, 10, 1):
        plt.scatter(x, y, color='blue', s=2, zorder=2, clip_on=False) 
for x in np.arange(0, 10, 1):
    for y in np.arange(0, 10, 1):
        value = np.random.randint(-3, 4)
        square[int(x), int(y)] = value
        plt.text(x-0.2, y-0.2, str(value), ha='center', va='center', fontsize=8, color='black')

r1 = 3
circle1 = Circle((0, 0), r1,  color="blue", alpha=0.5, ec='k', lw=1)
plt.gca().add_patch(circle1)

对应的矩阵：

[[ 3, -1,  1,  0, -1, -1, -3, -2, -2,  2],
 [ 0,  0,  3,  0,  0, -1,  2,  0, -2,  3],
 [ 2,  0,  3, -2,  3,  1,  2,  2,  1,  1],
 [-3,  0,  1,  0,  1,  2,  3,  1, -3, -1],
 [-3, -2,  1,  2,  1, -3, -2,  2, -2,  0],
 [-1, -3, -3,  1,  3, -2,  0,  2, -1,  1],
 [-2, -2, -1,  2, -2,  1, -1,  1,  3, -1],
 [ 1,  2, -1,  2,  0, -2, -1, -1,  2,  3],
 [-1, -2,  3, -1,  0,  0,  3, -3,  3, -2],
 [ 0, -3,  0, -1, -1,  0, -2, -3, -3, -1]]

当不可剔除任何区域时，最优半径为2√2=√8，圆形内元素和为11，这是最大值（半径增大后总和不会超过11）。计算该结果的代码：

import numpy as np
def make_data(N):
    np.random.seed(40)
    g = np.random.randint(-3, 4, (N, N))
    return g

def find_max(g):
    n = g.shape[0]
    sum_dist = np.zeros(2 * N * N, dtype=np.int32)
    for i in range(n):
        for j in range(n):
            dist = i**2 + j**2
            sum_dist[dist] += g[i, j]
    cusum = np.cumsum(sum_dist)
    return np.argmax(cusum), np.max(cusum)

N = 10
g = make_data(N)
g = g.T # 匹配示例图的矩阵方向
print(g)  

squared_dist, score = find_max(g)
print(np.sqrt(squared_dist), score)

带行剔除的示例

若剔除第0行以及第4-6行（共2个连续行块，k=2），最优半径增大至√58，对应元素和为24：
带行剔除的圆形覆盖网格

高效解决方案

核心思路

预处理事件点：将每行中元素的平方距离分界点（即当半径增大到某个值时，该行的贡献会新增元素）整理为事件，按距离升序排序，确保我们只在贡献变化的节点计算总和。
线性DP计算最优和：对每个距离节点对应的行贡献数组，用O(nk)的动态规划计算最多剔除k个连续行块的最大总和，利用k≤5的限制保证计算效率。
全局追踪最优值：遍历所有距离节点，记录最大总和对应的半径和剔除方案。

实现代码

import numpy as np
from itertools import chain

def precompute_event_points(g):
    """预处理所有行的贡献变化事件点，每个事件点是(平方距离d, 行号i, 增量)"""
    n = g.shape[0]
    events = []
    for i in range(n):
        current_sum = 0
        # 按j从小到大遍历，计算每个d对应的增量
        for j in range(n):
            d = i*i + j*j
            val = g[i][j]
            events.append((d, i, val))
    # 按d升序排序事件点
    events.sort(key=lambda x: x[0])
    # 合并相同d的事件，计算每行的总增量
    merged_events = []
    prev_d = None
    delta = np.zeros(n, dtype=np.int32)
    for d, i, val in events:
        if d != prev_d:
            if prev_d is not None:
                merged_events.append((prev_d, delta.copy()))
            prev_d = d
            delta[:] = 0
        delta[i] += val
    # 添加最后一个d的事件
    merged_events.append((prev_d, delta))
    return merged_events

def compute_max_sum_with_k_removals(row_sums, k_max):
    """给定每行的贡献数组row_sums，计算最多剔除k_max个连续行块的最大总和及对应剔除方案"""
    n = len(row_sums)
    # dp[j]表示处理到当前行，剔除j个块时的最大总和，以及前一个状态的信息（用于回溯方案）
    dp = [(-np.inf, None, None) for _ in range(k_max+1)]  # (max_sum, prev_j, is_removed)
    dp[0] = (0, None, False)
    # 维护每个j对应的max(dp[j-1][0] for l<=i)，用于快速计算剔除连续行的情况
    max_prev = [-np.inf]*(k_max+1)
    max_prev[0] = 0
    # 记录路径信息
    path = [[None for _ in range(k_max+1)] for __ in range(n)]
    
    for i in range(n):
        new_dp = [(-np.inf, None, None) for _ in range(k_max+1)]
        # 情况1：不剔除第i行
        for j in range(k_max+1):
            if dp[j][0] != -np.inf:
                candidate = dp[j][0] + row_sums[i]
                if candidate > new_dp[j][0]:
                    new_dp[j] = (candidate, j, False)
                    path[i][j] = (j, False)
        # 情况2：剔除第i行，作为第j个块的一部分
        for j in range(1, k_max+1):
            # 最大的dp[j-1][0] 到当前行之前的状态
            if max_prev[j-1] != -np.inf:
                candidate = max_prev[j-1]
                if candidate > new_dp[j][0]:
                    new_dp[j] = (candidate, j-1, True)
                    path[i][j] = (j-1, True)
        # 更新max_prev：对于每个j，取当前dp[j]和之前的max_prev[j]的最大值
        temp_max = max_prev.copy()
        for j in range(k_max+1):
            if new_dp[j][0] > temp_max[j]:
                temp_max[j] = new_dp[j][0]
        max_prev = temp_max
        dp = new_dp
    
    # 找到最大的总和对应的j
    max_total = -np.inf
    best_j = 0
    for j in range(k_max+1):
        if dp[j][0] > max_total:
            max_total = dp[j][0]
            best_j = j
    
    # 回溯得到剔除方案
    removals = []
    current_i = n-1
    current_j = best_j
    in_removal = False
    removal_start = -1
    while current_i >= 0:
        prev_j, is_removed = path[current_i][current_j]
        if is_removed:
            if not in_removal:
                in_removal = True
                removal_start = current_i
        else:
            if in_removal:
                removals.append((removal_start, current_i+1))
                in_removal = False
        current_j = prev_j
        current_i -= 1
    if in_removal:
        removals.append((removal_start, 0))
    # 反转得到正确的顺序
    removals = [(s, e) for s, e in reversed(removals)]
    return max_total, removals

def find_optimal_solution(g, k_max):
    """主函数：找到最优半径、最大元素和、最优剔除方案"""
    n = g.shape[0]
    events = precompute_event_points(g)
    
    max_total = -np.inf
    best_radius = 0
    best_removals = []
    current_row_sums = np.zeros(n, dtype=np.int32)
    
    for d, delta in events:
        current_row_sums += delta
        # 计算当前d下的最优和
        total, removals = compute_max_sum_with_k_removals(current_row_sums, k_max)
        if total > max_total:
            max_total = total
            best_radius = np.sqrt(d)
            best_removals = removals
    
    return best_radius, max_total, best_removals

# 测试示例
if __name__ == "__main__":
    np.random.seed(40)
    N = 10
    g = np.random.randint(-3, 4, (N, N)).T  # 匹配示例图的矩阵方向
    k_max = 2
    radius, total, removals = find_optimal_solution(g, k_max)
    print(f"最优半径: {radius:.4f}")
    print(f"最大元素和: {total}")
    print(f"最优剔除行块（格式：(起始行, 结束行)，左闭右开）: {removals}")

方案优势

时间效率：预处理阶段为O(n² log n)，每个事件点的DP计算为O(nk)，k≤5，整体复杂度可控，支持大规模矩阵。
空间效率：仅需O(nk)的空间存储DP数组和路径信息，占用内存小。
准确性：通过事件点精准捕捉贡献变化，结合动态规划确保找到全局最优的行块剔除方案。

内容的提问来源于stack exchange，提问作者Simd

如何寻找可剔除最多k个连续行区域的最优圆形区域

示例说明 #

基础示例（无行剔除） #

带行剔除的示例 #

高效解决方案 #

核心思路 #

实现代码 #

方案优势 #