使用 LaTeX 计算最长的 Collat​​z 序列

使用 LaTeX 计算最长的 Collat​​z 序列

最近,我开始尝试expl3编程语言。我发现它不仅能解决与排版相关的问题,还能解决其他编程语言可以解决的问题,非常强大。

为了练习,我选择了以下问题:欧拉计划在社区的帮助下tex.se,我设法解决了一些我曾经使用解决的问题Python

下面我引用了名为 的问题之一的解决方案代码Longest Collatz sequence。该代码并不完美,因为它耗费了大量的预估时间。也许有更优雅的算法可以使用 来解决这个问题expl3

在本节中,我想看到该问题的其他解决方案,不仅包括exp3,还包括可以使用的其他可能性LaTeX(例如LuaLaTeX等)


\documentclass[12pt]{article}

\usepackage{xparse, pgf}
\ExplSyntaxOn

\int_new:N \l_max_int 
\int_zero:N \l_max_int
\int_new:N \l_max_num_int 
\int_zero:N \l_max_num_int
\int_new:N \l_a_int 
\int_zero:N \l_a_int
\int_set:Nn \l_tmpa_int {1}

\NewDocumentCommand{\MaxNumberB}{ mmm }
{
    \int_while_do:nn {\l_tmpa_int <= #3}
    {
    \int_set:Nn \l_a_int { \l_tmpa_int }
    \int_set:Nn \l_tmpb_int { 1 }
        \int_while_do:nn{\l_a_int != 1}
        {
            \int_if_even:nTF {\l_a_int}
                {\int_set:Nn \l_a_int {\l_a_int/2} } 
                {\int_set:Nn \l_a_int {3*\l_a_int + 1} }
                \int_incr:N \l_tmpb_int
        }
        \int_compare:nTF{\l_tmpb_int > \l_max_int}
        {
        \int_set:Nn \l_max_int {\l_tmpb_int}
        \int_set:Nn \l_max_num_int {\l_tmpa_int}
        }{}

    \int_incr:N \l_tmpa_int
    } 

       \cs_new:Npx #1 { \int_eval:n { \l_max_int } }
       \cs_new:Npx #2 { \int_eval:n { \l_max_num_int } }

}

\ExplSyntaxOff

\begin{document}

\pgfmathsetmacro{\num}{1000000}
\MaxNumberB{\MaxNumEl}{\MaxNum}{\num}

For numbers from $0$ to $\num$, the number \MaxNum{} has the maximum number of elements in Collatz sequence~--- \MaxNumEl.


\end{document}

答案1

在 LuaTeX 中重新实现 Joseph 的答案。在 Lua 中进行计算的好处是它们始终是完全可扩展的。

\documentclass{article}
\usepackage{luacode}

\begin{luacode}
function collatz_next(n)
    if ( n \% 2 == 0 ) then
        return .5 * n
    else
        return 3 * n + 1
    end
end

function collatz_count(n, m)
    m = m or 1

    if ( n > 1 ) then
        return collatz_count( collatz_next(n), m + 1 )
    else
        return m
    end
end

function collatz_max_range(start, stop, count, len)
    count = count or 1
    len = len or start
    if start > stop then
        return len
    else
        if collatz_count(start) > count then
            return collatz_max_range(start+1, stop, collatz_count(start), start)
        else
            return collatz_max_range(start+1, stop, count, len)
        end
    end
end
\end{luacode}

\newcommand*\collatzcount[1]{\directlua{tex.sprint(collatz_count(#1))}}
\newcommand*\collatzmaxrange[2]{\directlua{tex.sprint(collatz_max_range(#1,#2))}}

\begin{document}

\edef\x{\collatzcount{13}}
\show\x

\edef\x{\collatzmaxrange{1}{10}}
\show\x

\end{document}

\collatzmaxrange{1}{1000000}性能相当不错英特尔® 酷睿™ i7-4790K

$ time lualatex --interaction=batchmode test.tex
This is LuaTeX, Version 1.0.4 (TeX Live 2017) 
 restricted system commands enabled.

luaotfload | main : initialization completed in 0.070 seconds
real    0m10.875s
user    0m10.860s
sys 0m0.012s

LuaJIT 比 Lua 快很多。

$ time luajittex -fmt=luajitlatex --interaction=batchmode test.tex
This is LuajitTeX, Version 1.0.4 (TeX Live 2017) 
 restricted system commands enabled.

luaotfload | main : initialization completed in 0.048 seconds
real    0m4.675s
user    0m4.644s
sys 0m0.008s

让我们实际启动 JIT。

$ time luajittex -jiton -fmt=luajitlatex --interaction=batchmode test.tex
This is LuajitTeX, Version 1.0.4 (TeX Live 2017) 
 restricted system commands enabled.

luaotfload | main : initialization completed in 0.048 seconds
real    0m0.807s
user    0m0.780s
sys 0m0.024s

好的,这已经是最快的了。

\documentclass{article}

\directlua{collatz = require"collatz"}

\newcommand*\collatzcount[1]{\directlua{tex.sprint(collatz.count(#1))}}
\newcommand*\collatzmaxrange[2]{\directlua{tex.sprint(collatz.max_range(#1,#2))}}

\begin{document}

\collatzmaxrange{1}{1000000}

\end{document}
$ time lualatex -interaction=batchmode test.tex
This is LuaTeX, Version 1.0.4 (TeX Live 2017) 
 restricted system commands enabled.

luaotfload | main : initialization completed in 0.069 seconds
real    0m0.367s
user    0m0.360s
sys 0m0.004s

编译为clang -shared -fpic -O3 -I/usr/include/lua5.2 test.c -o collatz.so

#include <lua.h>
#include <lauxlib.h>
#include <assert.h>

typedef unsigned long int_t;

inline int_t collatz_next(int_t n)
{
    if ( n % 2 == 0 )
      return n / 2;
    else
      return 3 * n + 1;
}

int_t collatz_count(int_t n, int_t m)
{
    if ( n > 1 )
      return collatz_count( collatz_next(n), m + 1 );
    else
      return m;
}

int_t collatz_max_range(int_t start, int_t stop, int_t count, int_t len)
{
  if ( start > stop )
    return len;
  else
  {
    if ( collatz_count(start,1) > count )
      return collatz_max_range(start+1, stop, collatz_count(start,1), start);
    else
      return collatz_max_range(start+1, stop, count, len);
  }
}

static int lcollatz_count(lua_State *L)
{
  assert( lua_gettop(L) == 1  && "Number of required arguments is 1" );
  assert( lua_isnumber(L, -1) && "First argument is not a number"    );

  int_t n = lua_tonumber(L, -1);

  int_t result = collatz_count(n,1);
  lua_pushnumber(L, result);

  return 1;
}

static int lcollatz_max_range(lua_State *L)
{
  assert( lua_gettop(L) == 2  && "Number of required arguments is 2" );
  assert( lua_isnumber(L, -2) && "First argument is not a number"    );
  assert( lua_isnumber(L, -1) && "Second argument is not a number"   );

  int_t start = lua_tonumber(L, -2);
  int_t stop  = lua_tonumber(L, -1);

  int_t result = collatz_max_range(start,stop,1,start);
  lua_pushnumber(L, result);

  return 1;
}

static const luaL_Reg lib_collatz[] = {
  {"count", lcollatz_count},
  {"max_range", lcollatz_max_range},
  {NULL, NULL} // sentinel
};

LUALIB_API int luaopen_collatz(lua_State *L)
{
  luaL_newlib(L, lib_collatz);
  return 1;
}

答案2

我会使用可扩展的实现来做到这一点:将所有内容保存在堆栈上很容易:

\documentclass{article}
\usepackage{expl3}
\usepackage{xparse}
\begin{document}
\ExplSyntaxOn
\cs_new:Npn \collatz_next:n #1
  {
    \int_if_even:nTF {#1}
      { \int_eval:n { #1 / 2 } }
      { \int_eval:n { 3 * #1 + 1 } }
  }
\cs_generate_variant:Nn  \collatz_next:n { f }
\cs_new:Npn \collatz_count:n #1
  { \__collatz_count:nn {#1} { 1 } }
\cs_new:Npn \__collatz_count:nn #1#2
  {
    \int_compare:nNnTF {#1} > { 1 }
      {
        \__collatz_count:ff
          { \collatz_next:n {#1} } { \int_eval:n { #2 + 1 } }
      }
      { #2 }
  }
\cs_generate_variant:Nn \__collatz_count:nn { ff }
\cs_new:Npn \collatz_max_range:nn #1#2
  {
    \__collatz_max_range:nnnn { 1 } {#1} {#1} {#2}
  }
\cs_new:Npn \__collatz_max_range:nnnn #1#2#3#4
  {
    \int_compare:nNnTF {#3} > {#4}
      {#2} % Could also dump length (#1) here
      {
        \__collatz_max_range:fnnnn { \collatz_count:n {#3} } 
          {#1} {#2} {#3} {#4}
      }
  }
\cs_generate_variant:Nn \__collatz_max_range:nnnn { nnf }
\cs_new:Npn \__collatz_max_range:nnnnn #1#2#3#4#5
  {
    \int_compare:nNnTF {#1} > {#2}
      { \__collatz_max_range:nnfn {#1} {#4} }
      { \__collatz_max_range:nnfn {#2} {#3} }
        { \int_eval:n { #4 + 1 } } {#5}
  }
\cs_generate_variant:Nn \__collatz_max_range:nnnnn { f }
\NewExpandableDocumentCommand \collatzcount { m }
  { \collatz_count:n {#1} }
\NewExpandableDocumentCommand \collatzmaxinrange { mm }
  {
    \collatz_max_range:nn {#1} {#2}
  }

\ExplSyntaxOff

\collatzcount{13}
\collatzmaxinrange{1}{13}

\end{document}

这里我利用了 -type 扩展本身可扩展的事实f:它停在第一个不可扩展的标记处,这在这里是理想的,因为它是我们想要的数字。

正如我在代码中指出的那样,人们可能想要安排“返回”两个值:这可以通过添加如下参数来实现

\cs_new:Npn \collatz_max_range:nnN #1#2#3
  {
    \__collatz_max_range:nnnnN { 1 } {#1} {#1} {#2} #3
  }
\cs_new:Npn \__collatz_max_range:nnnnN #1#2#3#4#5
  {
    \int_compare:nNnTF {#3} > {#4}
      { #5 {#1} {#2}} % Could also dump length (#1) here
      {
        \__collatz_max_range:fnnnnN { \collatz_count:n {#3} } 
          {#1} {#2} {#3} {#4}
      }
  }
\cs_generate_variant:Nn \__collatz_max_range:nnnnN { nnf }
\cs_new:Npn \__collatz_max_range:nnnnnN #1#2#3#4#5#6
  {
    \int_compare:nNnTF {#1} > {#2}
      { \__collatz_max_range:nnfnN {#1} {#4} }
      { \__collatz_max_range:nnfnN {#2} {#3} }
        { \int_eval:n { #4 + 1 } } {#5} #6
  }

第三个参数\collatz_max_range:nnN可能是简单的,\use_i:nn也可能是更复杂的(存储值,用填充文本排版,ETC。


也可以考虑使用“预先计算”方法来保存每个值的数据。例如,我们可以建立关于要选择的“下一个”值的数据:这可用于显示值链或链长度:

\prop_new:N \g__collatz_next_prop
\int_new:N \l__collatz__tmp_int
\cs_new_protected:Npn \collatz_precalc:n #1
  {
    \int_set:Nn \l__collatz__tmp_int { \collatz_next:n {#1} }
    \prop_gput:NnV \g__collatz_next_prop
      {#1} \l__collatz__tmp_int
    \prop_if_in:NVF \g__collatz_next_prop \l__collatz__tmp_int
      { \collatz_precalc:V  \l__collatz__tmp_int }
  }
\cs_generate_variant:Nn \collatz_precalc:n { V }
\int_step_inline:nnnn { 1 } { 1 } { 100 }
  { \collatz_precalc:n {#1} }

人们可能会考虑每个值使用一个tl(或int),但就 csname 使用而言,这会导致成本高昂(我们在int寄存器方面也受到一定限制)。根据实际用例,预先计算的值或预先计算的一组链可能会很有用:取决于数据是否被重复使用。

答案3

使用 Plain TeX(以及方便使用的 eTeX)

  1. 初始方法受 TeX 算法限制,最大起点为 N = 159486。

  2. 扩大方法进行大算术并能够掌握 N = 1,000,000

  3. 技术变体使用字体尺寸技术来存储所有 Collat​​z 长度,起始点从 1 到 N。最大 N 为 5,000,000(在我的 TeXLive 上,我可以达到略小于 8,000,000。超出需要扩大 TeX 内存)。最后\number\fontdimen<number>\cz存储了所有 Collat​​z 长度,<number>起始点从 1 到 N=5,000,000。

  4. 最后,一个稍微简单一些的算法将速度提高了近 2 倍。


編輯我最初的帖子根本没有执行处理已知事物在控制序列中的存储的那部分。因此它要慢得多,但另一方面,当我修正它时,我发现它很快就达到了 TeX 内存限制,因为在\begingroup\endgroup中包含扩展。但毕竟我们存储的是绝对数据,所以不需要组。而且,我决定只存储索引最大为 N 的数据,以便再次摆脱内存限制。

% pdftex (or etex)

\newcount\inputNmax

\newcount\intN
\newcount\intNa
\newcount\intNtop
\newcount\intL
\newcount\intLtop

\def\CollOne{%
    \advance\intNa 1
    \ifnum\intNa > \inputNmax
       \CollDone
    \else
      \intL 0 % will store number of steps starting at Na
      \intN = \intNa
      \expandafter\CollTwo
    \fi
}

% store in 
\def\CollTwo{%
    \ifcsname collatz\the\intN\endcsname
       \expandafter\CollThreeA
    \else
       \expandafter\CollThreeB
    \fi
}

\def\CollThreeA{%
    \advance\intL\csname collatz\the\intN\endcsname\relax
    \expandafter\edef\csname collatz\the\intNa\endcsname{\the\intL}%
    \ifnum\intL > \intLtop
          \intLtop = \intL
          \intNtop = \intNa
    \fi
    \intN = \intNa
    \CollUpdate
}

\def\CollUpdate{%
    \advance\intL -1
    \ifodd\intN
      \multiply\intN 3
      \advance\intN  1
    \else
      \divide\intN 2
    \fi
    \let\next\CollOne
    \ifnum\intN>\inputNmax
    \else
     \ifcsname collatz\the\intN\endcsname
     \else
        \expandafter\edef\csname collatz\the\intN\endcsname{\the\intL}%
        \let\next\CollUpdate
     \fi
    \fi
    \next
%% this variant seems to impact a bit negatively execution time
% \ifcsname collatz\the\intN\endcsname
%   \expandafter\CollOne
% \else
%   \ifnum\intN>\inputNmax
%   \else
%     \expandafter\edef\csname collatz\the\intN\endcsname{\the\intL}%
%   \fi
%   \expandafter\CollUpdate
% \fi
} 

\def\CollThreeB{%
    \advance\intL 1
    \ifodd\intN
      \multiply\intN 3
      \advance\intN  1
    \else
      \divide\intN 2
    \fi
    \CollTwo
}

\def\CollDone{%
    From 1 to \the\inputNmax, the longest sequence with smallest starting
    point was observed to start at \the\intNtop, and contained
    \the\intLtop\relax\
    elements.\par
}

\def\CollMax #1{% #1 integer at least 1
    %\begingroup
      \inputNmax=#1\relax
      \intNa   = 1
      \intNtop = 1
      \intLtop = 1
    % I would prefer counting steps to reach 1, so here 0
    % but it seems the question asks for number of elements, so here 1
      \expandafter\def\csname collatz1\endcsname{1}%
      \CollOne
    %\endgroup
}


\hsize10cm

\CollMax {10}

\CollMax {100}

\CollMax {1000}

\CollMax {10000}

\CollMax {100000}

\bye

输出:

在此处输入图片描述

执行时间处理时间:

$ time pdftex -interaction batchmode pcollatz.tex
This is pdfTeX, Version 3.14159265-2.6-1.40.18 (TeX Live 2017) (preloaded format=pdftex)
 restricted \write18 enabled.
entering extended mode

real    0m6.667s
user    0m6.636s
sys 0m0.027s

1000000 发生算术溢出。

我将使用添加变体信特在稍后的阶段。

扩展到大整数

实际上 159487 是导致 TeX 算术溢出的最小起点。(其第 56 次迭代是2265333694 >= 2**31

但是,在修改了上述代码以使用xintcore宏进行大运算之后,又出现了另一个问题,即“TeX 池”在 N=317893 的起点处被填满。我确实用 317892 运行了它,得到了

\CollMax {317892}
From 1 to 317892, the longest sequence with smallest starting point was
observed to start at 230631, and contained 443 elements.

然后我修改了代码,让它只将前 100000 个整数的序列长度存储在内存中,从而克服了“TeX 池”的限制。仍然通过以下方式处理大型算术新芯宏。

% pdftex (or etex)

\input xintcore.sty % to handle big integers.
% there is also bnumexpr but it is LaTeX interface only
% or xint, xintexpr which are more extensive than xintcore

\newcount\inputNmax

\newcount\intN
\newcount\intNa
\newcount\intNtop
\newcount\intL
\newcount\intLtop

\def\CollLoop{%
    \advance\intNa 1
% \immediate\write-1{\the\intNa}%
    \ifnum\intNa > \inputNmax
       \CollDone
    \else
      \intL 0 % will store number of steps starting at Na
      \intN = \intNa
      \expandafter\CollTwo
    \fi
}

\def\CollTwo{%
    \ifcsname collatz\the\intN\endcsname
       \expandafter\CollThreeA
    \else
       \expandafter\CollThreeB
    \fi
}

% \CollThreeA will be either \CollThreeAwithUpdate or \CollThreeAwithNoUpdate
\def\CollThreeAwithUpdate{%
    \advance\intL\csname collatz\the\intN\endcsname\relax
    \expandafter\edef\csname collatz\the\intNa\endcsname{\the\intL}%
    \ifnum\intL > \intLtop
          \intLtop = \intL
          \intNtop = \intNa
    \fi
    \intN = \intNa
    \CollUpdate
}

\def\CollUpdate{%
    \advance\intL -1
    \ifodd\intN
      \multiply\intN 3
      \advance\intN 1
    \else
      \divide\intN 2
    \fi
    \let\next\CollLoop
    \ifnum\intN>\inputNmax
    \else
     \ifcsname collatz\the\intN\endcsname
     \else
        \expandafter\edef\csname collatz\the\intN\endcsname{\the\intL}%
        \let\next\CollUpdate
     \fi
    \fi
    \next
} 

\def\CollThreeAwithNoUpdate{%
    \advance\intL\csname collatz\the\intN\endcsname\relax
    \ifnum\intL > \intLtop
          \intLtop = \intL
          \intNtop = \intNa
    \fi
    \CollLoop
}

\def\CollThreeB{%
    \let\next\CollTwo
    \ifodd\intN
      \ifnum\intN>\maxdimen
         % notice that necessarily this first happens with previous execution
         % had done (3x+1)/2, so the real antecedent was > "7FFFFFFF
         % and would have created arithmetic overflow if we had done
         % x->3x+1->(3x+1)/2
         \edef\bigintN{\the\intN}%
         \let\next\CollThreeBig
      \else
         \advance\intL 2
         \divide\intN 2
         \multiply\intN 3
         \advance\intN 2
      \fi
    \else
      \advance\intL 1
      \divide\intN 2
    \fi
    \next
}%

% \def\error{\immediate\write-1{\the\intNa, \the\intL}\csname end\endcsname}

\def\CollThreeBig{% 
% 159487 is the smallest starting integer which triggers this, as its
% 56th iterate 2265333694 exceeds 2**31, and the 57th is thus > \maxdimen
% \error
    \advance\intL 1
    % \xintLastItem does no expansion ...
    \ifodd\expandafter\xintLastItem\expandafter{\bigintN}
      \advance\intL 1
      \edef\bigintN{\xintHalf{\xintiiMul{\bigintN}3}}% Half truncates
      % possibly faster to use \xintDouble and an addition, not tested
    \else
      \edef\bigintN{\xintHalf{\bigintN}}%
    \fi
    % \xintLength does no expansion ...
    \ifnum\expandafter\xintLength\expandafter{\bigintN}>9
        \expandafter\CollThreeBig
    \else
        \intN = \bigintN\relax
        \expandafter\CollThreeB
    \fi
}%


\def\CollReport{%
    From 1 to \the\inputNmax, the longest sequence with smallest starting
    point was observed to start at \the\intNtop, and contained
    \the\intLtop\relax\
    elements.\par
}

\let\CollDone\CollReport

\def\CollMaxInitial {%
      \let\CollThreeA\CollThreeAwithUpdate % faster but uses macro storage
      \intNa   = 1
      \intNtop = 1
% I would prefer counting steps to reach 1, so here 0
% but it seems the question asks for number of elements, so 1
% (and not 4 although 1->4->2->1, as I consider that 1 is sequence in itself)
      \expandafter\def\csname collatz1\endcsname{1}%
      \intLtop = 1
      \CollLoop
}

\def\CollMax #1{% #1 integer at least 1
      \ifnum#1>100000
          \inputNmax 100000
          \let\CollDone\empty
          \CollMaxInitial
          \let\CollDone\CollReport
          \let\CollThreeA\CollThreeAwithNoUpdate
          \inputNmax=#1\relax
          \intNa = 100000
          \CollLoop
      \else
          \inputNmax = #1\relax
          \CollMaxInitial
      \fi
}


\hsize10cm

% \CollMax {10}

% \CollMax {100}

% \CollMax {1000}

% \CollMax {10000}

% \CollMax {100000}

% \CollMax {200000}
% From 1 to 200000, the longest sequence with smallest starting
% point was observed to start at 156159, and contained 383 elements.

\CollMax {1000000}
% From 1 to 1000000, the longest sequence with smallest start-
% ing point was observed to start at 837799, and contained 525
% elements.
\bye

这给出

在此处输入图片描述

执行时间为

$ time pdftex -interaction batchmode pcollatz-big.tex
This is pdfTeX, Version 3.14159265-2.6-1.40.18 (TeX Live 2017) (preloaded format=pdftex)
 restricted \write18 enabled.
entering extended mode

real    0m40.345s
user    0m40.247s
sys 0m0.087s

编辑:这是在一台通常比用于测试第一版代码(100000)的计算机快 15% 左右的计算机上进行的测试。但现在尝试使用 N=1000000 和较慢的计算机进行此代码,我观察到执行时间约为 1 分钟(因此 +50%...)。由于某些我不知道的原因,当 TeX 内存使用量很大时,pdftex 在我的笔记本电脑上确实变得更慢。我过去在测试xint数千位数字的计算时就观察到了这一点。

\fontdimen参数的使用

开始了

% pdftex (or etex)

\input xintcore.sty % to handle big integers.
% there is also bnumexpr but it is LaTeX interface only
% or xint, xintexpr which are more extensive than xintcore

% use fontdimen parameters to have an array where to store
% "Collatz lengths". At the end of the day
% \number\fontdimen<number>\cz gives the length of the
% sequence starting at <number> (for all numbers up to 5,000,000) and
% reaching 1 (included).

\newcount\inputNmax

\newcount\intN
\newcount\intNa
\newcount\intNtop
\newcount\intL
\newcount\intLtop

\newcount\czsize
\czsize 5000000
\font\cz=cmr10 at 1pt
\fontdimen \czsize\cz = 0sp % make room ...
% vz  texmf.cnf
% Words of font info for TeX (total size of all TFM files, approximately).
% Must be >= 20000 and <= 147483647 (without tex.ch changes).
% font_mem_size = 8000000


% make sure array entries are zero (only \fontdimen2 to 7 are populated for cmr10)
\intN 1
\loop
  \fontdimen\intN\cz = 0sp
  \advance\intN 1
\ifnum\intN < 8
\repeat

% do I need to do that for all, aren't they zero except first few ones ?
% no, it's ok (paranoide check here, done once)
% \intN 1
% \loop
%   \ifnum\fontdimen\intN\cz>0 \error\fi
% \ifnum\intN < \czsize
%   \advance\intN 1
% \repeat

\def\CollLoop{%
    \advance\intNa 1
    \ifnum\intNa > \inputNmax
       \CollDone
    \else
      \intL 0 % will store number of steps starting at Na
      \intN = \intNa
      \expandafter\CollTwo
    \fi
}

\def\CollTwo{%
  \let\next\CollThreeB
  \unless\ifnum\intN>\czsize
    \ifnum\fontdimen\intN\cz>0
       \let\next\CollThreeA
    \fi
  \fi
  \next
}

% \CollThreeA will be either \CollThreeAwithUpdate or \CollThreeAwithNoUpdate
% but we are going to use it always with update...
\def\CollThreeAwithUpdate{%
    \advance\intL\fontdimen\intN\cz
    \fontdimen\intNa\cz=\intL sp
    \ifnum\intL > \intLtop
          \intLtop = \intL
          \intNtop = \intNa
    \fi
    \intN = \intNa
    \CollUpdate
}

\def\CollUpdate{%
    \advance\intL -1
    \ifodd\intN
      \multiply\intN 3
      \advance\intN 1
    \else
      \divide\intN 2
    \fi
    \let\next\CollLoop
    \ifnum\intN>\inputNmax % always at most \cssize in this macro
    \else
     \ifnum\fontdimen\intN\cz=0
        \fontdimen\intN\cz=\intL sp
        \let\next\CollUpdate
     \fi
    \fi
    \next
} 

\def\CollThreeAwithNoUpdate{%
    \advance\intL\fontdimen\intN\cz
    \ifnum\intL > \intLtop
          \intLtop = \intL
          \intNtop = \intNa
    \fi
    \CollLoop
}

\def\CollThreeB{%
    \let\next\CollTwo
    \ifodd\intN
      \ifnum\intN>\maxdimen
         % notice that necessarily this first happens with previous execution
         % had done (3x+1)/2, so the real antecedent was > "7FFFFFFF
         % and would have created arithmetic overflow if we had done
         % x->3x+1->(3x+1)/2
         \edef\bigintN{\the\intN}%
         \let\next\CollThreeBig
      \else
         \advance\intL 2
         \divide\intN 2
         \multiply\intN 3
         \advance\intN 2
      \fi
    \else
      \advance\intL 1
      \divide\intN 2
    \fi
    \next
}%

% \def\error{\immediate\write-1{\the\intNa, \the\intL}\csname end\endcsname}

\def\CollThreeBig{% 
% 159487 is the smallest starting integer which triggers this, as its
% 56th iterate 2265333694 exceeds 2**31, and the 57th is thus > \maxdimen
% \error
    \advance\intL 1
    % \xintLastItem does no expansion ...
    \ifodd\expandafter\xintLastItem\expandafter{\bigintN}
      \advance\intL 1
      \edef\bigintN{\xintHalf{\xintiiMul{\bigintN}3}}% Half truncates
      % possibly faster to use \xintDouble and an addition, not tested
    \else
      \edef\bigintN{\xintHalf{\bigintN}}%
    \fi
    % \xintLength does no expansion ...
    \ifnum\expandafter\xintLength\expandafter{\bigintN}>9
        \expandafter\CollThreeBig
    \else
        \intN = \bigintN\relax
        \expandafter\CollThreeB
    \fi
}%


\def\CollReport{%
    From 1 to \the\inputNmax, the longest sequence with smallest starting
    point was observed to start at \the\intNtop, and contained
    \the\intLtop\relax\
    elements.\par
}

\let\CollDone\CollReport

\def\CollMaxInitial {%
      \let\CollThreeA\CollThreeAwithUpdate % limited by font array storage
      \intNa   = 1
      \intNtop = 1
% I would prefer counting steps to reach 1, so here 0
% but it seems the question asks for number of elements, so 1
% (and not 4 although 1->4->2->1, as I consider that 1 is sequence in itself)
      \fontdimen1 \cz= 1sp
      \intLtop = 1
      \CollLoop
}

\def\CollMax #1{% #1 integer at least 1
      \ifnum#1>\czsize
          \inputNmax \czsize
          \let\CollDone\empty
          \CollMaxInitial
          \let\CollDone\CollReport
          \let\CollThreeA\CollThreeAwithNoUpdate
          \inputNmax=#1\relax
          \intNa = \czsize
          \CollLoop
      \else
          \inputNmax = #1\relax
          \CollMaxInitial
      \fi
}


\hsize10cm

% check it does give same results as earlier...
% \CollMax {10}

% \CollMax {100}

% \CollMax {1000}

% \CollMax {10000}

%\CollMax {100000}

%\CollMax {1000000}

\CollMax {5000000}

% From 1 to 5000000, the longest sequence with smallest starting point was
% observed to start at 2929311, and contained 550 elements.
\bye

生成:

在此处输入图片描述

对于 1000000,它比使用 csname 的方法快 4 到 6 倍(很难说清楚,因为这取决于我使用哪台计算机进行测试)。也许我可以改进更多,我不太习惯那种技术。

这是用我的笔记本电脑(“慢速”计算机)计算 N=1,000,000 的结果:

$ time pdftex --interaction=batchmode pcollatz-array.tex
This is pdfTeX, Version 3.14159265-2.6-1.40.18 (TeX Live 2017) (preloaded format=pdftex)
 restricted \write18 enabled.
entering extended mode

real    0m10.408s
user    0m10.313s
sys 0m0.041s

对于 5,000,000 来说,所需的资金大约是 1,000,000 的五倍。


最终的重写带来了45%速度的提升。

\input xintcore.sty

\newcount\inputNmax
\newcount\inputNstop

\newcount\intN
\newcount\intNa
\newcount\intNb
\newcount\intNtop
\newcount\intL
\newcount\intLtop

\newcount\czsize
\czsize 5000000
\font\cz=cmr10 at 1pt
\fontdimen \czsize\cz = 0sp

\intN 1
\loop
  \fontdimen\intN\cz = 0sp
  \advance\intN 1
\ifnum\intN < 8
\repeat

\catcode`@ 11

\long\def\@gobble#1{}
\edef\sp{\expandafter\@gobble\string\sp}% (no space needed?)

\def\CollLoopA{%
    \advance\intNa 1
    \intNb \intNa
    \divide\intNb 2
    \intL\fontdimen\intNb\cz
    \advance\intL 1
    \fontdimen\intNa\cz=\intL\sp
    \advance\intNa 1
    \intNb \intNa
    \divide\intNb 3
    \multiply\intNb 2
    \advance\intNb 1
    \intL\fontdimen\intNb\cz
    \advance\intL -2
    \fontdimen\intNa\cz=\intL\sp
    \advance\intNa 1        
    \ifnum\intNa > \inputNstop
       \CollDone
    \else
      \let\CollLoopBack\CollLoopB
      \intL 0
      \intN = \intNa
      \expandafter\CollThreeB
    \fi
}

\def\CollLoopB{%
    \advance\intNa 1
    \ifnum\intNa > \inputNstop
       \CollDone
    \else
      \let\CollLoopBack\CollLoopC
      \intL 0
      \intN = \intNa
      \expandafter\CollThreeB
    \fi
}

\def\CollLoopC{%
    \advance\intNa 1
    \intNb \intNa
    \divide\intNb 2
    \intL\fontdimen\intNb\cz
    \advance\intL 1
    \fontdimen\intNa\cz=\intL\sp
    \advance\intNa 1
    \ifnum\intNa > \inputNstop
       \CollDone
    \else
      \let\CollLoopBack\CollLoopA
      \intL 0
      \intN = \intNa
      \expandafter\CollThreeB
    \fi
}

\def\CollTwo{%
  \ifnum\intN<\intNa
    \advance\intL\fontdimen\intN\cz
    \fontdimen\intNa\cz=\intL\sp
    \CollCheckTop
    \expandafter\CollLoopBack
  \else
    \expandafter\CollThreeB 
  \fi
}

\def\CollCheckTopyes {%
    \ifnum\intL > \intLtop
          \intLtop = \intL
          \intNtop = \intNa
    \fi
}%
\let\CollCheckTopno\empty
\let\CollCheckTop\CollCheckTopyes

\def\CollThreeB{%
    \let\next\CollTwo
    \ifodd\intN
      \ifnum\intN>\maxdimen
         \edef\bigintN{\the\intN}%
         \let\next\CollThreeBig
      \else
         \advance\intL 2
         \divide\intN 2
         \multiply\intN 3
         \advance\intN 2
      \fi
    \else
      \advance\intL 1
      \divide\intN 2
    \fi
    \next
}%

\def\CollThreeBig{% 
    \advance\intL 1
    % \xintLastItem does no expansion ...
    \ifodd\expandafter\xintLastItem\expandafter{\bigintN}
      \advance\intL 1
      \edef\bigintN{\xintHalf{\xintiiMul{\bigintN}3}}% Half truncates
      % possibly faster to use \xintDouble and an addition, not tested
    \else
      \edef\bigintN{\xintHalf{\bigintN}}%
    \fi
    % \xintLength does no expansion ...
    \ifnum\expandafter\xintLength\expandafter{\bigintN}>9
        \expandafter\CollThreeBig
    \else
        \intN = \bigintN\relax
        \expandafter\CollThreeB
    \fi
}%


\def\CollReport{%
    From 1 to \the\inputNmax, the longest sequence with smallest starting
    point was observed to start at \the\intNtop, and contained
    \the\intLtop\relax\
    elements.\par
}

\let\CollDone\empty

\def\CollMax #1{%
      \inputNmax #1\relax
      \inputNstop \inputNmax
      \divide\inputNstop 2
      \intNa = 1
      \intNtop = 1
      \fontdimen1 \cz= 1sp
      \intLtop = 1
      \let\CollCheckTop\CollCheckTopno
      \CollLoopC
      \inputNstop \inputNmax
      \let\CollCheckTop\CollCheckTopyes
      \ifcase \numexpr3+\intNa - 6*(\intNa/6)\relax
         \advance\intNa-2 \let\next\CollLoopC\or
         \error\or
         \error\or
         \advance\intNa-3 \let\next\CollLoopA\or
         \advance\intNa-1 \let\next\CollLoopB\else
         \error
      \fi 
      \next
      \CollReport
}
\catcode`@ 12

\hsize10cm

\CollMax {1000000}

\bye

与前者相比,time pdftex --interaction=batchmode pcollatz-arrayIII.tex我得到的用户时间为。在@HenriMenke的机器上,这种纯eTeX方法的执行时间应该在(对于)左右。0m5.780s0m10.313s2sN=1000000

相关内容