summaryrefslogtreecommitdiff
path: root/Master/texmf-dist/source/latex/l3kernel/l3fp-basics.dtx
diff options
context:
space:
mode:
Diffstat (limited to 'Master/texmf-dist/source/latex/l3kernel/l3fp-basics.dtx')
-rw-r--r--Master/texmf-dist/source/latex/l3kernel/l3fp-basics.dtx2135
1 files changed, 1072 insertions, 1063 deletions
diff --git a/Master/texmf-dist/source/latex/l3kernel/l3fp-basics.dtx b/Master/texmf-dist/source/latex/l3kernel/l3fp-basics.dtx
index 0b06d377232..591ccf5615e 100644
--- a/Master/texmf-dist/source/latex/l3kernel/l3fp-basics.dtx
+++ b/Master/texmf-dist/source/latex/l3kernel/l3fp-basics.dtx
@@ -35,8 +35,8 @@
%%
%
%<*driver>
-\RequirePackage{l3names}
-\GetIdInfo$Id: l3fp-basics.dtx 3986 2012-07-15 19:23:51Z joseph $
+\RequirePackage{l3bootstrap}
+\GetIdInfo$Id: l3fp-basics.dtx 4089 2012-08-14 04:52:20Z bruno $
{L3 Floating-point arithmetic}
\documentclass[full]{l3doc}
\begin{document}
@@ -77,63 +77,34 @@
%<@@=fp>
% \end{macrocode}
%
+% The \pkg{l3fp-basics} module implements addition, subtraction,
+% multiplication, and division of two floating points, and the absolute
+% value and sign-changing operations on one floating point.
% All operations implemented in this module yield the outcome of
% rounding the infinitely precise result of the operation to the
-% nearest representable number.
-%
-% ^^A begin[todo]: move
-% \section{Internal storage of floating points numbers}
-%
-% A floating point number \meta{X} is stored as
-% \begin{quote}
-% \cs{s_@@} \cs{@@_chk:w} \meta{case} \meta{sign} \meta{body} |;|
-% \end{quote}
-% Here, \meta{case} is 0 for $\pm 0$, 1 for normal numbers, 2 for $\pm
-% \infty$, and 3 for \texttt{nan}, and \meta{sign} is $0$ for positive
-% numbers, $1$ for \texttt{nan}s, and $2$ for negative numbers. The
-% \meta{body} of normal numbers is \Arg{exponent} \Arg{X_1} \Arg{X_2}
-% \Arg{X_3} \Arg{X_4}, with
-% \[
-% \meta{X} = (-1)^{\meta{sign}} 10^{-\meta{exponent}} \sum_i
-% \meta{X_i} 10^{-4i}.
-% \]
-% Calculations are done in base $10000$, \emph{i.e.} one myriad. The
-% \meta{exponent} lies between $\pm\cs{c_@@_max_exponent_int} = \pm
-% \the\csname\detokenize{c__fp_max_exponent_int}\endcsname$ inclusive.
-%
-% Additionally, positive and negative floating point numbers may only be
-% stored with $1000\leq\meta{X_1}<10000$. This requirement is necessary
-% in order to preserve accuracy and speed.
-%
-% ^^A end[todo]
-%
-% ^^A begin[todo]
+% nearest floating point.
%
% Some algorithms used below end up being quite similar to some
% described in \enquote{What Every Computer Scientist Should Know About
% Floating Point Arithmetic}, by David Goldberg, which can be found at
-% \texttt{http://cr.yp.to/2005-590/goldberg.pdf}. I need to compare them
-% very carefully.
-%
-% ^^A end[todo]
-%
-%^^A todo sanitize, pack.
+% \texttt{http://cr.yp.to/2005-590/goldberg.pdf}.
%
% \subsection{Common to several operations}
%
-% \begin{macro}[EXP]
+% \begin{macro}[aux, EXP]
% {
-% \@@_basics_pack_low:NNNNNw ,
-% \@@_basics_pack_high:NNNNNw ,
+% \@@_basics_pack_low:NNNNNw,
+% \@@_basics_pack_high:NNNNNw,
% \@@_basics_pack_high_carry:w
% }
% Addition and multiplication of mantissas are done in two steps:
-% first compute a (more or less) exact result,
-% then round and pack digits in the final (braced) form.
-% These functions take care of the packing, with special attention
-% given to the case where rounding has caused a carry.
-% In \cs{@@_basics_pack_high_carry:w}, |#1| should
-% always be $0000$.
+% first compute a (more or less) exact result, then round and pack
+% digits in the final (braced) form. These functions take care of the
+% packing, with special attention given to the case where rounding has
+% caused a carry. Since rounding can only shift the final digit by
+% $1$, a carry always produces an exact power of $10$. Thus,
+% \cs{@@_basics_pack_high_carry:w} is always followed by four times
+% |{0000}|.
% \begin{macrocode}
\cs_new:Npn \@@_basics_pack_low:NNNNNw #1 #2#3#4#5 #6;
{
@@ -154,284 +125,275 @@
% \end{macrocode}
% \end{macro}
%
-% \begin{macro}
+% \begin{macro}[aux, EXP]
% {
-% \@@_basics_return_nan_nan:NNww ,
-% \@@_basics_return_zero:NNww ,
-% \@@_basics_return_inf:NNww ,
-% \@@_basics_return_i:NNNNww ,
-% \@@_basics_return_ii:NNNNww ,
-% \@@_basics_return_nan:NNNNww
+% \@@_basics_pack_weird_low:NNNNw,
+% \@@_basics_pack_weird_high:NNNNNNNNw
% }
-% \begin{syntax}
-% \cs{@@_basics_return_...:NNww}
-% ~~\meta{sign_1} \meta{sign_2} \meta{body_1} |;| \meta{body_2} |;|
-% \end{syntax}
-% Used for binary operations, to return a value for some special
-% cases (common to several operations). All functions expand once
-% after their arguments.
-%
-% ^^A todo: redoc, changed.
-% The \texttt{nan_nan} function combines the \texttt{info} fields
-% of the two \texttt{nan}.\footnote{Bruno: check that messages are
-% kept.}
-% The \texttt{zero} and \texttt{inf} functions return $\pm 0$ or
-% $\pm\infty$ with a sign equal to the product of the two signs:
-% three \cs{exp_after:wN} are needed to escape out of the
-% conditional, and expand once after.
-% The \texttt{i} and \texttt{ii} functions return one of their
-% operands and expand after using \cs{@@_exp_after_o:w}.
-% In some cases, this could be optimized, since we know in advance
-% what case of number we have. However, it seems better to keep the
-% number of control sequences low: these functions are called only
-% in special cases anyways, so performance is not an issue.
+% I don't fully understand those functions, used for additions and
+% divisions. Hence the name.
% \begin{macrocode}
-\cs_new:Npn \@@_basics_return_nan_nan:NNww #1#2 #3; #4;
- { \@@_exp_after_o:w \s_@@ \@@_chk:w 3 1 #3 ; }
-\cs_new:Npn \@@_basics_return_zero:NNww #1#2 #3; #4;
- {
- \if_meaning:w #1 #2
- \exp_after:wN \exp_after:wN \exp_after:wN \c_zero_fp
- \else:
- \exp_after:wN \exp_after:wN \exp_after:wN \c_minus_zero_fp
- \fi:
- }
-\cs_new:Npn \@@_basics_return_inf:NNww #1#2 #3; #4;
+\cs_new:Npn \@@_basics_pack_weird_low:NNNNw #1 #2#3#4 #5;
{
- \if_meaning:w #1 #2
- \exp_after:wN \exp_after:wN \exp_after:wN \c_inf_fp
- \else:
- \exp_after:wN \exp_after:wN \exp_after:wN \c_minus_inf_fp
- \fi:
- }
-\cs_new:Npn \@@_basics_return_i:NNNNww #1#2 #3#4 #5; #6;
- { \@@_exp_after_o:w \s_@@ \@@_chk:w #1 #3 #5; }
-\cs_new:Npn \@@_basics_return_ii:NNNNww #1#2 #3#4 #5; #6;
- { \@@_exp_after_o:w \s_@@ \@@_chk:w #2 #4 #6; }
-\cs_new:Npn \@@_basics_return_nan:NNww #1#2
- {
- \if_meaning:w 1 #1
- \exp_after:wN \@@_basics_return_i:NNNNww
- \else:
- \exp_after:wN \@@_basics_return_ii:NNNNww
+ \if_meaning:w 2 #1
+ + \c_one
\fi:
- 3 3 #1 #2
+ \__int_eval_end:
+ #2#3#4; {#5} ;
}
+\cs_new:Npn \@@_basics_pack_weird_high:NNNNNNNNw
+ 1 #1#2#3#4 #5#6#7#8 #9; { ; {#1#2#3#4} {#5#6#7#8} {#9} }
% \end{macrocode}
% \end{macro}
%
% \subsection{Addition and subtraction}
%
-% \begin{macro}[int, EXP]{\@@_+_o:ww}
-% For addition, everything is easy. No need to grab the
-% \meta{body_2}.
-% \begin{macrocode}
-\cs_new:cpn { @@_+_o:ww }
- \s_@@ \@@_chk:w #1 #2 #3 ; \s_@@ \@@_chk:w #4 #5
- { \@@_add_cases:NN #1 #4 #2 #5 #3 ; }
-% \end{macrocode}
-% \end{macro}
+% We define here two functions, \cs{@@_-_o:ww} and \cs{@@_+_o:ww}, which
+% perform the subtraction and addition of their two floating point
+% operands, and expand the tokens following the result once.
+%
+% A more obscure function, \cs{@@_add_big_i_o:wNww}, is used in
+% \pkg{l3fp-expo}.
+%
+% The logic goes as follows:
+% \begin{itemize}
+% \item \cs{@@_-_o:ww} calls \cs{@@_+_o:ww} to do the work, with the
+% sign of the second operand flipped;
+% \item \cs{@@_+_o:ww} dispatches depending on the type of floating
+% point, calling specialized auxiliaries;
+% \item in all cases except summing two normal floating point numbers,
+% we return one or the other operands depending on the signs, or
+% detect an invalid operation in the case of $\infty - \infty$;
+% \item for normal floating point numbers, compare the signs;
+% \item to add two floating point numbers of the same sign or of
+% opposite signs, shift the mantissa of the smaller one to match the
+% bigger one, perform the addition or subtraction of mantissas,
+% check for a carry, round, and pack using the
+% |\__fp_basics_pack_...| functions.
+% \end{itemize}
+% The trickiest part is to round correctly when adding or subtracting
+% normal floating point numbers.
+%
+% \subsubsection{Sign, exponent, and special numbers}
%
% \begin{macro}[int, EXP]{\@@_-_o:ww}
-% Change the sign of the second argument.
+% A previous version of this function grabbed its two operands,
+% changed the sign of the second, and called \cs{@@_+_o:ww}. However,
+% for efficiency reasons, the operands were swapped in the process,
+% which means that error messages ended up wrong. Now, the
+% \cs{@@_+_o:ww} auxiliary has a hook: it takes one argument between
+% the first \cs{s_@@} and \cs{@@_chk:w}, which is applied to the sign
+% of the second operand. Positioning the hook there means that
+% \cs{@@_+_o:ww} can still check that it was followed by \cs{s_@@} and
+% not arbitrary junk.
% \begin{macrocode}
-\cs_new:cpn { @@_-_o:ww }
- \s_@@ \@@_chk:w #1 #2 #3 ; \s_@@ \@@_chk:w #4 #5
+\cs_new_nopar:cpx { @@_-_o:ww } \s_@@
{
- \exp_after:wN \@@_add_cases:NN
- \exp_after:wN #1
- \exp_after:wN #4
- \exp_after:wN #2
- \int_use:N \__int_eval:w \c_two - #5 \__int_eval_end:
- #3 ;
+ \exp_not:c { @@_+_o:ww }
+ \exp_not:n { \s_@@ \@@_neg_sign:N }
}
% \end{macrocode}
% \end{macro}
%
-% \subsubsection{Sign, exponent, and special numbers}
-%
-% \begin{macro}[EXP,aux]{\@@_add_cases:NN}
-% \begin{syntax}
-% \cs{@@_add_cases:NN} \meta{case_1} \meta{case_2}
-% ~~\meta{sign_1} \meta{sign_2} \meta{body_1} |;| \meta{body_2} |;|
-% \end{syntax}
-% This performs the addition. it also expands the following tokens
-% on the input stream once.
-%
-% Whenever \meta{case_1} is different from \meta{case_2}, the result
-% is simply the floating point number with the highest \meta{case}.
-% For instance, adding a normal number to a zero gives the normal
-% number, and adding a \texttt{nan} to any non-\texttt{nan} gives
-% that \texttt{nan}. Optimizing for addition of normal numbers,
-% we test for equality and then separate the \enquote{greater than}
-% and \enquote{less than} branches.
+% \begin{macro}[int, EXP]{\@@_+_o:ww}
+% This function is either called directly with an empty |#1| to
+% compute an addition, or it is called by \cs{@@_-_o:ww} with
+% \cs{@@_neg_sign:N} as |#1| to compute a subtraction (equivalent to
+% changing the \meta{sign_2} of the second operand). If the
+% \meta{types} |#2| and |#4| are the same, dispatch to case |#2| ($0$,
+% $1$, $2$, or $3$), where we call specialized functions: thanks to
+% \cs{__int_value:w}, those receive the tweaked \meta{sign_2}
+% (expansion of |#1#5|) as an argument. If the \meta{types} are
+% distinct, the result is simply the floating point number with the
+% highest \meta{type}. Since case $3$ (used for two \texttt{nan})
+% also picks the first operand, we can also use it when \meta{type_1}
+% is greater than \meta{type_2}. Also note that we don't need to
+% worry about \meta{sign_2} in that case since the second operand is
+% discarded.
% \begin{macrocode}
-\cs_new:Npn \@@_add_cases:NN #1 #2
+\cs_new:cpn { @@_+_o:ww }
+ \s_@@ #1 \@@_chk:w #2 #3 ; \s_@@ \@@_chk:w #4 #5
{
- \if_int_compare:w #1 = #2 \exp_stop_f:
- \exp_after:wN \@@_add_cases_eq:N
- \else:
- \if_int_compare:w #1 < #2 \exp_stop_f:
- \exp_after:wN \exp_after:wN
- \exp_after:wN \@@_basics_return_ii:NNNNww
+ \if_case:w
+ \if_meaning:w #2 #4
+ #2 \exp_stop_f:
\else:
- \exp_after:wN \exp_after:wN
- \exp_after:wN \@@_basics_return_i:NNNNww
+ \if_int_compare:w #2 > #4 \exp_stop_f:
+ \c_three
+ \else:
+ \c_minus_one
+ \fi:
\fi:
- \exp_after:wN #1
+ \exp_after:wN \@@_add_zeros_o:Nww \__int_value:w
+ \or: \exp_after:wN \@@_add_normal_o:Nww \__int_value:w
+ \or: \exp_after:wN \@@_add_inf_o:Nww \__int_value:w
+ \or: \@@_case_return_i_o:ww
+ \else: \exp_after:wN \@@_add_return_ii_o:Nww \__int_value:w
\fi:
- #2
+ #1 #5
+ \s_@@ \@@_chk:w #2 #3 ;
+ \s_@@ \@@_chk:w #4 #5
}
% \end{macrocode}
-% If the first \meta{case} is larger, then the first number remains
-% untouched, while the second number is ignored. On the other hand,
-% if the second \meta{case} is larger, the opposite happens: we retain
-% the second number. In both cases, there needs to be one step of
-% expansion after.
-% \begin{macrocode}
-% \end{macrocode}
-% We are then ready for the equality case: we split according
-% to the \meta{case}.
+% \end{macro}
+%
+% \begin{macro}[aux, EXP]{\@@_add_return_ii_o:Nww}
+% Ignore the first operand, and return the second, but using the sign
+% |#1| rather than |#4|. As usual, expand after the floating point.
% \begin{macrocode}
-\cs_new:Npn \@@_add_cases_eq:N #1
- {
- \if_case:w #1 \exp_stop_f:
- \exp_after:wN \@@_add_zeros:NNww
- \or: \exp_after:wN \@@_add_normal:NNww
- \or: \exp_after:wN \@@_add_inf:NNww
- \or: \exp_after:wN \@@_basics_return_nan_nan:NNww
- \fi:
- }
+\cs_new:Npn \@@_add_return_ii_o:Nww #1 #2 ; \s_@@ \@@_chk:w #3 #4
+ { \@@_exp_after_o:w \s_@@ \@@_chk:w #3 #1 }
% \end{macrocode}
-% Adding two zeros yields \cs{c_zero_fp}, except if both
-% zeros were $-0$.\footnote{Bruno: this should depend on the
-% rounding mode.}
+% \end{macro}
+%
+% \begin{macro}[aux, EXP]{\@@_add_zeros_o:Nww}
+% Adding two zeros yields \cs{c_zero_fp}, except if both zeros were
+% $-0$.
% \begin{macrocode}
-\cs_new:Npn \@@_add_zeros:NNww #1#2 #3;
+\cs_new:Npn \@@_add_zeros_o:Nww #1 \s_@@ \@@_chk:w 0 #2
{
- \if_int_compare:w #1 #2 = 02 \exp_stop_f:
- \@@_case_return_o:Nw \c_zero_fp
+ \if_int_compare:w #2 #1 = 20 \exp_stop_f:
+ \exp_after:wN \@@_add_return_ii_o:Nww
\else:
- \@@_case_return_same_o:w
+ \@@_case_return_i_o:ww
\fi:
+ #1
\s_@@ \@@_chk:w 0 #2
}
% \end{macrocode}
+% \end{macro}
+%
+% \begin{macro}[aux, EXP]{\@@_add_inf_o:Nww}
% If both infinities have the same sign, just return that infinity,
-% otherwise, it is an invalid operation.
+% otherwise, it is an invalid operation. We find out if that invalid
+% operation is an addition or a subtraction by testing whether the
+% tweaked \meta{sign_2} (|#1|) and the \meta{sign_2} (|#4|) are
+% identical.
% \begin{macrocode}
-\cs_new:Npn \@@_add_inf:NNww #1#2 #3;
+\cs_new:Npn \@@_add_inf_o:Nww
+ #1 \s_@@ \@@_chk:w 2 #2 #3; \s_@@ \@@_chk:w 2 #4
{
\if_meaning:w #1 #2
- \@@_case_return_same_o:w
+ \@@_case_return_i_o:ww
\else:
\@@_case_use:nw
{
- \@@_invalid_operation:Nnww \c_nan_fp { + }
- \s_@@ \@@_chk:w 2 #1 #3 ;
+ \if_meaning:w #1 #4
+ \exp_after:wN \@@_invalid_operation_o:Nww
+ \exp_after:wN +
+ \else:
+ \exp_after:wN \@@_invalid_operation_o:Nww
+ \exp_after:wN -
+ \fi:
}
\fi:
- \s_@@ \@@_chk:w 2 #2
+ \s_@@ \@@_chk:w 2 #2 #3;
+ \s_@@ \@@_chk:w 2 #4
}
% \end{macrocode}
% \end{macro}
%
-% \begin{macro}{\@@_add_normal:NNww}
-% \begin{syntax}
-% \cs{@@_add_normal:NNww} \meta{sign_1} \meta{sign_2}
-% ~~\Arg{exp_1} \meta{body_1} |;| \Arg{exp_2} \meta{body_2} |;|
-% \end{syntax}
+% \begin{macro}[aux, EXP]{\@@_add_normal_o:Nww}
+% \begin{quote}
+% \cs{@@_add_normal_o:Nww} \meta{sign_2}
+% \cs{s_@@} \cs{@@_chk:w} |1| \meta{sign_1}
+% \meta{exp_1} \meta{body_1} |;|
+% \cs{s_@@} \cs{@@_chk:w} |1| \meta{initial sign_2}
+% \meta{exp_2} \meta{body_2} |;|
+% \end{quote}
% We now have two normal numbers to add, and we have to check signs
% and exponents more carefully before performing the addition.
% \begin{macrocode}
-\cs_new:Npn \@@_add_normal:NNww #1#2
+\cs_new:Npn \@@_add_normal_o:Nww #1 \s_@@ \@@_chk:w 1 #2
{
\if_meaning:w #1#2
- \exp_after:wN \@@_add_npos:Nnwnw
+ \exp_after:wN \@@_add_npos_o:NnwNnw
\else:
- \exp_after:wN \@@_sub_npos:Nnwnw
+ \exp_after:wN \@@_sub_npos_o:NnwNnw
\fi:
- #1
+ #2
}
% \end{macrocode}
% \end{macro}
%
% \subsubsection{Absolute addition}
%
-% In this subsection, we perform the addition
-% of two positive normal numbers.
+% In this subsection, we perform the addition of two positive normal
+% numbers.
%
-% \begin{macro}[EXP]{\@@_add_npos:Nnwnw}
-% \begin{syntax}
-% \cs{@@_add_npos:Nnwnw} \meta{sign}
-% ~~\Arg{exp_1} \meta{body_1} |;| \Arg{exp_2} \meta{body_2} |;|
-% \end{syntax}
-% Since we are doing an addition, \meta{sign} will be the final sign.
-% The only special case which may arise is the case of an overflow.
-% This will be checked by \cs{@@_sanitize:Nw} at the end of
-% the calculation. We start an \cs{__int_eval:w}, responsible for
-% computing the exponent, which may receive a contribution of |+1|
-% in case of carry. The exponent should be stopped by |;| followed by
-% the overall \meta{sign} for the sanitizing to work properly.
-%
-% Grab and compare the exponents. The smaller number is decimated until
-% its exponent reaches that of the bigger number. We need to bring the
-% final sign down in the midst of the calculation to do the rounding
-% correctly.
+% \begin{macro}[aux, EXP]{\@@_add_npos_o:NnwNnw}
+% \begin{quote}
+% \cs{@@_add_npos_o:NnwNnw} \meta{sign_1} \meta{exp_1} \meta{body_1}
+% |;| \cs{s_@@} \cs{@@_chk:w} |1| \meta{initial sign_2} \meta{exp_2}
+% \meta{body_2} |;|
+% \end{quote}
+% Since we are doing an addition, the final sign is \meta{sign_1}.
+% Start an \cs{__int_eval:w}, responsible for computing the exponent:
+% the result, and the \meta{final sign} are then given to
+% \cs{@@_sanitize:Nw} which checks for overflow. The exponent is
+% computed as the largest exponent |#2| or |#5|, incremented if there
+% is a carry. To add the mantissas, we decimate the smaller number by
+% the difference between the exponents. This is done by
+% \cs{@@_add_big_i:wNww} or \cs{@@_add_big_ii:wNww}. We need to bring
+% the final sign with us in the midst of the calculation to round
+% properly at the end.
% \begin{macrocode}
-\cs_new:Npn \@@_add_npos:Nnwnw #1 #2#3; #4
+\cs_new:Npn \@@_add_npos_o:NnwNnw #1#2#3 ; \s_@@ \@@_chk:w 1 #4 #5
{
\exp_after:wN \@@_sanitize:Nw
\exp_after:wN #1
\int_use:N \__int_eval:w
- \if_int_compare:w #2 > #4 \exp_stop_f:
+ \if_int_compare:w #2 > #5 \exp_stop_f:
#2
- \exp_after:wN \@@_add_big_i:wNww \__int_value:w -
+ \exp_after:wN \@@_add_big_i_o:wNww \__int_value:w -
\else:
- #4
- \exp_after:wN \@@_add_big_ii:wNww \__int_value:w
+ #5
+ \exp_after:wN \@@_add_big_ii_o:wNww \__int_value:w
\fi:
- \__int_eval:w #4 - #2 ; #1 #3;
+ \__int_eval:w #5 - #2 ; #1 #3;
}
% \end{macrocode}
% \end{macro}
%
-% \begin{macro}{\@@_add_big_i:wNww,
-% \@@_add_big_ii:wNww}
-% \begin{syntax}
-% \cs{@@_add_big_i:wNww} \meta{shift} |;| \meta{sign}
-% ~~\meta{body_1} |;| \meta{body_2} |;|
-% \end{syntax}
-% Shift the mantissa of the small number, and then add with
-% \cs{@@_add_mantissa:NnnwnnnnN}.
+% \begin{macro}[int, rEXP]{\@@_add_big_i_o:wNww}
+% \begin{macro}[aux, rEXP]{\@@_add_big_ii_o:wNww}
+% \begin{quote}
+% \cs{@@_add_big_i_o:wNww} \meta{shift} |;| \meta{final sign}
+% \meta{body_1} |;| \meta{body_2} |;|
+% \end{quote}
+% Shift the mantissa of the small number, then add with
+% \cs{@@_add_mantissa_o:NnnwnnnnN}.
% \begin{macrocode}
-\cs_new:Npn \@@_add_big_i:wNww #1; #2 #3; #4;
+\cs_new:Npn \@@_add_big_i_o:wNww #1; #2 #3; #4;
{
\@@_decimate:nNnnnn {#1}
- \@@_add_mantissa:NnnwnnnnN
+ \@@_add_mantissa_o:NnnwnnnnN
#4
#3
#2
}
-\cs_new:Npn \@@_add_big_ii:wNww #1; #2 #3; #4;
+\cs_new:Npn \@@_add_big_ii_o:wNww #1; #2 #3; #4;
{
\@@_decimate:nNnnnn {#1}
- \@@_add_mantissa:NnnwnnnnN
+ \@@_add_mantissa_o:NnnwnnnnN
#3
#4
#2
}
% \end{macrocode}
% \end{macro}
+% \end{macro}
%
-% \begin{macro}{\@@_add_mantissa:NnnwnnnnN}
-% \begin{syntax}
-% \cs{@@_add_mantissa:NnnwnnnnN}
-% ~~\meta{rounding}
-% ~~\Arg{Y'_1} \Arg{Y'_2} \meta{extra-digits} |;|
-% ~~\Arg{X_1} \Arg{X_2} \Arg{X_3} \Arg{X_4}
-% ~~\meta{final sign}
-% \end{syntax}
+% \begin{macro}[aux, rEXP]{\@@_add_mantissa_o:NnnwnnnnN}
+% \begin{macro}[aux, rEXP]
+% {\@@_add_mantissa_pack:NNNNNNN, \@@_add_mantissa_test_o:N}
+% \begin{quote}
+% \cs{@@_add_mantissa_o:NnnwnnnnN}
+% \meta{rounding digit}
+% \Arg{Y'_1} \Arg{Y'_2} \meta{extra-digits} |;|
+% \Arg{X_1} \Arg{X_2} \Arg{X_3} \Arg{X_4}
+% \meta{final sign}
+% \end{quote}
% To round properly, we must know at which digit the rounding
% should occur. This requires to know whether the addition
% produces an overall carry or not. Thus, we do the computation
@@ -441,45 +403,43 @@
% give an exact power of $10$, for which it is easy to correct
% the result at the end.
% \begin{macrocode}
-\cs_new:Npn \@@_add_mantissa:NnnwnnnnN #1 #2#3 #4; #5#6#7#8
+\cs_new:Npn \@@_add_mantissa_o:NnnwnnnnN #1 #2#3 #4; #5#6#7#8
{
- \exp_after:wN \@@_add_mantissa_test:N
+ \exp_after:wN \@@_add_mantissa_test_o:N
\int_use:N \__int_eval:w 1#5#6 + #2
\exp_after:wN \@@_add_mantissa_pack:NNNNNNN
\int_use:N \__int_eval:w 1#7#8 + #3 ; #1
}
-% \end{macrocode}
-%
-% \begin{macrocode}
\cs_new:Npn \@@_add_mantissa_pack:NNNNNNN #1 #2#3#4#5#6#7
{
- \if:w 2 #1
+ \if_meaning:w 2 #1
+ \c_one
\fi:
; #2 #3 #4 #5 #6 #7 ;
}
-\cs_new:Npn \@@_add_mantissa_test:N #1
+\cs_new:Npn \@@_add_mantissa_test_o:N #1
{
- \if:w 2 #1
- \exp_after:wN \@@_add_mantissa_carry:wwNNNN
+ \if_meaning:w 2 #1
+ \exp_after:wN \@@_add_mantissa_carry_o:wwwNN
\else:
- \exp_after:wN \@@_add_mantissa_no_carry:wwNNNN
+ \exp_after:wN \@@_add_mantissa_no_carry_o:wwwNN
\fi:
}
% \end{macrocode}
+% \end{macro}
+% \end{macro}
%
+% \begin{macro}[aux, rEXP]{\@@_add_mantissa_no_carry_o:wwwNN}
% \begin{quote}
-% \cs{@@_add_mantissa_no_carry:wwNNNN}
-% ~~\meta{8d} |;| \meta{6d} |;| \meta{2d} |;|
-% ~~\meta{rounding} \meta{sign}
+% \cs{@@_add_mantissa_no_carry_o:wwwNN}
+% \meta{8d} |;| \meta{6d} |;| \meta{2d} |;|
+% \meta{rounding digit} \meta{sign}
% \end{quote}
-% If there's no carry, grab all the digits again, and just
-% set the rounding correctly.\footnote{Bruno: an optimization
-% would be to compute whether we need rounding or not,
-% and only grab digits if there is rounding.}
-%
+% If there's no carry, grab all the digits again and round. The
+% packing function \cs{@@_basics_pack_high:NNNNNw} takes care of the
+% case where rounding brings a carry.
% \begin{macrocode}
-\cs_new:Npn \@@_add_mantissa_no_carry:wwNNNN
+\cs_new:Npn \@@_add_mantissa_no_carry_o:wwwNN
#1; #2; #3#4 ; #5#6
{
\exp_after:wN \@@_basics_pack_high:NNNNNw
@@ -490,328 +450,407 @@
\exp_after:wN ;
}
% \end{macrocode}
+% \end{macro}
%
-% The case where there is a carry is very similar: rounding can even
-% raise the first digit from $1$ to $2$ (but we don't need to check that).
+% \begin{macro}[aux, rEXP]{\@@_add_mantissa_carry_o:wwwNN}
% \begin{quote}
-% \cs{@@_add_mantissa_carry:wwNNNN}
-% ~~\meta{8d} |;| \meta{6d} |;| \meta{2d} |;|
-% ~~\meta{rounding} \meta{sign}
+% \cs{@@_add_mantissa_carry_o:wwwNN}
+% \meta{8d} |;| \meta{6d} |;| \meta{2d} |;|
+% \meta{rounding digit} \meta{sign}
% \end{quote}
+% The case where there is a carry is very similar. Rounding can even
+% raise the first digit from $1$ to $2$, but we don't care.
% \begin{macrocode}
-\cs_new:Npn \@@_add_mantissa_carry:wwNNNN
+\cs_new:Npn \@@_add_mantissa_carry_o:wwwNN
#1; #2; #3#4; #5#6
{
+ \c_one
- \exp_after:wN \@@_add_mantissa_carry_pack:NNNNNNNNw
- \int_use:N \__int_eval:w 1 #1
- \exp_after:wN \@@_add_mantissa_carry_pack_ii:NNNNw
- \int_use:N \__int_eval:w 1 #2#3
- + \@@_round:NNNN #6 #3 #4 #5
+ \exp_after:wN \@@_basics_pack_weird_high:NNNNNNNNw
+ \int_use:N \__int_eval:w 1 1 #1
+ \exp_after:wN \@@_basics_pack_weird_low:NNNNw
+ \int_use:N \__int_eval:w 1 #2#3 +
+ \exp_after:wN \@@_round:NNN
+ \exp_after:wN #6
+ \exp_after:wN #3
+ \__int_value:w \@@_round_digit:Nw #4 #5 ;
\exp_after:wN ;
}
-\cs_new:Npn \@@_add_mantissa_carry_pack_ii:NNNNw #1 #2#3#4 #5;
- {
- \if:w 2 #1
- + \c_one
- \fi:
- \__int_eval_end:
- #2#3#4; {#5} ;
- }
-\cs_new:Npn \@@_add_mantissa_carry_pack:NNNNNNNNw
- #1#2#3#4 #5#6#7#8 #9; { ; {#1#2#3#4} {#5#6#7#8} {#9} }
% \end{macrocode}
% \end{macro}
%
% \subsubsection{Absolute subtraction}
%
-% \begin{macro}[EXP]{\@@_sub_npos:Nnwnw}
-% \begin{syntax}
-% \cs{@@_sub_npos:Nnwnw} \meta{sign}
-% ~~\Arg{exp_1} \meta{body_1} |;| \Arg{exp_2} \meta{body_2} |;|
-% \end{syntax}
-% Rounding properly in some modes requires to know what the sign
-% of the result will be. For addition, this was easy. Here, besides
-% comparing the exponents to know how to decimate, we need to
-% check carefully which number is bigger when they have the same
-% exponent.
+% \begin{macro}[aux, EXP]{\@@_sub_npos_o:NnwNnw}
+% \begin{macro}[aux, EXP]{\@@_sub_eq_o:Nnwnw, \@@_sub_npos_ii_o:Nnwnw}
+% \begin{quote}
+% \cs{@@_sub_npos_o:NnwNnw}
+% \meta{sign_1} \meta{exp_1} \meta{body_1} |;|
+% \cs{s_@@} \cs{@@_chk:w} |1|
+% \meta{initial sign_2} \meta{exp_2} \meta{body_2} |;|
+% \end{quote}
+% Rounding properly in some modes requires to know what the sign of
+% the result will be. Thus, we start by comparing the exponents and
+% mantissas. If the numbers coincide, return zero. If the second
+% number is larger, swap the numbers and call
+% \cs{@@_sub_npos_i_o:Nnwnw} with the opposite of \meta{sign_1}.
% \begin{macrocode}
-\cs_new:Npn \@@_sub_npos:Nnwnw #1 #2#3; #4 #5;
+\cs_new:Npn \@@_sub_npos_o:NnwNnw #1#2#3; \s_@@ \@@_chk:w 1 #4#5#6;
{
- \exp_after:wN \@@_sanitize:wN
- \int_use:N \__int_eval:w
- \if_int_compare:w #2 > #4 \exp_stop_f:
- #2
- \exp_after:wN \@@_sub_big_i:wNww \__int_value:w -
- \else:
- #4
- \if_int_compare:w #2 = #4 \exp_stop_f:
- \@@_sub_exponent_eq:nnnnnnnn #3 #5
- \else:
- \exp_after:wN \@@_sub_big_ii:wNww \__int_value:w
- \fi:
- \fi:
- \__int_eval:w #4 - #2 ; #1 #3; #5;
+ \if_case:w \@@_compare_npos:nwnw {#2} #3; {#5} #6; \exp_stop_f:
+ \exp_after:wN \@@_sub_eq_o:Nnwnw
+ \or:
+ \exp_after:wN \@@_sub_npos_i_o:Nnwnw
+ \else:
+ \exp_after:wN \@@_sub_npos_ii_o:Nnwnw
+ \fi:
+ #1 {#2} #3; {#5} #6;
+ }
+\cs_new:Npn \@@_sub_eq_o:Nnwnw #1#2; #3; { \exp_after:wN \c_zero_fp }
+\cs_new:Npn \@@_sub_npos_ii_o:Nnwnw #1 #2; #3;
+ {
+ \exp_after:wN \@@_sub_npos_i_o:Nnwnw
+ \int_use:N \__int_eval:w \c_two - #1 \__int_eval_end:
+ #3; #2;
}
% \end{macrocode}
% \end{macro}
+% \end{macro}
%
-% \begin{macro}{\@@_sub_exponent_eq:nnnnnnnn}
+% \begin{macro}[aux, EXP]{\@@_sub_npos_i_o:Nnwnw}
+% After the computation is done, \cs{@@_sanitize:Nw} checks for
+% overflow/underflow. It expects the \meta{final sign} and the
+% \meta{exponent} (delimited by |;|). Start an integer expression for
+% the exponent, which starts with the exponent of the largest number,
+% and may be decreased if the two numbers are very close. If the two
+% numbers have the same exponent, call the \texttt{near} auxiliary.
+% Otherwise, decimate $y$, then call the \texttt{far} auxiliary to
+% evaluate the difference between the two mantissas. Note that we
+% decimate by $1$ less than one could expect.
% \begin{macrocode}
-\cs_new:Npn \@@_sub_exponent_eq:nnnnnnnn #1#2#3#4 #5#6#7#8
+\cs_new:Npn \@@_sub_npos_i_o:Nnwnw #1 #2#3; #4#5;
{
- \if_int_compare:w #1#2 > #5#6 \exp_stop_f:
- \exp_after:wN \@@_sub_big_i:wNww \__int_value:w
- \else:
- \if_int_compare:w #1#2 < #5#6 \exp_stop_f:
- \exp_after:wN \@@_sub_big_ii:wNww \__int_value:w
+ \exp_after:wN \@@_sanitize:Nw
+ \exp_after:wN #1
+ \int_use:N \__int_eval:w
+ #2
+ \if_int_compare:w #2 = #4 \exp_stop_f:
+ \exp_after:wN \@@_sub_back_near_o:nnnnnnnnN
\else:
- \if_int_compare:w #3#4 > #7#8 \exp_stop_f:
- \exp_after:wN \@@_sub_big_i:wNww \__int_value:w
- \else:
- \if_int_compare:w #3#4 < #7#8 \exp_stop_f:
- \exp_after:wN \@@_sub_big_ii:wNww \__int_value:w
- \else:
- \exp_after:wN \@@_sub_eq:wNww \__int_value:w
- \fi:
- \fi:
+ \exp_after:wN \@@_decimate:nNnnnn \exp_after:wN
+ { \int_use:N \__int_eval:w #2 - #4 - \c_one \exp_after:wN }
+ \exp_after:wN \@@_sub_back_far_o:NnnwnnnnN
\fi:
- \fi:
+ #5
+ #3
+ #1
}
-\cs_new:Npn \@@_sub_eq:wNww #1; #2 #3; #4;
- { \exp_after:wN ; \exp_after:wN 1 \exp_after:wN ; }
% \end{macrocode}
% \end{macro}
%
-% \begin{macro}{\@@_sub_big_i:wNww,\@@_sub_big_ii:wNww}
-% \begin{syntax}
-% \cs{@@_sub_big_i:wNww} \meta{shift} |;| \meta{sign}
-% ~~\meta{body_1} |;| \meta{body_2} |;|
-% \end{syntax}
-% Shift the mantissa of the small number, and then subtract with
-% \cs{@@_sub_back_mantissa:NnnwNnnnn}.
+% \begin{macro}[aux, rEXP]{\@@_sub_back_near_o:nnnnnnnnN}
+% \begin{macro}[aux, rEXP]
+% {\@@_sub_back_near_pack:NNNNNNw, \@@_sub_back_near_after:wNNNNw}
+% \begin{quote}
+% \cs{@@_sub_back_near_o:nnnnnnnnN}
+% \Arg{Y_1} \Arg{Y_2} \Arg{Y_3} \Arg{Y_4}
+% \Arg{X_1} \Arg{X_2} \Arg{X_3} \Arg{X_4}
+% \meta{final sign}
+% \end{quote}
+% In this case, the subtraction is exact, so we discard the
+% \meta{final sign} |#9|. The very large shifts of $10^{9}$ and
+% $1.1\cdot10^{9}$ are unnecessary here, but allow the auxiliaries to
+% be reused later. Each integer expression produces a $10$ digit
+% result. If the resulting $16$ digits start with a $0$, then we need
+% to shift the group, padding with trailing zeros.
% \begin{macrocode}
-\cs_new:Npn \@@_sub_big_i:wNww #1; #2 #3; #4;
+\cs_new:Npn \@@_sub_back_near_o:nnnnnnnnN #1#2#3#4 #5#6#7#8 #9
{
- \@@_decimate:nNnnnn {#1}
- \@@_sub_back_mantissa:NnnwNnnnn
- #4
- #2
- #3
+ \exp_after:wN \@@_sub_back_near_after:wNNNNw
+ \int_use:N \__int_eval:w 10#5#6 - #1#2 - \c_eleven
+ \exp_after:wN \@@_sub_back_near_pack:NNNNNNw
+ \int_use:N \__int_eval:w 11#7#8 - #3#4 \exp_after:wN ;
}
-\cs_new:Npn \@@_sub_big_ii:wNww #1; #2 #3; #4;
+\cs_new:Npn \@@_sub_back_near_pack:NNNNNNw #1#2#3#4#5#6#7 ;
+ { + #1#2 ; {#3#4#5#6} {#7} ; }
+\cs_new:Npn \@@_sub_back_near_after:wNNNNw 10 #1#2#3#4 #5 ;
{
- \exp_after:wN \@@_sub_big_i:wNww
- \__int_value:w #1 \exp_after:wN ;
- \int_use:N \__int_eval:w 2 - #2 \__int_eval_end:
- #4; #3;
+ \if_meaning:w 0 #1
+ \exp_after:wN \@@_sub_back_shift:wnnnn
+ \fi:
+ ; {#1#2#3#4} {#5}
}
% \end{macrocode}
% \end{macro}
+% \end{macro}
%
-% \begin{macro}{\@@_sub_back_mantissa:NnnwNnnnn}
-% \begin{syntax}
-% \cs{@@_sub_back_mantissa:NnnwNnnnn}
-% ~~\meta{rounding} \Arg{Y'_1} \Arg{Y'_2} \meta{extra-digits}
-% ~~\meta{final sign}
-% ~~\Arg{X_1} \Arg{X_2} \Arg{X_3} \Arg{X_4}
-% \end{syntax}
-% At this stage, we know that \meta{Y} is less than \meta{X},
-% and we know the final sign.
+% \begin{macro}[aux, rEXP]{\@@_sub_back_shift:wnnnn}
+% \begin{macro}[aux, rEXP]
+% {
+% \@@_sub_back_shift_ii:ww,
+% \@@_sub_back_shift_iii:NNNNNNNNw,
+% \@@_sub_back_shift_iv:nnnnw
+% }
+% \begin{quote}
+% \cs{@@_sub_back_shift:wnnnn} |;|
+% \Arg{Z_1} \Arg{Z_2} \Arg{Z_3} \Arg{Z_4} |;|
+% \end{quote}
+% This function is called with $\meta{Z_1}\leq 999$. Act with
+% \tn{number} to trim leading zeros from \meta{Z_1} \meta{Z_2} (we
+% don't do all four blocks at once, since non-zero blocks would then
+% overflow \TeX{}'s integers). If the first two blocks are zero, the
+% auxiliary receives an empty |#1| and trims |#2#30| from leading
+% zeros, yielding a total shift between $7$ and~$16$ to the exponent.
+% Otherwise we get the shift from |#1| alone, yielding a result
+% between $1$ and~$6$. Once the exponent is taken care of, trim
+% leading zeros from |#1#2#3| (when |#1| is empty, the space before
+% |#2#3| is ignored), get four blocks of $4$~digits and finally clean
+% up. Trailing zeros are added so that digits can be grabbed safely.
% \begin{macrocode}
-\cs_new:Npn \@@_sub_back_mantissa:NnnwNnnnn #1 #2#3 #4; #5 #6#7#8#9
+\cs_new:Npn \@@_sub_back_shift:wnnnn ; #1#2
{
- \exp_after:wN \@@_sub_back_mantissa_i:NNwNNNNwN
- \exp_after:wN #1
- \exp_after:wN #5
- \int_use:N \__int_eval:w 2#6#7 - #2 - \c_two +
- \exp_after:wN \@@_sub_back_mantissa_round:wNN
- \int_use:N \__int_eval:w 2#8#9 - #3 ; #1 #5
+ \exp_after:wN \@@_sub_back_shift_ii:ww
+ \__int_value:w #1 #2 0 ;
}
-% \end{macrocode}
-% After the computation, we need to check whether the first digit of
-% the result is zero. This can only happen if the two numbers had the
-% same exponent, or exponents differing by $1$. In the latter case,
-% the \meta{rounding} digit is not quite enough to let us retrieve
-% the exact result (consider $\cdots25$ and $\cdots15$, both rounded
-% to $\cdots2$ in the usual mode), so we also move the result of
-% \cs{@@_round_neg:NNN} upstream as the digit $0$ or $1$.
-% \begin{macrocode}
-\cs_new:Npn \@@_sub_back_mantissa_round:wNN #1; #2 #3
+\cs_new:Npn \@@_sub_back_shift_ii:ww #1 0 ; #2#3 ;
{
- \exp_after:wN \@@_sub_back_mantissa_iii:N
- \__int_value:w
- \exp_after:wN \@@_round_neg:NNN
- \exp_after:wN #3
- \use_none:nnnnnnnn #1 #2
- + #1
+ \if_meaning:w @ #1 @
+ - \c_seven
+ - \exp_after:wN \use_i:nnn
+ \exp_after:wN \@@_sub_back_shift_iii:NNNNNNNNw
+ \__int_value:w #2#3 0 ~ 123456789;
+ \else:
+ - \@@_sub_back_shift_iii:NNNNNNNNw #1 123456789;
+ \fi:
+ \exp_after:wN \@@_pack_twice_four:wNNNNNNNN
+ \exp_after:wN \@@_pack_twice_four:wNNNNNNNN
+ \exp_after:wN \@@_sub_back_shift_iv:nnnnw
\exp_after:wN ;
+ \__int_value:w
+ #1 ~ #2#3 0 ~ 0000 0000 0000 000 ;
}
-\cs_new:Npn \@@_sub_back_mantissa_iii:N #1
- {
- \exp_after:wN \@@_sub_back_mantissa_ii:NNNNNNw
- \exp_after:wN #1
- \int_use:N \__int_eval:w
- - #1
- }
-\cs_new:Npn \@@_sub_back_mantissa_ii:NNNNNNw #1 #2 #3#4#5#6 #7;
- { #2 ; #1 {#3#4#5#6} {#7} ; }
+\cs_new:Npn \@@_sub_back_shift_iii:NNNNNNNNw #1#2#3#4#5#6#7#8#9; {#8}
+\cs_new:Npn \@@_sub_back_shift_iv:nnnnw #1 ; #2 ; { ; #1 ; }
% \end{macrocode}
% \end{macro}
+% \end{macro}
%
-% \begin{macro}{\@@_sub_back_mantissa_i:NNwNNNNwN}
-% Here, |#3| should always be $2$, but we have to take is
-% as a normal undelimited argument since that would break
-% if |#2| is $2$.
+% \begin{macro}[aux, rEXP]{\@@_sub_back_far_o:NnnwnnnnN}
+% \begin{quote}
+% \cs{@@_sub_back_far_o:NnnwnnnnN}
+% \meta{rounding} \Arg{Y'_1} \Arg{Y'_2} \meta{extra-digits} |;|
+% \Arg{X_1} \Arg{X_2} \Arg{X_3} \Arg{X_4}
+% \meta{final sign}
+% \end{quote}
+% If the difference is greater than $10^{\meta{expo_x}}$, call the
+% \texttt{very_far} auxiliary. If the result is less than
+% $10^{\meta{expo_x}}$, call the \texttt{not_far} auxiliary. If it is
+% too close a call to know yet, namely if $1 \meta{Y'_1} \meta{Y'_2} =
+% \meta{X_1} \meta{X_2} \meta{X_3} \meta{X_4} 0$, then call the
+% \texttt{quite_far} auxiliary. We use the odd combination of space
+% and semi-colon delimiters to allow the \texttt{not_far} auxiliary to
+% grab each piece individually, the \texttt{very_far} auxiliary to use
+% \cs{@@_pack_eight:wNNNNNNNN}, and the \texttt{quite_far} to ignore
+% the mantissas easily (using the |;| delimiter).
% \begin{macrocode}
-\cs_new:Npn \@@_sub_back_mantissa_i:NNwNNNNwN #1#2 #3 #4#5#6#7 #8; #9
+\cs_new:Npn \@@_sub_back_far_o:NnnwnnnnN #1 #2#3 #4; #5#6#7#8
{
- \if:w 0 #4
- \exp_after:wN \@@_sub_back_carry:NNwNnnnn
- \exp_after:wN #1
- \exp_after:wN #9
+ \if_case:w
+ \if_int_compare:w 1 #2 = #5#6 \use_i:nnnn #7 \exp_stop_f:
+ \if_int_compare:w #3 = \use_none:n #7#8 0 \exp_stop_f:
+ \c_zero
+ \else:
+ \if_int_compare:w #3 > \use_none:n #7#8 0 - \fi: \c_one
+ \fi:
+ \else:
+ \if_int_compare:w 1 #2 > #5#6 \use_i:nnnn #7 - \fi: \c_one
+ \fi:
+ \exp_after:wN \@@_sub_back_quite_far_o:wwNN
+ \or: \exp_after:wN \@@_sub_back_very_far_o:wwwwNN
+ \else: \exp_after:wN \@@_sub_back_not_far_o:wwwwNN
\fi:
- ; #2
- {#4#5#6#7} {#8}
+ #2 ~ #3 ; #5 #6 ~ #7 #8 ; #1
}
% \end{macrocode}
% \end{macro}
%
-% \begin{macro}{\@@_sub_back_carry:NNwNnnnn}
-% \begin{syntax}
-% \cs{@@_sub_back_carry:NNwNnnnn}
-% ~~\meta{rounding} \meta{0 or 1} |;| \meta{final sign}
-% ~~\Arg{Z_1} \Arg{Z_2} \Arg{Z_3} \Arg{Z_4} |;|
-% \end{syntax}
-% This function is called when $\meta{Z_1}\leq 999$. We revert
-% the carry, which is given by \meta{0 or 1}, and subtract the
-% \meta{rounding} digit as appropriate, then feed the result,
-% of the form \meta{$\leq$ 7d} |;| \meta{9d} |;| to
-% \cs{@@_sub_back_carry_i:wwN}. The result is always exact.
+% \begin{macro}[aux, EXP]{\@@_sub_back_quite_far_o:wwNN}
+% \begin{macro}[aux, EXP]{\@@_sub_back_quite_far_ii:NN}
+% The easiest case is when $x-y$ is extremely close to a power of
+% $10$, namely the first digit of $x$ is $1$, and all others vanish
+% when subtracting $y$. Then the \meta{rounding} |#3| and the
+% \meta{final sign} |#4| control whether we get $1$ or $0.9999 9999
+% 9999 9999$. In the usual round-to-nearest mode, we will get $1$
+% whenever the \meta{rounding} digit is less than or equal to $5$
+% (remember that the \meta{rounding} digit is only equal to $5$ if
+% there was no further non-zero digit).
% \begin{macrocode}
-\cs_new:Npn \@@_sub_back_carry:NNwNnnnn #1#2 ; #3 #4#5#6#7 ;
+\cs_new:Npn \@@_sub_back_quite_far_o:wwNN #1; #2; #3#4
{
- \exp_after:wN \@@_sub_back_carry_i:wwN
- \int_use:N \__int_eval:w #4 #5 - 1 + \exp_after:wN \@@_use_s:n
- \int_use:N \__int_eval:w 1 #6 #7 0 + #2 0 - #1 ; #3
+ \exp_after:wN \@@_sub_back_quite_far_ii:NN
+ \exp_after:wN #3
+ \exp_after:wN #4
}
-% \end{macrocode}
-% Unless the first block is zero, check how many digits is has,
-% and shift the exponent down by the corresponding amount. Then
-% pack digits into blocks of $4$ (there are between $10$ and $16$
-% digits in front of \cs{@@_sub_back_carry_large:NNNNNNNNw}).
-% \begin{macrocode}
-\cs_new:Npn \@@_sub_back_carry_i:wwN #1 ;
+\cs_new:Npn \@@_sub_back_quite_far_ii:NN #1#2
{
- \if:w 0 #1
- - 8
- \exp_after:wN \@@_sub_back_carry_small:wN \__int_value:w
+ \if_case:w \@@_round_neg:NNN #2 0 #1
+ \exp_after:wN \use_i:nn
\else:
- - \@@_sub_back_carry_ii:NNNNNNNNw #1 1234567;
- \exp_after:wN \@@_sub_back_carry_large:NNNNNNNNw
+ \exp_after:wN \use_ii:nn
\fi:
- #1
+ { ; {1000} {0000} {0000} {0000} ; }
+ { - \c_one ; {9999} {9999} {9999} {9999} ; }
}
% \end{macrocode}
-% The case where the number is non-zero is slightly easier.
+% \end{macro}
+% \end{macro}
+%
+% \begin{macro}[aux, rEXP]{\@@_sub_back_not_far_o:wwwwNN}
+% In the present case, $x$ and $y$ have different exponents, but
+% $y$~is large enough that $x-y$ has a smaller exponent than~$x$.
+% Decrement the exponent (with |- \c_one|). Then proceed in a way
+% similar to the \texttt{near} auxiliaries seen earlier, but
+% multiplying $x$ by~$10$ (|#30| and |#40| below), and with the added
+% quirk that the \meta{rounding} digit has to be taken into account.
+% Namely, we may have to decrease the result by one unit if
+% \cs{@@_round_neg:NNN} returns~$1$. This function expects the
+% \meta{final sign}~|#6|, the last digit of |1100000000+#40-#2|, and
+% the \meta{rounding} digit. Instead of redoing the computation for
+% the second argument, we note that \cs{@@_round_neg:NNN} only cares
+% about its parity, which is identical to that of the last digit
+% of~|#2|.
% \begin{macrocode}
-\cs_new:Npn \@@_sub_back_carry_ii:NNNNNNNNw #1#2#3#4#5#6#7#8#9; {#8}
-\cs_new:Npn \@@_sub_back_carry_large:NNNNNNNNw #1#2#3#4 #5#6#7#8 #9;
+\cs_new:Npn \@@_sub_back_not_far_o:wwwwNN #1 ~ #2; #3 ~ #4; #5#6
{
- \@@_sub_back_carry_large_ii:NNNNNNNNw
- #9 000000 ; {#1#2#3#4} {#5#6#7#8}
+ - \c_one
+ \exp_after:wN \@@_sub_back_near_after:wNNNNw
+ \int_use:N \__int_eval:w 1#30 - #1 - \c_eleven
+ \exp_after:wN \@@_sub_back_near_pack:NNNNNNw
+ \int_use:N \__int_eval:w 11 0000 0000 + #40 - #2
+ - \exp_after:wN \@@_round_neg:NNN
+ \exp_after:wN #6
+ \use_none:nnnnnnn #2 #5
+ \exp_after:wN ;
}
-\cs_new:Npn \@@_sub_back_carry_large_ii:NNNNNNNNw #1#2#3#4 #5#6#7#8 #9;
- { \@@_sub_back_carry_large_iii:nnnnN {#1#2#3#4} {#5#6#7#8} }
-\cs_new:Npn \@@_sub_back_carry_large_iii:nnnnN #1#2 #3#4 #5
- { ; #5 {#3}{#4} {#1}{#2} ; }
% \end{macrocode}
-% In the case of a \enquote{small} result, what comes after
-% \cs{@@_sub_back_carry_small:wN} has between $1$
-% and $9$ digits, and is not zero.
+% \end{macro}
+%
+% \begin{macro}[aux, EXP]{\@@_sub_back_very_far_o:wwwwNN}
+% \begin{macro}[aux, EXP]{\@@_sub_back_very_far_ii_o:nnNwwNN}
+% The case where $x-y$ and $x$ have the same exponent is a bit more
+% tricky, mostly because it cannot reuse the same auxiliaries. Shift
+% the $y$~mantissa by adding a leading~$0$. Then the logic is similar
+% to the \texttt{not_far} functions above. Rounding is a bit more
+% complicated: we have two \meta{rounding} digits |#3| and |#6| (from
+% the decimation, and from the new shift) to take into account, and
+% getting the parity of the main result requires a computation. The
+% first \cs{__int_value:w} triggers the second one because the number
+% is unfinished; we can thus not use $0$ in place of $2$ there.
% \begin{macrocode}
-\cs_new:Npn \@@_sub_back_carry_small:wN #1;
+\cs_new:Npn \@@_sub_back_very_far_o:wwwwNN #1#2#3#4#5#6#7
+ {
+ \@@_pack_eight:wNNNNNNNN
+ \@@_sub_back_very_far_ii_o:nnNwwNN
+ { 0 #1#2#3 #4#5#6#7 }
+ ;
+ }
+\cs_new:Npn \@@_sub_back_very_far_ii_o:nnNwwNN #1#2 ; #3 ; #4 ~ #5; #6#7
{
- - \exp_after:wN \@@_use_i_until_s:nw
- \use_none:nnnnnnnnn #1 012345678;
- \@@_sub_back_carry_small_ii:NNNNNNNN #1 00000000 ;
+ \exp_after:wN \@@_basics_pack_high:NNNNNw
+ \int_use:N \__int_eval:w 1#4 - #1 - \c_one
+ \exp_after:wN \@@_basics_pack_low:NNNNNw
+ \int_use:N \__int_eval:w 2#5 - #2
+ - \exp_after:wN \@@_round_neg:NNN
+ \exp_after:wN #7
+ \__int_value:w
+ \if_int_odd:w \__int_eval:w #5 - #2 \__int_eval_end:
+ 1 \else: 2 \fi:
+ \__int_value:w \@@_round_digit:Nw #3 #6 ;
+ \exp_after:wN ;
}
-\cs_new:Npn \@@_sub_back_carry_small_ii:NNNNNNNN #1#2#3#4 #5#6#7#8
- { \@@_sub_back_carry_small_iii:nnNwN {#1#2#3#4} {#5#6#7#8} }
-\cs_new:Npn \@@_sub_back_carry_small_iii:nnNwN #1 #2 #3 #4; #5
- { ; #5 {#1} {#2} {#3000} {0000} ; }
% \end{macrocode}
% \end{macro}
-%
+% \end{macro}
%
% \subsection{Multiplication}
%
-% \begin{macro}[int, EXP]{\@@_*_o:ww}
-% For multiplication, everything is easy. No need to grab the
-% \meta{body_2}.
-% \begin{macrocode}
-\cs_new:cpn { @@_*_o:ww }
- \s_@@ \@@_chk:w #1 #2 #3 ; \s_@@ \@@_chk:w #4 #5
- { \@@_mul_cases:NN #1 #4 #2 #5 #3 ; }
-% \end{macrocode}
-% \end{macro}
-%
% \subsubsection{Signs, and special numbers}
%
-% \begin{macro}[EXP,aux]{\@@_mul_cases:NN}
-% \begin{syntax}
-% \cs{@@_mul_cases:NN} \meta{case_1} \meta{case_2}
-% ~~\meta{sign_1} \meta{sign_2} \meta{body_1} |;| \meta{body_2} |;|
-% \end{syntax}
-% Expands the following tokens on the input stream once.
-% The special cases are coded at the start of this module,
-% and identical to the ones for division.\footnote{Bruno: \texttt{nan}
-% are not treated properly: $\infty\times 0$ should signal.}
-%^^A todo: use the faster \if_meaning:w.
+% \begin{macro}[int, EXP]{\@@_*_o:ww}
+% We go through an auxiliary, which is common with \cs{@@_/_o:ww}.
+% The first argument is the operation, used for the invalid operation
+% exception. The second is inserted in a formula to dispatch cases
+% slightly differently between multiplication and division. The third
+% is the operation for normal floating points. The fourth is there
+% for extra cases needed in \cs{@@_/_o:ww}.
% \begin{macrocode}
-\cs_new:Npn \@@_mul_cases:NN #1 #2
+\cs_new_nopar:cpn { @@_*_o:ww }
{
- \if_case:w \if_meaning:w 1 #1 #2 \else:
- \if_meaning:w 1 #2 #1 \else:
- \if_meaning:w #1#2 #1 \else:
- \if_int_compare:w \__int_eval:w #1 + #2 > \c_two
- 3 \else: 4 \fi: \fi: \fi: \fi:
- \exp_stop_f:
- \exp_after:wN \@@_basics_return_zero:NNww
- \or: \exp_after:wN \@@_mul_normal:NNww
- \or: \exp_after:wN \@@_basics_return_inf:NNww
- \or: \exp_after:wN \@@_basics_return_nan:NNww
- \or:
- \exp_after:wN \@@_mul_invalid:NNNNww
- \exp_after:wN #1
- \exp_after:wN #2
- \fi:
- }
-\cs_new:Npn \@@_mul_invalid:NNNNww #1#2#3#4#5; #6;
- {
- \@@_invalid_operation:Nnww \c_nan_fp { * }
- \s_@@ \@@_chk:w #1 #3 #5 ;
- \s_@@ \@@_chk:w #2 #4 #6 ;
+ \@@_mul_cases_o:NnNnww
+ *
+ { - \c_two + }
+ \@@_mul_npos_o:Nww
+ { }
}
% \end{macrocode}
% \end{macro}
%
-% \begin{macro}{\@@_mul_normal:NNww}
-% \begin{syntax}
-% \cs{@@_mul_normal:NNww} \meta{sign_1} \meta{sign_2}
-% ~~\Arg{exp_1} \meta{body_1} |;| \Arg{exp_2} \meta{body_2} |;|
-% \end{syntax}
-% We now have two normal numbers to multiply. Combine the signs.
+% \begin{macro}[int, EXP]{\@@_mul_cases_o:nNnnww}
+% Split into $10$ cases ($12$ for division).
+% If both numbers are normal, go to case $0$
+% (same sign) or case $1$ (opposite signs): in both cases, call
+% \cs{@@_mul_npos_o:Nww} to do the work. If the first operand is
+% \texttt{nan}, go to case $2$, in which the second operand is
+% discarded; if the second operand is \texttt{nan}, go to case $3$, in
+% which the first operand is discarded (note the weird interaction
+% with the final test on signs). Then we separate the case where the
+% first number is normal and the second is zero: this goes to cases
+% $4$ and $5$ for multiplication, $10$ and $11$ for division.
+% Otherwise, we do a computation which
+% dispatches the products $0\times 0 = 0\times 1 = 1\times 0 = 0$ to
+% case $4$ or $5$ depending on the combined sign, the products
+% $0\times\infty$ and $\infty\times0$ to case $6$ or $7$ (invalid
+% operation), and the products $1\times\infty = \infty\times1 =
+% \infty\times\infty = \infty$ to cases $8$ and $9$. Note that the
+% code for these two cases (which return $\pm\infty$) is inserted as
+% argument |#4|, because it differs in the case of divisions.
% \begin{macrocode}
-\cs_new:Npn \@@_mul_normal:NNww #1#2
+\cs_new:Npn \@@_mul_cases_o:NnNnww
+ #1#2#3#4 \s_@@ \@@_chk:w #5#6#7; \s_@@ \@@_chk:w #8#9
{
- \if:w #1#2
- \exp_after:wN \@@_mul_npos:Nnwnw
- \exp_after:wN 0
- \else:
- \exp_after:wN \@@_mul_npos:Nnwnw
- \exp_after:wN 2
+ \if_case:w \__int_eval:w
+ \if_int_compare:w #5 #8 = \c_eleven
+ \c_one
+ \else:
+ \if_meaning:w 3 #8
+ \c_three
+ \else:
+ \if_meaning:w 3 #5
+ \c_two
+ \else:
+ \if_int_compare:w #5 #8 = \c_ten
+ \c_nine #2 - \c_two
+ \else:
+ (#5 #2 #8) / \c_two * \c_two + \c_seven
+ \fi:
+ \fi:
+ \fi:
+ \fi:
+ \if_meaning:w #6 #9 - \c_one \fi:
+ \__int_eval_end:
+ \@@_case_use:nw { #3 0 }
+ \or: \@@_case_use:nw { #3 2 }
+ \or: \@@_case_return_i_o:ww
+ \or: \@@_case_return_ii_o:ww
+ \or: \@@_case_return_o:Nww \c_zero_fp
+ \or: \@@_case_return_o:Nww \c_minus_zero_fp
+ \or: \@@_case_use:nw { \@@_invalid_operation_o:Nww #1 }
+ \or: \@@_case_use:nw { \@@_invalid_operation_o:Nww #1 }
+ \or: \@@_case_return_o:Nww \c_inf_fp
+ \or: \@@_case_return_o:Nww \c_minus_inf_fp
+ #4
\fi:
+ \s_@@ \@@_chk:w #5 #6 #7;
+ \s_@@ \@@_chk:w #8 #9
}
% \end{macrocode}
% \end{macro}
@@ -821,39 +860,44 @@
% In this subsection, we perform the multiplication
% of two positive normal numbers.
%
-% \begin{macro}[EXP]{\@@_mul_npos:Nnwnw}
-% \begin{syntax}
-% \cs{@@_mul_npos:Nnwnw} \meta{sign}
-% ~~\Arg{exp_1} \meta{body_1} |;| \Arg{exp_2} \meta{body_2} |;|
-% \end{syntax}
-% As for addition, \meta{sign} is the final sign. After the computation,
-% \cs{@@_sanitize:Nw} checks for overflow or underflow.
-% As before, \cs{__int_eval:w} computes the exponent, catching any
-% shift coming from the computation in the mantissa. Again, the
-% \meta{sign} is needed for rounding to be done properly, so we move
-% it around with us. We setup the post-expansion here, triggered by
-% \cs{@@_mul_mantissa:nnnnNnnnn}.
+% \begin{macro}[int, EXP]{\@@_mul_npos_o:Nww}
+% \begin{quote}
+% \cs{@@_mul_npos_o:Nww} \meta{final sign}
+% \cs{s_@@} \cs{@@_chk:w} |1| \meta{sign_1} \Arg{exp_1} \meta{body_1} |;|
+% \cs{s_@@} \cs{@@_chk:w} |1| \meta{sign_2} \Arg{exp_2} \meta{body_2} |;|
+% \end{quote}
+% After the computation, \cs{@@_sanitize:Nw} checks for overflow or
+% underflow. As we did for addition, \cs{__int_eval:w} computes the
+% exponent, catching any shift coming from the computation in the
+% mantissa. The \meta{final sign} is needed to do the rounding
+% properly in the mantissa computation. We setup the post-expansion
+% here, triggered by \cs{@@_mul_mantissa_o:nnnnNnnnn}.
% \begin{macrocode}
-\cs_new:Npn \@@_mul_npos:Nnwnw #1 #2#3; #4 #5;
+\cs_new:Npn \@@_mul_npos_o:Nww
+ #1 \s_@@ \@@_chk:w #2 #3 #4 #5 ; \s_@@ \@@_chk:w #6 #7 #8 #9 ;
{
\exp_after:wN \@@_sanitize:Nw
\exp_after:wN #1
\int_use:N \__int_eval:w
- #2 + #4
- \@@_mul_mantissa:nnnnNnnnn #3 #1 #5
- \exp_after:wN ;
+ #4 + #8
+ \@@_mul_mantissa_o:nnnnNnnnn #5 #1 #9
}
% \end{macrocode}
% \end{macro}
%
-% \begin{macro}{\@@_mul_mantissa:nnnnNnnnn}
-% \begin{syntax}
-% \cs{@@_mul_mantissa:nnnnNnnnn}
-% ~~\Arg{X_1} \Arg{X_2} \Arg{X_3} \Arg{X_4} \meta{sign}
-% ~~\Arg{Y_1} \Arg{Y_2} \Arg{Y_3} \Arg{Y_4} \meta{;}
-% \end{syntax}
-% After one expansion, the token following \meta{Y_4} must be a semicolon
-% (represented by \meta{;}).
+% \begin{macro}[aux, rEXP]{\@@_mul_mantissa_o:nnnnNnnnn}
+% \begin{macro}[aux, EXP]
+% {\@@_mul_mantissa_drop:NNNNNw, \@@_mul_mantissa_keep:NNNNNw}
+% \begin{quote}
+% \cs{@@_mul_mantissa_o:nnnnNnnnn}
+% \Arg{X_1} \Arg{X_2} \Arg{X_3} \Arg{X_4} \meta{sign}
+% \Arg{Y_1} \Arg{Y_2} \Arg{Y_3} \Arg{Y_4}
+% \end{quote}
+% Note the three semicolons at the end of the definition. One is for
+% the last \cs{@@_mul_mantissa_drop:NNNNNw}; one is for
+% \cs{@@_round_digit:Nw} later on; and one, preceeded by
+% \cs{exp_after:wN}, which is correctly expanded (within an
+% \cs{__int_eval:w}), is used by \cs{@@_basics_pack_low:NNNNNw}.
%
% The product of two $16$ digit integers has $31$ or $32$ digits,
% but it is impossible to know which one before computing. The place
@@ -865,9 +909,9 @@
% known, and we can do the rounding within yet another set of
% \cs{__int_eval:w}.
% \begin{macrocode}
-\cs_new:Npn \@@_mul_mantissa:nnnnNnnnn #1#2#3#4 #5 #6#7#8#9
+\cs_new:Npn \@@_mul_mantissa_o:nnnnNnnnn #1#2#3#4 #5 #6#7#8#9
{
- \exp_after:wN \@@_mul_mantissa_after:NNN
+ \exp_after:wN \@@_mul_mantissa_test_f:NNN
\exp_after:wN #5
\int_use:N \__int_eval:w 99990000 + #1*#6 +
\exp_after:wN \@@_mul_mantissa_keep:NNNNNw
@@ -881,688 +925,641 @@
\exp_after:wN \@@_mul_mantissa_drop:NNNNNw
\int_use:N \__int_eval:w 99990000 + #3*#9 + #4*#8 +
\exp_after:wN \@@_mul_mantissa_drop:NNNNNw
- \int_use:N \__int_eval:w 100000000 + #4*#9 \exp_after:wN ;
+ \int_use:N \__int_eval:w 100000000 + #4*#9 ;
+ ; \exp_after:wN ;
}
\cs_new:Npn \@@_mul_mantissa_drop:NNNNNw #1#2#3#4#5 #6;
{ #1#2#3#4#5 ; + #6 }
\cs_new:Npn \@@_mul_mantissa_keep:NNNNNw #1#2#3#4#5 #6;
{ #1#2#3#4#5 ; #6 ; }
% \end{macrocode}
-% Once the first \cs{int_use:N} \cs{__int_eval:w}, and all the
-% \cs{@@_mul_mantissa_...:NNNNNw} have been expanded,
-% we get
+% \end{macro}
+% \end{macro}
+%
+% \begin{macro}[aux, rEXP]{\@@_mul_mantissa_test_f:NNN}
% \begin{quote}
-% \cs{@@_mul_mantissa_after:NNN} \meta{sign} |1|
-% ~~\meta{digits 1--8} |;| \meta{digits 9--12} |;| \meta{digits 13--16} |;|
-% ~~|+| \meta{digits 17--20} |+| \meta{digits 21--24}
-% ~~|+| \meta{digits 25--28} |+| \meta{digits 29--32} |;|
+% \cs{@@_mul_mantissa_test_f:NNN} \meta{sign} |1|
+% \meta{digits 1--8} |;| \meta{digits 9--12} |;| \meta{digits 13--16} |;|
+% |+| \meta{digits 17--20} |+| \meta{digits 21--24}
+% |+| \meta{digits 25--28} |+| \meta{digits 29--32} |;|
+% \cs{exp_after:wN} |;|
% \end{quote}
% If the \meta{digit 1} is non-zero, then for rounding we only care
-% about the digits $16$ and $17$, and whether all other digits are zero
+% about the digits $16$ and $17$, and whether further digits are zero
% or not (check for exact ties). On the other hand, if \meta{digit 1}
-% is zero, we care about digits $17$ and $18$, and whether all others are
-% zero.
+% is zero, we care about digits $17$ and $18$, and whether further
+% digits are zero.
% \begin{macrocode}
-\cs_new:Npn \@@_mul_mantissa_after:NNN #1 #2 #3
+\cs_new:Npn \@@_mul_mantissa_test_f:NNN #1 #2 #3
{
- \if:w 0 #3
- \exp_after:wN \@@_mul_mantissa_small:NNwwwN
+ \if_meaning:w 0 #3
+ \exp_after:wN \@@_mul_mantissa_small_f:NNwwwN
\else:
- \exp_after:wN \@@_mul_mantissa_large:NwwNNNN
+ \exp_after:wN \@@_mul_mantissa_large_f:NwwNNNN
\fi:
#1 #3
}
% \end{macrocode}
% \end{macro}
%
-% \begin{macro}{\@@_mul_mantissa_large:NwwNNNN}
+% \begin{macro}[aux, EXP]{\@@_mul_mantissa_large_f:NwwNNNN}
% In this branch, \meta{digit 1} is non-zero. The result is thus
% \meta{digits 1--16}, plus some rounding which depends on the digits
% $16$, $17$, and whether all subsequent digits are zero or not.
-% Here, \cs{@@_round_s:NNNw} takes the \meta{sign}, followed by
-% digits $16$, $17$, and an integer expression which is zero if and
-% only if all further digits are zero.
+% Here, \cs{@@_round_digit:Nw} takes digits $17$ and further (as an
+% integer expression), and replaces it by a \meta{rounding digit},
+% suitable for \cs{@@_round:NNN}.
% \begin{macrocode}
-\cs_new:Npn \@@_mul_mantissa_large:NwwNNNN #1 #2; #3; #4#5#6#7; +
+\cs_new:Npn \@@_mul_mantissa_large_f:NwwNNNN #1 #2; #3; #4#5#6#7; +
{
\exp_after:wN \@@_basics_pack_high:NNNNNw
\int_use:N \__int_eval:w 1#2
\exp_after:wN \@@_basics_pack_low:NNNNNw
- \int_use:N \__int_eval:w 1#3#4#5#6#7 + \@@_round_s:NNNw #1 #7
+ \int_use:N \__int_eval:w 1#3#4#5#6#7
+ + \exp_after:wN \@@_round:NNN
+ \exp_after:wN #1
+ \exp_after:wN #7
+ \__int_value:w \@@_round_digit:Nw
}
% \end{macrocode}
% \end{macro}
%
-% \begin{macro}{\@@_mul_mantissa_small:NNwwwN}
+% \begin{macro}[aux, rEXP]{\@@_mul_mantissa_small_f:NNwwwN}
% In this branch, \meta{digit 1} is zero. Our result will thus be
% \meta{digits 2--17}, plus some rounding which depends on the digits
% $17$, $18$, and whether all subsequent digits are zero or not.
% The $8$ digits |1#3| are followed, after expansion of the
% \texttt{small_pack} auxiliary, by the next digit, to form a $9$
-% digit number. Also, rounding may have caused a carry, which is
-% then converted to \cs{c_ten} rather than the usual \cs{c_one},
-% because of the shift.
+% digit number.
% \begin{macrocode}
-\cs_new:Npn \@@_mul_mantissa_small:NNwwwN #1 #2#3; #4; #5; + #6
+\cs_new:Npn \@@_mul_mantissa_small_f:NNwwwN #1 #2#3; #4#5; #6; + #7
{
- \c_one
\exp_after:wN \@@_basics_pack_high:NNNNNw
- \int_use:N \__int_eval:w 1#3
- \exp_after:wN \@@_mul_mantissa_small_pack:NNNNNNw
- \int_use:N \__int_eval:w 1#4#5#6 + \@@_round_s:NNNw #1 #6
- }
-\cs_new:Npn \@@_mul_mantissa_small_pack:NNNNNNw #1#2 #3#4#5#6 #7;
- {
- #2
- \if:w 2 #1
- + \c_ten
- \fi:
- ; {#3#4#5#6} {#7} ;
+ \int_use:N \__int_eval:w 1#3#4
+ \exp_after:wN \@@_basics_pack_low:NNNNNw
+ \int_use:N \__int_eval:w 1#5#6#7
+ + \exp_after:wN \@@_round:NNN
+ \exp_after:wN #1
+ \exp_after:wN #7
+ \__int_value:w \@@_round_digit:Nw
}
% \end{macrocode}
% \end{macro}
%
% \subsection{Division}
%
-% Time is now ripe to tackle the hardest of the four elementary
-% operations: division.
-%
-% \begin{macro}[EXP]{\@@_/_o:ww}
-% For division we swap the two floating point numbers.
-% \begin{macrocode}
-\cs_new:cpn { @@_/_o:ww }
- \s_@@ \@@_chk:w #1 #2 #3 ; \s_@@ \@@_chk:w #4 #5 #6 ;
- { \@@_div_back_cases:NN #4 #1 #5 #2 #6 ; #3 ; }
-% \end{macrocode}
-% \end{macro}
-%
% \subsubsection{Signs, and special numbers}
%
-% In the case of division, the order of the operands matters,
-% and it turns out to be slightly simpler if we internally
-% compute the \enquote{backwards} division.
-%
-% \begin{macro}[EXP,aux]{\@@_div_back_cases:NN}
-% \begin{syntax}
-% \cs{@@_div_back_cases:NN} \meta{case_2} \meta{case_1}
-% ~~\meta{sign_2} \meta{sign_1} \meta{body_2} |;| \meta{body_1} |;|
-% \end{syntax}
-% Expands the following tokens on the input stream once.
-% \begin{macrocode}
-\cs_new:Npn \@@_div_back_cases:NN #1 #2
- {
- \if_case:w \if_int_compare:w #1 = #2 \exp_stop_f:
- #1 \exp_stop_f:
- \else:
- \if_int_compare:w #1 < #2 \exp_stop_f:
- \if:w 3 #2 \c_four \else: \c_five \fi:
- \else:
- \if:w 3 #1 \c_six \else: \c_seven \fi:
- \fi:
- \fi:
- \exp_after:wN \@@_div_back_invalid:NNNww \exp_after:wN 0
- \or: \exp_after:wN \@@_div_back_normal:NNww
- \or: \exp_after:wN \@@_div_back_invalid:NNNww \exp_after:wN 2
- \or: \exp_after:wN \@@_basics_return_nan_nan:NNww
- \or:
- \exp_after:wN \@@_basics_return_ii:NNNNww
- \exp_after:wN #1
- \exp_after:wN #2
- \or: \exp_after:wN \@@_basics_return_inf:NNww
- \or:
- \exp_after:wN \@@_basics_return_i:NNNNww
- \exp_after:wN #1
- \exp_after:wN #2
- \or: \exp_after:wN \@@_basics_return_zero:NNww
- \fi:
- }
-% \end{macrocode}
-% Most of the special cases are common with some
-% previous operations. We only need to write the cases of
-% $0/0$ and $\infty/\infty$.
-% \begin{macrocode}
-\cs_new:Npn \@@_div_back_invalid:NNNww #1#2#3 #4; #5;
- {
- \@@_invalid_operation:Nnww \c_nan_fp { / }
- \s_@@ \@@_chk:w #1 #3 #5 ;
- \s_@@ \@@_chk:w #1 #2 #4 ;
- }
-% \end{macrocode}
-% \end{macro}
+% Time is now ripe to tackle the hardest of the four elementary
+% operations: division.
%
-% \begin{macro}{\@@_div_back_normal:NNww}
-% \begin{syntax}
-% \cs{@@_div_back_normal:NNww} \meta{sign_1} \meta{sign_2}
-% ~~\Arg{exp_1} \meta{body_1} |;| \Arg{exp_2} \meta{body_2} |;|
-% \end{syntax}
-% We now have two normal numbers to divide. Combine the signs.
+% \begin{macro}[int, EXP]{\@@_/_o:ww}
+% Filtering special floating point is very similar to what we did for
+% multiplications, with a few variations. Invalid operation
+% exceptions display |/| rather than |*|. In the formula for
+% dispatch, we replace |- \c_two +| by |-|. The case of normal
+% numbers is treated using \cs{@@_div_npos_o:Nww} rather than
+% \cs{@@_mul_npos_o:Nww}. There are two additionnal cases: if the
+% first operand is normal and the second is a zero, then the division
+% by zero exception is raised: cases $10$ and $11$ of the
+% \cs{if_case:w} construction in \cs{@@_mul_cases_o:NnNnww} are
+% provided as the fourth argument here.
% \begin{macrocode}
-\cs_new:Npn \@@_div_back_normal:NNww #1#2
+\cs_new_nopar:cpn { @@_/_o:ww }
{
- \if:w #1#2
- \exp_after:wN \@@_div_back_npos:Nnwnw
- \exp_after:wN 0
- \else:
- \exp_after:wN \@@_div_back_npos:Nnwnw
- \exp_after:wN 2
- \fi:
+ \@@_mul_cases_o:NnNnww
+ /
+ { - }
+ \@@_div_npos_o:Nww
+ {
+ \or:
+ \@@_case_use:nw
+ { \@@_division_by_zero_o:NNww \c_inf_fp / }
+ \or:
+ \@@_case_use:nw
+ { \@@_division_by_zero_o:NNww \c_minus_inf_fp / }
+ }
}
% \end{macrocode}
% \end{macro}
%
-% \subsubsection{Absolute (backwards) division}
-%
-% In this subsection, we perform the division
-% of two positive normal numbers.
-%
-% \begin{macro}[EXP]{\@@_div_back_npos:Nnwnw}
-% \begin{syntax}
-% \cs{@@_div_back_npos:Nnwnw} \meta{sign}
-% ~~\Arg{exp Z} \Arg{Z_1} \Arg{Z_2} \Arg{Z_3} \Arg{Z_4} |;|
-% ~~\Arg{exp A} \Arg{A_1} \Arg{A_2} \Arg{A_3} \Arg{A_4} |;|
-% \end{syntax}
-% We want to compute $A/Z$. As for addition and multiplication,
-% \meta{sign} is the final sign. Checking for underflow and
-% overflow is done using the same auxiliary as for multiplication.
-% As explained just below, we first compute $y$, which is
-% the $5$ first digits of $Z$, plus $1$, and then compute pieces
-% of the quotient roughly $4$ digits at a time. Here, |#1| is
-% a single digit, |#2| and |#7| are the exponents (integers),
-% |#8| is three brace groups, and all others are each $4$ digits.
+% \begin{macro}[aux, EXP]{\@@_div_npos_o:Nww}
+% \begin{quote}
+% \cs{@@_div_npos_o:Nww} \meta{final sign}
+% \cs{s_@@} \cs{@@_chk:w} |1| \meta{sign_A} \Arg{exp A}
+% \Arg{A_1} \Arg{A_2} \Arg{A_3} \Arg{A_4} |;|
+% \cs{s_@@} \cs{@@_chk:w} |1| \meta{sign_Z} \Arg{exp Z}
+% \Arg{Z_1} \Arg{Z_2} \Arg{Z_3} \Arg{Z_4} |;|
+% \end{quote}
+% We want to compute $A/Z$. As for multiplication,
+% \cs{@@_sanitize:Nw} checks for overflow or underflow; we provide it
+% with the \meta{final sign}, and an integer expression in which we
+% compute the exponent. We set up the arguments of
+% \cs{@@_div_mantissa_i_o:wnnw}, namely an integer \meta{y} obtained
+% by adding $1$ to the first $5$ digits of $Z$ (explanation given soon
+% below), then the four \Arg{A_{i}}, then the four \Arg{Z_{i}}, a
+% semi-colon, and the \meta{final sign}, used for rounding at the end.
% \begin{macrocode}
-\cs_new:Npn \@@_div_back_npos:Nnwnw #1 #2 #3#4#5#6; #7 #8;
+\cs_new:Npn \@@_div_npos_o:Nww
+ #1 \s_@@ \@@_chk:w 1 #2 #3 #4 ; \s_@@ \@@_chk:w 1 #5 #6 #7#8#9;
{
\exp_after:wN \@@_sanitize:Nw
\exp_after:wN #1
\int_use:N \__int_eval:w
- #7 - #2
- \@@_div_mantissa_i:wNwnn #3; #4;
- #8 {#3}{#4}{#5}{#6} #1
+ #3 - #6
+ \exp_after:wN \@@_div_mantissa_i_o:wnnw
+ \int_use:N \__int_eval:w #7 \use_i:nnnn #8 + \c_one ;
+ #4
+ {#7}{#8}#9 ;
+ #1
}
% \end{macrocode}
% \end{macro}
%
-% We are given two numbers, $A=0.A_{1}A_{2}A_{3}A_{4}$
-% and $Z=0.Z_{1}Z_{2}Z_{3}Z_{4}$, in blocks of $4$ digits,
-% and we know that the first digits of $A_{1}$ and of $Z_{1}$
-% are non-zero. To compute $A/Z$, we proceed as follows.
-% \begin{itemize}
-% \item Find an integer $Q_{1} \simeq 10^4 A / Z$.
-% \item Replace $A$ by $B = 10^4 A - Q_{1} Z$.
-% \item Find an integer $Q_{2} \simeq 10^4 B / Z$.
-% \item Replace $B$ by $C = 10^4 B - Q_{2} Z$.
-% \item Find an integer $Q_{3} \simeq 10^4 C / Z$.
-% \item Replace $C$ by $D = 10^4 C - Q_{3} Z$.
-% \item Find an integer $Q_{4} \simeq 10^4 D / Z$.
-% \item Consider $E = 10^4 D - Q_{4} Z$, and ensure
+% \subsubsection{Work plan}
+%
+% In this subsection, we explain how to avoid overflowing \TeX{}'s
+% integers when performing the division of two positive normal numbers.
+%
+% We are given two numbers, $A=0.A_{1}A_{2}A_{3}A_{4}$ and
+% $Z=0.Z_{1}Z_{2}Z_{3}Z_{4}$, in blocks of $4$ digits, and we know that
+% the first digits of $A_{1}$ and of $Z_{1}$ are non-zero. To compute
+% $A/Z$, we proceed as follows.
+% \begin{itemize}
+% \item Find an integer $Q_{A} \simeq 10^{4} A / Z$.
+% \item Replace $A$ by $B = 10^{4} A - Q_{A} Z$.
+% \item Find an integer $Q_{B} \simeq 10^{4} B / Z$.
+% \item Replace $B$ by $C = 10^{4} B - Q_{B} Z$.
+% \item Find an integer $Q_{C} \simeq 10^{4} C / Z$.
+% \item Replace $C$ by $D = 10^{4} C - Q_{C} Z$.
+% \item Find an integer $Q_{D} \simeq 10^{4} D / Z$.
+% \item Consider $E = 10^{4} D - Q_{D} Z$, and ensure
% correct rounding.
-% \end{itemize}
-% The calculations of $B$, $C$, $D$, and $E$ can be done
-% exactly with only $16$ (or $17$) digits.
-%
-% Unfortunately, things are not as easy as they seem.
-% Firstly, we make sure that all intermediate steps are positive,
-% since negative results would require extra calculations at the end.
-% This requires that $Q_{1} \leq 10^4 A / Z$ etc. A reasonable
-% attempt would be to define $Q_{1}$ as
-% \[
+% \end{itemize}
+% The result is then $Q = 10^{-4} Q_{A} + 10^{-8} Q_{B} + 10^{-12} Q_{C}
+% + 10^{-16} Q_{D} + \text{rounding}$. Since the $Q_{i}$ are integers,
+% $B$, $C$, $D$, and~$E$ are all exact multiples of $10^{-16}$, in other
+% words, computing with $16$ digits after the decimal separator yields
+% exact results. The problem will be overflow: in general $B$, $C$,
+% $D$, and $E$ may be greater than $1$.
+%
+% Unfortunately, things are not as easy as they seem. In particular, we
+% want all intermediate steps to be positive, since negative results
+% would require extra calculations at the end. This requires that
+% $Q_{A} \leq 10^{4} A / Z$ \emph{etc.} A reasonable attempt would be
+% to define $Q_{A}$ as
+% \begin{equation*}
% \cs{int_eval:n} \left\{
-% \frac{ A_{1} A_{2} }{ Z_{1} + 1 } - 1 \right\}.
-% \]
-% Subtracting $1$ at the end takes care of the fact that e\TeX{}'s
-% \cs{__int_eval:w} rounds instead of truncating. We add $1$ to $Z_{1}$
-% because $ Z_{1} \leq 10^4 Z < Z_{1}+1$ and we need $Q_{1}$
-% to be an underestimate. However, we are now underestimating
-% $Q_{1}$ too much: it can be wrong by up to $100$, for instance
-% when $Z = 0.1$ and $A \simeq 1$. Then $B$ could take values up to
-% $10$ (maybe more), and a few steps down the line, we would run into
-% arithmetic overflow, since \TeX{} can only handle integers less than
-% roughly $2\cdot 10^9$.
-%
-% A better formula is to take
-% \[
-% Q_{1} = \cs{int_eval:n} \left\{
+% \frac{ A_{1} A_{2} }{ Z_{1} + 1 } - 1 \right\}
+% \leq 10^{4} \frac{A}{Z}
+% \end{equation*}
+% Subtracting $1$ at the end takes care of the fact that \eTeX{}'s
+% \cs{__int_eval:w} rounds divisions instead of truncating (really,
+% $1/2$ would be sufficient, but we work with integers). We add $1$ to
+% $Z_{1}$ because $Z_{1} \leq 10^{4}Z < Z_{1}+1$ and we need $Q_{A}$ to
+% be an underestimate. However, we are now underestimating $Q_{A}$ too
+% much: it can be wrong by up to $100$, for instance when $Z = 0.1$ and
+% $A \simeq 1$. Then $B$ could take values up to $10$ (maybe more), and
+% a few steps down the line, we would run into arithmetic overflow,
+% since \TeX{} can only handle integers less than roughly $2\cdot
+% 10^{9}$.
+%
+% A better formula is to take
+% \begin{equation*}
+% Q_{A} = \cs{int_eval:n} \left\{
% \frac{ 10 \cdot A_{1} A_{2} }
% { \left\lfloor 10^{-3} \cdot Z_{1} Z_{2} \right\rfloor + 1 }
% - 1 \right\}.
-% \]
-% This is always less than $10^9 A / (10^5 Z)$, as we wanted.
-% In words, we take the $5$ first digits of $Z$ into account,
-% and the $8$ first digits of $A$, using $0$ as a $9$-th digit
-% rather than the true digit for efficiency reasons. We shall
-% prove that using this formula to define all the $Q_{i}$
-% leads to no overflow. For convenience, let us denote
-% \[
+% \end{equation*}
+% This is always less than $10^{9} A / (10^{5} Z)$, as we wanted. In
+% words, we take the $5$ first digits of $Z$ into account, and the $8$
+% first digits of $A$, using $0$ as a $9$-th digit rather than the true
+% digit for efficiency reasons. We shall prove that using this formula
+% to define all the $Q_{i}$ avoids any overflow. For convenience, let
+% us denote
+% \begin{equation*}
% y = \left\lfloor 10^{-3} \cdot Z_{1} Z_{2} \right\rfloor + 1,
-% \]
-% so that, taking into account the fact that e\TeX{} rounds ties
-% away from zero,
-% \[
-% Q_{1} = \left\lfloor A_{1}A_{2}0/y - 1/2 \right\rfloor.
-% \]
-% Note that $10^4<y\leq 10^5$, and $999 \leq Q_{1} \leq 99989$.
-% Also note that this formula does not cause an overflow as long as
-% $A<2.147\cdots$, since the numerator involves an integer slightly
-% smaller than $10^9A$.
-%
-% Let us bound $B$:
-% \begin{align*}
-% 10^5 B
+% \end{equation*}
+% so that, taking into account the fact that \eTeX{} rounds ties away
+% from zero,
+% \begin{align*}
+% Q_{A}
+% &= \left\lfloor \frac{A_{1}A_{2}0}{y} - \frac{1}{2} \right\rfloor
+% \\
+% &>\frac{A_{1}A_{2}0}{y} - \frac{3}{2}.
+% \end{align*}
+% Note that $10^{4}<y\leq 10^{5}$, and $999 \leq Q_{A} \leq 99989$.
+% Also note that this formula does not cause an overflow as long as $A <
+% (2^{31}-1) / 10^{9} \simeq 2.147\cdots$, since the numerator involves an
+% integer slightly smaller than $10^{9} A$.
+%
+% Let us bound $B$:
+% \begin{align*}
+% 10^{5} B
% &=
% A_{1}A_{2}0 + 10 \cdot 0.A_{3}A_{4}
-% - 10 \cdot Z_{1}.Z_{2}Z_{3}Z_{4}
-% \cdot \left\lfloor A_{1}A_{2}0/y - 1/2 \right\rfloor
+% - 10 \cdot Z_{1}.Z_{2}Z_{3}Z_{4} \cdot Q_{A}
% \\
% &<
% A_{1}A_{2}0
-% \cdot \left( 1 - 10 \frac{Z_{1}.Z_{2}Z_{3}Z_{4}}{y} \right)
+% \cdot \left( 1 - 10 \cdot \frac{Z_{1}.Z_{2}Z_{3}Z_{4}}{y} \right)
% + \frac{3}{2} \cdot 10 \cdot Z_{1}.Z_{2}Z_{3}Z_{4} + 10
% \\
% &\leq
-% \frac{A_{1}A_{2}0 \cdot (y - 10 Z_{1}.Z_{2}Z_{3}Z_{4})}{y}
+% \frac{A_{1}A_{2}0 \cdot (y - 10 \cdot Z_{1}.Z_{2}Z_{3}Z_{4})}{y}
% + \frac{3}{2} y + 10
% \\
% &\leq
-% \frac{A_{1}A_{2}0}{y} + \frac{3}{2} y + 10
+% \frac{A_{1}A_{2}0\cdot 1}{y} + \frac{3}{2} y + 10
% \leq
-% \frac{10^9 A}{y} + 1.6 y
-% \end{align*}
-% At the last step, we hide $10$ into the second term
-% for later convenience. The same reasoning yields\footnote{Bruno:
-% I need to find much better notations. These are not great.}
-% \begin{align*}
-% 10^5 B &< 10^9 A/y + 1.6 y, \\
-% 10^5 C &< 10^9 B/y + 1.6 y, \\
-% 10^5 D &< 10^9 C/y + 1.6 y, \\
-% 10^5 E &< 10^9 D/y + 1.6 y. \\
-% \end{align*}
-% The goal is now to prove that none of $B$, $C$, $D$, and $E$
-% can go beyond $2.147\cdots$. Simply bounding each term on the
-% right-hand side separately will not be tight enough: for instance,
-% we would get $10^5 B < 10^5 + 1.6\cdot 10^5 = 2.6 \cdot 10^5$,
-% which is too large.
-%
-% Combining the various inequalities together with $A<1$, we get
-% \begin{align*}
-% 10^5 B &< 10^9/y + 1.6 y, \\
-% 10^5 C &< 10^{13}/y^2 + 1.6 (y + 10^4), \\
-% 10^5 D &< 10^{17}/y^3 + 1.6 (y + 10^4 + 10^8/y), \\
-% 10^5 E &< 10^{21}/y^4 + 1.6 (y + 10^4 + 10^8/y + 10^{12}/y^2). \\
-% \end{align*}
-% All of those bounds are convex functions of $y$ (since every power
-% of $y$ involved is convex, and hte coefficients are positive), and
-% thus maximal at one of the end-points of the allowed range
-% $10^4<y\leq 10^5$. Thus,
-% \begin{align*}
-% 10^5 B &< \mathrm{max} ( 1.16\cdot 10^5, 1.7 \cdot 10^5), \\
-% 10^5 C &< \mathrm{max} ( 1.32\cdot 10^5, 1.77 \cdot 10^5), \\
-% 10^5 D &< \mathrm{max} ( 1.48\cdot 10^5, 1.777 \cdot 10^5), \\
-% 10^5 E &< \mathrm{max} ( 1.64\cdot 10^5, 1.7777 \cdot 10^5). \\
-% \end{align*}
-% All of those bounds are less than $2.147\cdot 10^5$, and
-% we are thus within \TeX{}'s bounds in all cases!\footnote{Bruno:
-% but I need to check this very carefully again.}
-%
-% We will later need to have a bound on the $Q_{i}$. Their
-% definitions imply that $Q_{1} < 10^9 A/y - 1/2 < 10^5 A$ and
-% similarly for the other $Q_{i}$. Thus each of them is at most
-% $177770$.
-%
-% The last step is to ensure correct rounding. We have
-% \[
-% A/Z = \sum_{i=1}^4 \left(10^{-4i} Q_{i}\right) + 10^{-16} E/Z
-% \]
-% exactly. Furthermore, we know that the result will be between
-% $0.1$ (inclusive) and $10$, so we only need to know the integer
-% part of $E/Z$, and a \enquote{rounding} digit encoding the rest
-% (see maybe addition for an explanation of why). Equivalently,
-% we need to find the integer part of $2E/Z$, and determine whether
-% it was an exact integer or not (this serves to detect ties). Since
-% \[
-% \frac{2E}{Z} = 2\frac{10^5 E}{10^5 Z}
-% \leq 2\frac{10^5 E}{10^4} < 36,
-% \]
-% this integer part is between $0$ and $35$ inclusive. We let
-% e\TeX{} round
-% \[
+% \frac{10^{9} A}{y} + 1.6\cdot y.
+% \end{align*}
+% At the last step, we hide $10$ into the second term for later
+% convenience. The same reasoning yields
+% \begin{align*}
+% 10^{5} B &< 10^{9} A/y + 1.6 y, \\
+% 10^{5} C &< 10^{9} B/y + 1.6 y, \\
+% 10^{5} D &< 10^{9} C/y + 1.6 y, \\
+% 10^{5} E &< 10^{9} D/y + 1.6 y. \\
+% \end{align*}
+% The goal is now to prove that none of $B$, $C$, $D$, and $E$ can go
+% beyond $(2^{31}-1) / 10^{9} = 2.147\cdots$.
+%
+% Combining the various inequalities together with $A<1$, we get
+% \begin{align*}
+% 10^{5} B &< 10^{9}/y + 1.6 y, \\
+% 10^{5} C &< 10^{13}/y^{2} + 1.6 (y + 10^{4}), \\
+% 10^{5} D &< 10^{17}/y^{3} + 1.6 (y + 10^{4} + 10^{8}/y), \\
+% 10^{5} E &< 10^{21}/y^{4} + 1.6 (y + 10^{4} + 10^{8}/y + 10^{12}/y^{2}). \\
+% \end{align*}
+% All of those bounds are convex functions of $y$ (since every power of
+% $y$ involved is convex, and the coefficients are positive), and thus
+% maximal at one of the end-points of the allowed range $10^{4} < y \leq
+% 10^{5}$. Thus,
+% \begin{align*}
+% 10^{5} B &< \mathrm{max} ( 1.16\cdot 10^{5}, 1.7 \cdot 10^{5}), \\
+% 10^{5} C &< \mathrm{max} ( 1.32\cdot 10^{5}, 1.77 \cdot 10^{5}), \\
+% 10^{5} D &< \mathrm{max} ( 1.48\cdot 10^{5}, 1.777 \cdot 10^{5}), \\
+% 10^{5} E &< \mathrm{max} ( 1.64\cdot 10^{5}, 1.7777 \cdot 10^{5}). \\
+% \end{align*}
+% All of those bounds are less than $2.147\cdot 10^{5}$, and we are thus
+% within \TeX{}'s bounds in all cases!
+%
+% We will later need to have a bound on the $Q_{i}$. Their definitions
+% imply that $Q_{A} < 10^{9} A/y - 1/2 < 10^{5} A$ and similarly for the
+% other $Q_{i}$. Thus, all of them are less than $177770$.
+%
+% The last step is to ensure correct rounding. We have
+% \begin{equation*}
+% A/Z = \sum_{i=1}^{4} \left(10^{-4i} Q_{i}\right) + 10^{-16} E/Z
+% \end{equation*}
+% exactly. Furthermore, we know that the result will be in $[0.1,10)$,
+% hence will be rounded to a multiple of $10^{-16}$ or of $10^{-15}$, so
+% we only need to know the integer part of $E/Z$, and a
+% \enquote{rounding} digit encoding the rest. Equivalently, we need to
+% find the integer part of $2E/Z$, and determine whether it was an
+% exact integer or not (this serves to detect ties). Since
+% \begin{equation*}
+% \frac{2E}{Z} = 2\frac{10^{5} E}{10^{5} Z}
+% \leq 2\frac{10^{5} E}{10^{4}} < 36,
+% \end{equation*}
+% this integer part is between $0$ and $35$ inclusive. We let \eTeX{}
+% round
+% \begin{equation*}
% P = \cs{int_eval:n} \left\{
-% 2 \frac{E_{1}E_{2}}{Z_{1}Z_{2}} \right\},
-% \]
-% which differs from $2E/Z$ by at most
-% \[
+% \frac{2\cdot E_{1}E_{2}}{Z_{1}Z_{2}} \right\},
+% \end{equation*}
+% which differs from $2E/Z$ by at most
+% \begin{equation*}
% \frac{1}{2}
% + 2 \left\lvert \frac{E}{Z} - \frac{E}{10^{-8} Z_{1}Z_{2}}\right\rvert
-% + 2 \left\lvert \frac{10^8 E - E_{1}E_{2}}{Z_{1}Z_{2}}\right\rvert
+% + 2 \left\lvert \frac{10^{8} E - E_{1}E_{2}}{Z_{1}Z_{2}}\right\rvert
% < 1,
-% \]
-% ($1/2$ comes from e\TeX{}'s rounding) because each absolute value
-% is less than $10^{-7}$. Thus $P$ is either the correct integer part,
-% or an overestimate by $1$ (impossible if $2E/Z$ is an integer). It
-% then suffices to compare $PZ$ with $2E$ to get the integer part of
-% $2E/Z$ and the information of whether it was an exact quotient or not.
-%
-% \begin{macro}[EXP]{\@@_div_mantissa_i:wNwnn}
-% \begin{syntax}
-% \cs{@@_div_mantissa_i:wNwnn}
-% ~~\meta{Z_1} |;| \meta{Z_2} |;|
-% ~~\Arg{A_1} \Arg{A_2} \Arg{A_3} \Arg{A_4}
-% ~~\Arg{Z_1} \Arg{Z_2} \Arg{Z_3} \Arg{Z_4} \meta{sign}
-% \end{syntax}
-% First compute $y$ from the first $5$ digits of $Z$, and
-% unbrace \meta{A_1} and \meta{A_2}.
-% \begin{macrocode}
-\cs_new:Npn \@@_div_mantissa_i:wNwnn #1; #2 #3; #4 #5
- {
- \exp_after:wN \@@_div_mantissa_ii:ww
- \int_use:N \__int_eval:w #1#2 + \c_one ;
- #4 #5 ;
- }
-% \end{macrocode}
+% \end{equation*}
+% ($1/2$ comes from \eTeX{}'s rounding) because each absolute value is
+% less than $10^{-7}$. Thus $P$ is either the correct integer part, or
+% is off by $1$; furthermore, if $2 E / Z$ is an integer, $P = 2 E / Z$.
+% We will check the sign of $2 E - P Z$. If it is negative, then $E / Z
+% \in \big((P - 1) / 2, P / 2\big)$. If it is zero, then $E / Z = P /
+% 2$. If it is positive, then $E / Z \in \big(P / 2, (P - 1) / 2\big)$.
+% In each case, we know how to round to an integer, depending on the
+% parity of $P$, and the rounding mode.
+%
+% \subsubsection{Implementing the mantissa division}
+%
+% \begin{macro}[aux, rEXP]{\@@_div_mantissa_i_o:wnnw}
% \begin{quote}
-% \cs{@@_div_mantissa_ii:ww}
-% ~~\meta{y} |;| \meta{A_1} \meta{A_2} |;| \Arg{A_3} \Arg{A_4}
-% ~~\Arg{Z_1} \Arg{Z_2} \Arg{Z_3} \Arg{Z_4} \meta{sign}
+% \cs{@@_div_mantissa_i_o:wnnw} \meta{y} |;|
+% \Arg{A_1} \Arg{A_2} \Arg{A_3} \Arg{A_4}
+% \Arg{Z_1} \Arg{Z_2} \Arg{Z_3} \Arg{Z_4} |;| \meta{sign}
% \end{quote}
-% Compute $Q_{1}$ by evaluating $\meta{A_1}\meta{A_2}0/y - 1$.
-% The result will be output to the left, in an \cs{__int_eval:w}
-% which we start now.
+% Compute $10^{6} + Q_{A}$ (a $7$~digit number thanks to the shift),
+% unbrace \meta{A_1} and \meta{A_2}, and prepare the
+% \meta{continuation} arguments for $4$ consecutive calls to
+% \cs{@@_div_mantissa_calc:wwnnnnnnn}. Each of these calls will need
+% \meta{y} (|#1|), and it turns out that we need post-expansion there,
+% hence the \cs{__int_value:w}. Here, |#4| is six brace groups, which
+% give the six first |n|-type arguments of the \texttt{calc} function.
% \begin{macrocode}
-\cs_new:Npn \@@_div_mantissa_ii:ww #1; #2;
+\cs_new:Npn \@@_div_mantissa_i_o:wnnw #1 ; #2#3 #4 ;
{
- \exp_after:wN \@@_div_mantissa_iii:www
- \__int_value:w #1 \exp_after:wN ;
- \__int_value:w
- \exp_after:wN \@@_div_mantissa_calc:Nwwnnnnnn
- \int_use:N \__int_eval:w #20/#1 + 999999 ; #2 ;
+ \exp_after:wN \@@_div_mantissa_test_o:w
+ \int_use:N \__int_eval:w
+ \exp_after:wN \@@_div_mantissa_calc:wwnnnnnnn
+ \int_use:N \__int_eval:w 999999 + #2 #3 0 / #1 ;
+ #2 #3 ;
+ #4
+ { \exp_after:wN \@@_div_mantissa_ii:wwn \__int_value:w #1 }
+ { \exp_after:wN \@@_div_mantissa_ii:wwn \__int_value:w #1 }
+ { \exp_after:wN \@@_div_mantissa_ii:wwn \__int_value:w #1 }
+ { \exp_after:wN \@@_div_mantissa_iii:wwnnnnn \__int_value:w #1 }
}
% \end{macrocode}
+% \end{macro}
+%
+% \begin{macro}[int, rEXP]{\@@_div_mantissa_calc:wwnnnnnnn}
+% \begin{macro}[aux, rEXP]
+% {
+% \@@_div_mantissa_calc_i:wwnnnnnnn,
+% \@@_div_mantissa_calc_ii:wwnnnnnnn,
+% }
% \begin{quote}
-% \cs{@@_div_mantissa_calc:Nwwnnnnnn} \meta{$10^6+{}$Q_1} |;|
-% ~~\meta{A_1} \meta{A_2} |;| \Arg{A_3} \Arg{A_4}
-% ~~\Arg{Z_1} \Arg{Z_2} \Arg{Z_3} \Arg{Z_4} \meta{sign}
+% \cs{@@_div_mantissa_calc:wwnnnnnnn} \meta{$10^{6}+{}$Q_{A}} |;|
+% \meta{A_1} \meta{A_2} |;| \Arg{A_3} \Arg{A_4}
+% \Arg{Z_1} \Arg{Z_2} \Arg{Z_3} \Arg{Z_4}
+% \Arg{continuation}
% \end{quote}
-% The goal here is to expand to
+% expands to
% \begin{quote}
-% \meta{$10^6+{}$Q_1} |;| \meta{B_1} \meta{B_2} |;| \Arg{B_3} \Arg{B_4}
-% ~~\Arg{Z_1} \Arg{Z_2} \Arg{Z_3} \Arg{Z_4} \meta{sign}
+% \meta{$10^{6}+{}$Q_{A}} \meta{continuation} |;|
+% \meta{B_1} \meta{B_2} |;| \Arg{B_3} \Arg{B_4}
+% \Arg{Z_1} \Arg{Z_2} \Arg{Z_3} \Arg{Z_4}
% \end{quote}
-% where $B = 10^4 A - Q_{1} \cdot Z$. More generally, this function
-% is used with $A\to B$, $B\to C$ and $Q_{1} \to Q_{2}$, etc.
-%
-% Computing the product $Q_{1} \cdot Z$ is almost simple, since
-% $Q_{1}$ is rather small, but not quite:
-% the product of $Q_{1}$ with each block of four digits $Z_{i}$
-% is within \TeX{}'s bounds, but we wouldn't be able to use the
-% usual trick of adding a large power of $10$ to ensure that the
-% number of digits is fixed (see other operations for many examples
-% of this). Instead, we split off the digit of $10^5$ in $Q_{1}$
-% (and more generally $Q_{i}$), and do something similar to the
-% case of the full multiplication.
-%
-% We know that $0<Q_{i}<1.8\cdot 10^5$, so $10^6+Q_{i}$ starts
-% with the digit $1$, followed by $\#1 = 1$ or $0$, then |#2|, which
-% is $5$ more digits. It would be somewhat simpler if we got |#1|
-% to be two digits, and |#2| four, but we are constrained by the $9$
-% arguments limit.
+% where $B = 10^{4} A - Q_{A} \cdot Z$. This function is also used to
+% compute $C$, $D$, $E$ (with the input shifted accordingly), and is
+% used in \pkg{l3fp-expo}.
+%
+% We know that $0<Q_{A}<1.8\cdot 10^{5}$, so the product of $Q_{A}$
+% with each $Z_{i}$ is within \TeX{}'s bounds. However, it is a
+% little bit too large for our purposes: we would not be able to use
+% the usual trick of adding a large power of $10$ to ensure that the
+% number of digits is fixed.
+%
+% The bound on $Q_{A}$, implies that $10^{6}+Q_{A}$ starts with the
+% digit $1$, followed by $0$ or $1$. We test, and call different
+% auxiliaries for the two cases. An earlier implementation did the
+% tests within the computation, but since we added a
+% \meta{continuation}, this is not possible because the macro has $9$
+% parameters.
%
% The result we want is then (the overall power of $10$ is arbitrary):
% \begin{align*}
-% &10^{-4} ( \#3 - \#2 \cdot \#6 - 10 \cdot \#1 \cdot \#6\#7 )
-% + 10^{-8} ( \#4 - \#2 \cdot \#7 - 10 \cdot \#1 \cdot \#8 ) \\
-% &+ 10^{-12}( \#5 - \#2 \cdot \#8 - 10 \cdot \#1 \cdot \#9 )
-% + 10^{-16}( - \#2 \cdot \#9 ).
+% &10^{-4} ( \#2 - \#1 \cdot \#5 - 10 \cdot \meta{i} \cdot \#5\#6 )
+% + 10^{-8} ( \#3 - \#1 \cdot \#6 - 10 \cdot \meta{i} \cdot \#7 ) \\
+% &+ 10^{-12}( \#4 - \#1 \cdot \#7 - 10 \cdot \meta{i} \cdot \#8 )
+% + 10^{-16}( - \#1 \cdot \#8 ),
% \end{align*}
-% The factors of $10$ come from the fact that
-% $Q_{i} = 10\cdot 10^4 \cdot \#1 + \#2$. As usual, to combine
-% all the terms, we need to choose some shifts which must ensure
-% that the number of digits of the second, third, and fourth terms
-% are each fixed. Here, a good choice is $2\cdot 10^9$.
-% We are flirting with \TeX{}'s limits once more.
-%
-% If $\#1=0$, then each term in parentheses (omitting the first)
-% is in the open interval $(-10^9, 10^4)$. Thus, adding
-% $2\cdot 10^9$ to it gives a $10$ digits number.\footnote{Bruno:
-% check that the carry from below does not screw that up. This
-% requires slightly tighter bounds.}
-%
-% If $\#1=1$, then $\#2 < 7.8 \cdot 10^4$, and each term
-% in parentheses (omitting the first) is in the interval
-% $(-8\cdot 10^8, 10^4)$, and we are even safer.
-%
-% We add the terms containing $\#1$ in a slightly tricky way
-% for efficiency reasons: if $\#1=0$, no need to do any computation,
-% while if $\#1=1$ we want $10$ times some number, simply obtained
-% by appending a $0$ digit.
+% where \meta{i} stands for the $10^{5}$ digit of $Q_{A}$, which is
+% $0$ or~$1$, and $\#1$, $\#2$, \emph{etc.\@} are the parameters of
+% either auxiliary. The factors of $10$ come from the fact that
+% $Q_{A} = 10\cdot 10^{4} \cdot \meta{i} + \#1$. As usual, to combine
+% all the terms, we need to choose some shifts which must ensure that
+% the number of digits of the second, third, and fourth terms are each
+% fixed. Here, the positive contributions are at most $10^{8}$ and
+% the negative contributions can go up to $10^{9}$. Indeed, for the
+% auxiliary with $\meta{i}=1$, |#1| is at most $80000$, leading to
+% contributions of at worse $-8\cdot 10^{8}4$, while the other
+% negative term is very small $<10^{6}$ (except in the first
+% expression, where we don't care about the number of digits); for the
+% auxiliary with $\meta{i}=0$, |#1| can go up to $99999$, but there is
+% no other negative term. Hence, a good choice is $2\cdot 10^{9}$,
+% which produces totals in the range $[10^{9}, 2.1\cdot 10^{9}]$. We
+% are flirting with \TeX{}'s limits once more.
% \begin{macrocode}
-\cs_new:Npn \@@_div_mantissa_calc:Nwwnnnnnn 1#1#2; #3;#4#5 #6#7#8#9
+\cs_new:Npn \@@_div_mantissa_calc:wwnnnnnnn 1#1
{
- 1 #1 #2 \exp_after:wN ;
- \int_use:N \__int_eval:w
- - 200000 + #3 - #2 * #6
- \if_meaning:w 1 #1
- - #6#70
- \fi:
- +
- \exp_after:wN \@@_div_mantissa_calc_last:NNNNNN
- \int_use:N \__int_eval:w
- 1999800000 + #4 - #2*#7
- \if_meaning:w 1 #1
- - #80
- \fi:
- +
- \exp_after:wN \@@_div_mantissa_calc_pack:NNNNNNw
- \int_use:N \__int_eval:w
- 1999800000 + #5 - #2*#8
- \if_meaning:w 1 #1
- - #90
- \fi:
- +
- \exp_after:wN \@@_div_mantissa_calc_pack:NNNNNNw
- \int_use:N \__int_eval:w 2000000000 - #2*#9 ;
- {#6}{#7}{#8}{#9}
+ \if_meaning:w 1 #1
+ \exp_after:wN \@@_div_mantissa_calc_i:wwnnnnnnn
+ \else:
+ \exp_after:wN \@@_div_mantissa_calc_ii:wwnnnnnnn
+ \fi:
}
-\cs_new:Npn \@@_div_mantissa_calc_pack:NNNNNNw #1#2#3#4#5#6 #7;
- { #1#2#3#4#5#6 ; {#7} }
-\cs_new:Npn \@@_div_mantissa_calc_last:NNNNNN #1#2#3#4#5#6
- { #1#2#3#4#5#6 \__int_eval_end: }
-% \end{macrocode}
-% \begin{quote}
-% \cs{@@_div_mantissa_iii:www} \meta{y} |;| \meta{$10^6+{}$Q_1} |;|
-% ~~\meta{B_1} \meta{B_2} |;| \Arg{B_3} \Arg{B_4}
-% ~~\Arg{Z_1} \Arg{Z_2} \Arg{Z_3} \Arg{Z_4} \meta{sign}
-% \end{quote}
-% \begin{macrocode}
-\cs_new:Npn \@@_div_mantissa_iii:www #1; #2; #3;
+\cs_new:Npn \@@_div_mantissa_calc_i:wwnnnnnnn #1; #2;#3#4 #5#6#7#8 #9
{
- \exp_after:wN \@@_div_mantissa_iii_after:w
- \int_use:N \__int_eval:w #2
- \exp_after:wN \@@_div_mantissa_iv:www
- \__int_value:w #1 \exp_after:wN ;
- \__int_value:w
- \exp_after:wN \@@_div_mantissa_calc:Nwwnnnnnn
- \int_use:N \__int_eval:w #30/#1 + 999999 ;
- #3 ;
+ 1 1 #1
+ #9 \exp_after:wN ;
+ \int_use:N \__int_eval:w \c_@@_Bigg_leading_shift_int
+ + #2 - #1 * #5 - #5#60
+ \exp_after:wN \@@_pack_Bigg:NNNNNNw
+ \int_use:N \__int_eval:w \c_@@_Bigg_middle_shift_int
+ + #3 - #1 * #6 - #70
+ \exp_after:wN \@@_pack_Bigg:NNNNNNw
+ \int_use:N \__int_eval:w \c_@@_Bigg_middle_shift_int
+ + #4 - #1 * #7 - #80
+ \exp_after:wN \@@_pack_Bigg:NNNNNNw
+ \int_use:N \__int_eval:w \c_@@_Bigg_trailing_shift_int
+ - #1 * #8 ;
+ {#5}{#6}{#7}{#8}
}
-% \end{macrocode}
-% \begin{quote}
-% \cs{@@_div_mantissa_iv:www} \meta{y} |;| \meta{$10^6+{}$Q_2} |;|
-% ~~\meta{C_1} \meta{C_2} |;| \Arg{C_3} \Arg{C_4}
-% ~~\Arg{Z_1} \Arg{Z_2} \Arg{Z_3} \Arg{Z_4} \meta{sign}
-% \end{quote}
-% \begin{macrocode}
-\cs_new:Npn \@@_div_mantissa_iv:www #1; #2; #3;
+\cs_new:Npn \@@_div_mantissa_calc_ii:wwnnnnnnn #1; #2;#3#4 #5#6#7#8 #9
{
- \exp_after:wN \@@_div_mantissa_pack:NNN
- \int_use:N \__int_eval:w #2
- \exp_after:wN \@@_div_mantissa_v:www
- \__int_value:w #1 \exp_after:wN ;
- \__int_value:w
- \exp_after:wN \@@_div_mantissa_calc:Nwwnnnnnn
- \int_use:N \__int_eval:w #30/#1 + 999999 ;
- #3 ;
+ 1 0 #1
+ #9 \exp_after:wN ;
+ \int_use:N \__int_eval:w \c_@@_Bigg_leading_shift_int
+ + #2 - #1 * #5
+ \exp_after:wN \@@_pack_Bigg:NNNNNNw
+ \int_use:N \__int_eval:w \c_@@_Bigg_middle_shift_int
+ + #3 - #1 * #6
+ \exp_after:wN \@@_pack_Bigg:NNNNNNw
+ \int_use:N \__int_eval:w \c_@@_Bigg_middle_shift_int
+ + #4 - #1 * #7
+ \exp_after:wN \@@_pack_Bigg:NNNNNNw
+ \int_use:N \__int_eval:w \c_@@_Bigg_trailing_shift_int
+ - #1 * #8 ;
+ {#5}{#6}{#7}{#8}
}
% \end{macrocode}
+% \end{macro}
+% \end{macro}
+%
+% \begin{macro}[aux, EXP]{\@@_div_mantissa_ii:wwn}
% \begin{quote}
-% \cs{@@_div_mantissa_v:www} \meta{y} |;| \meta{$10^6+{}$Q_3} |;|
-% ~~\meta{D_1} \meta{D_2} |;| \Arg{D_3} \Arg{D_4}
-% ~~\Arg{Z_1} \Arg{Z_2} \Arg{Z_3} \Arg{Z_4} \meta{sign}
+% \cs{@@_div_mantissa_ii:wwn} \meta{y} |;|
+% \meta{B_1} |;| \Arg{B_2} \Arg{B_3} \Arg{B_4}
+% \Arg{Z_1} \Arg{Z_2} \Arg{Z_3} \Arg{Z_4}
+% \meta{continuations} \meta{sign}
% \end{quote}
+% Compute $Q_{B}$ by evaluating $\meta{B_1}\meta{B_2}0 / y - 1$. The
+% result will be output to the left, in an \cs{__int_eval:w} which we
+% start now. Once that is evaluated (and the other $Q_{i}$ also,
+% since later expansions are triggered by this one), a packing
+% auxiliary takes care of placing the digits of $Q_{B}$ in an
+% appropriate way for the final addition to obtain $Q$. This
+% auxiliary is also used to compute $Q_{C}$ and $Q_{D}$ with the
+% inputs $C$ and $D$ instead of $B$.
% \begin{macrocode}
-\cs_new:Npn \@@_div_mantissa_v:www #1; #2; #3;
+\cs_new:Npn \@@_div_mantissa_ii:wwn #1; #2;#3
{
\exp_after:wN \@@_div_mantissa_pack:NNN
- \int_use:N \__int_eval:w #2
- \exp_after:wN \@@_div_mantissa_vi:wwnnnn
- \__int_value:w
- \exp_after:wN \@@_div_mantissa_calc:Nwwnnnnnn
- \int_use:N \__int_eval:w #30/#1 + 999999 ;
- #3 ;
+ \int_use:N \__int_eval:w
+ \exp_after:wN \@@_div_mantissa_calc:wwnnnnnnn
+ \int_use:N \__int_eval:w 999999 + #2 #3 0 / #1 ; #2 #3 ;
}
% \end{macrocode}
+% \end{macro}
+%
+% \begin{macro}[aux, rEXP]{\@@_div_mantissa_iii:wwnnnnn}
% \begin{quote}
-% \cs{@@_div_mantissa_vi:wwnnnn} \meta{$10^6+{}$Q_4} |;|
-% ~~\meta{E_1} \meta{E_2} |;| \Arg{E_3} \Arg{E_4}
-% ~~\Arg{Z_1} \Arg{Z_2} \Arg{Z_3} \Arg{Z_4} \meta{sign}
+% \cs{@@_div_mantissa_iii:wwnnnnn} \meta{y} |;|
+% \meta{E_1} |;| \Arg{E_2} \Arg{E_3} \Arg{E_4}
+% \Arg{Z_1} \Arg{Z_2} \Arg{Z_3} \Arg{Z_4} \meta{sign}
% \end{quote}
-% We compute $P$ by rounding $2 E_{1} E_{2}/Z_{1}Z_{2}$.
+% We compute $P \simeq 2E/Z$ by rounding $2 E_{1} E_{2}/Z_{1}Z_{2}$.
+% Note the first $0$, which multiplies $Q_{D}$ by $10$: we will later
+% add (roughly) $5\cdot P$, which amounts to adding $P/2 \simeq E/Z$
+% to $Q_{D}$, the appropriate correction from a hypothetical $Q_{E}$.
% \begin{macrocode}
-\cs_new:Npn \@@_div_mantissa_vi:wwnnnn #1; #2;#3#4 #5#6
+\cs_new:Npn \@@_div_mantissa_iii:wwnnnnn #1; #2;#3#4#5 #6#7
{
- \exp_after:wN \@@_div_mantissa_pack:NNN
- \int_use:N \__int_eval:w #10
- \exp_after:wN \@@_div_mantissa_vii:wwnnnnnn
- \int_use:N \__int_eval:w (\c_two*#2)/#5#6 ; % <- P
- #2;{#3}{#4}
- {#5}{#6}
+ 0
+ \exp_after:wN \@@_div_mantissa_iv:wwnnnnnnn
+ \int_use:N \__int_eval:w (\c_two * #2 #3) / #6 #7 ; % <- P
+ #2 ; {#3} {#4} {#5}
+ {#6} {#7}
}
% \end{macrocode}
-% Note that we used |#10| instead of |#2| which we had previously.
-% Two reasons: firstly, since we dropped $y$, the argument which
-% holds $Q_{i}$ has changed, and secondly, we will want the
-% fourth piece of the result to have $5$ digits, including the
-% \meta{rounding} digit, which we shall compute now from $P$.
+% \end{macro}
+%
+% \begin{macro}[aux, rEXP]
+% {
+% \@@_div_mantissa_iv:wwnnnnnnn,
+% \@@_div_mantissa_v:NNw,
+% \@@_div_mantissa_vi:Nw
+% }
% \begin{quote}
-% \cs{@@_div_mantissa_vii:wwnnnnnn} \meta{P} |;|
-% ~~\meta{E_1} \meta{E_2} |;| \Arg{E_3} \Arg{E_4}
-% ~~\Arg{Z_1} \Arg{Z_2} \Arg{Z_3} \Arg{Z_4} \meta{sign}
+% \cs{@@_div_mantissa_iv:wwnnnnnnn} \meta{P} |;|
+% \meta{E_1} |;| \Arg{E_2} \Arg{E_3} \Arg{E_4}
+% \Arg{Z_1} \Arg{Z_2} \Arg{Z_3} \Arg{Z_4} \meta{sign}
% \end{quote}
-% Then compute $2E-PZ$. Once more, we need to be careful and show
-% that the calculation $\#1\cdot\#5\#6$ below does not cause an
-% overflow: naively, $P$ can be up to $35$, and $\#5\#6$ up to
-% $10^8$, but both cannot happen simultaneously. To show that things
-% are fine, we split in two (non-disjoint) cases.
+% This adds to the current expression ($10^{7} + 10\cdot Q_{D}$) a
+% contribution of $5 \cdot P + \operatorname{sign}(T)$ with $T = 2 E -
+% P Z$. This amounts to adding $P / 2$ to $Q_{D}$, with an extra
+% \meta{rounding} digit. This \meta{rounding} digit is $0$ or $5$ if
+% $T$ does not contribute, \emph{i.e.,} if $0 = T = 2 E - P Z$, in
+% other words if $10^{16} A / Z$ is an integer or half-integer.
+% Otherwise it is in the appropriate range, $[1,4]$ or $[6,9]$. This
+% is precise enough for rounding purposes (in any mode).
+%
+% It seems an overkill to compute $T$ exactly as I do here, but I see
+% no faster way right now.
+%
+% Once more, we need to be careful and show that the calculation
+% $\#1\cdot\#6\#7$ below does not cause an overflow: naively, $P$ can
+% be up to $35$, and $\#6\#7$ up to $10^{8}$, but both cannot happen
+% simultaneously. To show that things are fine, we split in two
+% (non-disjoint) cases.
% \begin{itemize}
-% \item For small $P$, say, $P< 10$, the product obeys
-% $P\cdot\#5\#6 < 10^8 \cdot P < 10^9 $.
-% \item For large $P$, say, $P\geq 3$, the rounding error on $P$,
-% which is at most $1$, is less than a factor of $2$, hence
-% $P\leq 4E/Z$, and $P\cdot \#5\#6 \leq 4E\cdot 10^8 < 10^9$.
+% \item For $P < 10$, the product obeys $P\cdot\#6\#7 < 10^{8} \cdot P
+% < 10^{9} $.
+% \item For large $P\geq 3$, the rounding error on $P$, which is at
+% most $1$, is less than a factor of $2$, hence $P\leq 4E/Z$. Also,
+% $\#6\#7 \leq 10^{8} \cdot Z$, hence $P\cdot \#6\#7 \leq 4E\cdot
+% 10^{8} < 10^{9}$.
% \end{itemize}
% Both inequalities could be made tighter if needed.
%
-% Note however that $P\cdot \#7\#8$ may overflow,
-% since the two factors are now independent, and the result may reach
-% $3.5\cdot 10^9$.
-%
-% Also, we add $10\cdot P/2$ to the \enquote{fourth piece} of the result
-% as a first estimate of $10$ times $E/Z$. The goal is that the last digit
-% (for now $0$ or $5$) should be the \meta{rounding} digit. More precisely,
-% it will be corrected later by adding or subtracting $1$ depending on
-% whether $F$ was the correct integer part, or an overestimate (and nothing
-% is added when the quotient was exact). This does not give the
-% \enquote{correct} \meta{rounding} digit, but it always gives a digit
-% in the right \enquote{class} ($0$, $[1,4]$, $5$, or $[6-9]$), enough
-% for rounding purposes.
+% Note however that $P\cdot \#8\#9$ may overflow, since the two
+% factors are now independent, and the result may reach $3.5\cdot
+% 10^{9}$. Thus we compute the two lower levels separately. The rest
+% is standard, except that we use |+| as a separator (ending integer
+% expressions explicitly). $T$ is negative if the first character is
+% |-|, it is positive if the first character is neither |0| nor |-|.
+% It is also positive if the first character is |0| and second
+% argument of \cs{@@_div_mantissa_vi:Nw}, a sum of several terms, is
+% also zero. Otherwise, there was an exact agreement: $T = 0$.
% \begin{macrocode}
-\cs_new:Npn \@@_div_mantissa_vii:wwnnnnnn #1; #2;#3#4 #5#6#7#8
+\cs_new:Npn \@@_div_mantissa_iv:wwnnnnnnn #1; #2;#3#4#5 #6#7#8#9
{
+ \c_five * #1
- \exp_after:wN \@@_div_mantissa_ix:Nww
- \int_use:N \__int_eval:w -20 + 2*#2 - #1*#5#6 +
- \exp_after:wN \@@_div_mantissa_viii:NNw
- \int_use:N \__int_eval:w 199980 + 2*#3 - #1*#7 +
- \exp_after:wN \@@_div_mantissa_viii:NNw
- \int_use:N \__int_eval:w 200000 + 2*#4 - #1*#8 ; ;
+ \exp_after:wN \@@_div_mantissa_vi:Nw
+ \int_use:N \__int_eval:w -20 + 2*#2#3 - #1*#6#7 +
+ \exp_after:wN \@@_div_mantissa_v:NN
+ \int_use:N \__int_eval:w 199980 + 2*#4 - #1*#8 +
+ \exp_after:wN \@@_div_mantissa_v:NN
+ \int_use:N \__int_eval:w 200000 + 2*#5 - #1*#9 ;
}
-\cs_new:Npn \@@_div_mantissa_viii:NNw #1#2#3; { #1#2 ; + #3 }
-% \end{macrocode}
-% \begin{quote}
-% \cs{@@_div_mantissa_ix:Nww}
-% ~~\meta{F_1} \meta{F_2} |;| |+| \meta{F_3} |+| \meta{F_4} |;| \meta{sign}
-% \end{quote}
-% where $F=2E-PZ$. We only need to know whether it is positive,
-% negative, or exactly zero.
-% \begin{macrocode}
-\cs_new:Npn \@@_div_mantissa_ix:Nww #1#2;#3;
+\cs_new:Npn \@@_div_mantissa_v:NN #1#2 { #1#2 \__int_eval_end: + }
+\cs_new:Npn \@@_div_mantissa_vi:Nw #1#2;
{
\if_meaning:w 0 #1
- \exp_after:wN \@@_div_mantissa_x:w
- \int_use:N \__int_eval:w #3
+ \if_int_compare:w \__int_eval:w #2 > \c_zero + \c_one \fi:
\else:
- \if_meaning:w - #1
- -
- \else:
- +
- \fi:
- \c_one
- \fi:
- ;
- }
-\cs_new:Npn \@@_div_mantissa_x:w #1;
- {
- \if_int_compare:w #1 > \c_zero
- + \c_one
+ \if_meaning:w - #1 - \else: + \fi: \c_one
\fi:
;
}
% \end{macrocode}
-% We now obtain the following code, where \TeX{} is in the process
-% of expanding each of the integer expressions, and thus expands
-% the function at the bottom before the ones above it.
+% \end{macro}
+%
+% \begin{macro}[aux, EXP]{\@@_div_mantissa_pack:NNN}
+% At this stage, we are in the following situation: \TeX{} is in the
+% process of expanding several integer expressions, thus functions at
+% the bottom expand before those above.
% \begin{quote}
-% \cs{@@_div_mantissa_iii_after:w} $10^6 + Q_{1}$
-% \cs{@@_div_mantissa_pack:NNN} $10^6 + Q_{2}$
-% \cs{@@_div_mantissa_pack:NNN} $10^6 + Q_{3}$
+% \cs{@@_div_mantissa_test_o:w} $10^{6} + Q_{A}$
+% \cs{@@_div_mantissa_pack:NNN} $10^{6} + Q_{B}$
+% \cs{@@_div_mantissa_pack:NNN} $10^{6} + Q_{C}$
% \cs{@@_div_mantissa_pack:NNN}
-% $10^7 + 10\cdot Q_{4} + 5 \cdot P + \varepsilon$ |;| \meta{sign}
+% $10^{7} + 10\cdot Q_{D} + 5 \cdot P + \varepsilon$ |;| \meta{sign}
% \end{quote}
-% Here, $\varepsilon$ is $0$ in case $2E=PZ$ (\emph{i.e.}, $F=0$),
-% $1$ in case $2E>PZ$, which means that $P$ was the correct value,
-% but not with an exact quotient, and $-1$ if $2E<PZ$, \emph{i.e.},
-% $P$ was an overestimate.
+% Here, $\varepsilon = \operatorname{sign}(T)$ is $0$ in case $2E=PZ$,
+% $1$ in case $2E>PZ$, which means that $P$ was the correct value, but
+% not with an exact quotient, and $-1$ if $2E<PZ$, \emph{i.e.}, $P$
+% was an overestimate. The packing function we define now does
+% nothing special: it removes the $10^{6}$ and carries two digits (for
+% the $10^{5}$'s and the $10^{4}$'s).
% \begin{macrocode}
\cs_new:Npn \@@_div_mantissa_pack:NNN 1 #1 #2 { + #1 #2 ; }
% \end{macrocode}
-% Once those have been expanded, we get
+% \end{macro}
+%
+% \begin{macro}[aux, rEXP]{\@@_div_mantissa_test_o:w}
% \begin{quote}
-% \cs{@@_div_mantissa_iii_after:w} |1| |0| \meta{5d} |;|
+% \cs{@@_div_mantissa_test_o:w} |1| |0| \meta{5d} |;|
% ~~\meta{4d} |;| \meta{4d} |;| \meta{5d} |;| \meta{sign}
% \end{quote}
-% The reason we know that the first two digits are |1| and |0|
-% is that the final result is known to be between $0.1$ (inclusive)
-% and $10$, hence $\widetilde{Q_{1}}$ (the tilde denoting the
-% contribution from the other $Q_{i}$) is at most $99999$,
-% and $10^6+\widetilde{Q_{1}} = 10\cdots$.
-%
-% It is now time to round. This depends on how many digits
-% the final result will have.
+% The reason we know that the first two digits are |1| and |0| is that
+% the final result is known to be between $0.1$ (inclusive) and $10$,
+% hence $\widetilde{Q_{A}}$ (the tilde denoting the contribution from
+% the other $Q_{i}$) is at most $99999$, and $10^{6}+\widetilde{Q_{A}}
+% = 10\cdots$.
+%
+% It is now time to round. This depends on how many digits the final
+% result will have.
% \begin{macrocode}
-\cs_new:Npn \@@_div_mantissa_iii_after:w 10 #1
+\cs_new:Npn \@@_div_mantissa_test_o:w 10 #1
{
\if_meaning:w 0 #1
- \exp_after:wN \@@_div_mantissa_small:wwwNNNNwN
+ \exp_after:wN \@@_div_mantissa_small_o:wwwNNNNwN
\else:
- \exp_after:wN \@@_div_mantissa_large:wwwNNNNwN
+ \exp_after:wN \@@_div_mantissa_large_o:wwwNNNNwN
\fi:
#1
}
% \end{macrocode}
+% \end{macro}
+%
+% \begin{macro}[aux, EXP]{\@@_div_mantissa_small_o:wwwNNNNwN}
% \begin{quote}
-% \cs{@@_div_mantissa_small:wwwNNNNwN} |0| \meta{4d} |;|
-% ~~\meta{4d} |;| \meta{4d} |;| \meta{5d} |;| \meta{sign}
+% \cs{@@_div_mantissa_small_o:wwwNNNNwN} |0| \meta{4d} |;|
+% ~~\meta{4d} |;| \meta{4d} |;| \meta{5d} |;| \meta{final sign}
% \end{quote}
+% Standard use of \cs{@@_basics_pack_low:NNNNNw} and
+% \cs{@@_basics_pack_high:NNNNNw}. We finally get to use the
+% \meta{final sign} which has been sitting there for a while.
% \begin{macrocode}
-\cs_new:Npn \@@_div_mantissa_small:wwwNNNNwN
+\cs_new:Npn \@@_div_mantissa_small_o:wwwNNNNwN
0 #1; #2; #3; #4#5#6#7#8; #9
{
\exp_after:wN \@@_basics_pack_high:NNNNNw
@@ -1573,36 +1570,42 @@
\exp_after:wN ;
}
% \end{macrocode}
+% \end{macro}
+%
+% \begin{macro}[aux, rEXP]{\@@_div_mantissa_large_o:wwwNNNNwN}
% \begin{quote}
-% \cs{@@_div_mantissa_large:wwwNNNNwN} \meta{5d} |;|
+% \cs{@@_div_mantissa_large_o:wwwNNNNwN} \meta{5d} |;|
% ~~\meta{4d} |;| \meta{4d} |;| \meta{5d} |;| \meta{sign}
% \end{quote}
-% \footnote{Bruno: rename the \enquote{add mantissa carry pack} function.}
+% We know that the final result cannot reach $10$, hence |1#1#2|,
+% together with contributions from the level below, cannot reach
+% $2\cdot 10^{9}$. For rounding, we build the \meta{rounding digit}
+% from the last two of our $18$ digits.
% \begin{macrocode}
-\cs_new:Npn \@@_div_mantissa_large:wwwNNNNwN
+\cs_new:Npn \@@_div_mantissa_large_o:wwwNNNNwN
#1; #2; #3; #4#5#6#7#8; #9
{
+ \c_one
- \exp_after:wN \@@_div_mantissa_large_pack:NNNNNNNNw
- \int_use:N \__int_eval:w 1 #1 #2 %<- 1+9d
- \exp_after:wN \@@_add_mantissa_carry_pack_ii:NNNNw
- \int_use:N \__int_eval:w 1 #3 #4 #5 #6
- + \@@_round:NNNN #9 #6 #7 #8
- \exp_after:wN ;
+ \exp_after:wN \@@_basics_pack_weird_high:NNNNNNNNw
+ \int_use:N \__int_eval:w 1 #1 #2
+ \exp_after:wN \@@_basics_pack_weird_low:NNNNw
+ \int_use:N \__int_eval:w 1 #3 #4 #5 #6 +
+ \exp_after:wN \@@_round:NNN
+ \exp_after:wN #9
+ \exp_after:wN #6
+ \__int_value:w \@@_round_digit:Nw #7 #8 ;
+ \exp_after:wN ;
}
-\cs_new:Npn \@@_div_mantissa_large_pack:NNNNNNNNw
- 1 #1#2#3#4 #5#6#7#8 #9; { ; {#1#2#3#4} {#5#6#7#8} {#9} }
% \end{macrocode}
% \end{macro}
%
% \subsection{Unary operations}
%
-% \begin{macro}{\@@_neg:w}
-% This function flips the sign of the \meta{floating point}
-% and expands after it in the input stream, just like
-% \cs{@@_+_o:ww} etc.
+% \begin{macro}[int, EXP]{\@@_neg_o:w}
+% This function flips the sign of the \meta{floating point} and
+% expands after it in the input stream, just like \cs{@@_+_o:ww} etc.
% \begin{macrocode}
-\cs_new:Npn \@@_neg:w \s_@@ \@@_chk:w #1 #2
+\cs_new:Npn \@@_neg_o:w \s_@@ \@@_chk:w #1 #2
{
\exp_after:wN \@@_exp_after_o:w
\exp_after:wN \s_@@
@@ -1613,13 +1616,19 @@
% \end{macrocode}
% \end{macro}
%
-% \begin{macro}{\@@_abs:w}
+% \begin{macro}[int, EXP]{\@@_abs_o:w}
% This function sets the sign of the \meta{floating point} to be
% positive, and expands after itself in the input stream, just like
-% \cs{@@_neg:w}.
+% \cs{@@_neg_o:w}. We must leave the sign of \texttt{nan} invariant.
% \begin{macrocode}
-\cs_new:Npn \@@_abs:w \s_@@ \@@_chk:w #1 #2
- { \@@_exp_after_o:w \s_@@ \@@_chk:w #1 0 }
+\cs_new:Npn \@@_abs_o:w \s_@@ \@@_chk:w #1 #2
+ {
+ \exp_after:wN \@@_exp_after_o:w
+ \exp_after:wN \s_@@
+ \exp_after:wN \@@_chk:w
+ \exp_after:wN #1
+ \__int_value:w \if_meaning:w 1 #2 1 \else: 0 \fi: \exp_stop_f:
+ }
% \end{macrocode}
% \end{macro}
%
@@ -1631,4 +1640,4 @@
%
% \PrintChanges
%
-% \PrintIndex \ No newline at end of file
+% \PrintIndex