Nov
12
从高考完开始，每月保持至少一篇的节奏，终于还是在今年被打破了。幸好3月写得多，我还可以说今年平均每月还有一篇；但是明年恐怕就没这底气了。
其实并不是不想写博客了，只是因为无从下手敲键盘。工作上没什么值得写的事情。没有什么困难的东西，有点新意的也就是象征性地接触了一下以前想搞搞但是一直动力不足的hadoop什么的，安装说明什么的还是不写比较好...
八月以来，拾起三月看断的SICP，陆陆续续地看到了第四章，连自己都有点意外。虽然从书中学到了不少东西，但似乎也没有什么适合总结写下来的，也许是因为还没“悟道”吧，只能先把代码提交到 GitHub ，期望着哪天看完了能写点什么。
这几个月生活的主要亮点之一是买了Kindle Paperwhite以后，看了好多书。历史、哲学、金融、文学，都有涉猎，挺有意思的。最近把《量子物理史话》又看了一遍，然后开始看《明朝那些事儿》。不仅买了很多电子书，还莫名其妙心血来潮买了好多纸质书，其中有几块大部头，比如《30天自制操作系统》、《编译原理》、《什么是数学》……感觉有很多只能拿来垫显示器了，有点囧，尽量看吧。不管怎么说，买了kpw以后增加了自己的阅读量，挺好的，因此在这里顺便强烈推荐一下。
最后征集一下早餐的问题。在吃了两个月华夫饼、接着又泡了两个月永和豆浆粉以后，想换点什么新的、适合在办公室搞定的东西，求推荐……
p.s. 再顺便推荐几本好书吧，Kindle电子书 //亚马逊有Kindle for PC/Mac/iOS/Android
《上帝掷骰子吗:量子物理史话》 ￥2.99
《时间的形状:相对论史话》 ￥5.00
哲学家们都干了些什么? ￥0.99
搞懂金融的第一本书 ￥3.99
重说中国近代史 ￥1.60
明朝那些事儿(1-7套装) ￥15.00
其实并不是不想写博客了，只是因为无从下手敲键盘。工作上没什么值得写的事情。没有什么困难的东西，有点新意的也就是象征性地接触了一下以前想搞搞但是一直动力不足的hadoop什么的，安装说明什么的还是不写比较好...
八月以来，拾起三月看断的SICP，陆陆续续地看到了第四章，连自己都有点意外。虽然从书中学到了不少东西，但似乎也没有什么适合总结写下来的，也许是因为还没“悟道”吧，只能先把代码提交到 GitHub ，期望着哪天看完了能写点什么。
这几个月生活的主要亮点之一是买了Kindle Paperwhite以后，看了好多书。历史、哲学、金融、文学，都有涉猎，挺有意思的。最近把《量子物理史话》又看了一遍，然后开始看《明朝那些事儿》。不仅买了很多电子书，还莫名其妙心血来潮买了好多纸质书，其中有几块大部头，比如《30天自制操作系统》、《编译原理》、《什么是数学》……感觉有很多只能拿来垫显示器了，有点囧，尽量看吧。不管怎么说，买了kpw以后增加了自己的阅读量，挺好的，因此在这里顺便强烈推荐一下。
最后征集一下早餐的问题。在吃了两个月华夫饼、接着又泡了两个月永和豆浆粉以后，想换点什么新的、适合在办公室搞定的东西，求推荐……
p.s. 再顺便推荐几本好书吧，Kindle电子书 //亚马逊有Kindle for PC/Mac/iOS/Android
《上帝掷骰子吗:量子物理史话》 ￥2.99
《时间的形状:相对论史话》 ￥5.00
哲学家们都干了些什么? ￥0.99
搞懂金融的第一本书 ￥3.99
重说中国近代史 ￥1.60
明朝那些事儿(1-7套装) ￥15.00
Aug
22
本篇来自这个问题：python中创建父类对象的问题
也许可以认为这是Python设计的缺陷导致的，因为Python 3.0只需要用 `super().some_attr` 就行了。至于为什么需要两个参数，可以大概分析一下，如果有错，欢迎指正。
先介绍下背景知识：super是从Python 2.2开始随着 `new-style class` 一起引入的，也只能应用于所谓的 `new-style class` (即直接或间接继承于 `object` 的class)，可以在一定程度上解决所谓`钻石继承`的问题。2.2起所有python内置class都是new-style class。
那么super是什么呢？实际上它是一个builtin的type，你调用 `super(Base, self)` 会返回一个 `__class__ = super` 的对象，这个对象可以作为一个代理，通过BFS的顺序（实际上列表已经保存在`Base.__mro__`里了）去访问第一个**直接定义**了目标属性的基类里的属性。`super(type, obj_or_subtype).attr` 基本上可以认为是 `find_in_mro_provide_by(obj_or_subtype, "attr")` ，有点晦涩。
**注意**，super()返回的不是“父类对象”，而是一个super类型的对象；实例化的时候，第一个参数是一个类型(type)，第二个参数可以是type的instance，也可以是type的subclass。
介绍完基础知识，可以先说说为什么要有第二个参数。理由很简单 —— 我们知道 self 代表当前对象，每个对象的方法在定义的时候都需要显示地把self作为第一个参数，你本来应该写成这样
some_class.some_method(obj, *args, **kwargs)
但是因为Python语法允许 `obj.some_metod(*args, **kwargs)` (本质上是个语法糖) ，所以你可以写得简单点（不必显式地给出方法的第一个参数）。而super对象则不同，它没有语法上的直接支持，所以在内部invoke some_method的时候必须指定某个对象，而这个对象的得你自己塞给它，也就是说把这个super对象绑定(bound)到第二个参数上。所以实际上并不是非要用self作为super的第二个参数，甚至super并不是必须在class内部才能调用：
至于为什么要有第一个参数，如果你看了 `help(super)` 就会发现它居然提供了个单参数的版本：
super(type) -> unbound super object
这个理解起来比较困难点。先说这个`unbound`，没绑定，就是说这个super对象没有实际绑定到某个对象上，它并不是可以直接用的。要怎么用呢？那就得注意到`help(super)`里面有一个 `__get__` 方法，也就是说，super类还是一个Descriptor类！哎，这个Descriptor类展开又是好大一段，简单地说就是它可以把自己和某个object的属性绑定，使得访问这个属性的时候，实际上是在调用它的`__get__/__set__`方法：
回到unbound super, 为了使用它，就要通过 __get__ 方法把它再绑定到一个对象上。由于Descriptor的属性，特别适合这么用：
这是它的一个可能用法。我不确定是否还有其他更合适的用途，但是在这里实际上也并不是很好，尤其是遇到staticmethod的时候还会出错，再加上绕了这么大一个弯，实在不是很推荐使用。
据说unbound super的使用非常少，不知道现实意义有多少，也不知道当初为什么设计成这个样子，总之由于Python3.0已经改了（虽然仍然保留了unbound super），所以基本上还是可以认为这是设计的历史遗留问题。
参考文献（unbound super主要参考了这个系列）：Things to Know About Python Super
[1] http://www.artima.com/weblogs/viewpost.jsp?thread=236275
[2] http://www.artima.com/weblogs/viewpost.jsp?thread=236278
[3] http://www.artima.com/weblogs/viewpost.jsp?thread=237121
也许可以认为这是Python设计的缺陷导致的，因为Python 3.0只需要用 `super().some_attr` 就行了。至于为什么需要两个参数，可以大概分析一下，如果有错，欢迎指正。
先介绍下背景知识：super是从Python 2.2开始随着 `new-style class` 一起引入的，也只能应用于所谓的 `new-style class` (即直接或间接继承于 `object` 的class)，可以在一定程度上解决所谓`钻石继承`的问题。2.2起所有python内置class都是new-style class。
那么super是什么呢？实际上它是一个builtin的type，你调用 `super(Base, self)` 会返回一个 `__class__ = super` 的对象，这个对象可以作为一个代理，通过BFS的顺序（实际上列表已经保存在`Base.__mro__`里了）去访问第一个**直接定义**了目标属性的基类里的属性。`super(type, obj_or_subtype).attr` 基本上可以认为是 `find_in_mro_provide_by(obj_or_subtype, "attr")` ，有点晦涩。
**注意**，super()返回的不是“父类对象”，而是一个super类型的对象；实例化的时候，第一个参数是一个类型(type)，第二个参数可以是type的instance，也可以是type的subclass。
介绍完基础知识，可以先说说为什么要有第二个参数。理由很简单 —— 我们知道 self 代表当前对象，每个对象的方法在定义的时候都需要显示地把self作为第一个参数，你本来应该写成这样
some_class.some_method(obj, *args, **kwargs)
但是因为Python语法允许 `obj.some_metod(*args, **kwargs)` (本质上是个语法糖) ，所以你可以写得简单点（不必显式地给出方法的第一个参数）。而super对象则不同，它没有语法上的直接支持，所以在内部invoke some_method的时候必须指定某个对象，而这个对象的得你自己塞给它，也就是说把这个super对象绑定(bound)到第二个参数上。所以实际上并不是非要用self作为super的第二个参数，甚至super并不是必须在class内部才能调用：
class A(object):
def foo(self):
print self.name
class B(A):
def __init__(self, name):
self.name = name
b1 = B('b1')
super(B, b1).foo() #produces 'b1'
def foo(self):
print self.name
class B(A):
def __init__(self, name):
self.name = name
b1 = B('b1')
super(B, b1).foo() #produces 'b1'
至于为什么要有第一个参数，如果你看了 `help(super)` 就会发现它居然提供了个单参数的版本：
super(type) -> unbound super object
这个理解起来比较困难点。先说这个`unbound`，没绑定，就是说这个super对象没有实际绑定到某个对象上，它并不是可以直接用的。要怎么用呢？那就得注意到`help(super)`里面有一个 `__get__` 方法，也就是说，super类还是一个Descriptor类！哎，这个Descriptor类展开又是好大一段，简单地说就是它可以把自己和某个object的属性绑定，使得访问这个属性的时候，实际上是在调用它的`__get__/__set__`方法：
class Descr(object):
def __get__(self, obj, tp): #最后一句传进来的 obj=x, tp=X
return 'get!'
class X(object):
t = Descr()
x = X()
print x.t #this produces 'get!'
def __get__(self, obj, tp): #最后一句传进来的 obj=x, tp=X
return 'get!'
class X(object):
t = Descr()
x = X()
print x.t #this produces 'get!'
回到unbound super, 为了使用它，就要通过 __get__ 方法把它再绑定到一个对象上。由于Descriptor的属性，特别适合这么用：
class A(object):
def foo(self):
print 'foo'
class B(A):
def foo(self):
self._super.foo()
B._super = super(B) #这还不能直接写在B的定义里，多蛋疼啊。
b = B()
b.foo() #produces 'foo'
def foo(self):
print 'foo'
class B(A):
def foo(self):
self._super.foo()
B._super = super(B) #这还不能直接写在B的定义里，多蛋疼啊。
b = B()
b.foo() #produces 'foo'
这是它的一个可能用法。我不确定是否还有其他更合适的用途，但是在这里实际上也并不是很好，尤其是遇到staticmethod的时候还会出错，再加上绕了这么大一个弯，实在不是很推荐使用。
据说unbound super的使用非常少，不知道现实意义有多少，也不知道当初为什么设计成这个样子，总之由于Python3.0已经改了（虽然仍然保留了unbound super），所以基本上还是可以认为这是设计的历史遗留问题。
参考文献（unbound super主要参考了这个系列）：Things to Know About Python Super
[1] http://www.artima.com/weblogs/viewpost.jsp?thread=236275
[2] http://www.artima.com/weblogs/viewpost.jsp?thread=236278
[3] http://www.artima.com/weblogs/viewpost.jsp?thread=237121
Aug
16
因为不满ssh tunnel的使用效果，所以2012年12月某天（大概是17号）心血来潮写了这个小东西，由于 [socks5协议]( http://www.openssh.com/txt/rfc1928.txt ) 本身很简单、加上gevent/greenlet使得异步开发跟同步似的，所以200行就搞定了。但是性能上问题很大——主要是加密有问题。尽管加密就是最简单的xor，但是因为python不适合处理大量的小对象，所以当时写了一个python扩展，性能上就没问题了，但是又多了一项麻烦的依赖。后来发现已经有更成熟的shadowsocks，于是就弃坑了，也一直没有发布。
今天[2013.08.16]心血来潮，用ctypes来实现同样的功能，似乎也挺合适的；不过跟shadowsocks比起来有两个地方做得不好，一是没有更“高级”的加密方式（他家用了M2Crypto，代码看起来很复杂），另一个是shadowsocks在本地先回应socks5请求，只把必要的host:port信息发送给server，减少了一个来回，而我原先的实现则是在server端实现完整的socks5(现在把step1搬到client了，因为改动很小)。
总之好歹也是个凑合能用的东西了，发布出来晾着吧，也许哪天有人就用上了呢。
项目地址： https://github.com/felix021/mixo
今天[2013.08.16]心血来潮，用ctypes来实现同样的功能，似乎也挺合适的；不过跟shadowsocks比起来有两个地方做得不好，一是没有更“高级”的加密方式（他家用了M2Crypto，代码看起来很复杂），另一个是shadowsocks在本地先回应socks5请求，只把必要的host:port信息发送给server，减少了一个来回，而我原先的实现则是在server端实现完整的socks5(现在把step1搬到client了，因为改动很小)。
总之好歹也是个凑合能用的东西了，发布出来晾着吧，也许哪天有人就用上了呢。
项目地址： https://github.com/felix021/mixo
Aug
16
为了扩展Python，我们可以[用C/C++编写模块](http://docs.python.org/2/extending/ )，但是这要求对Python的底层有足够的了解，包括Python对象模型、常用模块、引用计数等，门槛较高，且不方便利用现有的C库。而 [ctypes](http://docs.python.org/2/library/ctypes.html ) 则另辟蹊径，通过封装`dlopen/dlsym`之类的函数，并提供对C中数据结构的包装/解包，让Python能够加载动态库、导出其中的函数直接加以利用。
快速上手
-------
最简单的，我们可以从 libc 开始：
由于ctypes库能够直接处理 None, integers, longs, byte strings 和 unicode strings 类型的转换，因此在这里不需要任何额外的操作就能完成任务。（注意，最后一行的25是printf的返回值，标识输出了25个字符）
自己动手
-------
[这里](http://wolfprojects.altervista.org/articles/dll-in-c-for-python/)有一个最简单的例子——实现一个 `int sum(int a, int b)`：
编译之得到 foo.so:
$ gcc -fPIC -shared -o foo.so foo.c
然后在Python里头：
真是超级简单。
丰衣足食
-------
下面来个稍微复杂点的例子：反转一个字符串。尽管ctypes可以直接处理string的转换，但是你不能直接修改string里的内容，所以这里需要变通一下。可能的解决方案有两个：
1. 在foo.c里面malloc一点空间然后返回。这样可以不修改原string，但是却要考虑释放该空间的问题，不太好处理。
2. 让Python分配好空间，将原字符串拷贝进去，再让foo里面的函数来修改它。ctypes为这种方案提供了`create_string_buffer`这个函数，正适合。
那就让我们用方案2来实现它：
然后在Python里头：
这样就OK啦。
补充说明
-------
以上的操作都是在Linux下完成的，在Windows下基本上一致，除了windows不能直接load libc.so.6，而是应该找 msvcrt :
>>> libc = ctypes.cdll.msvcrt
或者直接导入文件
>>> libc = ctypes.CDLL("msvcr90.dll")
小结
-------
前面演示了ctypes库的简单使用，已经能够完成许多任务了；但是它可以完成更多的任务，包括操作更复杂的例如结构体、数组等C数据结构，具体的这里不细说了，可以参考详细的[官方文档](http://docs.python.org/2/library/ctypes.html)。
好了，就到这里。
快速上手
-------
最简单的，我们可以从 libc 开始：
felix021@ubserver:~/code$ python
Python 2.7.3 (default, Jul 5 2013, 08:39:51)
[GCC 4.6.3] on linux2
Type "help", "copyright", "credits" or "license" for more information.
>>> import ctypes
>>> libc = ctypes.CDLL("libc.so.6")
>>> libc.printf("hello, %s! %d + %d = %d\n", "world", 1, 2, 3)
hello, world! 1 + 2 = 3
25
Python 2.7.3 (default, Jul 5 2013, 08:39:51)
[GCC 4.6.3] on linux2
Type "help", "copyright", "credits" or "license" for more information.
>>> import ctypes
>>> libc = ctypes.CDLL("libc.so.6")
>>> libc.printf("hello, %s! %d + %d = %d\n", "world", 1, 2, 3)
hello, world! 1 + 2 = 3
25
由于ctypes库能够直接处理 None, integers, longs, byte strings 和 unicode strings 类型的转换，因此在这里不需要任何额外的操作就能完成任务。（注意，最后一行的25是printf的返回值，标识输出了25个字符）
自己动手
-------
[这里](http://wolfprojects.altervista.org/articles/dll-in-c-for-python/)有一个最简单的例子——实现一个 `int sum(int a, int b)`：
//file foo.c
#ifdef __WIN32
#define DLLEXPORT __declspec(dllexport)
#else
#define DLLEXPORT
#endif
DLLEXPORT int sum(int a, int b)
{
return a + b;
}
#ifdef __WIN32
#define DLLEXPORT __declspec(dllexport)
#else
#define DLLEXPORT
#endif
DLLEXPORT int sum(int a, int b)
{
return a + b;
}
编译之得到 foo.so:
$ gcc -fPIC -shared -o foo.so foo.c
然后在Python里头：
>>> foo = ctypes.CDLL("./foo.so")
>>> foo.sum(5, 3)
8
>>> foo.sum(5, 3)
8
真是超级简单。
丰衣足食
-------
下面来个稍微复杂点的例子：反转一个字符串。尽管ctypes可以直接处理string的转换，但是你不能直接修改string里的内容，所以这里需要变通一下。可能的解决方案有两个：
1. 在foo.c里面malloc一点空间然后返回。这样可以不修改原string，但是却要考虑释放该空间的问题，不太好处理。
2. 让Python分配好空间，将原字符串拷贝进去，再让foo里面的函数来修改它。ctypes为这种方案提供了`create_string_buffer`这个函数，正适合。
那就让我们用方案2来实现它：
//file foo.c
DLLEXPORT int reverse(char *str)
{
int i, j, t;
for (i = 0, j = strlen(str) - 1; i < j; i++, j--)
t = str[i], str[i] = str[j], str[j] = t;
return 0;
}
DLLEXPORT int reverse(char *str)
{
int i, j, t;
for (i = 0, j = strlen(str) - 1; i < j; i++, j--)
t = str[i], str[i] = str[j], str[j] = t;
return 0;
}
然后在Python里头：
>>> s = ctypes.create_string_buffer("hello world!")
>>> foo.reverse(s)
0
>>> print s.value
!dlrow olleh
>>> foo.reverse(s)
0
>>> print s.value
!dlrow olleh
这样就OK啦。
补充说明
-------
以上的操作都是在Linux下完成的，在Windows下基本上一致，除了windows不能直接load libc.so.6，而是应该找 msvcrt :
>>> libc = ctypes.cdll.msvcrt
或者直接导入文件
>>> libc = ctypes.CDLL("msvcr90.dll")
小结
-------
前面演示了ctypes库的简单使用，已经能够完成许多任务了；但是它可以完成更多的任务，包括操作更复杂的例如结构体、数组等C数据结构，具体的这里不细说了，可以参考详细的[官方文档](http://docs.python.org/2/library/ctypes.html)。
好了，就到这里。
Jul
30
简单地说，MPI 就是个并行计算框架，模型也很直接——就是多进程。和hadoop不同，它不提供计算任务的map和reduce，只提供了一套通信接口，需要程序员来完成这些任务；它也不提供冗余容错等机制，完全依赖于其下层的可靠性。但是因为把控制权几乎完全交给了程序员，所以有很大的灵活性，可以最大限度地榨取硬件性能。超级计算机上的运算任务，基本上都是使用MPI来开发的。
～ 下载编译安装：
现在貌似一般都用MPICH，开源的MPI库，可以从这里获取： http://www.mpich.org/ ，现在的最新版本是3.0.4，编译安装过程可以参考安装包里的README的说明，基本步骤如下（万恶的configure）：
$ wget http://www.mpich.org/static/downloads/3.0.4/mpich-3.0.4.tar.gz
$ tar zxf mpich-3.0.4.tar.gz
$ cd mpich-3.0.4
$ mkdir ~/mpich
$ ./configure --prefix=$HOME/mpich --disable-f77 --disable-fc 2>&1 | tee c.txt #我禁用了fortran的支持
$ make -j4 2>&1 | tee m.txt
$ make install 2>&1 | tee i.txt
$ echo 'export PATH=$PATH:~/mpich/bin' >> ~/.bashrc
下面给出三个例子，参考教程：http://wenku.baidu.com/view/ee8bf3390912a216147929f3.html （注：22页有BUG，它把 MPI_Comm_XXX 错写成了 MPI_Common_xxx //包括全大写版本，共四处），给出了MPI框架中最常用、最基础的6个API的例子。更复杂的API可以参考mpich的手册。这些例子只是简单地演示了MPI框架的使用；实际上在使用MPI开发并行计算的软件时，还需要考虑到很多方面的问题，这里就不展开说了（其实真相是我也不会-.-，有兴趣的话可以请教 @momodi 和 @dumbear 两位）。
1. 最简单的：Hello world
代码如下： hello.c
编译：
$ mpicc -o hello hello.c
运行：
$ mpiexec -n 4 ./hello
Hello world!
Hello world!
Hello world!
Hello world!
可以看到这里启动了4个进程。注意 -n 和 4 之间一定要有空格，否则会报错。
2. 进程间通信
MPI最基本的通信接口是 MPI_Send/MPI_Recv:
编译运行：
$ mpicc comm.c
$ mpiexec -n 4 ./a.out
I'm 0 of 4
from 1: hello from 1
from 2: hello from 2
I'm 1 of 4
I'm 2 of 4
from 3: hello from 3
I'm 3 of 4
3. 来个复杂点的：数数前1亿个自然数里有几个 雷劈数
代码后附，答案是97（真少），不过这不是重点，重点是MPI对硬件的利用率是怎样 :D
测试机器是 16核 AMD Opteron 6128HE @2GHz，32G内存
单进程（无MPI版本）：56.9s
4进程：15.3s
8进程：7.85s
12进程：5.35s
考虑到16核跑满可能会受到其他进程的影响（性能不稳定，4.2～4.9s），这个数据就不列进来比较了。
可以看出来，在这个例子里，因为通信、同步只有在计算完之后才有那么一点点，所以在SMP架构下，耗费的时间基本上是跟进程数成反比的，说明MPI框架对硬件性能的利用率还是相当高的。
具体代码如下：
～ 下载编译安装：
现在貌似一般都用MPICH，开源的MPI库，可以从这里获取： http://www.mpich.org/ ，现在的最新版本是3.0.4，编译安装过程可以参考安装包里的README的说明，基本步骤如下（万恶的configure）：
引用
$ wget http://www.mpich.org/static/downloads/3.0.4/mpich-3.0.4.tar.gz
$ tar zxf mpich-3.0.4.tar.gz
$ cd mpich-3.0.4
$ mkdir ~/mpich
$ ./configure --prefix=$HOME/mpich --disable-f77 --disable-fc 2>&1 | tee c.txt #我禁用了fortran的支持
$ make -j4 2>&1 | tee m.txt
$ make install 2>&1 | tee i.txt
$ echo 'export PATH=$PATH:~/mpich/bin' >> ~/.bashrc
下面给出三个例子，参考教程：http://wenku.baidu.com/view/ee8bf3390912a216147929f3.html （注：22页有BUG，它把 MPI_Comm_XXX 错写成了 MPI_Common_xxx //包括全大写版本，共四处），给出了MPI框架中最常用、最基础的6个API的例子。更复杂的API可以参考mpich的手册。这些例子只是简单地演示了MPI框架的使用；实际上在使用MPI开发并行计算的软件时，还需要考虑到很多方面的问题，这里就不展开说了（其实真相是我也不会-.-，有兴趣的话可以请教 @momodi 和 @dumbear 两位）。
1. 最简单的：Hello world
代码如下： hello.c
#include <stdio.h>
#include <mpi.h>
int main(int argc, char *argv[])
{
MPI_Init(&argc, &argv); //初始化MPI环境
printf("Hello world!\n");
MPI_Finalize(); //结束MPI环境
return 0;
}
#include <mpi.h>
int main(int argc, char *argv[])
{
MPI_Init(&argc, &argv); //初始化MPI环境
printf("Hello world!\n");
MPI_Finalize(); //结束MPI环境
return 0;
}
编译：
$ mpicc -o hello hello.c
运行：
$ mpiexec -n 4 ./hello
Hello world!
Hello world!
Hello world!
Hello world!
可以看到这里启动了4个进程。注意 -n 和 4 之间一定要有空格，否则会报错。
2. 进程间通信
MPI最基本的通信接口是 MPI_Send/MPI_Recv:
#include <stdio.h>
#include <mpi.h>
int main(int argc, char *argv[])
{
int myid, numprocs, source, msg_tag = 0;
char msg[100];
MPI_Status status;
MPI_Init(&argc, &argv);
MPI_Comm_size(MPI_COMM_WORLD, &numprocs); //共启动几个进程
MPI_Comm_rank(MPI_COMM_WORLD, &myid); //当前进程的编号（0～n-1）
printf("I'm %d of %d\n", myid, numprocs);
if (myid != 0)
{
int len = sprintf(msg, "hello from %d", myid);
MPI_Send(msg, len, MPI_CHAR, 0, msg_tag, MPI_COMM_WORLD); //向id=0的进程发送信息
}
else
{
for (source = 1; source < numprocs; source++)
{
//从id=source的进程接受消息
MPI_Recv(msg, 100, MPI_CHAR, source, msg_tag, MPI_COMM_WORLD, &status);
printf("from %d: %s\n", source, msg);
}
}
MPI_Finalize();
return 0;
}
#include <mpi.h>
int main(int argc, char *argv[])
{
int myid, numprocs, source, msg_tag = 0;
char msg[100];
MPI_Status status;
MPI_Init(&argc, &argv);
MPI_Comm_size(MPI_COMM_WORLD, &numprocs); //共启动几个进程
MPI_Comm_rank(MPI_COMM_WORLD, &myid); //当前进程的编号（0～n-1）
printf("I'm %d of %d\n", myid, numprocs);
if (myid != 0)
{
int len = sprintf(msg, "hello from %d", myid);
MPI_Send(msg, len, MPI_CHAR, 0, msg_tag, MPI_COMM_WORLD); //向id=0的进程发送信息
}
else
{
for (source = 1; source < numprocs; source++)
{
//从id=source的进程接受消息
MPI_Recv(msg, 100, MPI_CHAR, source, msg_tag, MPI_COMM_WORLD, &status);
printf("from %d: %s\n", source, msg);
}
}
MPI_Finalize();
return 0;
}
编译运行：
$ mpicc comm.c
$ mpiexec -n 4 ./a.out
I'm 0 of 4
from 1: hello from 1
from 2: hello from 2
I'm 1 of 4
I'm 2 of 4
from 3: hello from 3
I'm 3 of 4
3. 来个复杂点的：数数前1亿个自然数里有几个 雷劈数
代码后附，答案是97（真少），不过这不是重点，重点是MPI对硬件的利用率是怎样 :D
测试机器是 16核 AMD Opteron 6128HE @2GHz，32G内存
单进程（无MPI版本）：56.9s
4进程：15.3s
8进程：7.85s
12进程：5.35s
考虑到16核跑满可能会受到其他进程的影响（性能不稳定，4.2～4.9s），这个数据就不列进来比较了。
可以看出来，在这个例子里，因为通信、同步只有在计算完之后才有那么一点点，所以在SMP架构下，耗费的时间基本上是跟进程数成反比的，说明MPI框架对硬件性能的利用率还是相当高的。
具体代码如下：
#include <stdio.h>
#include <mpi.h>
int is_lp(long long x)
{
long long t = x * x, i = 10;
while (i < t)
{
long long l = t / i, r = t % i;
if (l + r == x)
return 1;
i *= 10;
}
return 0;
}
int main(int argc, char *argv[])
{
int myid, numprocs, source;
const int N = 100000000;
MPI_Status status;
MPI_Init(&argc, &argv);
MPI_Comm_rank(MPI_COMM_WORLD, &myid);
MPI_Comm_size(MPI_COMM_WORLD, &numprocs);
printf("I'm %d of %d\n", myid, numprocs);
int start = myid * (N / numprocs), stop = (myid + 1) * (N / numprocs);
if (myid == numprocs - 1)
stop = N;
printf("start from %d to %d\n", start, stop);
int ans = 0, i;
for (i = start; i < stop; i++)
if (is_lp(i))
ans += 1;
printf("%d finished calculation with %d numbers\n", myid, ans);
if (myid != 0)
{
MPI_Send(&ans, 1, MPI_INT, 0, 0, MPI_COMM_WORLD);
}
else
{
int tmp;
for (source = 1; source < numprocs; source++)
{
MPI_Recv(&tmp, 1, MPI_INT, source, 0, MPI_COMM_WORLD, &status);
printf("from %d: %d\n", source, tmp);
ans += tmp;
}
printf("final ans: %d\n", ans);
}
MPI_Finalize();
return 0;
}
#include <mpi.h>
int is_lp(long long x)
{
long long t = x * x, i = 10;
while (i < t)
{
long long l = t / i, r = t % i;
if (l + r == x)
return 1;
i *= 10;
}
return 0;
}
int main(int argc, char *argv[])
{
int myid, numprocs, source;
const int N = 100000000;
MPI_Status status;
MPI_Init(&argc, &argv);
MPI_Comm_rank(MPI_COMM_WORLD, &myid);
MPI_Comm_size(MPI_COMM_WORLD, &numprocs);
printf("I'm %d of %d\n", myid, numprocs);
int start = myid * (N / numprocs), stop = (myid + 1) * (N / numprocs);
if (myid == numprocs - 1)
stop = N;
printf("start from %d to %d\n", start, stop);
int ans = 0, i;
for (i = start; i < stop; i++)
if (is_lp(i))
ans += 1;
printf("%d finished calculation with %d numbers\n", myid, ans);
if (myid != 0)
{
MPI_Send(&ans, 1, MPI_INT, 0, 0, MPI_COMM_WORLD);
}
else
{
int tmp;
for (source = 1; source < numprocs; source++)
{
MPI_Recv(&tmp, 1, MPI_INT, source, 0, MPI_COMM_WORLD, &status);
printf("from %d: %d\n", source, tmp);
ans += tmp;
}
printf("final ans: %d\n", ans);
}
MPI_Finalize();
return 0;
}
Jul
29
nth_element是一个用烂了的面试题，09年的时候我也曾经跑过一点数据（回头一看好像有点不太对，那时候CPU有那么慢吗？应该是当时没有把读取数据的时间给去掉吧），昨天看到 从一道笔试题谈算法优化（上） 和 从一道笔试题谈算法优化（下） 这个系列，里面提到了从简单到复杂的6种算法，根据作者的说法，经过各种优化以后solution6达到了相当的效率。受到作者启发，我也想到了一种算法，于是花了些时间，把常用的几种算法实现了进行对比（包括对比stl里的nth_element）。
回到 nth_element 的问题：从 n 个数里头，找出其中的 top k （top可以是找最大 也可以是找最小）。简单起见，后面都以找最小为例。
那篇文章里的后面的3种算法和我的算法都是基于它的solution3进行优化的，所以先介绍一下solution3：
1. 将 n 个数的前 k 个拷贝出来
2. 找出这 k 个数的最大值 m
3. 对于每一个剩下的 n - k 个数 i ：如果 i 小于 m ，将 m 替换为 i ，然后再找出新的 k 个数中的最大值 m
4. 返回过滤后得到的 k 个数
虽然第3步的判断条件成立时，这一步是 O(k) 的复杂度，但是在大部分情况下，这个判断条件都不成立，所以它需要执行的次数很少（这种效果越往后越明显，因为 m 的值越来越小，可以滤掉更多的数）。但是对于极端情况（或者接近极端的情况）——如果数组完全是降序排列的话，那么这个条件每次都会成立，会导致算法退化成 O(n * m) ，性能就不可接受了。
solution4~6的优化我觉得读起来有点晦涩，有兴趣的同学可以自己去看，这里不展开说了。它们都存在算法退化的情况。
我的算法可能理解起来要简单些，基本思路是偷懒：把缓冲区设置为 2 * k ，在扫数组的时候，如果这个数比当前保留的最大值还要小，就把它塞进缓冲区，直到缓冲区塞满了，再排序、取最大值、删多余的数。
1. 开辟一个 2 × k 的临时数组 r
2. 将 N[1..k] 拷贝到 r[1..k]，并记录 r 的末尾 e = k + 1
3. 找出 r[1..k] 中的最大值 m
4. 对于每一个剩下的 n - k 个数 i：
4.1 如果 i 小于 m，将 m 塞到 r 的末尾 r[e]；e = e + 1
4.1.1 如果 e == 2×k（r被塞满了），对 r 进行排序，得出前k个数中的最大值；抛弃后面的k个数，e = k + 1
5. 对 r[1..e] 进行排序
6. 返回 r[1..k]
这个算法的实际效果很好，在 n = 1亿、k = 1万 的情况下，对于随机的n个数，4.1条件成立的次数通常只有十几次，几乎可以忽略，因此大约只需要扫描整个数组时间的2～3倍即可。
但是这个算法同样存在退化的问题：如果全都是倒序排列的话，也变成O(n*m)了。幸而解决方案也很简单：对这个数组进行 n / 100 次的随机化（考虑到随机化耗时也比较大，100这个数字是试了几次以后发现的合适值，不一定最优，但是都差不多了）：
p.s. 很不幸的是随机化这个方法对于solution3来说虽然提升也很明显，但是远远不够。solution4～6也不太行。
这里给出随机情况下和完全倒序情况下的性能对比吧（ubuntu 12.04 x86_64，i5 2400@3.1G）：
补充说明一下，nth_element算法在面试的时候可能给出的 n 会更大，以至于在内存中存不下，这时候通常会认为用堆来实现最好——因为只要O(k)的空间就可以搞定，而且最差时间复杂度是O(logk * n)，但是要注意logk的常数实际上是很大的，所以你可以看到前面的数据里 堆算法 的效率并不是最好的。事实上这里给出的所有算法（不考虑随机化的话）也都只需要O(k)的空间；如果需要随机化的话也很简单，分段处理，只要保证每段能在内存中保存下来就行了。
最后上代码存档（为了方便测试用了些全局变量，看起来可能有点挫）：
回到 nth_element 的问题：从 n 个数里头，找出其中的 top k （top可以是找最大 也可以是找最小）。简单起见，后面都以找最小为例。
那篇文章里的后面的3种算法和我的算法都是基于它的solution3进行优化的，所以先介绍一下solution3：
1. 将 n 个数的前 k 个拷贝出来
2. 找出这 k 个数的最大值 m
3. 对于每一个剩下的 n - k 个数 i ：如果 i 小于 m ，将 m 替换为 i ，然后再找出新的 k 个数中的最大值 m
4. 返回过滤后得到的 k 个数
虽然第3步的判断条件成立时，这一步是 O(k) 的复杂度，但是在大部分情况下，这个判断条件都不成立，所以它需要执行的次数很少（这种效果越往后越明显，因为 m 的值越来越小，可以滤掉更多的数）。但是对于极端情况（或者接近极端的情况）——如果数组完全是降序排列的话，那么这个条件每次都会成立，会导致算法退化成 O(n * m) ，性能就不可接受了。
solution4~6的优化我觉得读起来有点晦涩，有兴趣的同学可以自己去看，这里不展开说了。它们都存在算法退化的情况。
我的算法可能理解起来要简单些，基本思路是偷懒：把缓冲区设置为 2 * k ，在扫数组的时候，如果这个数比当前保留的最大值还要小，就把它塞进缓冲区，直到缓冲区塞满了，再排序、取最大值、删多余的数。
1. 开辟一个 2 × k 的临时数组 r
2. 将 N[1..k] 拷贝到 r[1..k]，并记录 r 的末尾 e = k + 1
3. 找出 r[1..k] 中的最大值 m
4. 对于每一个剩下的 n - k 个数 i：
4.1 如果 i 小于 m，将 m 塞到 r 的末尾 r[e]；e = e + 1
4.1.1 如果 e == 2×k（r被塞满了），对 r 进行排序，得出前k个数中的最大值；抛弃后面的k个数，e = k + 1
5. 对 r[1..e] 进行排序
6. 返回 r[1..k]
这个算法的实际效果很好，在 n = 1亿、k = 1万 的情况下，对于随机的n个数，4.1条件成立的次数通常只有十几次，几乎可以忽略，因此大约只需要扫描整个数组时间的2～3倍即可。
但是这个算法同样存在退化的问题：如果全都是倒序排列的话，也变成O(n*m)了。幸而解决方案也很简单：对这个数组进行 n / 100 次的随机化（考虑到随机化耗时也比较大，100这个数字是试了几次以后发现的合适值，不一定最优，但是都差不多了）：
void rand_factor()
{
//randomization
const int factor = 100;
int r, t, i;
for (i = n / factor; i > 0; i--)
{
r = rand() % n;
t = a[0], a[0] = a[r], a[r] = t;
}
}
{
//randomization
const int factor = 100;
int r, t, i;
for (i = n / factor; i > 0; i--)
{
r = rand() % n;
t = a[0], a[0] = a[r], a[r] = t;
}
}
p.s. 很不幸的是随机化这个方法对于solution3来说虽然提升也很明显，但是远远不够。solution4～6也不太行。
这里给出随机情况下和完全倒序情况下的性能对比吧（ubuntu 12.04 x86_64，i5 2400@3.1G）：
引用
随机数据
========
随机化生成数据： 7.087s
遍历耗时： 0.053s
1/100随机化耗时： 0.064s
STL的nth_element： 0.841s
最大堆+随机化： 0.162s
Solution3+随机化： 3.121s
solution4+随机化： 2.825s
solution6+随机化： 0.249s
我的算法+随机化： 0.163s
倒序数据
========
无随机化生成数据： 0.150s
遍历耗时： 0.052s
1/100随机化耗时： 0.067s
STL的nth_element： 1.420s
最大堆+随机化： 0.301s
Solution3+随机化： 67.552s
solution4+随机化： 66.405s
solution6+随机化： 2.388s
我的算法+随机化： 0.170s
========
随机化生成数据： 7.087s
遍历耗时： 0.053s
1/100随机化耗时： 0.064s
STL的nth_element： 0.841s
最大堆+随机化： 0.162s
Solution3+随机化： 3.121s
solution4+随机化： 2.825s
solution6+随机化： 0.249s
我的算法+随机化： 0.163s
倒序数据
========
无随机化生成数据： 0.150s
遍历耗时： 0.052s
1/100随机化耗时： 0.067s
STL的nth_element： 1.420s
最大堆+随机化： 0.301s
Solution3+随机化： 67.552s
solution4+随机化： 66.405s
solution6+随机化： 2.388s
我的算法+随机化： 0.170s
补充说明一下，nth_element算法在面试的时候可能给出的 n 会更大，以至于在内存中存不下，这时候通常会认为用堆来实现最好——因为只要O(k)的空间就可以搞定，而且最差时间复杂度是O(logk * n)，但是要注意logk的常数实际上是很大的，所以你可以看到前面的数据里 堆算法 的效率并不是最好的。事实上这里给出的所有算法（不考虑随机化的话）也都只需要O(k)的空间；如果需要随机化的话也很简单，分段处理，只要保证每段能在内存中保存下来就行了。
最后上代码存档（为了方便测试用了些全局变量，看起来可能有点挫）：
Jul
11
在模板里引入其他模板应该是很常见的一种需求，但是webpy默认的template居然没有提供这种机制，挺神奇的。
官方的解决办法是把某个模板的输出给另一个模板，看起来和用起来都超级不爽；知乎网友给出的解决方案也很不爽：用render初始化时指定的layout，但是这个跟include差很多，不灵活。
实际上，由于render渲染模板后的输出本身是一个字符串，所以如果能在模板里头直接调用render渲染其他模板就最好了：而且这是可以实现的，只是略带tricky。
由于webpy给render的默认globals是空的，所以模板里只能用基本的python语法，默认连builtin的东西都用不了（比如zip、str），但是可以通过初始化render时指定globals的方式来引入：
所以我们只需要把render自己也加到这个globals里头去，就可以在模板里引用它了：
第三句看起来虽然tricky，但是由于python的对象传的是引用，所以可以达到预期的效果。
这样只要在模板里这样写就行了：
官方的解决办法是把某个模板的输出给另一个模板，看起来和用起来都超级不爽；知乎网友给出的解决方案也很不爽：用render初始化时指定的layout，但是这个跟include差很多，不灵活。
实际上，由于render渲染模板后的输出本身是一个字符串，所以如果能在模板里头直接调用render渲染其他模板就最好了：而且这是可以实现的，只是略带tricky。
由于webpy给render的默认globals是空的，所以模板里只能用基本的python语法，默认连builtin的东西都用不了（比如zip、str），但是可以通过初始化render时指定globals的方式来引入：
render = web.template.render('view/', cache=False, globals = __builtins__.__dict__)
所以我们只需要把render自己也加到这个globals里头去，就可以在模板里引用它了：
render_globals = {}
render = web.template.render('view/', cache=False, globals = render_globals)
render_globals['render'] = render
render = web.template.render('view/', cache=False, globals = render_globals)
render_globals['render'] = render
第三句看起来虽然tricky，但是由于python的对象传的是引用，所以可以达到预期的效果。
这样只要在模板里这样写就行了：
$:render.header('首页') #那个冒号别漏掉:)
Jul
3
先吐槽一下，顺便说下起因。
去年9月买的这个ML-G3000，青轴全无冲（回想起来有点2，青轴又不是拿来玩游戏的，买什么全无冲。。还贵了50块……），419入手的。它的键位基本上是标准的，但是R-Win变成了“禁用Win”功能键（直接在键盘电路逻辑里实现的，OS里检测不到keyCode的）、右键菜单（属性）键变成了功能键（配合Fx实现多媒体/音量控制什么的）。禁用Win对于游戏控来说是挺不错了，但是对于我这种几乎不怎么玩PC游戏的人来说就太废了点，而且还经常把这个当成Win键误按，例如我要Win+L锁屏，经常就变成了禁用Win，然后再换"左Win"+L，可是这时候Win已经被禁用了……各种抓狂。于是就用拔键器拔了一下，没想到键帽和轴结合太紧，直接把轴给拔坏了，里面弹簧蹦了出来。然后心里很不爽，就拿CapsLock开刀，果然又拔坏了……于是就成了这个样子：
去年9月买的这个ML-G3000，青轴全无冲（回想起来有点2，青轴又不是拿来玩游戏的，买什么全无冲。。还贵了50块……），419入手的。它的键位基本上是标准的，但是R-Win变成了“禁用Win”功能键（直接在键盘电路逻辑里实现的，OS里检测不到keyCode的）、右键菜单（属性）键变成了功能键（配合Fx实现多媒体/音量控制什么的）。禁用Win对于游戏控来说是挺不错了，但是对于我这种几乎不怎么玩PC游戏的人来说就太废了点，而且还经常把这个当成Win键误按，例如我要Win+L锁屏，经常就变成了禁用Win，然后再换"左Win"+L，可是这时候Win已经被禁用了……各种抓狂。于是就用拔键器拔了一下，没想到键帽和轴结合太紧，直接把轴给拔坏了，里面弹簧蹦了出来。然后心里很不爽，就拿CapsLock开刀，果然又拔坏了……于是就成了这个样子：
