旧 C 代码的 TeX 翻译 - 在 foreach 循环*内部*评估数组元素

旧 C 代码的 TeX 翻译 - 在 foreach 循环*内部*评估数组元素

我正在尝试将算法中概念上简单的部分从 -engine 的原始源代码R(用 编写C)转换为 TeX。这应该只需要一些逻辑树和 for/while 循环。

根据@percusse 的评论,我将解释完整的用例以了解我想要做什么。

我希望模仿的代码如下:

int attribute_hidden chebyshev_init(double *dos, int nos, double eta)
{
    int i, ii;
    double err;

    if (nos < 1)
        return 0;

    err = 0.0;
    i = 0;  /* just to avoid compiler warnings */
    for (ii=1; ii<=nos; ii++) {
        i = nos - ii;
        err += fabs(dos[i]);
        if (err > eta) {
        return i;
        }
    }
    return i;
}

来源:https://github.com/wch/r-source/blob/776708efe6003e36f02587ad47b2eaaaa19e2f69/src/nmath/chebyshev.c#L48(第 48-66 行)

实际用例可以在以下位置找到:https://github.com/wch/r-source/blob/0b493ddbefbc3cc06f675841c4672b3197dbfbbb/src/nmath/gamma.c#L103. 论据表明:

ngam = chebyshev_init(gamcs, 42, DBL_EPSILON/20)

Gamcs定义为:

const static double gamcs[42] = {
+.8571195590989331421920062399942e-2,
+.4415381324841006757191315771652e-2,
+.5685043681599363378632664588789e-1,
-.4219835396418560501012500186624e-2,
+.1326808181212460220584006796352e-2,
-.1893024529798880432523947023886e-3,
+.3606925327441245256578082217225e-4,
-.6056761904460864218485548290365e-5,
+.1055829546302283344731823509093e-5,
-.1811967365542384048291855891166e-6,
+.3117724964715322277790254593169e-7,
-.5354219639019687140874081024347e-8,
+.9193275519859588946887786825940e-9,
-.1577941280288339761767423273953e-9,
+.2707980622934954543266540433089e-10,
-.4646818653825730144081661058933e-11,
+.7973350192007419656460767175359e-12,
-.1368078209830916025799499172309e-12,
+.2347319486563800657233471771688e-13,
-.4027432614949066932766570534699e-14,
+.6910051747372100912138336975257e-15,
-.1185584500221992907052387126192e-15,
+.2034148542496373955201026051932e-16,
-.3490054341717405849274012949108e-17,
+.5987993856485305567135051066026e-18,
-.1027378057872228074490069778431e-18,
+.1762702816060529824942759660748e-19,
-.3024320653735306260958772112042e-20,
+.5188914660218397839717833550506e-21,
-.8902770842456576692449251601066e-22,
+.1527474068493342602274596891306e-22,
-.2620731256187362900257328332799e-23,
+.4496464047830538670331046570666e-24,
-.7714712731336877911703901525333e-25,
+.1323635453126044036486572714666e-25,
-.2270999412942928816702313813333e-26,
+.3896418998003991449320816639999e-27,
-.6685198115125953327792127999999e-28,
+.1146998663140024384347613866666e-28,
-.1967938586345134677295103999999e-29,
+.3376448816585338090334890666666e-30,
-.5793070335782135784625493333333e-31
    }

DBL_EPSILON定义为DBL_EPSILON = 2^(-52)

要知道我想要的结果,我们可以模拟的整个计算chebyshev_init。为此,您可以在中运行以下代码R

DBL_EPSILON = 2^(-52)
err = 0.0
i = 0
dos = c(+.8571195590989331421920062399942e-2,
+.4415381324841006757191315771652e-2,
+.5685043681599363378632664588789e-1,
-.4219835396418560501012500186624e-2,
+.1326808181212460220584006796352e-2,
-.1893024529798880432523947023886e-3,
+.3606925327441245256578082217225e-4,
-.6056761904460864218485548290365e-5,
+.1055829546302283344731823509093e-5,
-.1811967365542384048291855891166e-6,
+.3117724964715322277790254593169e-7,
-.5354219639019687140874081024347e-8,
+.9193275519859588946887786825940e-9,
-.1577941280288339761767423273953e-9,
+.2707980622934954543266540433089e-10,
-.4646818653825730144081661058933e-11,
+.7973350192007419656460767175359e-12,
-.1368078209830916025799499172309e-12,
+.2347319486563800657233471771688e-13,
-.4027432614949066932766570534699e-14,
+.6910051747372100912138336975257e-15,
-.1185584500221992907052387126192e-15,
+.2034148542496373955201026051932e-16,
-.3490054341717405849274012949108e-17,
+.5987993856485305567135051066026e-18,
-.1027378057872228074490069778431e-18,
+.1762702816060529824942759660748e-19,
-.3024320653735306260958772112042e-20,
+.5188914660218397839717833550506e-21,
-.8902770842456576692449251601066e-22,
+.1527474068493342602274596891306e-22,
-.2620731256187362900257328332799e-23,
+.4496464047830538670331046570666e-24,
-.7714712731336877911703901525333e-25,
+.1323635453126044036486572714666e-25,
-.2270999412942928816702313813333e-26,
+.3896418998003991449320816639999e-27,
-.6685198115125953327792127999999e-28,
+.1146998663140024384347613866666e-28,
-.1967938586345134677295103999999e-29,
+.3376448816585338090334890666666e-30,
-.5793070335782135784625493333333e-31)
nos = 42
eta = DBL_EPSILON/20
for (ii in 1:nos) {
i = nos - ii
err = err + abs(dos[i])
if (err>eta)
{
break
}
}

的结果值为i23,如果 ,情况确实如此。也就是说,我们通过将其设置为大于( )而不是大于或等于( ) err=>eta,将预期值提高了一。预期结果还显示>=>https://github.com/wch/r-source/blob/0b493ddbefbc3cc06f675841c4672b3197dbfbbb/src/nmath/gamma.c#L115。是的,我知道知道结果会让这个练习变得毫无用处,但是这个案例是有解释力的,而且算法不仅仅在这种情况下是必需的。

所以R实际上确实给出了正确的结果。我制作了这个算法的自己的版本,理想情况下看起来有点像这样:

\documentclass{article}
\usepackage{tikz}
\usetikzlibrary{fpu}
\usepackage{xintexpr}
\def\gamcs{{
+.8571195590989331421920062399942e-2,
+.4415381324841006757191315771652e-2,
+.5685043681599363378632664588789e-1,
-.4219835396418560501012500186624e-2,
+.1326808181212460220584006796352e-2,
-.1893024529798880432523947023886e-3,
+.3606925327441245256578082217225e-4,
-.6056761904460864218485548290365e-5,
+.1055829546302283344731823509093e-5,
-.1811967365542384048291855891166e-6,
+.3117724964715322277790254593169e-7,
-.5354219639019687140874081024347e-8,
+.9193275519859588946887786825940e-9,
-.1577941280288339761767423273953e-9,
+.2707980622934954543266540433089e-10,
-.4646818653825730144081661058933e-11,
+.7973350192007419656460767175359e-12,
-.1368078209830916025799499172309e-12,
+.2347319486563800657233471771688e-13,
-.4027432614949066932766570534699e-14,
+.6910051747372100912138336975257e-15,
-.1185584500221992907052387126192e-15,
+.2034148542496373955201026051932e-16,
-.3490054341717405849274012949108e-17,
+.5987993856485305567135051066026e-18,
-.1027378057872228074490069778431e-18,
+.1762702816060529824942759660748e-19,
-.3024320653735306260958772112042e-20,
+.5188914660218397839717833550506e-21,
-.8902770842456576692449251601066e-22,
+.1527474068493342602274596891306e-22,
-.2620731256187362900257328332799e-23,
+.4496464047830538670331046570666e-24,
-.7714712731336877911703901525333e-25,
+.1323635453126044036486572714666e-25,
-.2270999412942928816702313813333e-26,
+.3896418998003991449320816639999e-27,
-.6685198115125953327792127999999e-28,
+.1146998663140024384347613866666e-28,
-.1967938586345134677295103999999e-29,
+.3376448816585338090334890666666e-30,
-.5793070335782135784625493333333e-31
    }}
\xintNewExpr\DBLEPSILON{\xintPow{2}{-52}}
\ExplSyntaxOn
\cs_new:Npn \breakloop #1 #2
  {
    \fp_compare:nNnT { #1 } > { #2 }
      { \breakforeach }
  }
\ExplSyntaxOff
\makeatletter
\pgfmathdeclarefunction{chebyshevinit}{3}{%
  \begingroup
    \pgfkeys{/pgf/fpu,/pgf/fpu/output format=fixed}
      \pgfmathparse{#1}%
      \edef\dos{\pgfmathresult}%
      \pgfmathparse{#2-1}% This will allow to refer to the first element as 1 instead of 0
      \edef\nos{\pgfmathresult}%
      \pgfmathparse{#3}%
      \edef\eta{\pgfmathresult}%
        \pgfmathparse{0}%
        \edef\err{\pgfmathresult}%
        \edef\i{\pgfmathresult}%
          \foreach \ii in {1,...,\nos}
          {
            \pgfmathparse{\nos-\ii}%
            \xdef\i{\pgfmathresult}%
            \pgfmathparse{\err}%
            \xdef\err@{\pgfmathresult}%
            \pgfmathparse{\err@ + abs(\dos[\i])}%
            \xdef\err{\pgfmathresult}%
            \breakloop{\err}{\eta}%
          }
        \pgfmathparse{\i}%
    \pgfmath@smuggleone\pgfmathresult%
  \endgroup
}
\makeatother
\begin{document}
\pgfmathparse{chebyshevinit(\gamcs,42,\DBLEPSILON/20)}\pgfmathresult    %
\end{document}

显然上面的例子不起作用:数组不是一个实际的“值”,而是包含以各种方式扩展的多个值,这使得它在设置中使用起来相当困难\pgfmathparse,在实际函数中使用起来更困难(Carlisle,2014)。

我曾经遇到过 Andrew Swann 解决方案的一个更简单的案例PGF 循环 - 如何在 pgfmathdeclare 函数中使用数组?,但这会导致下溢。该fpu库无法处理较小的数字!

我在之前的问题中也提出过这样的想法:将循环的某些部分移到循环之外,然后首先进行评估。

但是,这现在还不能回答我的具体问题。如果我必须评估一系列元素(现在就是这种情况),那么事先评估值就行不通了。

那么我如何在 中进行这样的评估pgfmath? 外面的解决方案pgfmath也可以,只要我可以在、 等里面重用最终结果pgfmathtikzpictures只要最终结果至少与R输出的某些数字相似即可。(10^{-4}?此时我会选择任何东西。:-/


编辑我稍微把这个问题带回到了本质。这是另一个有同样问题的例子:如何在循环内评估数组?考虑以下代码:

\documentclass{article}
\usepackage{tikz}
\usetikzlibrary{fpu}    
\gdef\myarray{{1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,16,17,18,19,20,21,22,23,24,25,26,27,28,29,30,31,32,33,34,35,36,37,38,39,40,41,42}}

\makeatletter
\pgfmathdeclarefunction{blah}{2}{%
  \begingroup
    \pgfkeys{/pgf/fpu,/pgf/fpu/output format=fixed}
    \pgfmathparse{#1}%
    \edef\x{\pgfmathresult}%
    \pgfmathparse{#2}%
    \edef\n{\pgfmathresult}%
    \pgfmathparse{\x * 2}%
    \edef\twox{\pgfmathresult}%
    \pgfmathparse{0}%
    \edef\b{\pgfmathresult}%
    \edef\b@{\pgfmathresult}%
    \edef\b@@{\pgfmathresult}%
      \foreach \i in {1,...,\n}
      {
      \pgfmathparse{\b@}%
      \xdef\b@@{\pgfmathresult}%
      \pgfmathparse{\b}%
      \xdef\b@{\pgfmathresult}%
      \pgfmathparse{\myarray[\n-\i]}%
      \xdef\arrayvalue{\pgfmathresult}%
      \pgfmathparse{\twox * \b@ - \b@@ + \arrayvalue}%
      \xdef\b{\pgfmathresult}%
      }
    \pgfmathparse{(\b-\b@@) * .5}%
    \pgfmath@smuggleone\pgfmathresult
  \endgroup
}
\makeatother

\begin{document}

\pgfmathparse{blah(5,22)}\pgfmathresult %

\end{document}

因此,如您所见,Andrew Swann 的技巧在这里不适用,我正在努力寻找解决方法。如何myarray在我的 中使用foreach

答案1

这里有一个 luatex 实现。它非常简单:

\documentclass{article}
% Define the function in lua, and the constant DBL_EPSILON
\directlua{
DBL_EPSILON = math.pow(2,-52)
function chebyshev_init(dos, nos, eta)
    local i, ii, err;
    if nos<1 then
        return 0
    end

    err = 0.0
    i = 0
    for ii=0,nos do
        i = nos - ii
        err = err + math.abs(dos[i])
        if err > eta then
            return i
        end
    end
    return i
end
}

% Define a tex macro to call the lua function
% and "print" the result
\def\chebyshevInit#1#2#3{%
   \directlua{tex.print(chebyshev_init(#1,#2,#3))}
}

% Define example array to test if it works
\xdef\gmcs{{
     .8571195590989331421920062399942e-2,
     .4415381324841006757191315771652e-2,
     .5685043681599363378632664588789e-1,
    -.4219835396418560501012500186624e-2,
     .1326808181212460220584006796352e-2,
    -.1893024529798880432523947023886e-3,
     .3606925327441245256578082217225e-4,
    -.6056761904460864218485548290365e-5,
     .1055829546302283344731823509093e-5,
    -.1811967365542384048291855891166e-6,
     .3117724964715322277790254593169e-7,
    -.5354219639019687140874081024347e-8,
     .9193275519859588946887786825940e-9,
    -.1577941280288339761767423273953e-9,
     .2707980622934954543266540433089e-10,
    -.4646818653825730144081661058933e-11,
     .7973350192007419656460767175359e-12,
    -.1368078209830916025799499172309e-12,
     .2347319486563800657233471771688e-13,
    -.4027432614949066932766570534699e-14,
     .6910051747372100912138336975257e-15,
    -.1185584500221992907052387126192e-15,
     .2034148542496373955201026051932e-16,
    -.3490054341717405849274012949108e-17,
     .5987993856485305567135051066026e-18,
    -.1027378057872228074490069778431e-18,
     .1762702816060529824942759660748e-19,
    -.3024320653735306260958772112042e-20,
     .5188914660218397839717833550506e-21,
    -.8902770842456576692449251601066e-22,
     .1527474068493342602274596891306e-22,
    -.2620731256187362900257328332799e-23,
     .4496464047830538670331046570666e-24,
    -.7714712731336877911703901525333e-25,
     .1323635453126044036486572714666e-25,
    -.2270999412942928816702313813333e-26,
     .3896418998003991449320816639999e-27,
    -.6685198115125953327792127999999e-28,
     .1146998663140024384347613866666e-28,
    -.1967938586345134677295103999999e-29,
     .3376448816585338090334890666666e-30,
    -.5793070335782135784625493333333e-31
}}

\begin{document}
The result is: \chebyshevInit{\gmcs}{42}{DBL_EPSILON/20}
\end{document}

用 编译它lualatex,你会得到以下结果(抱歉,粘贴的图像在视觉上很单调,但希望对 OP 来说令人兴奋:-)

结果

相关内容