Saddleback search

[Xinyu LIU]

看了一篇Dijkstra前辈手写的论文：
http://www.cs.utexas.edu/users/EWD/ewd09xx/EWD934.PDF

经典！

Larry, LIU Xinyu
https://sites.google.com/site/algoxy/
https://github.com/liuxinyu95/AlgoXY

[Xinyu LIU]

简单地介绍一下Saddleback search。

Binary search 是非常popular的一种divide and conquer的search 方法。Knuth说：它的思想极其朴素，但是它的实现惊人的trick。

Jon Bentley在Programming pearls里用了一整章讲binary search。他说70%的seasoned programer写出的binary search有bug。
他的学生写信告诉他，10本教科书里，至少有4本的binary search有bug。然后Bentley给出的binary search在20年后发现里面也
藏了个bug。

所以binary search的实现，是一个不错的课外作业。

一个特别自然的想法是，binary search 能从1维扩展到2维上来么？

Richard Bird，用这道题目作为Oxford大学的入学面试：

给一个表格，m行n列，然后里面的数字在行列两个方向都是strict increasing的，然后他要求面试者用一个好的方法找出某个
数字z在表格里的位置。

例如：
    1 & 2 & 3 & 4 & ... \\
    2 & 4 & 5 & 6 & ... \\
    3 & 5 & 7 & 8 & ... \\
    4 & 6 & 8 & 9 & ...
    ... \\

要求找出数字5出现在哪些位置。

结果Bird发现一个有趣的现象，如果学生以前了接触过编程，就会不自觉的使用binary search的思想，上来就去表格的中间找，然后和z比较，但接下来就不知道怎么做了。

这个故事告诉我们把binary search扩展到2维空间并不容易。

在讲解Dijkstra的saddleback search前，先看看brute-force search如何解这道题目。首先，我们把此题generalize然后再formalize一下：

给一个strict increasing的函数f(x, y)，（例如f(x, y) = 2^x + 3^y)，和一个指定的整数z，我们希望找出所有符合f(x, y) = z的非负整数解。

上述的表格其实可以特化为这一定义的特例：

f(x, y) = a[x][y] if 0 <= x < n, 0 <= y < m, where m, n是表格的行数和列数。

brute-force的解法：

function Solve(f, z)
    A←Φ
    for x ∈ {0,1,2,...,z} do
        for y ∈ {0,1,2,...,z} do
            if f(x,y) = z then A ← A ∪ {(x, y)}
  return A

或者用公式表示为：
    solve(f,z) = {(x,y)|x ∈ {0,1,...,z},y ∈ {0,1,...,z},f(x,y) = z}

明显brute-force的解法调用了f (z+1)^2次。

Dijkstra 指出，这道题的窍门在于不要从原点（0,0)开始search,而要从左上角(0, z)开始。
对于任何candidate (p, q)，我们和z做比较：
若f(p, q) < z，由于f strict increasing，所以所有的0 <= y < q都有f(p, y) < z因此我们可以把(p, q)这点向下划一条线段，把上面的点都丢弃掉；
若f(p, q) > z，由于f strict increasing， 所以所有的p < x <= z都有f(x, q)>z 因此可以把(p, q)向右划一条线段，把上面的点都丢弃掉；
否则f(p, q) = z，我们找到一个解，我们可以把(p, q)所在的水平和垂直方向上线内的点都丢弃掉。

这带来了一个很直观的解：

function Solve(f, z)
  p ← 0, q ← z
  S←Φ
   while p ≤ z ∧ q ≥ 0 do
     z′ ←f(p,q)
     if z′ < z then
        p←p+1
     else if z′ > z then
        q←q−1
    else
        S ← S ∪ {(p, q)}
    p ← p + 1, q ← q − 1

或者写成函数式的：
solve f z = search 0 z where
  search p q | p > z || q < 0 = []
             | f p q < z = search (p + 1) q
             | f p q > z = search p (q - 1)
             | otherwise = (p, q) : search (p + 1) (q - 1)

这一算法的复杂度是多少呢？首先看3个best cases

1. p一直增加，直到超过z，共调用f, z+1次；

2. q一直减小，知道变为负数，共调用f, z+1次；
3. 每次都找到f(x, x) = z，沿着对角线一直到达(z, z)，共调用f, z+1次。

再看worst case：

从左上角(0, z)曲折向下，最终到达(z, z)，把水平方向的线段投射到x轴，垂直方向的投射到y轴，就会发现共调用了f，2(z+1)次。

对比上面的brute-force算法，这一saddleback search的basic version仅仅调用f,O(z)次，从平方提高到了线性。

今天先到这里，稍后我讲如何把binary search结合到saddleback search中去。

有兴趣，不怕英文、会翻墙的朋友可以看这里：
https://sites.google.com/site/algoxy/search
页面中的pdf文件。

Larry, LIU Xinyu
https://sites.google.com/site/algoxy/
https://github.com/liuxinyu95/AlgoXY

[Tian Guo]

这个就是杨式表， 是很经典的数据结构。

但说其是二维二分，感觉有点牵强。多维数据搜索， 基于二分思想的主要是kd-tree 和range tree.

--
 
---
您收到此邮件是因为您订阅了 Google 网上论坛的“TopLanguage”论坛。
要退订此论坛并停止接收此论坛的电子邮件，请发送电子邮件到 pongba+un...@googlegroups.com。
要查看更多选项，请访问 https://groups.google.com/groups/opt_out。
 
 

[Xinyu LIU]

首先这个不是Young-tableau.

Young-tableau可以理解为binary heap with implicit array扩展到2D上的情形，其extract-min操作(对应siftup term in some materials)使用infinity来填充。
可以参见CLRS《算法导论》Problem 6-3. 请注意区别。

Saddleback search见我前面引述的
Edsger W. Dijkstra. ``The saddleback search''. EWD-934. 1985. http://www.cs.utexas.edu/users/EWD/index09xx.html.

其次，我前面说了，二分查找还没有用上呢，这个是basic version的saddleback search，我会在后继讲如何从saddle back search演进到2D binary search。
如果你等不及，请参考论文：
Bird, R. (2006). Improving saddleback search: a lesson in algorithm design. Mathematics of Program Construction, LNCS 4014, pp. 82–9.

最后，我前4年在TL陆陆续续post过很多相关的内容，这里不再重头写了。

Larry, LIU Xinyu
https://sites.google.com/site/algoxy/
https://github.com/liuxinyu95/AlgoXY

2013/5/3 Tian Guo <tian...@epfl.ch>

[Tian Guo]

哦 看到了 不好意思  刚才没仔细看， 多谢指点！

期待您后面的讲解.  

[qiaojie]

我尝试写了个二维版的binary search，请指教

int Mat[M][N];

void search(int l, int t, int r, int b, int z)

{

if(Mat[l][t] == z)

output(l,t,z);

else if(Mat[r][b] == z)

output(r,b,z);

else if(Mat[l][t] < z && Mat[r][b] > z)

{

int m0 = (l + r) / 2;

int m1 = (t + b) / 2;

search(l, t, m0, m1, z);

search(m0 + 1, t, r, m1, z);

search(l, m1 + 1, m0, b, z);

search(m0 + 1, m1 + 1, r, b, z);

}

}

search(Mat, 0,0, M-1,N -1, z);

算法复杂度应该是O(ln(m*n))

[lor fal]

我不明白 function Solve( f, z)中为什么要把z赋值给q？？

我找了一篇文章：http://www.cs.geneseo.edu/~baldwin/math-thinking/saddleback.html
是从(0,n-1)开始搜索的。我觉得这才比较合理吧。

2013/5/3 Xinyu LIU <liuxi...@gmail.com>

--
 
---
您收到此邮件是因为您订阅了 Google 网上论坛的“TopLanguage”论坛。
要退订此论坛并停止接收此论坛的电子邮件，请发送电子邮件到 pongba+un...@googlegroups.com。
要查看更多选项，请访问 https://groups.google.com/groups/opt_out。
 
 

--
Best regards,
LORFAL

[Tian Guo]

这怎么会是log 的复杂度？ 按照那个代码，递推方程是 f(n)=4f(n/4)+O(1)  是线性的。

在 2013年5月3日下午11:05，Leo <enpe...@gmail.com>写道：

用个bool做返回值好些，如果找到了就不需要再找了。
复杂度，最糟糕的情况应该是log2(M*M*N) 假设 M<=N

[Keren Zhou]

Tian Guo

根据你的递推公式来看，是没有大问题的。这里的f(n)你应该指的是n=M*N吧。

因此：

f(n)=4f(n/4)+O(1)

时间复杂度是O(n)显而易见。根据master method分析即可。

试想一下，一维数组的binary search和这个是一个原理，它的搜索复杂度是O(log(n))。

因为其本身就是建立在已经排序的基础上的。

binary search：

f(n)=f(n/2)+O(1)

因此，正如前面Leo说的，在这里做一个返回值会好一些。

但是这样的话，怎样分析average cost还没有想好，不过worst case还是O(n)无疑。

还请各位指教！

[qiaojie]

你没看懂代码，前两个if条件满足是不会进入后面那个else if的，所以找到了就不会继续搜索了。

在 2013年5月4日上午5:05，Leo <enpe...@gmail.com>写道：

用个bool做返回值好些，如果找到了就不需要再找了。
复杂度，最糟糕的情况应该是log2(M*M*N) 假设 M<=N

On Friday, May 3, 2013 12:34:38 PM UTC-4, qiaojie wrote:

[qiaojie]

不知道你这个递推公式怎么算的，因为数组是严格递增的，所以不存在最坏情况让所有的块都满足搜索条件，大部分块会不满足条件而停止搜索，如何能得出f(n)=4f(n/4)+O(1) ？

在 2013年5月4日上午6:27，Tian Guo <tian...@epfl.ch>写道：

[qiaojie]

不过说log复杂度应该是不对的，至少是跟z的大小线性相关的

[Xinyu LIU]

Hi, 各位，

找到了不能立即返回，这样找的的解不全。这是2D与1D的一个不同之处。

虽然f是strict increasing的，但是某个z会有多个对应值，请参考我前面给出的那个表格。

例如4，还有5，都出现了多次。

也就是说，虽然在某一行，或者某一列是严格单调增的（无重复），但是不能保证在表格内没有重复元素。

各位可以用函数f(x, y) = x^2 + y^2生成表格看看，很有趣。

to Lor Fal:

将q初始化为z，就是要从top-left corner: (0, z)开始search, 请参考下面链接的python代码：

https://github.com/liuxinyu95/AlgoXY/blob/algoxy/search/binary-search/src/saddleback.py

我粘贴一下其中的一个片段：

def saddleback(f, z):
    (p, q) = (0, z)
    res = []
    while p <= z and q >= 0:
        z1 = f(p, q)
        if z1 < z:
            p = p + 1
        elif z1 > z:
            q = q - 1
        else:
            res.append((p, q))
            (p, q) = (p + 1, q - 1)
    return res

<<<

请特别注意上述链接中的invariant测试用例。确保结果和brute-force是一致的。

另外，先在中心(z/2, z/2)检查，比较f(z/2, z/2)和z的大小后，关键点在于4个矩形中，可以丢弃的只有一个，搜索区域变成了一个L型。

（要么丢弃top-left矩形，要么丢弃bottom-right矩形）。

对于这剩下的L型，如何高效的进行divide and conquer的search是关键。这需要仔细的设计解法，我后继会给出。

Larry, LIU Xinyu
https://sites.google.com/site/algoxy/
https://github.com/liuxinyu95/AlgoXY

2013/5/4 qiaojie <qia...@gmail.com>

[qiaojie]

在 2013年5月4日下午9:46，Xinyu LIU <liuxi...@gmail.com>写道：

Hi, 各位，

找到了不能立即返回，这样找的的解不全。这是2D与1D的一个不同之处。

虽然f是strict increasing的，但是某个z会有多个对应值，请参考我前面给出的那个表格。

例如4，还有5，都出现了多次。

我的解法是找全部解的，请问有问题吗？

 

 

另外，先在中心(z/2, z/2)检查，比较f(z/2, z/2)和z的大小后，关键点在于4个矩形中，可以丢弃的只有一个，搜索区域变成了一个L型。

（要么丢弃top-left矩形，要么丢弃bottom-right矩形）。

对于这剩下的L型，如何高效的进行divide and conquer的search是关键。这需要仔细的设计解法，我后继会给出。

L型是最坏情况，也可能只有一块或者2块矩形满足条件的。

[Tian Guo]

哦 对， 可能是代码习惯不一样， 这么写也对。

不过感觉楼主提到的 saddleback search 不是这种简单的情况，还没来得及看dij 的手稿。 

期待楼主下文。

[Ray Song]

Richard S. Bird 的 Improving Saddleback Search: A Lesson in Algorithm Design 被選作 Pearls of Functional Algorithm Design

的一節了。主要思想有這麽幾個：

z 對應值在 Young tableau 裏的出現情況有 C(m+n, m) 種，探測一次 f(x,y) 返回三種結果(<,=,>)，

所以下界是 Ω(log C(m+n,m)) = Ω(m log (n/m)) (設m > n)

假設m<n，在n這一維度binary search，可以得到T(m,n) = log n + 2T(m/2, n/2)，達到下界

另外一種實現簡單的方式(也適合人腦使用)是不用二分，從右上到左下……O(m+n)

--

Ray

http://maskray.me

[Xinyu LIU]

Hi,

2年前读Richard Bird的Pearls of functional programming就觉得这章很好。但是这次我给出一个不同角度的解释（或者说如何进化到最终的那个算法）。
首先，前面时间比较紧，有些特别容易混乱的地方我再次强调一下：

1. 对于表格和矩阵，我们通常规定左上角top-left是(0, 0)，右下角bottom-right是(m, n)

2. 对于一个二元函数的作用域，这个结果恰好和1不一样，top-left是(0, m), 而bottom-right是(n, 0)
前面的saddleback search是从(0, m)开始，向(n, 0)进行搜索。

现在我们回到2D。并且看看使用中点的判定结果。在中点(p, q)划两条直线（水平+垂直），作用域矩形被划分为4个区域：

A, B

C, D

如果f(p, q) < z，我们可以丢弃C这个小矩形，如果f(p, q) > z我们可以丢弃B这个小矩形；但是接下来，我们要在剩余的3个矩形组成的L形状内搜索。

但是还有一个情况，f(p, q) == z，这种情况，由于f是单调增的。所以(p, q)所在的水平、垂直直线上的所有点都可以丢弃。
并且，我们可以同时丢弃C和B两个矩形，而递归地在A, D这两个矩形内搜索。显然这种情况是最理想的。因为我们把作用域一下子砍掉1半。

divide conquer search的一个要点就是，丢掉的不可能区域越大越快，搜索的就越快。但是使用中点划分，我们大多数情况，只能丢弃一个
小矩形，只有相等的时候，才能丢弃两个矩形（砍掉一半）。

有没有什么方法，使得我们每次都丢弃一半的矩形呢？有！只要我们不执着于使用矩形的终点。我们可以找到矩形内任何一点(p, q)
如果f(p, q) == z，那么使用这点划分，就保证能丢弃一半。

我们看看下面的方法：
假设要搜索矩形(a, b) - (c, d)，我们先比较矩形的长和宽，如果矩形的长>宽，我们就划一条垂直平分线；否则划一条水平平分线；
然后我们在这条中线上进行二分查找。（由于我们事先比较了长和宽，所以我们保证二分查找的范围尽量小）。
二分查找的目标是是否存在一个直线上的点(p, q)，使得f(p, q) == z，找到了我们就可以砍掉矩形的一半，然后递归搜索了。

但是这个方法有个问题：存在这样的可能，在这条中线上，不存在一个点(p, q)使得f(p, q) == z。
解决方案是，我们放宽二分查找的条件，查找这样一个点:
如果是水平线，我们找f(p, q) <= z < f(p+1, q)
如果是垂直线，我们找f(p, q) <= z < f(p, q+1)
这样的f(p, q)同样可以让我们砍掉矩形的一半（请读者思考为什么？）。

最后再考虑一个edge case：如果中线上所有的点求得f值都< z或者都> z怎么半？只要我们的二分查找能保证在这种条件下返回lower bound或者upper bound就可以。这样我们就可以把整个中线的一边全部砍掉。

要完整实现这个算法，还需要考虑一个细节，找到某个点(p, q)后，如果f(p, q) ！= z，那么我们要砍掉的矩形是否包含经过(p, q)的中线线段？这个留给读者思考。

那么这个算法有多快呢？显然，由于我们每次砍掉一半的矩形面积，所以我们需要log (z^2) = 2 log z 次砍的动作；
但是每次我们需要使用binary search寻找合适的(p, q)，binary search的复杂度是log(min(m, n))，其中m, n是矩形的长和宽。
我们列一下递归公式，假设m > n，有

T(m, n) = log n + 2 T(m/2, n/2)
解方程得 T(m, n) = O(m log (n/m))
而这是一个复杂度的一个上限！见上面Ray Song的post。

原理就是如此，我现在自己的实现还有一些bug，先不贴了。

最后再重复一下这个算法的设计关键：扔得越快越好。

Larry, LIU Xinyu
https://sites.google.com/site/algoxy/
https://github.com/liuxinyu95/AlgoXY

2013/5/5 Ray Song <emac...@gmail.com>

[Xinyu LIU]

还有一个事实特别值得我们思考：
什么时候使用上述的2D二分查找比basic ver. 的Saddleback search快？ 快多少？
如果我们把Saddleback search稍作改动，不从(0, z), 搜索到(z, 0)，而是从(0, m)搜索到(n, 0)，其中m, n是沿着x轴和y轴进行二分查找找到的边界。

我们发现，大多数情况下，saddleback search表现的已经非常好了。一个简单而朴素的策略从1D扩展到2D，结果完全没有按照人们的预期去变化发展，实在是非常有趣的现象。

Larry, LIU Xinyu
https://sites.google.com/site/algoxy/
https://github.com/liuxinyu95/AlgoXY

2013/5/7 Xinyu LIU <liuxi...@gmail.com>

[qiaojie]

为什么你只用矩阵中的中心点去和z比较？根据矩阵单调递增的特性，左上角元素是最小值，右下角是

最大值，所以应该以最小值和最大值去和z比较，只有z值落在最小值和最大值之间才需要递归搜索，否

则可以直接抛弃或者输出结果。

一般来说z值在矩阵中的分布规律接近于一条线，经过几次递归之后，需要搜索的子矩阵会收敛到这条

线附近，因此这个搜索算法的平均复杂度是线性的。

[Xinyu LIU]

Hi,

z值的分布和函数有关，你所说的一条直线的观察是基于这个生成函数的：

f(x, y) = x + y
可以看一下：[[x + y | x<-[0..10]] | y<-[0..10]]

[[0,1,2,3,4,5,6,7,8,9,10],
[1,2,3,4,5,6,7,8,9,10,11],
[2,3,4,5,6,7,8,9,10,11,12],
[3,4,5,6,7,8,9,10,11,12,13],
[4,5,6,7,8,9,10,11,12,13,14],
[5,6,7,8,9,10,11,12,13,14,15],
[6,7,8,9,10,11,12,13,14,15,16],
[7,8,9,10,11,12,13,14,15,16,17],
[8,9,10,11,12,13,14,15,16,17,18],
[9,10,11,12,13,14,15,16,17,18,19],
[10,11,12,13,14,15,16,17,18,19,20]]

但是如果不是这样的对称函数，则并没有直线这一结论：
例如f(x, y) = x + y^2
可以看一下：[[x + y^2 | x<-[0..10]] | y<-[0..10]]
[[0,1,2,3,4,5,6,7,8,9,10],
[1,2,3,4,5,6,7,8,9,10,11],
[4,5,6,7,8,9,10,11,12,13,14],
[9,10,11,12,13,14,15,16,17,18,19],
[16,17,18,19,20,21,22,23,24,25,26],
[25,26,27,28,29,30,31,32,33,34,35],
[36,37,38,39,40,41,42,43,44,45,46],
[49,50,51,52,53,54,55,56,57,58,59],
[64,65,66,67,68,69,70,71,72,73,74],
[81,82,83,84,85,86,87,88,89,90,91],
[100,101,102,103,104,105,106,107,108,109,110]]

然后再给一个对称的函数：

f(x, y) = x^2 + y^2
[[x^2 + y^2 | x<-[0..10]] | y<-[0..10]]

[[0,1,4,9,16,25,36,49,64,81,100],
[1,2,5,10,17,26,37,50,65,82,101],
[4,5,8,13,20,29,40,53,68,85,104],
[9,10,13,18,25,34,45,58,73,90,109],
[16,17,20,25,32,41,52,65,80,97,116],
[25,26,29,34,41,50,61,74,89,106,125],
[36,37,40,45,52,61,72,85,100,117,136],
[49,50,53,58,65,74,85,98,113,130,149],
[64,65,68,73,80,89,100,113,128,145,164],
[81,82,85,90,97,106,117,130,145,162,181],
[100,101,104,109,116,125,136,149,164,181,200]]

Larry, LIU Xinyu
https://sites.google.com/site/algoxy/
https://github.com/liuxinyu95/AlgoXY

2013/5/8 qiaojie <qia...@gmail.com>

[Xinyu LIU]

给一个Improved Saddleback search with binary search的实现，具体解释见我昨天的post.

solve'' f z = search f z (0, m) (n, 0) where
  m = bsearch (f 0) z (0, z)
  n = bsearch (\x -> f x 0) z (0, z)
 
search f z (a, b) (c, d) | c < a || b < d = []
                         | c - a < b - d = let q = (b + d) `div` 2 in csearch (bsearch (\x -> f x q) z (a, c), q)
                         | otherwise = let p = (a + c) `div` 2 in rsearch (p, bsearch (f p) z (d, b))
  where 
    csearch (p, q) | z < f p q = search f z (p, q - 1) (c, d)
                   | f p q == z = search f z (a, b) (p - 1, q + 1) ++ (p, q) : search f z (p + 1, q - 1) (c, d)
                   | otherwise = search f z (a, b) (p, q + 1) ++ search f z (p + 1, q - 1) (c, d)
    rsearch (p, q) | z < f p q = search f z (a, b) (p - 1, q)
                   | f p q == z = search f z (a, b) (p - 1, q + 1) ++ (p, q) : search f z (p + 1, q - 1) (c, d)
                   | otherwise = search f z (a, b) (p - 1, q + 1) ++ search f z (p + 1, q) (c, d)


大致说明一下。首先我们使用改进的1D binary search，沿x轴和y轴找到更加精确的一个范围。
这个改进的binary search 是针对(l, u)的范围和一个指定值z，搜索这样的一个x 使得 f(x) <= z < f(x+1)

bsearch f y (l, u) | u <= l = l
                   | f m <= y = if f (m + 1) <= y then bsearch f y (m + 1, u) else m
                   | otherwise = bsearch f y (l, m-1)
  where m = (l + u) `div` 2

注意，如果最小的值f(l)都比z大，则返回l。

然后，我们在矩形(0, m)和(n, 0)内搜索f(x, y) = z
如果矩形大小是空，则自然结果为空；
如果矩形是瘦高的，我们使用csearch沿着矩形的水平中线找到一个(p, q)，可以有效地把矩形丢弃掉一半
如果矩形是矮胖的，我们使用rsearch沿着矩形的垂直中线找到一个(p, q)，可以有效地把矩形丢弃一半

该代码使用随机生成的100条case测试通过。所搜寻结果和brute force search一致。

Larry, LIU Xinyu
https://sites.google.com/site/algoxy/
https://github.com/liuxinyu95/AlgoXY

2013/5/7 Xinyu LIU <liuxi...@gmail.com>

[Xinyu LIU]

Hi，

一直没有对qiaojie的算法进行复杂度计算。我刚才算了一下。

先重复一遍这个算法：
初始搜索矩形为(0, 0) - (z, z)

搜索任何矩形(a, b) - (c, d)时
- 若f(a, b) < z 或 f(c, d) > z，则返回空；
- 否则，取矩形中点将其一分为4，然后递归搜索4个矩形。

先看best case:
每次矩形中点都满足f(p, q) = z, 其中p = (a + b)/2, q = (c+d)/2。由于f单调增，故而四个矩形中，有两个矩形的搜索立即返回空，另外两个继续递归。
故而算法共调用f函数次数如下：

T(z) = 8 + 2* T(z/2)
解此递归方程得：

T(z) = 8(2^(i+1) - 1)，其中z = 2^i

= 8（2 z - 1) = 16 z - 8 次

再看worst case:
每次矩形中点要么大于，要么小于z, 由于f单调增，故而四个矩形中，有一个矩形搜索立即返回空，另外3个继续递归。
故而算法共调用f函数次数如下：

T(z) = 8 + 3 T(z/2)
解此递归方程得：

T(z) = 8(3^(i+1) -1)/2 = 8 (3 * z^(log 3) -1) = 24 z^(log 3) - 8 次

如有计算错误，请及时指出。

basic ver 的 saddleback search的复杂度为：
  best case计算f共计z次，worst case共计2z次。 

Larry, LIU Xinyu
https://sites.google.com/site/algoxy/
https://github.com/liuxinyu95/AlgoXY

2013/5/8 Xinyu LIU <liuxi...@gmail.com>

[Xinyu LIU]

Hi,

抱歉，计算有错误，寅吃牟粮把下次递归的计算次数算入本次了。

best case:

T(z) = 2 + 2 * T(z/2)

worst case:

T(z) = 2 + 3 * T(z/2)

不是8，是2

[qiaojie]

我是说一条线，没有说是直线。如果你把矩阵看成一个高度场，那么因为其单调性，任意高度的等高线都是唯一的。所以你后面的复杂度分析是错的，最好情况是O(1)，比如找左上角或右下角元素。每次递归时最坏情况只能扔掉一个子矩阵，但是继续往下递归就会发现不可能每次都是最坏情况，最后子矩阵范围会收敛到等于z值的等高线附近。所以定性分析的话这个算法是和z相关的线性复杂度算法。定量计算的话还没思路怎么算

在 2013-5-8 下午1:55，"Xinyu LIU" <liuxi...@gmail.com>写道：

[Leo]

抛开块递归，可以用另外一个方式思考下。

先考虑个简单情况，假设m=n，整个矩阵元素可以分解成3个不相交的单调集合：

1. 对角线

2. 对角线的右上区域，每个对角线上的点，都对应有一条先水平再竖直的折线

3. 对角线的左下区域，类似...

最简单的就是对每个单调数列做2分查找，当然这意义不大。

可以做些优化：

1. 先对对角线做2分查找，结果是可以找到2个点， 假设是p0(x, x), p1(x+1, x+1), f(p0)<=z<f(p1)

2. 考虑和对角线相邻的右上部分折线，根据前一个搜索结果，我们可以缩小当前折线的搜索范围，从p0点向上做投影，从p1点向右做投影，即q0 (x, x-1), q1(x+1,x)，显然z肯定是在这个范围内，在这个范围内应用2分搜索，然后记录下找到的位置，再应用到下一个折线...

3. 左下也是类似

m != n的情况，其实也是一样，只不过对角线是一个折线，其他的也是大同小异。

有兴趣的可以推演下复杂度. 

[Earthson Lu]

哈，LS这个我也正想说呢:)

每次对对角线（取折线，长度是m+n）进行二分查找，找到两个点。这样能排除两个块。令x=m+n，即矩阵两个维度的和，那么最坏情况下有如下复杂度

f(x) = log(x) + f(a) + f(b)，其中a+b = x

~~~~~~~~~~

这其实相当于对Leo提到的第一个算法的“优化”。原本二分是随着对角线而进行的，假设右上和左下区域我们先不急着递归下去，而是等对角线的区域先递归（只找对角线上的）。那么我们再利用这个递归的结果重新选择两块区域（把砍掉的区域重组到右上和左下，并且去掉原本就不需要的块）。实际上我们产生的区域将趋向于更加扁平，我觉得越扁平，就越倾向于有序，也许会更快一些？

2013/5/9 Leo <enpe...@gmail.com>

--
 
---
您收到此邮件是因为您订阅了 Google 网上论坛的“TopLanguage”论坛。
要退订此论坛并停止接收此论坛的电子邮件，请发送电子邮件到 pongba+un...@googlegroups.com。
要查看更多选项，请访问 https://groups.google.com/groups/opt_out。
 
 

--

~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

Perfection is achieved 

not when there is nothing more to add

 but when there is nothing left to take away

[Earthson Lu]

另外，a/b ~ m/n。这个结构和Heap初始化很像，应该是O(m+n)的吧。

2013/5/9 Earthson Lu <earth...@gmail.com>

[Earthson Lu]

囧了，貌似看走眼了。Leo上面提的那个和我上面说的这个不太一样。说实话，那个优化的第二步我没看懂的说。@Leo

我上面提到的这个其实是这么个过程。

1. 取对角折线。二分找到两个点a，b（在折线上），其中a <= z <= b

2. 取a的左下投影为需要递归的矩阵A。取b的右下投影为矩阵B

3. 对矩阵A，B做适当裁剪后递归，这里需要分挺多情况的，就不展开了（当然，我觉得不裁也不太会影响复杂度）。

2013/5/9 Earthson Lu <earth...@gmail.com>

[Earthson Lu]

上面笔误，第二步是取b的右上投影。

算法很容易扩展到非严格递增的情况（允许相等），虽然结果最多可能有O(m*n)个，但是描述复杂度依然是O(m+n)的应该，因为如果左上和右下元素相等，那么整个矩阵就只有“一个”元素。

2013/5/9 Earthson Lu <earth...@gmail.com>

[Earthson Lu]

@qiaojie有提到那条“线”。实际上我觉得lz给的saddleback search就是先找到线，然后沿着线走。很容易证明线在行和列的方向上都是单调的。不知道我有没有“看走眼。。。

2013/5/9 Earthson Lu <earth...@gmail.com>

[Xinyu LIU]

Hi,

这就是Saddleback这个名字的由来（Dijkstra说这个名字是Gries给起的）。在三维空间上绘制f(x, y) =z， 会得到一个类似saddleback的形状。这个search的要点就是，从saddleback的一端开始搜索，沿着saddleback走到另一端。

@Qiaojie: 我的worst-case计算可能不对，但是best case不会出错： T(z, z) = 2 + 2 T(z/2, a/2)。你说best case是O(1)无法令人理解。

最信服的方法是，你在每次调用f (在你最初的源代码是二维数组access)时将一个计数器加1，然后你看看结果是否是4z。

Larry, LIU Xinyu
https://sites.google.com/site/algoxy/
https://github.com/liuxinyu95/AlgoXY

2013/5/9 Earthson Lu <earth...@gmail.com>

[Xinyu LIU]

OK，我猜到你这个O(1)是怎么来的了。这一特殊情况就是矩阵，你上来比较发现矩阵的M[m][n]比z还小，所以就结束了。
我考虑的是general的单调函数f(x, y)，所以不管你给出什么样的z (z是non negative integer)，都不会出现上述情况。因为f(z, z) >= z

本质上将，矩阵并非一个在N上的单调函数，（因为我们可以定义f(x, y) = M[y][x] if x in [1..n], y in [1..m] else -1）。

Larry, LIU Xinyu
https://sites.google.com/site/algoxy/
https://github.com/liuxinyu95/AlgoXY

2013/5/10 Xinyu LIU <liuxi...@gmail.com>

[qiaojie]

O(1)的情况有很多种啊，比如z等于矩阵的左上角或者右下角元素，这个递归函数就立刻结束了。或者递归几次以后找到z元素，同时其他子矩阵的范围超出z。

另外为了验证我的复杂度猜想，我写了个暴力程序，分别选取3种函数及256*256、1024*1024两种矩阵规模，穷举每一个矩阵数值范围内的z值所需要的递归搜索次数，记录搜索次数最多的时候所对应的z值，列表如下：

f(x,y)=x+y

矩阵大小：256*256      最大搜索次数：1497  对应的z：253 

矩阵大小：1024*1024  最大搜索次数：6097  对应的z：1021

f(x,y)=x*x+y

矩阵大小：256*256      最大搜索次数：1021  对应的z：16639

矩阵大小：1024*1024  最大搜索次数：4093  对应的z：263167

f(x,y)=x*x+y*y

矩阵大小：256*256      最大搜索次数：1441  对应的z：63809

矩阵大小：1024*1024  最大搜索次数：5685  对应的z：1029641

可见最坏情况下的算法复杂度也是线性的，当然O(1)出现的情况就很多了。

[Xinyu LIU]

Hi,

Saddleback这个名字的由来可以运行下面的python代码看画出的图像。

import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
import numpy as np

fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')

xs, ys = np.meshgrid(range(100), range(100))
zs = xs**2 + ys**2
ax.plot_wireframe(xs, ys, zs)
ax.set_xlabel('X')
ax.set_ylabel('Y')
ax.set_zlabel('Z')

plt.show()

Larry, LIU Xinyu
https://sites.google.com/site/algoxy/
https://github.com/liuxinyu95/AlgoXY

2013/5/11 qiaojie <qia...@gmail.com>