Fractran是这样一种自动机：

该自动机由一个只包含正分数（即$a/b$，且$a,b\in\mathbb N$）的有序且有限的列表，和一个状态值s组成。状态值总是正整数，其初值（也叫种子值）在初始化时指定。

每次迭代时，找到列表中中第一个与s相乘后能结果为整数的分数$a/b$，并将状态更新为$as/b$。

例如， J.H.Conway曾给出如下Fractran自动机，它的初始状态为2，内部的列表为： $$\frac{17}{91}, \; \frac{78}{85}, \; \frac{19}{51}, \; \frac{23}{38}, \; \frac{29}{33}, \; \frac{77}{29}, \; \frac{95}{23}, \; \frac{77}{19}, \; \frac{1}{17}, \; \frac{11}{13}, \; \frac{13}{11}, \; \frac{15}{2}, \; \frac{1}{7}, \; \frac{55}{1}.$$ 该自动机迭代产生的序列为： $$15, 825, 725, 1925, 2275, 425, \ldots, 68, \text{\textbf{4}}, 30, \ldots, 136, \text{\textbf{8}}, 60, \ldots, 544, \text{\textbf{32}}, 240, \ldots$$ 该序列中依次出现$2^2, 2^3, 2^5, \ldots$等2的幂，其指数均为素数，且恰好按素数本身的顺序而出现。

若用上述自动机来解决欧拉计划的第7题（找出第10001个素数$p_{10001}$），则需要迭代多少步才能生成$2^{p_{10001}}$？

本题难度:

解答

本题相当繁琐，但索性文献资料比较齐全，可以参考OEIS A007542及其页面上的资料。

观察这些分数，可以发现状态$s$总是可以写成$2^n3^a5^b7^cx$的形式，其中x是素因子在$\{11, 13, 17, 19, 23, 29\}$之中的无平方因子数（或1），n,a,b,c都是非负整数。

可以把n,a,b,c理解为寄存器，x理解为控制器，把乘以这14个分数的操作依次记作 A, B, ..., N，这就得到了本题中Fractran的抽象计算模型，寄存器中值的变化可以用以下伪代码来描述：

while (true) {
  while (n > 0) { n--; a++: b++: }           		/* (L)n    */
  while (c > 0) { c--; }                    		/* (M)c    */
  b++;                                       		/*  N      */
  while (a > 0) { a--; c++; }               		/* (EF)a K */
    
  while (b > 0) {                       	        /*  b-max(x: 0 < x < b, x|b)-1 次  */
    while (b >= c) {                                    /*  floor(b/c) 次, c = b-1, ..., d */
      while (c > 0) { n++; a++; b--; c--; }  		/* (AB)c J */
      while (a > 0) { a--; c++; }           		/* (EF)a K */
    }
    
    if (b > 0) {
      while (b > 0) { n++; a++; b--; c--; }  		/* (AB)b   */
      a--; c--;                              		/*  A C    */
      while (n > 0) { n--; b++; }            		/* (DG)n   */
      c++;                                  		/*  H      */
      while (a > 0) { a--; c++; }                       /* (EF)a K */
    } else {
      c--;                                              /*  A I    */
    }
  }
}

其中形如(AB)c的注释记号表示每轮依次执行操作A和操作B，一共执行c轮。

仔细上述审查状态变化和对应的操作次数可以归纳出计算公式（公式和对应的计算方法不唯一），最终结果是$1539669807660924$。

注：以下代码改编自官网论坛，代码中打印了进度信息。

n=10001
i=2
r=0
while n>0:
    j=i-1
    r+=j
    while i%j>0:
        r+=2+6*i+(i/j)*2
        j-=1
    r+=2+6*i+(i/j)*2+j-1
    if j==1:
        n-=1
        if n%100==0:
            print 100-n/100,"percent completed"
    i+=1

print r