您的数字足够小，可以使用智能蛮力方法。
我使用两种方法来量化和最小化与“干净”等距值 (linspace(low, high, 7)) 的偏差。 "abserr" 表示平方差，"relerr" 表示平方误差除以平方净值。最后我也检查了corrcoefs，但我从未见过低于99.8%的东西
下面的代码首先找到具有最小错误的 clean 值的随机播放。这只需几秒钟，因为我们使用了以下技巧：
- 将 4 列分成两对
- 每对有 7 个！ 相对排列，即使平方也是一个可控的数字（每对一个因素）
- 计算这些 (7!)^2 次随机播放并对 对 求和
- 不必遍历配对之间的所有相对洗牌，我们观察到，如果两组配对和以相反的顺序排列，则总误差最小化，
"abserr" 和 "relerr" 是这样
最后，这些值被更正以使行总和为 100。这里我们再次使用这样一个事实，即当均匀分布时，求和误差最小化。
下面的代码包含两个变体，一个是旧版本solve，它在最小化relerr 和一个更正版本improved_solve 时包含一个小错误。他们经常找到不同的解决方案，但在 100 多个随机问题中，只有一个导致 improved_solve 的错误非常小。
几个例子的答案：
OP 的例子：
((75, 98), (6, 15), (2, 8), (0.05, 0.5))
solve relerr improved_solve relerr
table: table:
76.14213 15.22843 8.12183 0.50761 76.14213 15.22843 8.12183 0.50761
79.02431 13.53270 7.01696 0.42603 79.02431 13.53270 7.01696 0.42603
81.83468 11.87923 5.93961 0.34648 81.83468 11.87923 5.93961 0.34648
84.57590 10.26644 4.88878 0.26888 84.57590 10.26644 4.88878 0.26888
87.25048 8.69285 3.86349 0.19317 87.25048 8.69285 3.86349 0.19317
89.86083 7.15706 2.86282 0.11928 89.86083 7.15706 2.86282 0.11928
92.40924 5.65771 1.88590 0.04715 92.40924 5.65771 1.88590 0.04715
avgerr: avgerr:
0.03239 0.03239
corrcoefs: corrcoefs:
0.99977 0.99977 0.99977 0.99977 0.99977 0.99977 0.99977 0.99977
对某些列进行升序排序不是最佳的示例：
((11, 41), (4, 34), (37, 49), (0.01, 23.99))
请注意，求解器找到不同的解决方案，但错误是相同的。
solve relerr improved_solve relerr
table: table:
10.89217 18.81374 46.53926 23.75483 11.00037 24.00080 49.00163 15.99720
26.00087 9.00030 49.00163 15.99720 16.00107 19.00127 45.00300 19.99467
31.00207 4.00027 45.00300 19.99467 25.74512 13.86276 36.63729 23.75483
16.00000 29.00000 43.00000 12.00000 35.99880 8.99970 46.99843 8.00307
20.99860 33.99773 40.99727 4.00640 41.00000 4.00000 43.00000 12.00000
40.99863 13.99953 36.99877 8.00307 20.99860 33.99773 40.99727 4.00640
36.35996 24.23998 39.38996 0.01010 31.30997 29.28997 39.38996 0.01010
avgerr: avgerr:
0.00529 0.00529
corrcoefs: corrcoefs:
0.99993 0.99994 0.99876 0.99997 0.99989 0.99994 0.99877 0.99997
这就是improved_solve 实际击败旧版solve 的问题：
((36.787862883725872, 43.967159949544317),
(40.522239654303483, 47.625869880574164),
(19.760537036548321, 49.183056694462799),
(45.701873101046154, 48.051424087501672))
solve relerr improved_solve relerr
table: table:
21.36407 23.53276 28.56241 26.54076 20.25226 26.21874 27.07599 26.45301
22.33545 24.52391 26.03695 27.10370 21.53733 26.33278 25.10656 27.02333
23.33149 25.54022 23.44736 27.68093 22.90176 26.45386 23.01550 27.62888
24.35314 26.58266 20.79119 28.27301 24.35314 26.58266 20.79119 28.27301
25.40141 27.65226 18.06583 28.88050 25.90005 26.71994 18.42047 28.95953
26.47734 28.75009 15.26854 29.50403 27.55225 26.86656 15.88840 29.69279
27.58205 29.87728 12.39644 30.14424 29.32086 27.02351 13.17793 30.47771
avgerr: avgerr:
0.39677 0.39630
corrcoefs: corrcoefs:
0.99975 0.99975 0.99975 0.99975 0.99847 0.99847 0.99847 0.99847
代码：
import numpy as np
import itertools
import math
N_CHUNKS = 3
def improved_solve(LH, errtype='relerr'):
N = math.factorial(7)
# accept anything that looks like a 2d array
LH = np.asanyarray(LH)
# build equidistant columns
C = np.array([np.linspace(l, h, 7) for l, h in LH])
# subtract offset; it's cheaper now than later
c0, c1, c2, c3 = C - 25
# list all permutiations of a single column
p = np.array(list(itertools.permutations(range(7))))
# split into left and right halves, compute all relative permutiations
# and sort them by their sums of corresponding elements.
# Left pairs in ascending, right pairs in descending order.
L = np.sort(c0 + c1[p], axis=1)
R = np.sort(c2 + c3[p], axis=1)[:, ::-1]
# For each pair of permutations l in L, r in R compute the smallest
# possible error (sum of squared deviations.)
if errtype == 'relerr':
err = np.empty((N, N))
split = np.linspace(0, N, N_CHUNKS+1, dtype=int)[1:-1]
for LCH, ECH in zip(np.split(L, split, axis=0),
np.split(err, split, axis=0)):
dev = LCH[:, None] + R[None, :]
((dev / (100+dev))**2).sum(axis=-1, out=ECH)
del dev
elif errtype == 'abserr':
err = (np.add.outer(np.einsum('ij,ij->i', L, L),
np.einsum('ij,ij->i', R, R))
+ np.einsum('ik, jk->ij', 2*L, R))
else:
raise ValueError
# find pair of pairs with smallest error
i = np.argmin(err.ravel())
i1, i3 = np.unravel_index(i, (N, N))
# recreate shuffled table
c0, c1, c2, c3 = C
lidx = np.argsort(c0 + c1[p[i1]])
ridx = np.argsort(c2 + c3[p[i3]])[::-1]
C = np.array([c0[lidx], c1[p[i1]][lidx], c2[ridx], c3[p[i3]][ridx]])
# correct rowsums, calculate error and corrcoef and return
if errtype == 'relerr':
result = C * (100.0 / C.sum(axis=0, keepdims=True))
err = math.sqrt((((result-C)/C)**2).mean())
else:
result = C + (25 - C.mean(axis=0, keepdims=True))
err = math.sqrt(((result-C)**2).mean())
rs = np.sort(result, axis=1)
cc = tuple(np.corrcoef(ri, range(7))[0, 1] for ri in rs)
return dict(table=result.T, avgerr=err, corrcoefs=cc)
def solve(LH, errtype='relerr'):
LH = np.asanyarray(LH)
if errtype=='relerr':
err1 = 200 / LH.sum()
diff = np.diff(LH * err1, axis=1).ravel()
elif errtype=='abserr':
err1 = 25 - LH.mean()
diff = np.diff(LH, axis=1).ravel()
else:
raise ValueError
C = np.array([np.linspace(-d/2, d/2, 7) for d in diff])
c0, c1, c2, c3 = C
p = np.array(list(itertools.permutations(range(7))))
L = np.sort(c0 + c1[p], axis=1)
R = np.sort(c2 + c3[p], axis=1)[:, ::-1]
err = (np.add.outer(np.einsum('ij,ij->i', L, L),
np.einsum('ij,ij->i', R, R))
+ np.einsum('ik, jk->ij', 2*L, R)).ravel()
i = np.argmin(err)
i1, i3 = np.unravel_index(i, (math.factorial(7), math.factorial(7)))
L = np.argsort(c0 + c1[p[i1]])
R = np.argsort(c2 + c3[p[i3]])[::-1]
ref = [np.linspace(l, h, 7) for l, h in LH]
if errtype=='relerr':
c0, c1, c2, c3 = [np.linspace(l, h, 7) for l, h in LH * err1]
C = np.array([c0[L], c1[p[i1]][L], c2[R], c3[p[i3]][R]])
err2 = 100 / np.sum(C, axis=0)
C *= err2
cs = list(map(sorted, C))
err = math.sqrt(sum((c/r-1)**2 for ci, ri in zip(cs, ref) for c, r in zip(ci, ri)) / 28)
elif errtype=='abserr':
c0, c1, c2, c3 = [np.linspace(l, h, 7) for l, h in LH + err1]
C = np.array([c0[L], c1[p[i1]][L], c2[R], c3[p[i3]][R]])
err2 = 25 - np.mean(C, axis=0)
C += err2
cs = list(map(sorted, C))
err = math.sqrt(sum((c-r)**2 for ci, ri in zip(cs, ref) for c, r in zip(ci, ri)) / 28)
else:
raise ValueError
cc = tuple(np.corrcoef(ci, range(7))[0, 1] for ci in cs)
return dict(table=C.T, avgerr=err, corrcoefs=cc)
for problem in [((75, 98), (6, 15), (2, 8), (0.05, 0.5)),
((11, 41), (4, 34), (37, 49), (0.01, 23.99)),
((80, 94), (5, 14), (0.5, 5), (0.05, 0.5)),
((36.787862883725872, 43.967159949544317),
(40.522239654303483, 47.625869880574164),
(19.760537036548321, 49.183056694462799),
(45.701873101046154, 48.051424087501672))]:
for errtype in ('relerr', 'abserr'):
print()
columns = []
for solver in (solve, improved_solve):
sol = solver(problem, errtype)
column = [[' '.join((solver.__name__, errtype))]] + \
[[k + ':'] + [' '.join([f'{e:8.5f}' for e in r])
for r in np.atleast_2d(v)]
for k, v in sol.items()]
column = (line for block in column for line in block)
columns.append(column)
for l, r in zip(*columns):
print(f"{l:39s} {r:39s}")
problems = []
for i in range(0):
problem = np.sort(np.random.random((4, 2)), axis=1) * 50
for errtype in ('relerr', 'abserr'):
sol0 = solve(problem, errtype)
sol1 = improved_solve(problem, errtype)
if not np.allclose(sol0['table'], sol1['table']):
print(i, end= " ")
if np.abs((sol0['avgerr']-sol1['avgerr'])
/(sol0['avgerr']+sol1['avgerr']))>1e-6:
print(problem)
problems.append(problem)
columns = []
for sol, name in [(sol0, 'old '), (sol1, 'improved ')]:
column = [[name + errtype]] + \
[[k + ':'] + [' '.join([f'{e:8.5f}' for e in r])
for r in np.atleast_2d(v)]
for k, v in sol.items()]
column = (line for block in column for line in block)
columns.append(column)
for l, r in zip(*columns):
print(f"{l:39s} {r:39s}")