spigot: an exact real calculatorspigot?spigot
spigot's expression language
spigot's command line options
spigotspigot
This manual documents spigot, a command-line exact real calculator.
spigot?
spigot is a calculating program. It supports the usual arithmetic operations, square and cube roots, trigonometric and exponential functions, and a few other special functions such as erf.
spigot differs from the average calculating program in that it is an exact real calculator. This means that it does not suffer from rounding errors; in principle, it can keep generating more and more digits of the number you asked for until it runs out of memory.
In particular, if you ask for a complex expression such as sin(sqrt(pi)), then most calculating systems would compute first pi, then sqrt(pi) and finally sin(sqrt(pi)), accumulating a rounding error at each step, so that the final result had a build-up of error and you would have to do some additional error analysis to decide how much of the output you could trust.
spigot, on the other hand, does not output any digit until it is sure that digit is correct, so if you ask for (say) 100 digits of sin(sqrt(pi)) then you can be sure they are the right 100 digits.
(The downside of doing this is that spigot is slow, compared to the more usual methods of computer arithmetic. You wouldn't want to use even its individual arithmetic operations and functions for any kind of bulk computation. And because it evaluates a complicated expression as a whole, it's especially unsuited to any iterative algorithm which involves looking at the first result you get and then deciding what further arithmetic to do on it: normal arithmetic systems can do the next step of the algorithm with only the cost of the extra operations, but spigot would have to re-evaluate the entire thing from the start every time. So spigot can be useful if you really do need a lot of digits, or if you're doing a computation prone to numerical error, or as a cross-check on other calculating systems, but it wouldn't be usable outside that kind of specialist niche.)
spigot can evaluate expressions starting from numbers it can construct internally (integers, rationals, π, etc), and it can also read a number from a file or a pipe. In the latter mode, it operates ‘on-line’, i.e. it writes output as it goes along, and once it has read enough digits of an input number to know some of the output, it will print it.
spigot
The simplest thing spigot can do is to compute well-known mathematical constants such as π. For example, try running this command:
$ spigot pi
3.1415926535897932384626433832795028841971693993751058209749445923078164
062862089986280348253421170679821480865132823066470938446095505822317253
594081284811174502841027019385211055596446229489549303819644288109756659
(and so on)
In this default mode, spigot generates unbounded output: it will continue writing digits of π to the terminal until you interrupt it (e.g. by pressing Ctrl-C).
As spigot generates more digits, it will slow gradually down and consume more and more memory to keep track of where it's got to, and one of time and memory will ultimately be the limit to how much data it can generate. But in principle, given unlimited time and memory, it could just keep on going.
spigot can also evaluate more complicated mathematical expressions. (Some will be much slower than others.) You could try some of these, for example:
$ spigot pi*e
$ spigot pi/e
$ spigot pi^e
$ spigot -- '-pi^2'
$ spigot 'sin(sqrt(pi))'
$ spigot 'exp(pi*sqrt(163))'
(The above command lines assume that you are invoking spigot from a POSIX shell, so that expressions containing parentheses need to be single-quoted in order for spigot to receive them unmodified without the shell interfering. On other shells or operating systems, quoting conventions may vary.)
Since spigot supports a variety of command-line options beginning with -, if your expression also begins with a - (as one of the examples above does), then you need to precede it with the special argument word --, which instructs spigot to treat further command-line arguments as an expression rather than an option.
See chapter 3 for a full description of the range of expressions you can give to spigot.
In its default mode, spigot keeps generating digits until you interrupt it, or until it determines that it has output the exactly correct answer and needs no more digits at all.
You can limit spigot to a finite number of digits by using the -d option. For example:
$ spigot -d40 pi
3.1415926535897932384626433832795028841971
(The argument to -d counts decimal places, rather than significant figures. So here, there are 40 digits after the decimal point.)
By default, spigot does no rounding, i.e. it simply truncates the full decimal expansion of the number you asked for. (This is equivalent to rounding towards zero in all cases.)
If you want to change that, there's a variety of rounding-mode options, detailed in section 4.2. For example, you can ask spigot to round the final digit towards whichever value is nearer to the true result using --rn:
$ spigot --rn -d40 pi
3.1415926535897932384626433832795028841972
(This example rounds up, because the next digit following the 1 is a 6.)
By default, spigot's output is in base 10. You can use the -b option to output in another base instead:
$ spigot -b2 -d40 pi
11.0010010000111111011010101000100010000101
$ spigot -b3 -d40 pi
10.0102110122220102110021111102212222201112
$ spigot -b16 -d40 pi
3.243f6a8885a308d313198a2e03707344a4093822
$ spigot -b36 -d40 pi
3.53i5ab8p5fsa5jhk72i8asc47wwzlacljj9zn98l
For bases larger than 10, letters of the alphabet are used as additional digits. If you use the -B option in place of -b, those letters will be displayed in upper case:
$ spigot -B16 -d40 pi
3.243F6A8885A308D313198A2E03707344A4093822
As well as generating numbers in ordinary base notation, spigot can also generate them in the form of continued fractions. You can generate the continued fraction terms of a number, one per line, using the -c option:
$ spigot -c pi
3
7
15
1
292
1
(and so on)
(For those unfamiliar with continued fractions, this is equivalent to saying that π is equal to 3 + 1 / (7 + 1 / (15 + 1 / ...)).)
This one-per-line format is useful for applications which are consuming the numbers automatically, but for readability, you might prefer to add -l to replace the newlines with commas, so that you can see many more terms at once:
$ spigot -c -l pi
3;7,15,1,292,1,1,1,2,1,3,1,14,2,1,1,2,2,2,2,1,84,2,1,1,15,3,13,1,4,2,6,6
,99,1,2,2,6,3,5,1,1,6,8,1,7,1,2,3,7,1,2,1,1,12,1,1,1,3,1,1,8,1,1,2,1,6,1
,1,5,2,2,3,1,2,4,4,16,1,161,45,1,22,1,2,2,1,4,1,2,24,1,2,1,3,1,2,1,1,10,
(and so on)
(The semicolon emphasises that the first number is the integer part.)
Finally, you can use -C to output the continued fraction convergents (that is, the rational numbers obtained by evaluating successive truncations of the continued fraction, which are the best approximations to the target number by a particular metric):
$ spigot -C pi
3/1
22/7
333/106
355/113
103993/33102
104348/33215
(and so on)
In any of these modes, the -d option allows you to limit the number of continued fraction terms or convergents that spigot outputs. Again, the integer part is not counted, so that for example -d3 gives you three terms or convergents as well as the integer one:
$ spigot -c -d3 pi
3
7
15
1
$ spigot -C -d3 pi
3/1
22/7
333/106
355/113
As well as applying its collection of mathematical operators and functions to real numbers that it has derived itself from first principles, spigot can also apply the same operations to numbers received as input from elsewhere, which makes it at least theoretically able to operate on any real number.
To begin with, you can write a number into a file in ordinary decimal notation, and use a special syntax in the spigot expression language to read from that file:
$ echo 1.12358132134558914423337761098715972584418167651094617711 > z
$ spigot -d30 'sin( base10file:z )'
0.901654985409730168388244848164
You might wonder why spigot needs a feature like this: since the file contains only a finite amount of data, you surely could just have pasted the same number on to spigot's command line directly. The difference is that if you write a decimal number on spigot's command line, spigot will assume it is exact, i.e. that the precise number you meant is the rational number with that terminating decimal expansion. But if you read from a file like this, spigot will treat it as a prefix of some longer decimal expansion, and if it reaches the end of the file then it will report a special error saying that it can't compute any more output digits because it ran out of input. So this permits you to compute with imperfectly known inputs, without the risk of generating any incorrect digits of output due to the input imprecision.
You can also tell spigot to read from a file in number bases other than 10, or in continued fraction format. See section 3.5 for full details of all the available options.
But spigot can also read numeric data from its own standard input (or, if it's compiled on Unix, a numbered file descriptor of your choice), which can point at an unbounded data source such as a pipe. So if you invoke spigot with its standard input pointing at such a pipe, with a program at the other end of the pipe prepared to generate as many digits of the number as requested, then you need only use another special piece of expression syntax to tell spigot to read a real number from it, as far as is necessary.
For example, let's generate the Champernowne constant, which you get by concatenating all the positive integers in order immediately after the decimal point, with no leading zeroes. Here's a piece of Perl which generates that:
$ perl -e 'print"0.";print while++$_'
0.1234567891011121314151617181920212223242526272829303132333435363738394
041424344454647484950515253545556575859606162636465666768697071727374757
677787980818283848586878889909192939495969798991001011021031041051061071
(and so on)
Now we can use spigot to compute with this number (or any other input real), by piping that script's output into spigot, and telling spigot to expect to read a base-10 number from its standard input. As a really simple example, let's ask for spigot to output the same number, but in continued fraction form:
$ perl -e 'print"0.";print while++$_' | spigot -c -l base10stdin
0;8,9,1,149083,1,1,1,4,1,1,1,3,4,1,1,1,15,457540111391031076483646628242
956118599603939710457555000662004393090262659256314937953207747128656313
8641209375503552094607183089984575801469863148833592141783010987,6,1,1,2
1,1,9,1,1,2,3,1,7,2,1,83,1,156,4,58,8,54,4457353800911178833959067671634
293788437292958096324947188556700067877659324583930837874799958333344441
(and so on)
(Apart from it being very easy to generate, another reason I chose this example constant is because it has interestingly large numbers in its continued fraction – the huge number spanning the first three lines occurs because there's an extremely good rational approximation to the constant at the point where numbers change from 2 to 3 digits. The number starting on the fourth line is the corresponding term at the 3/4 digit boundary, and goes on for far longer!)
You can also use the expression base10stdin as part of a larger expression, so that you could compute mathematical functions of an input number, or add it to things (perhaps another input pipe, if you execute spigot in such a way as to assign multiple input fds).
spigot's expression language
This section gives a full specification of the language spigot uses for mathematical expressions.
spigot's expression language is semantically simpler than most such languages. In deference to the limitations of efficient computing, most expression languages have to carefully distinguish a collection of numeric types representing different subsets of the real numbers with different tradeoffs between range, precision and performance. But spigot, because its whole purpose is to deal in exact real numbers regardless of performance, does not have to do this. Every expression and subexpression in spigot has the same type: real.
So you don't have to remember to avoid writing 1/2 for a half (which in, say, gnuplot, would be interpreted as integer division and yield a rounded-down quotient of zero), or include any explicit syntax to change between formats as appropriate. Just write what looks like the obvious thing, and it will probably mean what you wanted.
spigot's expression language provides the following mathematical operators:
+ b (for two subexpressions a and b).
- b.
* b. You can also write this by juxtaposition, i.e. just write a b. For example, (pi+2)(pi+3) is a valid spigot expression, synonymous with (pi+2)*(pi+3).
/ b. As mentioned in the previous section, this is the real-number division operation, with no implicit rounding to integers in any circumstances. It is an error to use this operator with b = 0.
% b or a mod b (the two syntaxes are synonymous). This returns a value of the form a - nb, with n chosen to be whichever integer puts the result in the range [0,b) (if b is positive) or (b,0] (if b is negative).
rem b. This also returns a value of the form a - nb, but with slightly different semantics: this operation chooses n to be whichever integer puts the result in the range [-b/2,+b/2]. In case of a tie (if the result is precisely half way between two integer multiples of b, so that two different values of n meet that constraint), the answer is chosen to also make n an even number. (This is the same semantics as the IEEE 754 remainder operation.)
^ b or a ** b (the two syntaxes are synonymous).
- a.
!. This is implemented using the gamma function, so it can apply to non-integer values. It is an error to use this operator with a being a negative integer.
+ a is supported for symmetry with - a, but does nothing (just returns a unchanged).
The priorities and associativity of these operators, from lowest to highest, are as follows.
+ and - have lowest priority, and associate left to right. (E.g. an expression such as a-b+c-d is treated as if it were (((a-b)+c)-d).)
*, /, % and its synonym mod, and rem have the next higher priority, and also associate left to right. Multiplication by juxtaposition is included in this: when two numbers or parenthesised expressions appear adjacent in an expression, the parser imagines a * between them, and treats that * with the same priority as if it had been explicit.
- and + have the next higher priority. Since these are unary operators, there is no question of which way they associate: --a has only one possible parse in any case, namely -(-a).
^ and its synonym ** have the next higher priority. These operators associate right to left, i.e. a^b^c parses as if it had been a^(b^c). Note that since these operators have higher priority than the unary operators, -a^b parses as -(a^b), rather than the usually less helpful (-a)^b.
! has the highest priority.
As well as infix and prefix operators, spigot also supports a range of mathematical functions, all written in the style function(argument) or function(argument1,argument2). The following functions are provided:
sqrt(x) returns the square root of x. It is an error to use this function with x < 0.
cbrt(x) returns the cube root of x.
hypot(x,y) returns the hypotenuse of a right triangle with legs x and y, or equivalently, the square root of x^2 + y^2.
hypot(x0,…,xn) generalises the two-argument hypot function to return the square root of the sum of the squares of all its arguments.
sin(x) returns the sine of x, interpreted in radians.
cos(x) returns the cosine of x, interpreted in radians.
tan(x) returns the tangent of x, interpreted in radians. In theory, it would be an error to use this function with x being an odd integer multiple of π/2; in practice, spigot will never actually notice. (See chapter 5 for some discussion of this.)
asin(x) returns an angle θ in radians, in the range [-π/2,+π/2], such that sin(θ) = x. It is an error to use this function if x is not in the range [-1,+1].
acos(x) returns an angle θ in radians, in the range [0,+π], such that cos(θ) = x. It is an error to use this function if x is not in the range [-1,+1].
atan(x) returns an angle θ in radians, in the range [-π/2,+π/2], such that tan(θ) = x.
atan2(y,x) returns an angle θ in radians, in the range (-π,+π], such that the vector (sin(θ),cos(θ)) is a positive real multiple of the input vector (y,x). It is an error to use this function with y = x = 0.
sind(x) returns the sine of x, interpreted in degrees.
cosd(x) returns the cosine of x, interpreted in degrees.
tand(x) returns the tangent of x, interpreted in degrees. It is an error to use this function with x being an odd integer multiple of 90.
asind(x) returns an angle θ in degrees, in the range [-90,+90], such that sind(θ) = x. It is an error to use this function if x is not in the range [-1,+1].
acosd(x) returns an angle θ in degrees, in the range [0,+180], such that cosd(θ) = x. It is an error to use this function if x is not in the range [-1,+1].
atand(x) returns an angle θ in degrees, in the range [-90,+90], such that tand(θ) = x.
atan2d(y,x) returns an angle θ in degrees, in the range (-180,+180], such that the vector (sind(θ),cosd(θ)) is a positive real multiple of the input vector (y,x). It is an error to use this function with y = x = 0.
sinc(x) returns sin(x)/x (sometimes known as the ‘cardinal sine’ function). sinc(0) returns 1 (which makes the function continuous).
sincn(x) is the ‘normalised’ sinc function, equivalent to sinc(π x).
exp(x) returns e to the power of x.
exp2(x) returns 2 to the power of x.
exp10(x) returns 10 to the power of x.
log(x) returns the natural logarithm of x. It is an error to use this function if x <= 0.
log(x,b) returns the logarithm of x to the base b. It is an error to use this form of the log function if x <= 0, or if b <= 0, or if b = 1.
log2(x) returns the logarithm of x to the base 2. It is an error to use this function if x <= 0.
log10(x) returns the logarithm of x to the base 10. It is an error to use this function if x <= 0.
expm1(x) returns the same as exp(x) - 1.
log1p(x) returns the same as log(1 + x). It is an error to use this function if x <= -1.
pow(x,y) returns x raised to the power y. This is synonymous with the ^ and ** operators described in section 3.1, and included for convenience if pasting expressions out of languages such as C.
sinh(x) returns the hyperbolic sine of x.
cosh(x) returns the hyperbolic cosine of x.
tanh(x) returns the hyperbolic tangent of x.
asinh(x) returns the value y such that sinh(y) = x.
acosh(x) returns the non-negative value y such that cosh(y) = x. It is an error to use this function if x < 1.
atanh(x) returns the value y such that tanh(y) = x. It is an error to use this function if x is not in the range [-1,+1].
gamma(x) and tgamma(x) return the gamma function of x. (The two function names are synonymous – tgamma is provided for convenience if pasting expressions out of C.) It is an error to use this function if x is a negative integer or zero.
lgamma(x) returns the logarithm of gamma(x), or the logarithm of -gamma(x) if gamma(x) < 0. It is an error to use this function if x is a negative integer or zero.
factorial(x) returns the factorial of x, i.e. the same thing as gamma(x+1). It is an error to use this function if x is a negative integer.
erf(x) returns the error function of x.
erfc(x) returns the same as 1 - erf(x).
Phi(x) and norm(x) return the cumulative normal distribution function of x. (The two names are synonymous.)
erfinv(x) and inverf(x) both return the value y such that erf(y) = x. It is an error to use this function if x is not in the range (-1,+1).
erfcinv(x) and inverfc(x) both return the value y such that erfc(y) = x. It is an error to use this function if x is not in the range (0,2).
Phiinv(x), norminv(x), invPhi(x), invnorm(x) and probit(x) all return the value y such that Phi(y) = x. It is an error to use this function if x is not in the range (0,1).
W(x) and Wn(x) are the two branches of the Lambert W function. Each of them returns a value y such that y exp(y) = x. W(x) returns y >= -1, and it is an error to use it if x < -1/e; Wn(x) returns y <= -1, and it is an error to use it if x is not in the range [-1/e,0).
Ei(x) is the indefinite integral of exp(x)/x, defined in the negative domain so that its limit at negative infinity is zero, and defined in the positive domain (to deal with the pole at x = 0) so that the limit of Ei(-ε)-Ei(+ε) is zero as ε → 0. It is an error to use this function with x = 0.
En(n,x) is defined for non-negative integer n, and is the n-times-iterated indefinite integral of exp(-x)/x. That is, En(0,x) is just equal to exp(-x)/x; and for all larger n, En(n,x) is the indefinite integral of -En(n-1,x), defined so that its limit at infinity is zero. (The sign is flipped each time so that the function is positive for all n.) It is an error to use this function with x < 0, or with x = 0 and n equal to either 0 or 1.
E1(x) is a shorthand for En(1,x), and is also equal to -Ei(-x). It is an error to use this function with x <= 0.
Ein(x) is the indefinite integral of (1-exp(-x))/x, defined so that Ein(0) = 0. (This is the only one of the exponential integral family which is actually defined everywhere.)
li(x) is the indefinite integral of 1/log(x), defined so that li(0) = 0, and dealing with the pole at x = 1 by defining the limit of li(1-ε)-li(1+ε) to be zero as ε → 0 (the same trick as Ei above). It is an error to use this function with x < 0, or with x = 1.
Li(x) is also the indefinite integral of 1/log(x), but this time, defined so that Li(2) = 0. (That is, Li and li only differ by a constant.) Again, it is an error to use this function with x < 0 or with x = 1.
Li2(x) is the indefinite integral of log(1-x)/x, defined so that Li2(0) = 0. It is an error to use this function with x > 1.
Si(x) is the indefinite integral of sin(x)/x, defined so that Si(0) = 0.
si(x) is also the indefinite integral of sin(x)/x, defined so that its limit at positive infinity is zero.
Ci(x) is the indefinite integral of cos(x)/x, defined so that its limit at positive infinity is zero.
Cin(x) is the indefinite integral of (1-cos(x))/x, defined so that Cin(0) = 0.
FresnelS(x) and FresnelC(x) are the normalised Fresnel integrals: the indefinite integral of sin(π x^2/2) and cos(π x^2/2) respectively, both defined to be zero at zero.
UFresnelS(x) and UFresnelC(x) are the ‘unnormalised’ Fresnel integrals, i.e. the integrals of sin(x^2) and cos(x^2), also defined to be zero at zero.
BesselJ(a,x) is the Bessel function of the first kind of order a, which can be defined as the definite integral from t=0 to t=π of cos(at - x sin(t)) / π. It is an error to call this function with a not obviously an integer.
BesselI(a,x) is the modified Bessel function of the first kind of order a, which is obtained from the corresponding J function by extending it to the complex plane, giving it an imaginary input value, and rotating the output value back to the positive real axis. It is an error to call this function with a not obviously an integer.
zeta(s) is the Riemann zeta function for real arguments. That is, for s > 1, zeta(s) is the infinite sum of n^-s over all natural numbers n. The definition can be extended to s > 0 by rearranging it into the alternating sum of n^-s divided by 1-2^(1-s), and to negative s by Riemann's functional equation relating zeta(1-s) to zeta(s) (namely that their ratio is 2*gamma(s)*cos(pi*s/2)/((2*pi)^s)). It is an error to use this function with s = 1.
agm(a,b) is the arithmetic-geometric mean. This is the limiting value obtained by replacing the pair of numbers a and b with their arithmetic mean and their geometric mean, then repeating that operation so that both numbers converge inwards to a common limit. It is an error to use this function with a <= 0.
Hg(numerator-factors;denominator-factors;x) is the generalised hypergeometric function, defined by a power series in x in which the ratio between successive coefficients is given by a rational function with a term on top for each numerator factor, one on the bottom for each denominator factor, and an extra denominator factor consisting of factorials. Only rational factors are supported, but x can be irrational. It is an error to use this function with the number of numerator factors exceeding the number of denominator factors by more than 1, with any denominator factor being zero or a negative integer, or with the primary input value x being outside the function's domain of convergence (which varies with the parameters).
abs(x) returns the absolute value of x, i.e. either x or -x, whichever is non-negative.
sign(x) returns +1 if x is positive, -1 if x is negative, or zero if x is zero. (However, be wary of exactness hazards; many instances of the zero case may not manage to return any value at all See chapter 5.)
ceil(x) returns the smallest integer n such that n >= x.
floor(x) returns the largest integer n such that n <= x.
frac(x) returns the fractional part of x. This is always positive, even if x is negative; i.e. this is the same as x - floor(x).
fmod(x,y) returns the remainder of x by y, as the standard C function fmod would compute it. Unlike the mod and rem operators, the answer (if non-zero) is chosen to have the same sign as x.
round_mode(x) returns x rounded to an integer, by the rounding mode specified by the mode suffix on the function name. mode can be any of the same rounding-mode names that are available as command-line options, described in section 4.2. For example, round_rne rounds to the nearest integer, breaking a tie in favour of to whichever of the two equally near integers is even; round_ru rounds upwards (equivalent to ceil), and so on.
fracpart_mode(x) returns x - round_mode(x), i.e. the ‘fractional part’ of x after the integer part (determined by any rounding mode you choose) has been subtracted.
remainder_mode(x,y) returns the value x - qy, where q is an integer obtained by rounding the quotient x / y using the specified rounding mode. In particular, remainder_rne is equivalent to the rem operator; remainder_rd is equivalent to the mod operator; remainder_rz is equivalent to the fmod function.
algebraic(lo,hi,a0,a1,…,an) returns a root of the polynomial a0 + a1 x + … + an x^n which lies within the interval (lo,hi). It is an error (which spigot will not reliably detect) to use this function if the specified polynomial does not have a unique real root within that interval, and also an error if either of lo or hi or any polynomial coefficient is not obviously rational.
Some of these functions are considerably slower than others. The inverses of erf and Phi and the W functions, in particular, are implemented by laborious interval-bisection and are very slow indeed.
You can write actual numbers in spigot expressions, of course. By default a string of digits, or a string of digits with a decimal point somewhere in it, will be interpreted in decimal.
spigot supports C-style scientific notation, in which a decimal number (with or without a decimal point) is suffixed with e followed by an optional + or - and then another decimal integer; this denotes the first number times 10 to the power of the second. For example, 1.2 means 1.2, but 1.2e10 (or 1.2e+10) means 12000000000, and 1.2e-10 means 0.00000000012.
spigot also supports hex numbers, by prefixing 0x or 0X to a number, e.g. 0xabc means 2748. Hex numbers can be suffixed with an exponent part similar to the decimal one described above, only using p in place of e as the separator character; this means the hex number should be multiplied by the specified power of two. For example, 0x1.2 means 1.125 (one and an eighth); then 0x1.2p1 (or 0x1.2p+1) means twice that, i.e. 2.25, and 0x1.2p-1 means half of it, i.e. 0.5625.
If you want to input numbers in bases other than 10 and 16, you can prefix a numeric literal with base2:, base3:, …, base36:. For example, you could write base2:11.011, or base36:ZYX.WVU. No exponent suffix is recognised in this mode, because there's no solid convention for what it should look like or even what should be raised to the specified power (since the usual exponent-bases for bases 10 and 16 are 10 and 2 respectively). If you want to specify a number in scientific notation and an arbitrary base, write the exponent as an explicit power of something, e.g. base7:1.234 * 7^13.
A final way you can specify literal numbers to spigot is by specifying them as a hex number interpreted according to IEEE 754. You do this by writing ieee: followed by a number of hex digits, which must be exactly 4, 8, 16 or 32. These lengths correspond to the following formats:
An IEEE hex bit pattern can be followed by a . and further hex digits, in which case those digits are treated as an extension to the mantissa field. For example, ieee:3f800000 represents 1 (in IEEE single precision); ieee:3f800001 is the next representable number, namely 1+2^-23; and spigot will interpret ieee:3f800000.8 as the number half way between those two, i.e. 1+2^-24.
Infinities and NaNs are not permitted as IEEE format input.
spigot also supports a few mathematical constants under built-in names.
pi means π, the ratio between a circle's radius and half its circumference.
tau is a shorthand for 2*pi: the ratio between a circle's radius and its circumference. (See the Tau Manifesto.)
e means e, the base of natural logarithms.
phi means the golden ratio, i.e. (1+sqrt(5))/2.
eulergamma means the Euler-Mascheroni constant, i.e. the limit of the difference between log(n) and the sum of the reciprocals of the first n integers, as n tends to infinity.
apery means Apéry's constant, i.e. the sum of the reciprocals of the cubes of the natural numbers.
catalan means Catalan's constant, i.e. the alternating sum of the reciprocals of the squares of the odd natural numbers.
gauss means Gauss's constant, i.e. 1/agm(1,sqrt(2)).
(Many of these numbers can be generated by other methods, such as 4*atan(1) or exp(1), but it's more convenient to have them available as predefined constants.)
spigot can read a number from a file in various formats, and use it as input to its range of mathematical operations and functions (or just output it directly in a different format).
To read a number from a file in ordinary decimal, write base10file: followed by the file name. Any sequence of non-space characters following the : will be assumed to be the file name – even if they're punctuation or delimiters, e.g. if you write sin(base10file:myfile) then spigot will interpret the ) as part of the file name, and (even if a file with that strange name does exist) will complain that the expression is incomplete.
If you need to read a file that does have a space in its name, there's a quoting syntax to permit it. If the first character after the : is either ' or ", then spigot will look for the next occurrence of the same character, and treat all characters in between as the file name. A doubled quote character will be treated as a literal quote. For example:
base10file:'my file' will load the file ‘my file’.
base10file:"my file" will do the same.
base10file:'this isn''t sensible' will load the file ‘this isn't sensible’ (the doubled ' turns into a single one and does not terminate the quoted name).
base10file:"this isn't sensible" is a neater way to specify the same name (since the ' is not special in a filename quoted with ").
If your file contains the number in a base other than 10, you can prefix its name with base2file:, base3file:, …, base36file:. Bases larger than 36 are not supported, because after that the letters of the alphabet run out.
When spigot reads a file in base notation, it ignores newlines and white space interspersed between the digits.
You can also prefix a file name with cfracfile:, in which case spigot will expect it to contain a sequence of continued fraction coefficients (written in decimal). The coefficients can be separated by any non-digit characters you like (including newlines, like spigot's -c output mode, or ; and , like spigot's -c -l output mode); a minus sign is permitted before the very first coefficient to indicate a negative number, but ignored thereafter.
As mentioned in section 2.4, spigot will treat numbers read from files using any of the above keywords as if they are not an exact representation of a rational number, but a prefix of an infinitely long expansion of an irrational. So spigot will compute the output value as far as the input precision permits, and if it reaches the point where it can't generate any further output without more input precision, it will stop, print an error message indicating which input file ran out first, and return a special exit status 2.
If you really did want the contents of a file to be interpreted as the exact terminating expansion of a rational number, you can use the keywords base10xfile: (or likewise with bases other than 10) and cfracxfile, where the x stands for ‘exact’. The advantage of doing this rather than just specifying your expansion on spigot's command line directly is performance: if your input file is extremely large, then spigot will only have to read enough of it to generate the amount of output you asked for, whereas if you put the same number directly on the command line then spigot would have to parse all of it before starting.
spigot can also read from its standard input in the same way that it would read from a file, if you write an expression containing the keyword base10stdin (or another base, as above) or cfracstdin. For example, if you write base10stdin, then spigot will expect to read a number in base 10 from its own standard input (so you might pipe the output of another program into it).
If spigot has been compiled with the feature enabled (which it typically will have been on Unix systems; use spigot --version to check), you can also read from a numbered file descriptor of your choice in place of standard input, by writing base10fd:n (or another base, as above) or cfracfd:n. If you do this, spigot will expect its own file descriptor numbered n to already be open and readable, and will expect to read a number from there in the appropriate format. For example, base10fd:0 is another way of writing base10stdin; but, more interestingly, you could invoke spigot with a shell redirection operator such as 3<filename (where filename, in turn, could point at something interesting such as a named pipe, not necessarily a regular file), and then write (say) base10fd:3 to refer to that number.
Numbers read from stdin or numbered file descriptors are always treated as exact if the file descriptor signals EOF. (The expectation is that the whole point of using a file descriptor is so that you can connect it to a program that just keeps generating data.)
(If spigot is told to read from one or more file descriptors which refer to Unix terminal devices, you can use the -T option to set those terminals into raw mode. See section 4.3 for more detail.)
If you want spigot to run in a mode where it will refuse to read from files or file descriptors at all (for example, as a mild safety measure if input expressions are coming from an untrusted source), you can use the --safe option.
Sometimes, it's convenient to define your own variables or functions in a spigot expression. For example, suppose you want to evaluate a polynomial at some particular value. You could write something like this:
$ spigot 'sin(1.1)^3 - 2*sin(1.1)^2 + 5*sin(1.1) - 7'
but you probably got annoyed at the repetitiveness just reading that. So instead you could define a variable to have the repeated value:
$ spigot 'let x=sin(1.1) in x^3 - 2*x^2 + 5*x - 7'
or alternatively define a function to represent the polynomial:
$ spigot 'let f(x) = x^3 - 2*x^2 + 5*x - 7 in f(sin(1.1))'
In this example, it's more or less a matter of taste which of those you prefer, but the function syntax becomes more useful if you want to evaluate the function multiple times, e.g.
$ spigot 'let f(x) = x^3 - 2*x^2 + 5*x - 7 in f(f(sin(1.1)))'
The full syntax of the let statement is as follows:
let definitions in expression
where definitions is a comma-separated list of definitions, each in one of the following forms:
= expression
(parameters) = expression
In the latter case, parameters is in turn a comma-separated list of identifiers, and expression is allowed to refer to those identifiers as if they were variables.
Each definition comes into scope as soon as it is complete, but not before. So a sequence of definitions can refer back to each other:
$ spigot 'let x=pi/2, y=x^3, z=exp(y) in sin(z)'
$ spigot 'let f(x)=x+1, g(x)=sin(f(x)) in g(3)'
but a definition cannot refer to itself – in particular, functions cannot be defined recursively. (This is a fundamental limitation of spigot's evaluation system.)
spigot's command line options
spigot supports the following options to control the format in which numbers are output:
-b and -B, followed by a number between 2 and 36 inclusive, tell spigot to output in that number base. (Bases above 36 are not supported because the alphabet runs out.) In bases above 10, the Latin alphabet is used for extra digits, in lower case if -b was used or upper case if -B was used. The default is -b 10.
-w, followed by a positive number, tells spigot to output at least that many digits of the number's integer part, by writing leading zeroes if necessary.
-s, followed by a positive integer e, tells spigot to output a form of ‘scientific notation’, in which the number is expressed as a power of e times a number in the range [1,e). The power of e in this output mode is written first, unlike the more typical form of scientific notation which writes it at the end (because in spigot's case there might not be a far end of the number to write it at). You can use the special string ‘b’ in place of a positive integer, which means to use the same base for the exponent as the base the rest of the number is written in; for example, -sb on its own will produce output such as ‘10^8 * 4.85165’, with a power of 10 followed by a base-10 number.
-c tells spigot to output the number as a list of continued fraction convergents, separated by newlines.
-l is only effective in -c mode, and tells spigot to output the continued fraction convergents all on one long line, separated by ,, except that the initial integer part is separated from the rest by ;.
-C tells spigot to output the continued fraction convergents of the number. The convergents are separated by newlines, and each one is in the form of two decimal integers with a / between them.
-S, -D, -Q and -H tell spigot to output the number in the form of its IEEE 754 bit pattern, translated into hex. The four options specify single, double, quad and half precision respectively; see section 3.3 for more information on those formats. If the number is not exactly representable in the specified format, then the normal-length IEEE bit pattern will be followed by a . and further hex digits extending the mantissa field. If the number is too small, it will underflow to denormals or to zero; if it is too big, spigot will output the IEEE bit pattern of the appropriate sign of infinity.
-d followed by a number tells spigot to limit the number of digits or continued fraction terms it outputs. In the ordinary base output modes, the number counts digits after the integer part (i.e. decimal places); in continued fraction mode, it counts coefficients or convergents not including the initial one (which is conceptually the integer part, again); in IEEE mode, it counts extra bits of mantissa generated after the ‘.’. -d is only an upper limit: if the number is exactly representable with fewer digits than specified, then spigot will not output extra trailing zeroes. The argument to -d can be negative, in which case (in base or IEEE mode) rounding will occur at the specified distance before the point.
-R tells spigot to output the number's value as a rational. If the number is not rational, spigot will compute for ever without working that out, because some numbers it handles do turn out to be rational after some computation.
--printf, followed by a C-style printf format string, tells spigot to format the number in the way the C printf function would. The format string must consist of a single floating-point formatting directive (i.e. it must begin with % and end with the subsequent e, f, g or a conversion specifier). Most printf formatting modes also limit the number of digits, except for %a with no precision specification, which will continue to output digits until the number has been exactly represented – and if that is never (i.e. if the number does not have a terminating representation in binary) then spigot will continue to print mantissa digits for ever and never print the trailing exponent.
--printf mode with the a or A conversion specifier, --nibble adjusts the choice of exponent. spigot's default behaviour in hex printf mode is to choose the largest possible exponent, so that the leading hex digit of the output is always 1 (unless the number is exactly zero). --nibble changes the behaviour so that the exponent is always a multiple of 4, so that any leading digit is possible, and so that the hex digits always align to the digits that would be output in ordinary -b16 mode.
When spigot is given a digit limit via -d, and is not in continued fraction output mode, it defaults to printing only digits from the real expansion of the number, i.e. it always truncates toward zero and never rounds up.
You can make spigot round the output at the last digit position in various modes, using one of the following options:
--rz tells spigot to always round towards zero. This is the default.
--ri tells spigot to always round away from zero, i.e. the final digit is always incremented (unless the number was exact).
--ru tells spigot to always round up, i.e. towards positive infinity. This behaves like --rz for negative numbers, and like --ri for positive numbers.
--rd tells spigot to always round down, i.e. towards negative infinity. This behaves like --ri for negative numbers, and like --rz for positive numbers.
--rn tells spigot to round to nearest, i.e. the final digit will be incremented or not depending on which of the output numbers is closer in value to the true mathematical result. If there is an exact tie (the output is exactly half way between two representable numbers, e.g. asking for 0.25 to one decimal place) then spigot will break the tie by choosing whichever of the two outputs is an even multiple of the place-value of the final digit. Another name for this option is --rne (for ‘round to nearest, tie-breaking to even’).
(In an even base, this rounding rule is the same as the IEEE 754 round-to-nearest rule, which specifies breaking ties by choosing the output in which the final digit itself is even. In an odd base, the IEEE rule as written doesn't quite make sense, because in the case where rounding up causes a carry, both of the candidate values have an even final digit. But if you reinterpret the IEEE wording as meaning an even multiple of the final digit's place value, which is equivalent in an even base, then that version of the rule does generalise to odd bases and give an unambiguous answer in all cases.)
--rno tells spigot to round to nearest, but this time, tie-break to an odd multiple of the place value of the final digit.
--rnz, --rni, --rnu and --rnd tell spigot to round to nearest, and break ties by rounding them according to the corresponding directed-rounding option: --rz, --ri, --ru, or --rd respectively.
If spigot is told to print a limited number of digits via --printf rather than -d, then it defaults to rounding to nearest and tiebreaking to even (i.e. equivalent to --rn), because that's the default behaviour of the C printf function which spigot is imitating. This is done by making the --printf option behave as if --rn had been specified alongside it. So if you want to use a different rounding mode in printf-format output, you will need to specify the rounding-mode option second, or else the implicit one in the --printf will override it. For example, --printf --ru will work, but --ru --printf will revert to the --rn implied by the --printf option.
Here are some options which didn't fit nicely into the above categories.
--safe tells spigot to disallow all the keywords in the expression syntax that ask to read from files, file descriptors or standard input, as described in section 2.4 and section 3.5.
-T tells spigot to set Unix terminal devices into raw mode, if it is told to read from any by the baseNstdin, cfracstdin, baseNfd: or cfracfd: literal syntaxes (see section 3.5). Specifically, ‘raw mode’ means clearing the ICANON and ECHO bits. For example, you could run ‘spigot -T -b2 base10fd:0’, and then type a decimal number at the keyboard, and spigot will output the same number in binary, printing each bit as soon as it has enough information to know it; the use of -T means that your keystrokes are received instantly (without waiting for Return), and are not echoed to the screen to confuse the output.
-n tells spigot not to print a trailing newline, in base output mode (-b or -B), --printf mode, or one-line continued fraction mode (-c -l). Normally, if the output terminates, spigot will print a newline character before exiting.
-o, followed by a file name, tells spigot to send its output to that file instead of to its ordinary standard output channel.
--tentative=off, --tentative=on and --tentative=auto control the printing of ‘tentative output’. Tentative output consists of digits that spigot is not completely certain of yet, but is displaying temporarily in case it never manages to compute something better. The default is ‘auto’, meaning that tentative output is printed only when spigot's standard output points to a terminal device. For more detail, see chapter 5.
--version tells spigot to report its version number, and also any other build options, such as whether file descriptor input is supported, and which library is being used to provide the big-integer computation that spigot depends on.
--help displays an abbreviated list of spigot's command line options.
--licence displays spigot's copyright notice and licence text.
So far, this manual has more or less portrayed spigot as a magic device that can compute anything it likes. Perhaps some expressions are evaluated more slowly than others, but it'll get there in the end.
Well … not quite, sorry.
spigot's core algorithm works by finding an interval of rationals bracketing the target number, and gradually narrowing that interval further and further as more information is received or computed. Its various output modes work by waiting until the interval has narrowed to the point that the next digit is uniquely determined, and then writing that digit out.
The problem is: suppose that in order to work out the right output value, spigot has to decide whether a number falls on one side or the other of some boundary value. (This can happen in several different kinds of situation, which I categorise below.) It will do this by narrowing its interval until the boundary value doesn't fall inside it any more, and then it knows which side the number is on.
So the closer to the boundary the number is, the longer spigot will take before its interval narrows enough to work out which side it's on. And here's the problem: if the number it's trying to compute is exactly on the boundary value, then spigot may very well not ever notice – it will just keep narrowing its interval further and further, and the boundary value will still be inside the interval no matter how much it does that, so it will never manage to make a decision.
It's not standard terminology in exact real computation, but I've tended to call this sort of problem an ‘exactness hazard’. The general rule of thumb is that spigot is good at generating difficult output – complicated and fiddly irrational numbers – but can get confused if the answer to any computation is too easy, in particular if it's zero, or an integer, or has a finite number of digits in the output representation. So if you're computing any kind of long and fiddly expression, you need to watch out for too-simple exact numbers cropping up anywhere in it, because if you're unlucky, they can cause problems that prevent the expression as a whole from being computed.
In the following subsections, I list some possible effects of this kind of problem, and go into more detail about what numbers can trigger it.
One obvious situation in which spigot has to determine which side of a boundary a number falls is during final output, if the entire result of the computation has a terminating representation in the output base or number system. For example, suppose you run this command:
$ spigot 'sin(asin(0.12345))'
Obviously, we can see that the true answer is 0.12345 exactly. But spigot, with no symbolic algebra system, can't see that. So what can it do? It will manage to print ‘0.1234’ easily enough, but then its interval of possible output values will narrow for ever without allowing it to be certain of the next digit – no matter how much the interval narrows, it will still have one end looking like ‘0.123449999…’ and the other like ‘0.123450000…’, and any further work will just add more 9s to the former and more 0s to the latter, so there will never be a point at which it can be sure of that fifth digit.
spigot mitigates this situation by means of its tentative output feature. (See section 4.3 for options to control this.) What will actually happen if you run the above command, at least on a default Unix terminal, is that spigot will print ‘0.1234’ in the normal way, then follow that by printing the next digit ‘5’ in red. The red text indicates that the output is tentative, i.e. it might be retracted if more information comes along. So you might see a display along these lines:
$ spigot 'sin(asin(0.12345))'
0.12345 (10^-205)
This should be read as spigot saying: ‘The number definitely starts with 0.1234, and it looks as if it's exactly 0.12345, but I've only looked as far as 205 digits beyond that point, so there might still be surprises as I go further.’ As spigot computes more and more digits without finding a divergence from 0.12345, the power of 10 in the suffix will keep growing. And the red tentative text can be retracted if further information is discovered: if you had instead asked spigot to compute a number that was not exactly equal to 0.12345 but merely very very close to it, then eventually spigot would find that out, delete the red text, and print normal text in its place once it knew what it should be.
So if you see this kind of tentative output from spigot, that's your cue to have a closer look at the expression you asked it to compute, and see if you can see some mathematical reason why the answer is exactly what spigot is suggesting it might be. In this example case above, of course, that reason is pretty obvious – that's what asin means.
One exceptional case in which this kind of thing will work is that spigot makes a special case for numbers it knows from first principles to be rational. So while the above expression with answer 0.12345 ran into trouble, the following much simpler thing does work sensibly, and just prints the right answer immediately:
$ spigot 12345/100000
0.12345
because this time, spigot can remember that the number being output is the same as the known-rational number it got as input, so it can spot that it's an exact digit boundary and output it correctly.
(The special case for rationals is the only reason why section 4.2 needs to have all those different options for methods of tie-breaking in round-to-nearest mode – if no terminating representation could ever be successfully generated, then no tie-breaking options would be required, because spigot could never recognise a tie anyway.)
A more difficult case arises when spigot has the same problem of locating a number on one side or another of a boundary value, but in an intermediate result in a complex expression. For example, suppose you run one of these commands:
$ spigot 'tan(pi/2)'
$ spigot 'floor(sin(pi))'
In each of these situations, spigot will hang completely, and never manage to print any output at all. This is because, in each case, the value of a subexpression (respectively pi/2 and sin(pi)) is exactly on a point of discontinuity of the function it's being fed to next. In order to even start outputting the value of tan(x), spigot first needs to know whether it's positive or negative – and it can't find that out until it knows whether x is on one side or the other of π/2, which in this case it will never manage to decide. Similarly with floor(sin(pi)): to compute floor(x), spigot needs to know whether x < 0 or x >= 0, and again, it can't work that out in the case where x is exactly 0.
In this kind of situation, there really is nothing spigot can do but hang; even tentative output can't mitigate the situation.
As in the previous section, if the number is obviously a rational, then spigot can do better. If you actually ask it for floor(0), then it will handle that fine. But in the example above, the input to ‘floor’ is non-obviously zero, in the sense that in order to know it spigot would have to actually think about maths rather than just computing.
Therefore, if this happens to you, you can sometimes debug it by pulling out subexpressions and trying to evaluate them on their own; when you find one with a simple rational answer (probably signalled in turn by some red tentative output), see if you can prove that to be exactly the number you wanted, and substitute it in. For example:
$ spigot 'floor(sin(pi))'
(spigot hangs)
$ spigot 'sin(pi)'
0 (10^-951)
$ spigot 'floor(0)'
0
Here, the user observes spigot hanging, and suspects an exactness hazard somewhere in their expression. The only interesting subexpression is sin(pi), so the user evaluates that on its own, and gets back tentative output suggesting that the answer is exactly zero. The user thinks about it, realises that of course the answer is exactly zero, and simplifies the original expression by substituting a literal zero in place of the needlessly complicated sin(pi).
(Of course, in this simple example case, once the user had spotted the exact zero it hardly needed to ask spigot for floor of it! But if something more complicated were being done after the difficult point, it might well be easiest to substitute in the simple answer and try again with spigot.)
However, this debugging technique can only work in cases where the troublesome intermediate result is rational, because it's only rationals for which spigot has a special-case handler. There would be no equivalent way to debug the tan(pi/2) example above – spigot actually cannot recognise an input value to ‘tan’ as being an odd multiple of pi/2.
(This is the reason, mentioned in section 3.2, why the tan function can never return an error, even if you try to pass it an invalid input value. Because all the values where tan has no finite answer are irrational, so spigot will always run into this exactness hazard which prevents it attempting to compute an infinity.)
A third place where this same problem can occur, potentially, is internal to spigot itself. spigot's internal implementation of its various functions will often need to locate the input value on one or other side of some boundary value, and if spigot does that without sufficient caution, then it might hang forever when passed exactly the boundary value.
This should never happen. Exactness hazards internal to spigot itself are bugs. I've found and fixed all the ones I know of, but in case any more turn up, this section demonstrates how to recognise one, so that you can report it as a bug.
Here's an example of an exactness hazard that used to exist. spigot computes the power function a^b by computing exp(b*log(a)), in most cases. But if a < 0, then the right answer will have magnitude exp(b*log(abs(a))), but could be positive, negative or undefined depending on what kind of number b is. So the implementation of pow, having first ruled out a variety of special cases that don't require taking logs at all, would test the sign of a, and run straight into a hazard if a was non-obviously zero. For example, this could happen:
$ spigot '0^pi'
0
$ spigot 'sin(pi)'
0 (10^-951)
$ spigot 'sin(pi)^pi'
(spigot used to just hang)
If spigot knew that a was zero, then it could easily decide that zero to any power is zero. And it's at least capable of realising that sin(pi) is arbitrarily close to zero – but it wasn't able to produce even tentative output for sin(pi)^pi, because that internal sign test of a hung forever trying to decide which side of zero sin(pi) fell on.
Now the bug is fixed, and the last of those commands produces perfectly good tentative output:
$ spigot 'sin(pi)^pi'
0 (10^-807)
So if you suspect you're in a situation like this, where some spigot computation hangs (without even tentative output) for no reason you can see, try this diagnostic procedure:
f() is the outermost one in that subexpression, check that spigot can evaluate each argument to f at least as far as producing tentative output. (f might be a literal function, written with brackets, such as sin(x) or atan2(y,x), or it might be a mathematical operator such as addition or negation; the same procedure applies no matter how f is written.)
f is continuous at the point in question. If not, then this is an unavoidable hang as described in section 5.2, not an internal exactness hazard.
f is continuous at the specified point, and that spigot is producing at least tentative output for all the inputs to f, then please report it as a bug, including the exact spigot commands you ran to demonstrate the failing subexpression and the successful evaluation of the outermost function's arguments.
For example, you might have reported the above example as follows:
The following
spigotcommand hangs without even tentative output:$ spigot 'sin(pi)^pi'I think this is an internal exactness hazard, because
spigotcan evaluate the two operands to the power function safely (the former produces tentative output and the latter produces definite output):$ spigot 'sin(pi)' $ spigot 'pi'And the power function is mathematically continuous at the point
(0,π), so I thinkspigotshould not hang in this case.
That example report contains everything needed to check the facts, reproduce the failure, and decide whether it is indeed a bug.
One last class of hang in spigot is due not to a bug in the implementation, but a lack of generality in the fundamental design, which couldn't really be fixed without writing a totally different program from scratch.
As mentioned above, spigot works by narrowing an interval of rationals bracketing the number. But sometimes you know that the number is somewhere within one of two intervals of rationals, and can keep narrowing both of those intervals without ever working out which one contains the number. And sometimes, depending on what you want to do with the number next, that ought in principle to be good enough.
For example, it's possible to imagine a system that could cope with this:
$ spigot 'abs(atan2(sin(pi),-1))'
because the more you compute the input values to atan2, the closer you know the output is to one of +π and -π, even though you still don't know which – and since both of those values come out the same once they've been through the subsequent abs function, you could imagine the system being clever enough to cope with the expression as a whole. But spigot's design fundamentally is not.
spigot
Whatever output format you have chosen (see section 4.1), spigot will write data in the specified format to its standard output channel, i.e. Unix file descriptor 1.
If spigot successfully generates all the precision you asked for, it will print a terminating newline (if the output is any of the one-line types and you did not specify -n; see section 4.3), and then it will terminate with exit status 0 (the usual signal for success).
If spigot cannot generate enough output precision because it was reading a number from an input file and that input file ran out, it will still print its terminating newline (subject to the same conditions as above) on standard output; it will also write an error message to its standard error channel (Unix file descriptor 2) indicating which input file ran out first (in case there was more than one), and terminate with exit status 2.
If spigot encounters any other failure to evaluate the expression, it will print an error message to its standard error, and terminate with exit status 1.
spigot's central algorithm is derived from the paper ‘An unbounded spigot algorithm for the digits of π’, by Jeremy Gibbons. (American Mathematical Monthly, 113(4):318-328, 2006.) At the time of writing this, Jeremy Gibbons's web site has a PDF copy of the paper available.
Gibbons's algorithm is readily adapted to produce a continued fraction representation as output in place of base notation, and to accept numbers other than π as input. It forms the core of everything spigot does, and hence it seemed appropriate to name the entire program after it.
(Perhaps ironically, one part of Gibbons's algorithm that spigot does not use any more is his representation of π! It turned out that once spigot had grown the ability to compute square roots and do arithmetic, a representation based on Chudnovsky's formula was much faster.)
My algorithm for basic arithmetic is a sort of hybrid of Gibbons's spigot algorithm and William Gosper's algorithm for doing basic arithmetic on continued fractions (which I formerly used unmodified, but had to change it because using a continued fraction representation introduces exactness hazards). Gosper's algorithm comes originally from ‘HAKMEM’: Memo 239, Artificial Intelligence Laboratory, Massachusetts Institute of Technology, Cambridge, Mass., 1972. HAKMEM can also be downloaded at the time of writing this. Also, I haven't checked with great care, but the algorithm I constructed by combining both ideas may very well be the same as the one used by Imperial College's exact real arithmetic system (see the link below), in which case, they definitely had the same idea before I did.
Other algorithms used in this program came from my own head (though I'm not aware that any of them came from my head first).
spigot is not the only exact real calculator around, and was not really written with the aim of filling any particular gap in the available software. (I mostly wrote it for the fun of playing with Gibbons's elegant algorithm, and only realised some time later that it had become useful enough to be worth publishing.) Hence, if spigot is of interest to you, other implementations might also be of interest. Here are some that I managed to find:
spigot is copyright 2007-2019 Simon Tatham. All rights reserved.
Permission is hereby granted, free of charge, to any person obtaining a copy of this software and associated documentation files (the "Software"), to deal in the Software without restriction, including without limitation the rights to use, copy, modify, merge, publish, distribute, sublicense, and/or sell copies of the Software, and to permit persons to whom the Software is furnished to do so, subject to the following conditions:
The above copyright notice and this permission notice shall be included in all copies or substantial portions of the Software.
THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
spigot
spigot – command-line exact real calculator
spigot [ options ] expression
spigot is an exact real calculator: that is, you give it a mathematical expression to evaluate, and it computes it to any desired precision, by default simply printing digits to standard output until it is interrupted.
spigot provides command-line options to control the format of the output, restrict it to a specified number of digits, and apply rounding at the end of those digits. It can produce output in any base between 2 and 36 (after that it runs out of digit characters), or as a continued fraction, and it can read input numbers from files in any of those formats as well.
This man page gives only a brief summary of spigot's functionality. For full detail, you should read the main manual spigot.html; if that is not installed on your system, you can find it on the web at
https://www.chiark.greenend.org.uk/~sgtatham/spigot/spigot.html
The following options control spigot's basic output format:
-b base, -B base
-b and -B respectively. The default is -b 10.
-c
-C
/ between them.
-R
/ between them, or just one decimal integer if the number is a rational. If spigot does not know the number to be rational immediately, it will start evaluating it to see if it turns out rational later, so if it is not rational then spigot will compute for ever.
-S, -D, -Q, -H
--printf format, --printf=format
printf(3) would, given the formatting directive format. format must begin with a % and end with the associated conversion specifier, which must be a floating-point one (one of efgaEFGA).
The following options modify the details of those output formats:
-d limit
-b mode, no more than limit digits after the decimal point are printed. In -c or -C mode, no more than limit continued fraction coefficients or convergents are printed, not counting the initial one representing the number's integer part. In the IEEE 754 output modes, no more than limit additional bits of precision are generated after the end of the official mantissa. limit may be negative.
-l
-c mode, output continued fraction terms all on one line, separated by a ; after the first term and , after each subsequent term.
-w min-int-digits
-b mode, output at least min-int-digits of the number's integer part, by printing leading zeroes if necessary.
-s exponent-base (or -s b or -sb)
-b mode, output a prefix which is the largest power of exponent-base less than or equal to the absolute value of the number, and then print the number after dividing by that power. The special string ‘b’ for exponent-base means to use the same base for the exponent as the digit base used for the main output.
--nibble
--printf mode with the ‘a’ or ‘A’ conversion specifier, choose the output exponent to always be a multiple of 4, instead of the default behaviour of choosing it as large as possible.
-n
spigot prints output on a single line, suppress the usual trailing newline if spigot's output terminates.
The following options control rounding, when spigot's output is limited by the -d option. (Rounding does not occur in continued fraction modes.)
--rz
--ri
--ru
--rd
--rn, --rne
--rno
--rnz, --rni, --rnu, --rnd
--rz, --ri, --ru or --rd respectively.
Miscellaneous options:
-o output-file
spigot sends its output. The default is standard output; this option can be used to send it to a file instead.
--tentative=state
spigot does not know for sure what the next digit of the number is because it's starting to look as if it's exactly on a digit boundary. Tentative output is in red, and followed by an indication of about how many digits spigot has examined beyond that point (i.e. how close to exact that digit is known to be); spigot will retract it later if it finds out something definite.
state can be ‘on’, ‘off’ or ‘auto’. ‘auto’ is the default, and means that spigot should only print tentative output if its output is directed to a terminal device.
--safe
spigot to read from arbitrary files or from its own file descriptors. Might be useful if an expression is coming from an untrusted source (although, in that situation, you should still beware of other risks such as the expression author forcing spigot into an endless loop).
-T
ICANON and ECHO modes) while doing so.
spigot's expression language supports the following options, in order of priority from lowest to highest:
+ and -
*, /, %, mod, rem
% and mod are synonyms, which both return a remainder between 0 and the denominator; rem returns a remainder of either sign, with absolute value at most half that of the denominator, and ties broken by rounding to even in IEEE 754 style.
- and +
^, **
!
You can define variables and functions of your own in subexpressions using the let expression, as follows:
let var=value in expression
spigot input.
let fn(params)=defn in expression
(args) to refer to the expression defn with the arguments substituted in for the parameters. params must be a comma-separated list of identifiers; args is a comma-separated list of expressions.
A let expression can contain multiple definitions, separated by commas, e.g. ‘let x=1,y=2 in x+y’. Each definition is in scope for subsequent definitions, so you can write ‘let x=1,y=x+1 in’ expr. But definitions are not in scope for themselves; in particular, functions may not be recursive.
spigot also provides the following built-in functions:
sqrt, cbrt
hypot, atan2 (two arguments, or more for hypot)
hypot can also take a number of arguments other than two.
sin, cos, tan, asin, acos, atan
sind, cosd, tand, asind, acosd, atand, atan2d
d’ on the end except that angles are measured in degrees.
sinc, sincn
sinc is the ‘unnormalised’ form, i.e. just sin(x)/x; sincn is the ‘normalised’ form equal to sinc(pi*x).
sinh, cosh, tanh, asinh, acosh, atanh
exp, exp2, exp10, log, log2, log10
log.
expm1, log1p
exp(x)-1 and log(1+x).
pow (two arguments)
^ operator.
gamma, tgamma, lgamma
gamma and tgamma are synonyms for this), and the log of the absolute value of the gamma function.
factorial
! operator.
erf, erfc, Phi, norm
Phi and norm are synonyms for the cumulative normal distribution function.
erfinv, erfcinv, Phiinv, norminv
W, Wn
x exp(x). W is the branch with value at least -1, and Wn is the branch with value at most -1.
Ei, En (two arguments), E1, Ein
exp(x)/x. Ei(x) is the indefinite integral of exp(x)/x itself; En(n,x) (for non-negative integer n) is the result of integrating exp(-x)/x n times, flipping the sign each time; E1(x) is shorthand for En(1,x); and Ein(x) is the integral of (1-exp(-x))/x.
Li, li
1/log(x). Li(x) and li(x) are both the indefinite integral of 1/log(x); only their constants differ, in that Li(2) and li(0) are each defined to be zero.
Li2
-log(1-x)/x.
Si, si, Ci, Cin
sin(x)/x and cos(x)/x. Si(x) and si(x) are both the indefinite integral of sin(x)/x, differing only in the constant: Si(0)=0, but si(x) has limit 0 as x tends to positive infinity. Ci(x) is the indefinite integral of cos(x)/x, also with limit 0 at positive infinity; Cin(x) is the indefinite integral of (1-cos(x))/x, with Cin(0)=0.
UFresnelS, UFresnelC, FresnelS, FresnelC
UFresnelS and UFresnelC are the indefinite integrals of sin(x^2) and cos(x^2); FresnelS and FresnelC are the ‘normalised’ versions, i.e. integrals of sin(π x^2/2) and cos(π x^2/2). All are zero at the origin.
BesselJ, BesselI (two arguments, the first one an integer order)
zeta
agm (two arguments)
Hg (two lists of parameters and one primary input value, with semicolons separating the three kinds of argument)
abs
sign
-1 for a negative input, +1 for a positive input, or 0 for a zero input (provided spigot can determine that it is zero without an exactness hazard).
ceil, floor
frac
x - floor(x).
fmod
fmod(3), chosen to have the same sign as the numerator.
round_rz, round_rd, round_rne, ...
fracpart_rz, fracpart_rd, fracpart_rne, ...
round_xx function.
remainder_rz, remainder_rd, remainder_rne, ...
round_xx function.
algebraic (variable number of arguments)
spigot supports the following names for built-in constants:
pi, tau
e
phi
(1+sqrt(5))/2.
eulergamma
1/n.
apery
catalan
gauss
1/agm(1,sqrt(2)).
Numbers can be input in the following formats:
e+exponent or e-exponent for scientific notation
0x, and an optional C99-style p+exponent or p-exponent representing a power of 2 multiplier
baseN:, e.g. base7:0.123456
ieee:, followed by optional decimal point and extra mantissa digits
baseNfile: followed by a filename, e.g. base10file:pi.txt. The filename is taken to be the maximal sequence of non-space characters following the prefix, unless it starts with ' or ", in which case it is taken to be everything up to a matching closing quote, with doubled quote marks in between representing a literal quote character.
cfracfile: followed by a filename.
file: replaced by xfile: to indicate that end of file should be taken as the number being exactly represented rather than running out of precision.
baseNfd: or cfracfd: followed by an fd number, e.g. base10fd:0 to read from standard input. Also, baseNstdin and cfracstdin are available as synonyms for baseNfd:0 and cfracfd:0.
spigot returns 0 if its output terminates (because the result is exact, or because it reached the specified -d limit) with no problems.
In case of a parse error, or an invalid operand to a function, or any other kind of fatal error, spigot returns 1.
If spigot is unable to generate output to the desired precision because more precision was needed from a number read from an input file using baseNfile: or cfracfile:, then spigot returns 2, and prints an error message indicating which input file (in case there was more than one) ran out first.
Due to inherent limitations of its exact real arithmetic strategy, spigot is generally unable to recognise when a number it is computing is exactly equal to a specific boundary value.
One effect of this is that spigot will not behave as you'd like if the output number has a terminating representation in the selected base. For example, asking for sin(asin(0.12345)) will not be able to print 0.12345 and exit. Instead, spigot will get as far as printing ‘0.1234’, and then print tentative output (mentioned above) to indicate that it thinks the next digit might be exactly 5, but it will never reach a point where it's sure of that.
Another effect is that if you ask spigot to evaluate an expression in which an intermediate result is precisely on a point of discontinuity of the function it is passed to, then it may never manage to even start producing output. For example, spigot will hang completely if you ask it for floor(sin(pi)), since sin(pi) = 0 is a point of discontinuity of the floor function, and spigot will never be able to work out that the value of the input to floor is exactly zero, only that it seems to be closer and closer to zero the more it computes.
(An exception is numbers that spigot knows from first principles to be rational. For example, if you ask spigot to evaluate the simpler expressions ‘0.12345’ or ‘floor(0)’, it will print the complete output and terminate successfully, in both cases.)
spigot is free software, distributed under the MIT licence. Type ‘spigot --licence’ to see the full licence text.