遥远的街市

这个 blog 用到的一些算法

henix — Fri, 11 May 2012 17:34:35 +0800

　　开发个人博客除了写写 HTML/CSS/Javascript ，也是会遇到不少算法问题的。下面是我对一些问题的思考，以及相关的方案：

1. 热门度 ranking

　　即如何确定热门文章的排序。用于度量文章热门程度的主要有两项指标：点击数和评论数。我看到不少博客程序的办法是，有两个列表，一个是按点击量的 top 10 ，一个是按评论数的 top 10 。或者有些网站是只有评论数 top 10 。

　　我这个 blog 采用的是一个很简单的加权：热门度 = 点击数 + 7 * 评论数。

　　但这样也有一些问题，比如，“热门文章”一栏很少改变，排名第一的文章是历史上获得点击最多的文章，它在那里已经呆了将近一年了，后来的文章很难反超。于是可以考虑的改进是：引入时间因素。可以按时间加权，离现在越近的权值越高。

　　其实本来我以前也有一个“本月热门”，本月热门是只用当月数据排序出来的结果，反映了当前的热门。但后来发现这个结果基本上跟首页是一样的，读者还不如直接看首页。另一方面，博客跟新闻类网站不同，不是实时的，发文频率也很低，其实没必要这么搞。

　　目前的加权方式相当于所有时间点上的数据的权都是一样的，比较照顾历史数据。但是这样也有好处：有些东西的价值是不会随着时间的流逝而减少的。这样的加权照顾了经典文章。从最近的情况看，真正的有价值的文章也是有机会从默默无闻冲进前 10 的（从搜索引擎过来的流量帮了不少忙）。

　　所以各种加权方式也是各有优劣，究竟选用哪一种还需要仔细比较。

　　阮一峰前不久有几篇介绍 ranking 的文章可供参考：基于用户投票的排名算法（六）：贝叶斯平均。

2. 相关文章

　　两篇文章的相关度由它们的 tag 决定。直观上，两篇文章共有的 tag 越多，它们就越相关。所以一种直接的做法是：相关度 = 公共 tag 数。

　　我采用的算法是经典的余弦相似度：

$$ 相似度(A,B) = \frac{AB 共同拥有的 tag 数量}{\sqrt{A的tag数}\sqrt{B的tag数}} $$

　　这个算法不但考虑了公共的 tag 数，也考虑到另外一个因素：AB 单独拥有的 tag 越多，则它们越不相关。

　　效率上，由于要同时计算所有文章与其他每篇文章的相似度，时间是 n² 的，也有一些算法来改善这个时间，但需要用 Jaccard 相似度，以后有时间研究一下。

　　其他的玩法还可以：对每篇文章进行中文分词，然后用 TF-IDF 提取文章的关键词，这可以用于 auto tagging（现在我的 tag 都是手动添加的）。不过只是一个个人博客做到这种程度。。我还没有那么蛋疼...

　　当然，实践上，也可以用无觅的推荐插件，不过我觉得自己做的可以随时调整，准确率更高。

3. HTML 摘要

　　我的首页并没有显示每篇文章的全文，只是显示了开头的一段。这就需要设计一种 HTML 摘要算法。

　　第一种想法是，先去掉所有 tag ，再按字符串长度截整。比如 CSDN 博客就是这样做的。这样算法简单，但问题是，文章的格式信息完全被删除了，文章的内容全部都挤在一段里，看着觉得排版很乱。

　　另一种办法是，对 HTML 进行解析，取出前面的几个标签，直到总长度超过一个阈值。这种方案编程复杂，需要考虑各种特殊情况。

　　我的方案是：考虑对第一个算法进行改进。实际上读者并不需要标签跟原文一模一样，只需要换行跟原文一样就行了。于是我设计了这样一个算法：

将 p 的结束符 </p> 替换成一个不会在其他地方出现的特殊字符，如 0x7f
将 <br /> 替换成 0x7f
删去所有标签
将 0x7f 替换成 <br />

　　这样，标签跟原文不一样，但原文换了行的地方，摘要里也换了行。读者看的时候就能得到跟原文差不多的效果了。目前我的博客首页的摘要就是由这个算法生成的。

　　在具体实现的时候，还遇到过一些其他问题。比如，如果最后一个字符恰好是 HTML 转义符，如 & ，截断的位置不能出现在 HTML 实体的中间。所以，如果遇到了 & 就需要把整个 HTML 实体完整地“吞下”。Lua 代码：

function digestHTML(str, length)
	local ret = string.gsub(str, "</p>", string.char(127))
	ret = string.gsub(ret, "<br />", string.char(127))
	ret = string.gsub(ret, "<.->", "")

	-- truncate UTF-8 string
	local i = 1
	local curlen = 0
	local len = string.len(ret)
	while curlen < length and i <= len do
		local t = string.byte(ret, i)
		if t > 127 then
			i = i + 3
		else
			if t == string.byte('&') then
				while string.byte(ret, i) ~= string.byte(';') and i < len do
					i = i + 1
				end
			end
			i = i + 1
		end
		curlen = curlen + 1
	end

	return string.gsub(string.sub(ret, 1, i - 1), string.char(127), "<br />")
end

[bookmarklet]人人网批量已读脚本

henix — Sun, 15 Apr 2012 18:20:18 +0800

　　这个脚本可以将当前页面上已加载的所有新鲜事置为“已读”。

　　注意这个跟人人网自带的最下面的“全部标记为已读”是不同的。“全部标记为已读”是将你的全部新鲜事——不管是否已经在当前页面加载，标记为已读。我这个是只清除当前已加载的新鲜事。

人人-全部已读

用法

　　将上面的链接拖到书签栏。然后在人人网的页面点击书签栏的按钮，会提示你是否确定。如果确定，最后会显示已清除的新鲜事的数量。

　　对于 IE 用户：此脚本只支持 IE8 及以上，而且由于 Windows 上 IE 不支持 UTF-8 编码的字符串，会出现乱码。请使用这个英文版：人人-全部已读。Windows 上的其他浏览器（Chrome / Firefox / Opera）都可以正常工作。

源代码

　　实际上就是模拟点击每条新鲜事右边的“X”。

(function () {
	function fireEvent(element, eventType) {
		var evt;
		if (document.createEvent) {
			evt = document.createEvent("HTMLEvents");
			evt.initEvent(eventType, true, true); // type, bubbling, cancelable
			return !element.dispatchEvent(evt);
		} else {
			evt = document.createEventObject();
			return element.fireEvent('on' + eventType, evt);
		}
	}
	var goon = confirm('确定？');
	if (goon) {
		var feedlist = document.querySelector('div.feed-list');
		var deletes = feedlist.querySelectorAll('a.delete');
		var len = deletes.length;
		for (var i = 0; i < len; i++) {
			fireEvent(deletes[i], 'click');
		}
		alert('已清除 ' + len + ' 条新鲜事');
	}
})();

一些我使用和开发的油猴脚本

henix — Sun, 15 Apr 2012 17:08:38 +0800

Google Real Link: 觉得 Google 搜索的重定向很烦？这个脚本可以让 Google 的搜索结果不包含重定向，点击链接直接到对应网站。
豆沙绿#C7EDCC: 我在之前的介绍 awesome 的文章提到，浏览器的背景太白，感觉很刺眼。这个脚本可以把所有网页的白色部分改成豆沙绿。
t.co Bypasser[modified]: 让 twitter 上的短地址不再经过 t.co 。

　　下面是我开发的 userscripts（当然也都在用）：

Douban Timeline Marker: 为豆瓣的友邻广播增加一个书签，下次就很容易找到上次看到哪儿了。
Weibo Bookmark: 为微博增加一个书签，功能同上。

什么是油猴脚本，我为什么要用它？

　　扩展浏览器功能的方式，有这么几种：

各个浏览器的扩展（Extensions）或称插件（Plugins/Add-ons）：功能很强大，浏览器间不能通用
油猴脚本，或称 userscripts：就是一段 javascript ，当访问指定的网站时由浏览器加载，浏览器间可以通用
bookmarklets：一段 javascript 直接保存在书签栏里，要用时手动点

　　可见油猴脚本具有很高的通用性，并且适用于那些需要跟目标网站一起加载的任务。

　　如果一个任务可以用 JavaScript 搞定，就可以写成油猴脚本或者 bookmarklets 而不需要浏览器插件。

如何安装油猴脚本？

chrome 用户直接点。实际上，当打开一个以 .user.js 结尾的链接的时候就会提示安装，chrome 是把它当扩展安装了。
firefox 用户需要安装扩展 GreaseMonkey 。
opera 用户打开 Settings -> Advanced -> Content -> JavaScript Options -> User JavaScript Folder 然后选择一个文件夹并把 userjs 放在那个文件夹下面。

树形 dp 初探

henix — Wed, 04 Apr 2012 00:06:51 +0800

　　树形 dp 就是在一颗树上做 dp ，其基本思想是由子节点的信息推出父节点的信息。

　　由于树都是一般树，可能不只二叉，所以要用到一般树的“左儿子右兄弟”表示法（详见代码中的 first_child 和 next_sibling）。

hdu 1520 Anniversary party

　　最基本的树形 dp 。题目大意是，有一群人之间有上下级关系，在一个 party 中，有直接的上下级关系（即树中的父子关系）的人不能同时出席，每个人都有个 rating ，闻如何选择出席的人，使得所有人的 rating 之和最大。

　　每个节点有两个值，表示这个节点及其子节点所能取得的最大 rating 。max_with 是选择此节点的情况，max_without 不选择此节点。

　　有状态转移方程：

选择此节点的最大值 = 所有子节点不选择的最大值之和

不选择此节点的最大值 = 每个子节点选择或者不选择的最大值之和

　　说得太绕了...还是直接看代码吧

　　另外用了一个技巧，由于题目中没说谁是跟，用了一种很 tricky 的方法得到根节点的 id ，详见 root_id

　　另外我也尝试过非递归的，用了发现非递归反而比递归更耗时，所以就用递归的写法就好。

g++ 31ms

#include <stdio.h>
#include <string.h>
#include <stdlib.h>

#define N 6000
#define MAX(a, b) ((a)>(b)?(a):(b))

struct node_t {
	struct node_t *first_child;
	struct node_t *next_sibling;
	int max_with;
	int max_without;
};

struct node_t nodes[N];

void dp(struct node_t *node) {
	struct node_t *child_node = node->first_child;
	// node->max_with = node->rating;
	// node->max_without = 0;
	while (child_node != NULL) {
		dp(child_node);
		node->max_with += child_node->max_without;
		node->max_without += MAX(child_node->max_with, child_node->max_without);
		child_node = child_node->next_sibling;
	}
}

int getint(char end) {
	int s = 0;
	int ch;
	ch = getchar();
	while (ch != end && ch != EOF) {
		s = s * 10 + ch - '0';
		ch = getchar();
	}
	return s;
}

int main(int argc, char const *argv[])
{
	int n;
	while (scanf("%d", &n) == 1) {
		getchar(); // ignore a \n
		memset(nodes, 0, sizeof(struct node_t) * n);
		int i;
		for (i = 0; i < n; ++i) {
			// scanf("%d", &nodes[i].max_with);
			nodes[i].max_with = getint('\n');
		}
		int l, k;
		int root_id = (n - 1) * n / 2;
		while(1) {
			// scanf("%d%d", &l, &k);
			l = getint(' ');
			k = getint('\n');
			if (l == 0 && k == 0) {
				break;
			}
			l--;
			k--;
			/* add the parent-child relation: k - parent, l - child */
			nodes[l].next_sibling = nodes[k].first_child;
			nodes[k].first_child = nodes + l;
			root_id -= l;
		}
		dp(nodes + root_id);
		printf("%d\n", MAX(nodes[root_id].max_with, nodes[root_id].max_without));
	}
	return 0;
}

poj 1463 Strategic game

　　跟上题基本上是对偶的，不过这题的输入格式比较变态。

gcc 16ms

#include <stdio.h>
#include <string.h>
#include <stdlib.h>

#define N 1500
#define MIN(a, b) ((a)<(b)?(a):(b))

struct node_t {
	struct node_t *first_child;
	struct node_t *next_sibling;
	int min_with;
	int min_without;
};

struct node_t nodes[N];

void dp(struct node_t *node) {
	struct node_t *child_node = node->first_child;
	node->min_with = 1;
	while (child_node != NULL) {
		dp(child_node);
		node->min_without += child_node->min_with;
		node->min_with += MIN(child_node->min_with, child_node->min_without);
		child_node = child_node->next_sibling;
	}
}

int getint(char end) {
	// static char buf[20];
	// int i = 0;
	int s = 0;
	int ch;
	ch = getchar();
	while (ch != end && ch != EOF) {
		// buf[i++] = ch;
		s = s * 10 + ch - '0';
		ch = getchar();
	}
	// buf[i] = '\0';
	return s;
}

int getint2(int *out) {
	// static char buf[20];
	// int i = 0;
	register int s = 0;
	register int ch;
	ch = getchar();
	while (ch < '0' || ch > '9') {
		ch = getchar();
	}
	while (ch >= '0' && ch <= '9') {
		// buf[i++] = ch;
		s = s * 10 + ch - '0';
		ch = getchar();
	}
	// buf[i] = '\0';
	*out = s; // atoi(buf);
	return ch;
}

int main(int argc, char const *argv[])
{
	int n;
	while (scanf("%d", &n) == 1) {
		getchar(); // ignore a \n
		memset(nodes, 0, sizeof(struct node_t) * n);
		int i;
		int root_id;
		for (i = 0; i < n; ++i) {
			int parent_id, num;
			// scanf("%d:(%d)", &parent_id, &num);
			// parent_id = getint(':');
			getint2(&parent_id);
			if (i == 0) {
				root_id = parent_id;
			}
			getchar(); // (
			// num = getint(')');
			getint2(&num);
			// getchar(); // space
			int ch;
			for (num--; num >= 0; num--) {
				int node_id;
				ch = getint2(&node_id);
				// node_id = getint2(' ', '\n');
				// scanf("%d", &node_id);
				nodes[node_id].next_sibling = nodes[parent_id].first_child;
				nodes[parent_id].first_child = nodes + node_id;
			}
			while (ch != '\n') {
				ch = getchar();
			}
		}
		dp(nodes + root_id);
		printf("%d\n", MIN(nodes[root_id].min_with, nodes[root_id].min_without));
	}
	return 0;
}

poj 1947 Rebuilding Roads

　　这题的递推方法跟前面的有所不同，需要用子树及其兄弟的信息递推。也就是要设 F[i][j] 为：“从节点 i 的子节点及其右边的兄弟节点中选择 j 个节点所需的最小的边切割的数量。”具体的还是看 Beyond the Void 的解题报告吧，他讲得比我好得多了。

　　这题各种边界条件，各种特殊情况，很恶心。

gcc 16ms

#include <stdio.h>
#include <string.h>
#include <limits.h>

#define N 150

#define MIN(a,b) ((a)<(b)?(a):(b))

struct node_t {
	struct node_t *first_child;
	struct node_t *next_sibling;
	int capacity;
	int min_edge_cut[N+1];
} nodes[N];

int infadd(int a, int b) {
	if (a == INT_MAX || b == INT_MAX) {
		return INT_MAX;
	}
	return a + b;
}

int dp(struct node_t *node, int node_num) {
	if (node->min_edge_cut[node_num] != -1) {
		return node->min_edge_cut[node_num]; // calculated
	}

	// case1: not choose node
	int case1_min = INT_MAX;
	if (node->next_sibling) {
		int ret = dp(node->next_sibling, node_num);
		case1_min = infadd(ret, 1);
	} else if (node_num == 0) {
		case1_min = 1;
	}

	// case2: choose node
	int case2_min = INT_MAX;
	if (node_num == 0) {
		// impossible, pass
	} else {
		if (node->first_child == NULL && node->next_sibling == NULL) {
			if (node_num == 1) {
				case2_min = 0;
			}
		} else if (node->first_child == NULL) {
			case2_min = dp(node->next_sibling, node_num - 1);
		} else if (node->next_sibling == NULL) {
			case2_min = dp(node->first_child, node_num - 1);
		} else {
			int k;
			for (k = 0; k < node_num; ++k) {
				int ret = infadd(dp(node->first_child, k), dp(node->next_sibling, node_num - 1 - k));
				if (ret < case2_min) {
					case2_min = ret;
				}
			}
		}
	}

	int tmp = MIN(case1_min, case2_min);
	node->min_edge_cut[node_num] = tmp;
	// printf("node %d provide %d nodes: %d\n", node - nodes + 1, node_num, tmp);
	return tmp;
}

/**
 * init capacity
 */
void init_node(struct node_t *node) {
	struct node_t *child_node = node->first_child;
	node->capacity = 1;
	memset(node->min_edge_cut, -1, sizeof(node->min_edge_cut));
	while (child_node != NULL) {
		init_node(child_node);
		node->capacity += child_node->capacity;
		child_node = child_node->next_sibling;
	}
}

int main(int argc, char const *argv[])
{
	int n, p;
	scanf("%d%d", &n, &p);
	int i;
	int root_id = (n - 1) * n / 2;
	for (i = 0; i < n-1; ++i) {
		int p, c;
		scanf("%d%d", &p, &c);
		p--; c--;
		nodes[c].next_sibling = nodes[p].first_child;
		nodes[p].first_child = nodes + c;
		root_id -= c;
	}
	init_node(nodes + root_id);
	int min_cut = INT_MAX;
	for (i = 0; i < n; ++i) {
		int ret = (i != root_id);
		if (nodes[i].first_child) {
			ret = infadd(ret, dp(nodes[i].first_child, p - 1));
		} else {
			ret = infadd(ret, p == 1 ? 0 : INT_MAX);
		}
		if (ret < min_cut) {
			min_cut = ret;
		}
	}
	printf("%d\n", min_cut);
	return 0;
}

其他树形 dp 题目：

hdu 1011 Starship Troopers
poj 3345 Bribing FIPA
hdu 2196 Computer

为 bash 转义文件名

henix — Sat, 31 Mar 2012 12:23:40 +0800

　　最近遇到的一个问题：程序中有一个文件名，需要把这个文件名放在 shell 中执行，但文件名中可能包含特殊字符，所以需要转义。

　　比如，如果文件名是：

[SumiSora&CASO&HKG][Tears_to_Tiara][02][GB].rmvb

　　这个文件名肯定不能直接放到 bash 中的，因为“&”和 [ 、] 等都是 bash 的特殊字符。

　　bash 的自动补全默认采用反斜线转义：

\[SumiSora\&CASO\&HKG\]\[Tears_to_Tiara\]\[02\]\[GB\].rmvb

　　或者用单引号转义：

'[SumiSora&CASO&HKG][Tears_to_Tiara][02][GB].rmvb'

　　所以问题是，如何正确地实现转义？

　　经过一些搜索：

　　找到了两个东西可以实现这个功能：

python 3.3a 的 shlex.quote
bash 的内置命令 printf "%q" str（这货可不是 coreutils 的 printf ！）

　　但既然要实现 bash 的文件名转义，没有什么比 bash 本身的代码更权威的了。于是下载了 bash-4.2 的源代码来看。先花了很多时间定位，最终定位到 builtins/printf.def 这个文件，在大约 500 多行 case 'q' 的部分调用了以下函数：

ansic_shouldquote
ansic_quote
sh_backslash_quote

　　前两个函数在 lib/sh/strtrans.c 中，后一个函数在 lib/sh/shquote.c 中。所以最后终于定位到 shquote.c 这个文件。

如果要使用单引号转义，那么使用 sh_single_quote 的算法
想用反斜线转义，那么使用 sh_blackslash_quote 的算法

　　这两个函数的代码如下：

　　单引号转义：

/* **************************************************************** */
/*								    */
/*	 Functions for quoting strings to be re-read as input	    */
/*								    */
/* **************************************************************** */

/* Return a new string which is the single-quoted version of STRING.
   Used by alias and trap, among others. */
char *
sh_single_quote (string)
     const char *string;
{
  register int c;
  char *result, *r;
  const char *s;

  result = (char *)xmalloc (3 + (4 * strlen (string)));
  r = result;
  *r++ = '\'';

  for (s = string; s && (c = *s); s++)
    {
      *r++ = c;

      if (c == '\'')
	{
	  *r++ = '\\';	/* insert escaped single quote */
	  *r++ = '\'';
	  *r++ = '\'';	/* start new quoted string */
	}
    }

  *r++ = '\'';
  *r = '\0';

  return (result);
}

　　用 lua 实现如下：

function shquote(s)
	return "'"..string.gsub("'", "'\''").."'"
end

　　由于在单引号里面再用 \' 转义单引号也是非法的（我想这是因为单引号里面连 \ 也不是特殊字符），所以对于文件名里面出现的单引号，必须先结束上一个串，插入单引号，再开始下一个串。

　　反斜线转义：

/* Quote special characters in STRING using backslashes.  Return a new
   string.  NOTE:  if the string is to be further expanded, we need a
   way to protect the CTLESC and CTLNUL characters.  As I write this,
   the current callers will never cause the string to be expanded without
   going through the shell parser, which will protect the internal
   quoting characters. */
char *
sh_backslash_quote (string)
     char *string;
{
  int c;
  char *result, *r, *s;

  result = (char *)xmalloc (2 * strlen (string) + 1);

  for (r = result, s = string; s && (c = *s); s++)
    {
      switch (c)
	{
	case ' ': case '\t': case '\n':		/* IFS white space */
	case '\'': case '"': case '\\':		/* quoting chars */
	case '|': case '&': case ';':		/* shell metacharacters */
	case '(': case ')': case '<': case '>':
	case '!': case '{': case '}':		/* reserved words */
	case '*': case '[': case '?': case ']':	/* globbing chars */
	case '^':
	case '$': case '`':			/* expansion chars */
	case ',':				/* brace expansion */
	  *r++ = '\\';
	  *r++ = c;
	  break;
#if 0
	case '~':				/* tilde expansion */
	  if (s == string || s[-1] == '=' || s[-1] == ':')
	    *r++ = '\\';
	  *r++ = c;
	  break;

	case CTLESC: case CTLNUL:		/* internal quoting characters */
	  *r++ = CTLESC;			/* could be '\\'? */
	  *r++ = c;
	  break;
#endif

	case '#':				/* comment char */
	  if (s == string)
	    *r++ = '\\';
	  /* FALLTHROUGH */
	default:
	  *r++ = c;
	  break;
	}
    }

  *r = '\0';
  return (result);
}

最大流之最短增广路算法

henix — Mon, 19 Mar 2012 10:39:57 +0800

　　最短增广路算法在黑书上有介绍，其时间复杂度跟 Edmonds-Karp 是一样的 O(VE²) ，而且编程还复杂些。那为什么我还要介绍这个算法呢？因为这个算法引入了“距离标号”的概念，是后面的预推流算法的基础。

距离标号

　　节点 i 的距离标号 d[i] 首先被初始化成 i 到汇 sink 的实际距离，即最短路径的长度。

　　然后，在增广的时候要遵循一条规则：增广路上的每一条弧，其起始点的距离标号必须比结束点的距离标号大 1 ，这样的弧称为允许弧（admissible arc）。

　　如果我们把距离标号看成高度（实际上在预推流里面就叫做 height），把图的流看成水流，那么上述规则可以理解成：水流只能从高的地方流向低的地方，且每次只能下降一级。这样一比喻，整个算法就容易理解了。

距离标号的修改

　　关键的问题是：为什么要修改距离标号呢？

　　我们先看看最短增广路算法的框架：

　　这个算法框架主要参考了《网络流：算法、理论与应用》一书的 7.4 节，要注意跟黑书上讲的最短增广路算法略有不同。

初始化所有点的距离标号
i = source
while d(s) < n do
	if i 有允许弧 e(i,j) then
		-- 将 (i,j) 这条弧记录进当前的增广路
		pred[j] = i
		i = j
		if i == sink then -- 如果已经到达汇点
			增广
			i = s
		end
	else
		-- 从 i 出发没有允许弧，则修改 i 的距离标号
		d(i) = min(d(j) + 1 | 弧 e(i,j) 在残余网络中且可增广容量 c(i,j) > 0)
		if i != s then
			i = pred[i] -- 若不是源，则退回一步
		end
	end
end

　　黑书上的是如果没有允许弧则还要判断是否有弧，而且如果修改了也不会回退。

　　从 i 出发没有允许弧的时候，有两种情况：

从 i 出发有可增广弧，但 i 的 d 值太低，所以修改距离标号，然后就可以增广了
从 i 出发根本就没有可增广路，即满足上面的条件的 j 不存在，此时置 d(i) = n ，然后退一步。在我的代码中的体现是，min 的默认值是 n ，所以如果 j 不存在的话，就会把 n 赋给 d(i) 。

算法运行过程

　　下面以 poj 1273 为例，数据是：

　　0 是源，3 是汇。每个节点上标注了该节点的 d 值。

　　Step 1. 现在除汇外所有节点的 d 都等于 1 ，找到一条增广路：

　　Step 2. 灰色节点表示此时必须修改距离标号。我们看到从 0 出发有 0 -> 1 这条路还有容量，但由于 0 的距离标号太小，所以这不是允许弧：

　　Step 3. 修改距离标号后，找到一条增广路：

　　Step 4. 然后 1 又找不到允许弧了，于是修改它的距离标号，并回退到 0 ：

　　Step 5. 从 0 开始，发现也没有允许弧，于是修改距离标号：

　　Step 6. 0 的距离标号修改成 3 后，找到一条允许路：

　　Step 7. 注意这张图的布局跟前面的已经不同，大家只要记住 0 是源，3 是汇就行了。

　　Step 8.

　　Step 9. 0 号节点的距离标号终于被修改成 n = 4 ，算法结束。

源代码

　　我的实现，在 oj 上测试，gcc 0ms 。

　　其中 init_distance() 用来获得初始距离。本来以为要用 Dijkstra 单源点最短路径算法的，但实际上，由于这里相当于费用全部等于 1 的特殊情形，所以用一个 BFS 遍历就可以了。

#include <stdio.h>
#include <string.h>
#include <stdbool.h>
#include <limits.h>
#include <assert.h>

struct edge_t {
	int capacity;
};

int source, sink;

int maxflow(int n, struct edge_t edges[n][n], int source, int sink) {
	int pred[n];
	int d[n]; // distance to sink
	int sum = 0;

	void augment_path() {
		int delta = INT_MAX;
		int nv = sink;
		int v;
		while (nv != source) {
			v = pred[nv];
			int t = edges[v][nv].capacity;
			if (t < delta) {
				delta = t;
			}
			nv = v;
		}
		nv = sink;
		while (nv != source) {
			v = pred[nv];
			edges[v][nv].capacity -= delta;
			edges[nv][v].capacity += delta;
			nv = v;
		}
		sum += delta;
	}

	void init_distance() {
		int queue[n];
		int qhead = 0;
		int qend = 0;
		queue[qend++] = sink;

		memset(d, 0, sizeof(d));

		while (qhead < qend) {
			int v = queue[qhead++];
			int i;
			for (i = 0; i < n; ++i) {
				if (edges[i][v].capacity > 0 && d[i] == 0) {
					d[i] = d[v] + 1;
					queue[qend++] = i;
				}
			}
		}
	}

	init_distance();

	int i = source;
	while (d[source] < n) {
		// find an admissible arc of i
		int j;
		bool found = false;
		int min = n;
		for (j = 0; j < n; ++j) {
			if (edges[i][j].capacity > 0) {
				if (d[i] == d[j] + 1) {
					found = true;
					break;
				}
				if (d[j] + 1 < min) {
					min = d[j] + 1; // btw, find the min
				}
			}
		}
		if (found) {
			pred[j] = i;
			i = j;
			if (i == sink) {
				augment_path();
				i = source;
			}
		} else {
			d[i] = min;
			if (i != source) {
				i = pred[i];
			}
		}
	}

	return sum;
}

int main(int argc, char const *argv[])
{
	int n, m;
	while (scanf("%d%d", &m, &n) == 2) {

		struct edge_t edges[n][n];
		memset(edges, 0, sizeof(edges));

		for (; m > 0; m--) {
			int s, e, c;
			scanf("%d%d%d", &s, &e, &c);
			s--; // my index start from 0
			e--;
			edges[s][e].capacity += c;
		}
		source = 0;
		sink = n - 1;
		int max = maxflow(n, edges, source, sink);
		printf("%d\n", max);
	}
	return 0;
}

与 Dinic 算法的比较

　　Dinic 算法跟最短增广路算法很类似，但还是不一样的。Dinic 算法是建立了距离标号后在上面不断增广，不能增广后重新计算整个图的距离标号。详见参考资料。

参考资料

最大流之 Edmonds-Karp

henix — Sun, 18 Mar 2012 19:27:07 +0800

　　网络流问题定义什么的我就不介绍了，网上到处都找得到。

　　我第一次听说“最大流有 4 种做法”，大约是在这篇帖子里：http://tieba.baidu.com/p/1359675218 ，现在知道最大流的求法不止 4 种。

寻找增广路派：

原始的 Ford-Fulkerson 算法，用 DFS 寻找增广路
Edmonds-Karp 算法，用 BFS 寻找增广路
最短增广路算法，引入距离标号、允许弧，每次在允许弧上增广
Dinic 算法，引入层次网络

预推流派：

原始预推流
FIFO 实现预推流
最高标号预推流

　　还有 R.E.Tarjan 大神的 link/cut tree ，使得时间复杂度降到 log 级。

　　Edmonds-Karp 就是用 BFS 在残余网络上寻找增广路，时间复杂度 O(VE²) 。以 poj 1273 这题为例，下面是我第一次写的代码：

#include <stdio.h>
#include <string.h>
#include <stdbool.h>
#include <limits.h>
#include <assert.h>

struct edge_t {
	int capacity;
	int flow;
};

int source, sink;

int step = 0;

int maxflow(int n, struct edge_t edges[n][n], int source, int sink) {
	int parent[n];

	void augment_path() {
		int delta = INT_MAX;
		int nv = sink;
		int v;
		while (nv != source) {
			v = parent[nv];
			int t = edges[v][nv].capacity - edges[v][nv].flow;
			if (t < delta) {
				delta = t;
			}
			nv = v;
		}
		nv = sink;
		while (nv != source) {
			v = parent[nv];
			edges[v][nv].flow += delta;
			edges[nv][v].flow -= delta;
			nv = v;
		}
	}

	bool find_path() {
		int queue[n];
		int qhead = 0;
		int qend = 0;
		queue[qend++] = source;

		bool occured[n];
		memset(occured, 0, sizeof(occured));
		occured[source] = true;

		while (qhead < qend) {
			int v = queue[qhead++];
			if (edges[v][sink].flow < edges[v][sink].capacity) {
				parent[sink] = v;
				return true;
			}
			int i;
			for (i = 0; i < n; ++i) {
				if (!occured[i]) {
					if (edges[v][i].flow < edges[v][i].capacity) {
						occured[i] = true;
						parent[i] = v;
						queue[qend++] = i;
					}
				}
			}
		}
		return false;
	}

	while (find_path()) {
		augment_path();
		// dump(sink);
	}

	int sum = 0;
	int i;
	for (i = 0; i < n; ++i) {
		sum += edges[source][i].flow;
	}
	return sum;
}

int main(int argc, char const *argv[])
{
	int n, m;
	while (scanf("%d%d", &m, &n) == 2) {

		struct edge_t edges[n][n];
		memset(edges, 0, sizeof(edges));

		bool soccured[n];
		bool eoccured[n];
		memset(soccured, 0, sizeof(soccured));
		memset(eoccured, 0, sizeof(eoccured));

		for (; m > 0; m--) {
			int s, e, c;
			scanf("%d%d%d", &s, &e, &c);
			s--; // my index start from 0
			e--;
			edges[s][e].capacity += c;
			edges[s][e].flow = 0;
			soccured[s] = true;
			eoccured[e] = true;
		}
		source = 0;
		sink = n - 1;
		// puts("digraph G {");
		// puts("rankdir=LR");
		int max = maxflow(n, edges, source, sink);
		printf("%d\n", max);
		// puts("}");
	}
	return 0;
}

　　主要参考了 wiki 。但后来看了别人写的代码，发现其实 flow 可以不用保存，第二版：

#include <stdio.h>
#include <string.h>
#include <stdbool.h>
#include <limits.h>
#include <assert.h>

struct edge_t {
	int capacity;
};

int source, sink;

int step = 0;

int maxflow(int n, struct edge_t edges[n][n], int source, int sink) {
	int parent[n];
	int sum = 0;

	void augment_path() {
		int delta = INT_MAX;
		int nv = sink;
		int v;
		while (nv != source) {
			v = parent[nv];
			int t = edges[v][nv].capacity;
			if (t < delta) {
				delta = t;
			}
			nv = v;
		}
		nv = sink;
		while (nv != source) {
			v = parent[nv];
			edges[v][nv].capacity -= delta;
			edges[nv][v].capacity += delta;
			nv = v;
		}
		sum += delta;
	}

	bool find_path() {
		int queue[n];
		int qhead = 0;
		int qend = 0;
		queue[qend++] = source;

		bool occured[n];
		memset(occured, 0, sizeof(occured));
		occured[source] = true;

		while (qhead < qend) {
			int v = queue[qhead++];
			if (edges[v][sink].capacity > 0) {
				parent[sink] = v;
				return true;
			}
			int i;
			for (i = 0; i < n; ++i) {
				if (!occured[i]) {
					if (edges[v][i].capacity > 0) {
						occured[i] = true;
						parent[i] = v;
						queue[qend++] = i;
					}
				}
			}
		}
		return false;
	}

	while (find_path()) {
		augment_path();
		// dump(sink);
	}

	return sum;
}

int main(int argc, char const *argv[])
{
	int n, m;
	while (scanf("%d%d", &m, &n) == 2) {

		struct edge_t edges[n][n];
		memset(edges, 0, sizeof(edges));

		for (; m > 0; m--) {
			int s, e, c;
			scanf("%d%d%d", &s, &e, &c);
			s--; // my index start from 0
			e--;
			edges[s][e].capacity += c;
		}
		source = 0;
		sink = n - 1;
		// puts("digraph G {");
		// puts("rankdir=LR");
		int max = maxflow(n, edges, source, sink);
		printf("%d\n", max);
		// puts("}");
	}
	return 0;
}

　　capacity 可以看成“这条边上可以增加的流的大小”，当前找到的增广路保存在 parent 数组里。

实测效率

　　在 OJ 上实测，如果用原始的 Ford-Fulkerson 算法，则 TLE ，因为原始算法的时间复杂度跟最大流的大小有关。而用 Edmonds-Karp 就 0ms 过了。

参考资料

Edmonds–Karp algorithm - Wikipedia, the free encyclopedia

archlinux 上 awesome 的安装与配置

henix — Wed, 29 Feb 2012 09:46:14 +0800

　　先来张图：

　　我以前用的是 wmii ，但后来转到 awesome ，主要是因为 awesome 的脚本化能力很强大。

　　使用 linux 的过程中我们经常面对的东西有两样：terminal 和浏览器。前者漆黑，而后者很白。在两者之前来回切换造成的后果是让你的眼睛低亮度的环境还没适应过来，一下子就面对高亮度，眼睛很容易疲劳。

　　所以我首先想到的解决方法是：为不同的程序设置不同的壁纸。对浏览器，设置一个全黑的壁纸；对终端，设置一个比较亮的壁纸。然后让窗口半透明，这样壁纸就可以与窗口本身的亮度中和。

　　当然，实际上我们是做不到对不同的程序设置不同的壁纸的，于是我改成：在不同的桌面上设置不同的壁纸。比如 1 号桌面是黑的，2-6 是白的，7-9 都是黑的。然后只在白的桌面上开终端，只在黑的桌面上开浏览器。

awesome 的安装

　　先编译 cairo-xcb ，再编译 awesome 。我在编译 awesome-3.4.11-2 的时候遇到一个编译错误，大意是找不到 libiconv.a 。我发现系统里本来就没有这个东西，而且也没有一个叫做 libiconv 的库。最后我发现 libiconv 是包含在 gcc 里面的，在链接的时候直接加上 -liconv 就可以了。所以，需要修改 PKGBUILD ，下面是 diff ：

--- PKGBUILD.orig	2011-12-09 16:43:16.517210395 +0800
+++ PKGBUILD	2011-12-09 16:44:28.363881864 +0800
@@ -26,7 +26,7 @@
   cd "$srcdir/$pkgname-$pkgver"
 
   make CMAKE_ARGS=" -DPREFIX=/usr -DSYSCONFDIR=/etc \
-	-DCMAKE_BUILD_TYPE=RELEASE"
+	-DCMAKE_BUILD_TYPE=RELEASE -DCMAKE_EXE_LINKER_FLAGS=-liconv"
 }
 
 package() {

　　就是在链接器选项里加一个 -liconv ，我现在已经不知道当时是如何查阅那庞大的 cmake 手册的了。。。

awesome 的配置

1. 前面说的为不同的桌面设置不同的壁纸的代码如下：

-- {{{ Tags
-- Define a tag table which hold all screen tags.
last_desk = 'dark'
tags = {}
for s = 1, screen.count() do
    -- Each screen has its own tag table.
    tags[s] = awful.tag({ 1, 2, 3, 4, 5, 6, 7, 8, 9 }, s, layouts[2])
	awful.tag.attached_add_signal(s, "property::selected", function(t)
		-- local selected = awful.tag.selected(s)
		if last_desk == 'light' and (t.name == "1" or t.name == '7' or t.name == "8" or t.name == "9") then
			awful.util.spawn('feh --bg-fill /home/henix/black.png')
			last_desk = 'dark'
		elseif last_desk == 'dark' and (t.name ~= "1" and t.name ~= '7' and t.name ~= "8" and t.name ~= "9") then
			awful.util.spawn('feh --bg-center /home/henix/57a2ef86e99ba83166096e62.jpg')
			last_desk = 'light'
		end
	end)
end
-- }}}

　　通过 feh 设置壁纸，且仅在当前桌面和要切换过去的桌面的壁纸不同的时候在调用 feh ，减少 fork 的开销。

2. 一些自定义快捷键

-- Win + p 打开 dmenu ：
awful.key({modkey}, "p", function ()
    awful.util.spawn('dmenu_run')
end),

-- Ctrl + Atl + l 锁屏：
-- 因为 Windows 中是 Win + L ，但在 awesome 中 Win + L 已经被用来干其他事情了
awful.key({"Mod1", "Control"}, "l", function ()
    awful.util.spawn('xscreensaver-command -lock')
end),

-- 取消 Win + Shift + q 的退出，
-- 退出的时候用菜单，因为我老是不小心按了...
-- awful.key({modkey, "Shift"}, "q", awesome.quit),

　　当然这里又引发了我的另一个思考，如何设计用户界面，我称之为用户界面设计的 henix's law ：

用户越需要频繁使用的功能应该越容易被调用；相反，越不需要频繁使用的功能应该越难被找到。

　　容易调用的例子：快捷键、工具栏；难找到的例子：庞大的菜单、一个套一个的对话框、配置文件。

　　比如退出这种事情，实际上不需要快捷键，直接点菜单就行了，因为开一天机可能就用到一两次退出。但锁屏就需要快捷键，因为我离开电脑就需要，使用得更频繁。

3. 窗口背景透明：

　　在 manage 的 signal 里加入：

    if c.class ~= "MPlayer" then
        c.opacity = 0.72
    end

　　将除了 mplayer 外所有窗口的透明度设置为 0.72 。

4. 添加 vicious 控件：

-- vicious

batwidget = widget({ type = "textbox" })
vicious.register(batwidget, vicious.widgets.bat, "$1 $2% $3 | ", 5, "BAT0")

uptimewidget = widget({ type = "textbox" })
vicious.register(uptimewidget, vicious.widgets.uptime, "$4 $5 $6 | ", 7)

cpuwidget = widget({ type = "textbox" })
vicious.register(cpuwidget, vicious.widgets.cpu, "$1% | ", 3)

memwidget = widget({ type = "textbox" })
vicious.register(memwidget, vicious.widgets.mem, "$1% ($2MB/$3MB) | ", 13)

mystatusbar = awful.wibox({ position = "bottom", name = "mystatusbar" })
mystatusbar.widgets = {
    mytextclock,
    memwidget,
    cpuwidget,
    uptimewidget,
    batwidget,
    layout = awful.widget.layout.horizontal.rightleft
}
mystatusbar.screen = 1

　　在底部状态栏显示电池、uptime 、CPU 和内存。

IE CSS Hacks

henix — Wed, 22 Feb 2012 09:57:08 +0800

　　最近做右边“热门文章”的彩虹表效果，为了让它在 IE 上能优雅降级，研究了下针对 IE 的 CSS Hacks ：

selector{ 
  property:value; /* all browsers */ 
  property:value\9; /* IE 6 7 8 9 */ 
  property:value\0; /* IE 8 9 Opera */ 
  property:value\0/; /* IE 8 9 */
  property:value\9\0; /* IE 9 */ 
  *property:value; /* IE6 IE7 */
  +property:value; /* IE7 */ 
  _property:value; /* IE6 */ 
}
* html selector {
  property:value; /* IE 6 */
}
@media \0screen {
  selector {
    property:value; /* IE 8 */
  }
}
:root selector {
  property:value\9; /* IE 9 */
}

　　关键是 \0 这个 hack ：很多地方都说是 IE 8 9 的 hack ，但经我测试后发现最新版的 Opera 也识别（Opera 是 SB），恰好我平时用 Opera 用得比较多……最后终于找到了 \0/ 这个 hack ，算是圆满解决。

　　所以，CSS Hacks 这种东西，一定要自己测试，要注意文章的实效性。

References:

POJ 2778 AC 自动机上的 dp

henix — Fri, 10 Feb 2012 17:09:41 +0800

　　poj 2778 DNA Sequence 题目大意：给你一些疾病 DNA 片段，求长度为 n 的 DNA 串中不包含这些片段的串的数量。

　　一开始不知道怎么做，直到看到网上的一句话，“长度为 n 的字符串可以由长度为 n - 1 的字符串后面加一个字符构成”，由此大彻大悟，原来这就是传说中的自动机 dp ！

　　举个例子：{AG, CG} ，首先构造 AC 自动机：

　　给每个节点加了编号，只有 0 1 2 是合法的，其他状态到达就直接死了，不用处理。

　　然后将上面的 AC 自动机转化成下面这张表：

	A	C	T	G
0	1	2	0	0
1	1	2	0	3
2	1	2	0	4

　　这张表表示当前状态遇到一个字母后跳转到什么状态。我称之为“跳转表”，对应于我代码中的 jump 域。这个貌似跟传说中的“trie 图”类似，都是把自动机确定化。

　　然后数一下每一行有多少个 0 ，多少个 1 ……将上表转化成下面的递推关系：

$f_0(n) = 2 * f_0(n-1) + 1 * f_1(n-1) + 1 * f_2(n-1)$
$f_1(n) = 1 * f_0(n-1) + 1 * f_1(n-1) + 1 * f_2(n-1)$
$f_2(n) = 1 * f_0(n-1) + 1 * f_1(n-1) + 1 * f_2(n-1)$

　　其中 $f_i(n)$ 表示到第 n 步，有多少个字符串处于状态 i 。$f_0(n-1)$ 前面的系数 2 是 0 那一行包含的 0 的数量，表示如果上一步有一个处于状态 0 ，那么下一步就会有 2 个处于状态 0（分别对应遇到 T 和 G 的情况）。其他类似。

　　将上面的系数矩阵提出来，立即得到：

$$ \left[ \begin{array}{c} f_0(n)\\f_1(n)\\f_2(n) \end{array} \right] = \left[ \begin{array}{ccc} 2&1&1\\1&1&1\\1&1&1 \end{array} \right]^n \left[ \begin{array}{c} 1\\0\\0 \end{array} \right] $$

　　这就是所谓“常系数线性递推式”，详见黑书练习题 2.1.8 。比如对于 n = 3 ，先计算系数矩阵的 3 次方：

$$ \left[ \begin{array}{ccc} 2&1&1\\1&1&1\\1&1&1 \end{array} \right]^3 = \left[ \begin{array}{ccc} 20&14&14\\14&10&10\\14&10&10 \end{array} \right] $$

　　然后把结果矩阵的第一列加起来即得到答案 48 。

　　题目中 n 可能非常大，这时就轮到快速幂运算出场了。使用快速幂运算，计算 $A^n$ 的时间复杂度为 O(log n) 次乘法（此题是矩阵乘法）。我的代码中这部分用了两个指针，每次计算都交换，这是为了避免矩阵复制。

　　最后，还要考虑一种情况，即有些疾病 DNA 片段互相包含的情况，以 {ACG, C} 为例：

　　其中 C 是死节点，那么 AC（红色的）是不是呢？表面上看 AC 不是疾病片段，但实际上，匹配了 AC 就相当于匹配了 C ，所以 AC 也是死节点！我称这种情况为“疾病的传染”，见代码中用 spread 注释的部分。

　　代码（gcc 63MS）：

#include <stdio.h>
#include <string.h>
#include <stdlib.h>

int dnas[26] = {0, 0, 1, 0, 0, 0, 2, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 3};

int dna_seq(char c) {
	return dnas[c - 'A'];
}

struct trie_node {
	char isend;
	struct trie_node *fail;
	struct trie_node *child[4];
	struct trie_node *jump[4];
} nodes[101];

int endm;

void reset() {
	endm = 1;
	memset(nodes, 0, sizeof(struct trie_node));
}

void insert(const char *str) {
	struct trie_node *index = nodes;
	const char *p = str;
	for (; *p != '\0'; p++) {
		int i = dna_seq(*p);
		if (index->child[i] == NULL) {
			index->child[i] = nodes + endm;
			memset(nodes + endm, 0, sizeof(struct trie_node));
			endm++;
		}
		index = index->child[i];
	}
	index->isend = 1;
}

struct trie_node *queue[100];
int qhead, qend;

void add_fail() {
	qhead = -1;
	qend = 0;
	/* add root to the queue */
	queue[qend++] = nodes;
	while (qhead + 1 < qend) {
		struct trie_node *x = queue[++qhead];
		int i;
		for (i = 0; i < 4; i++) {
			struct trie_node *child = x->child[i];
			struct trie_node *t = x->fail;
			while ((t != NULL) && (t->child[i] == NULL)) {
				t = t->fail;
			}
			if (child != NULL) {
				queue[qend++] = child;
				child->fail = t ? t->child[i] : nodes;
				if (child->fail->isend) {
					child->isend = 1; // spread
				}
				x->jump[i] = child;
			} else {
				x->jump[i] = t ? t->child[i] : nodes;
			}
		}
	}
}

/**
 * c = a * b
 */
void matrix_mul(int n, int a[n][n], int b[n][n], int c[n][n]) {
	int i, j, k;
	for (i = 0; i < n; i++) {
		for (j = 0; j < n; j++) {
			long long sum = 0;
			for (k = 0; k < n; k++) {
				sum += ((long long)a[i][k] * b[k][j]);
			}
			c[i][j] = sum % 100000;
		}
	}
}

int main(int argc, const char *argv[])
{
	int m, n;
	scanf("%d%d", &m, &n);
	getchar(); // ignore a \n
	char buf[16];
	reset();
	int i;
	for (i = 0; i < m; i++) {
		gets(buf);
		insert(buf);
	}
	add_fail();
	// construct the matrix
	int len = endm - m; // max len, but may not correct
	// printf("len = %d, endm = %d\n", len, endm);
	int matrix[len][len];
	int j, k;
	int id = 0;
	for (i = 0; i < endm; i++) {
		if (! nodes[i].isend) {
			// collect how many i occured in others
			int id2 = 0;
			for (j = 0; j < endm; j++) {
				if (! nodes[j].isend) {
					int count = 0;
					for (k = 0; k < 4; k++) {
						if (nodes[j].jump[k] == nodes + i) {
							count++;
						}
					}
					matrix[id][id2] = count;
					id2++;
				}
			}
			id++;
		}
	}
	len = id; // true len
	/*for (i = 0; i < len; i++) {
		for (j = 0; j < len; j++) {
			printf("%d ", matrix[i][j]);
		}
		putchar('\n');
	}*/
	int *d1 = malloc(len * len * sizeof(int));
	int *d2 = malloc(len * len * sizeof(int));
	int *c1 = malloc(len * len * sizeof(int));
	int *c2 = malloc(len * len * sizeof(int));
	int *result;
	if (n == 1) {
		result = matrix;
	} else {
		memset(c1, 0, len * len * sizeof(int));
		for (i = 0; i < len; i++) {
			c1[i * len + i] = 1;
		}
		for (i = 0; i < len; i++) {
			for (j = 0; j < len; j++) {
				d1[i * len + j] = matrix[i][j];
			}
		}
		while (n > 0) {
			if (n & 1) {
				matrix_mul(len, c1, d1, c2);
				int *tmp = c1;
				c1 = c2;
				c2 = tmp;
			}
			n >>= 1;
			matrix_mul(len, d1, d1, d2);
			int *tmp = d1;
			d1 = d2;
			d2 = tmp;
		}
		result = c1;
	}
	/*for (i = 0; i < len; i++) {
		for (j = 0; j < len; j++) {
			printf("%d ", result[i * len + j]);
		}
		putchar('\n');
	}*/
	int sum = 0;
	for (i = 0; i < len; i++) {
		sum = (sum + result[i * len]) % 100000;
	}
	printf("%d\n", sum);
	free(d1);
	free(d2);
	free(c1);
	free(c2);
	return 0;
}

　　其他相关题目：

poj 3691 DNA repair
poj 1625 Censored!
hdu 2825 Wireless Password

　　贴两个解题报告：