我正在尝试将算法中概念上简单的部分从 -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
}
}
的结果值为i
23,如果 ,情况确实如此。也就是说,我们通过将其设置为大于( )而不是大于或等于( ) 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
也可以,只要我可以在、 等里面重用最终结果pgfmath
,tikzpictures
和只要最终结果至少与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 来说令人兴奋:-)