NARS2000 - User contributions [en]

Ball Arithmetic

2018-10-26T10:17:05Z

Drab:

==Introduction==

Ball<ref name="BallArithmetic">NARS2000 Wiki, [http://arblib.org Ball Arithmetic]</ref> and Interval Arithmetic<ref name="IntervalArithmetic">Wikipedia, [https://en.wikipedia.org/wiki/Interval_arithmetic Interval Arithmetic]</ref> both refer to the same idea of calculating with arbitrary precision inexact numbers — that is, with a range of floating point values instead of a single representative. This allows such numbers to reflect their inexactness explicitly.

"The most common use is to keep track of and handle rounding errors directly during the calculation and of uncertainties in the knowledge of the exact values of physical and technical parameters. The latter often arise from measurement errors and tolerances for components or '''due to limits on computational accuracy'''. Interval arithmetic also helps find reliable and guaranteed solutions to equations and optimization problems."<ref name="IntervalArithmetic" />

The difference between the two is in how each approach represents its values: Interval Arithmetic uses the two ends of the interval <apll>[a, b]</apll> whereas Ball Arithmetic uses a Ball with a Midpoint and Radius. Each has its advantages. Interval Arithmetic is implemented by several libraries such as MPFI<ref name="MPFI">[http://perso.ens-lyon.fr/nathalie.revol/software.html MPFI]</ref> and Ball Arithmetic is implemented by ARB<ref name="BallArithmetic" />.

This feature uses the ARB<ref name="BallArithmetic" /> library to implement Ball Arithmetic.

==Environment==

In the examples below, unless otherwise specified, the precision of multi-precision floating-point numbers is set to <apll>53</apll> bits via <apll>⎕FPC←53</apll> which is the same precision as standard double-precision numbers. The Printing Precision is set to <apll>17</apll> digits via <apll>⎕PP←17</apll> so as to display the maximum number of significant digits at <apll>53</apll> bits of precision. The Comparison Tolerance is set via <apll>⎕CT←1<_E/>¯10</apll>.

==Calculating With Ball Arithmetic==

If you have ever entered an expression such as

<apll><pre>
÷3
0.3333333333333333
</pre></apll>

undoubtedly, you are aware that the result is inexact. Moreover, as you incorporate into your code the many other functions that produce inexact numbers, you should be wondering "How does the accumulation of inexactness affect my results?".

The answer to that question lies in the province of [https://en.wikipedia.org/wiki/Numerical_analysis Numerical Analysis] where it is described as "the study of algorithms that use numerical approximation". In general, applying the techniques of Numerical Analysis often involves time-consuming analysis. Fortunately, Ball Arithmetic is an easy to use alternative to that. That is, instead of analyzing the algorithm, run the algorithm with all calculations done with Ball Arithmetic.

Here is the same expression in Ball Arithmetic (use Alt-‘q’ or Ctrl-‘q’, depending on your keyboard layout, to enter the <apll><_pom/></apll> symbol):

<apll><pre>
÷3<_pom/>
0.3333333333333333<_pom/>5.6<_E/>¯17
</pre></apll>

Essentially, each Ball Arithmetic number carries with it an explicit Radius that describes the interval which contains the true value. Because Ball Arithmetic data propagates throughout a calculation and is never demoted to a lower datatype, its advantage is that after each operation the resulting Ball reliably includes the true answer. Across an entire algorithm, the final result is a Ball that indicates explicitly how precise it is, thus eliminating the guesswork as to the effect on the algorithm of the inexactness of the input and any other inexact constants used.

For example,

<apll><pre>
sqrt←{⍺←= ⋄ {0.5×⍵+⍺÷⍵}⍣⍺⍨⍵} ⍝ Calculate square root using Newton's method
√2<_pom/>
1.414213562373095<_pom/>2.2<_E/>¯16
sqrt 2<_pom/>
1.414213562373095<_pom/>1.1<_E/>¯15
10 sqrt 2<_pom/>
1.414213562373095<_pom/>2.8<_E/>¯15
100 sqrt 2<_pom/>
1.414213562373095<_pom/>3.3<_E/>¯14
1000 sqrt 2<_pom/>
1.414213562373095<_pom/>3.3<_E/>¯13
10000 sqrt 2<_pom/>
1.414213562373095<_pom/>3.5<_E/>¯12
</pre></apll>

As you can see, after a few iterations, the accuracy of the result doesn't improve. However, because it is an iterative algorithm, the precision worsens as the number of iterations increases. Only by using Ball Arithmetic can you see just how the precision of the result is affected.

To demonstrate with an extreme example, consider the following function:

<apll><pre>
f←{⍪⍵-3√(⍵*3)-⍵*2}
f 10*3×⍳5
0.33344450621439137
0.3333334452472627
0.33333444595336914
0.3348388671875
2.25
</pre></apll>

Without knowing much about this ([https://en.wikipedia.org/wiki/Condition_number ill-conditioned]) algorithm, you have no reason to suspect these results. As it turns out, there's a mathematical proof that as the argument increases to infinity, the algorithm converges to <apll>÷3</apll>. Only when you calculate with Ball Arithmetic do you see that there is a problem:

<apll><pre>
f 10*3×⍳5<_pom/>
0.3334445062153009<_pom/>2.2<_E/>¯12
0.3333334461785853<_pom/>4.3<_E/>¯9
0.333337664604187<_pom/>0.0000062
0.3385009765625<_pom/>0.0085
1.25<_pom/>5
</pre></apll>

As you can see, the precision worsens quickly to the point where the result is meaningless. The last line even says that the answer lies somewhere between <apll>¯3.75</apll> and <apll>6.25</apll>. In other words, that answer is so '''imprecise''' as to be '''inaccurate'''! At higher powers of ten, the answers only worsen.

The good news with Ball Arithmetic is that if you need greater precision, just increase the value of <apll>⎕FPC</apll>. For example,

<apll><pre>
⎕FPC←128 ⋄ ⎕PP←40
sqrt 2<_pom/>
1.41421356237309504880168872420969807857<_pom/>5.9<_E/>¯39
</pre></apll>

In the case of an ill-conditioned algorithm, we just push the problem off a bit, but not very far — the problem is that the algorithm requires infinite precision:

<apll><pre>
f 10*3×⍳10<_pom/>
0.333444506214021971364848444464298690637<_pom/>3.4<_E/>¯35
0.333333444444506172880658466392407138179<_pom/>1.1<_E/>¯31
0.333333333444444444506172839632668381652<_pom/>1.6<_E/>¯28
0.333333333333444444444444547458328123363<_pom/>2.2<_E/>¯25
0.333333333333333444444619474517508476812<_pom/>3.2<_E/>¯22
0.333333333333333333552432522356445687706<_pom/>1.3<_E/>¯19
0.333333333333333491771410805881714622956<_pom/>4.5<_E/>¯16
0.333333333333417414223731611855328083038<_pom/>4.5<_E/>¯13
0.33333333363771089352667331695556640625<_pom/>4.5<_E/>¯10
0.3333334587514400482177734375<_pom/>2.6<_E/>¯7
</pre></apll>

While your floating point code might not be so ill-conditioned, the only way to be sure of its results is one of

<ul>
<li>Hire a Numerical Analyst</li>
<li>Become a Numerical Analyst</li>
<li>Use Ball Arithmetic</li>
</ul>

==Ball Arithmetic Numbers==

A Ball Arithmetic number consist of two parts: a Midpoint and Radius. The Midpoint value represents the middle of the Ball and the Radius represents the distance around the Ball of possible values. There is no way to tell which value in the Ball is the correct answer — you should assume that all values are potential candidates.

The precision at which the Midpoint and Radius are stored depends upon the current value of <apll>⎕FPC</apll>. For example, the following display

<apll><pre>
⎕FPC←53
÷3<_pom/>
0.3333333333333333<_pom/>5.6<_E/>¯17
</pre></apll>

indicates that the number <apll>÷3<_pom/></apll> has a Midpoint of <apll>0.3333333333333333</apll> at <apll>53</apll>-bit precision and a Radius of <apll>5.6<_E/>¯17</apll>. The Midpoint is stored as a variable precision number, and the Radius is stored as a low-precision number.

Alternatively,

<apll><pre>
⎕FPC←128 ⋄ ⎕PP←40
÷3<_pom/>
0.333333333333333333333333333333333333332<_pom/>1.5<_E/>¯39
</pre></apll>

shows the corresponding value at a higher precision and correspondingly narrower Radius.


The approximate number of digits of precision for a given value of <apll>⎕FPC</apll> is about <apll>10⍟2*⎕FPC</apll>. The value of <apll>⎕FPC</apll> for a given required precision is about <apll>⌈2⍟10*digits</apll>.



The current <apll>⎕FPC</apll> precision
should be sufficient to cover both the midpoint and the radius
of all ball numbers being used.
Otherwise,
some results might be misleading.



For example,

<apll>
2v*⎕fpc←128
3.402823669E38
⍟ ⎕← * 1±1e¯1000
2.718281828±1.2E¯38
1±7.3E¯39
</apll>


What happened? -- <apll>⎕FPC</apll> was too small!


<apll>
2v*⎕fpc←4096
1.044388881E1233
⍟ ⎕← * 1±1e¯1000
2.718281828±2.7E¯1000
1±1E¯1000
</apll>

 

A Ball describes a range of '''variable-precision floating point numbers''', and so the expression <apll>23<_pom/>1</apll> represents the set of Real numbers between <apll>22</apll> and <apll>24</apll>, inclusive, '''not''' the three integers <apll>22</apll>, <apll>23</apll>, <apll>24</apll>.

In general, suffix notation on Ball Arithmetic constants should be sufficient for '''almost all purposes'''. That is, the Radius is already correctly calculated for constants such as <apll>23<_pom/></apll> whose Radius is zero, and <apll>1.2<_pom/></apll> whose Radius at <apll>53</apll> bits of precision is on the order of <apll>2*¯53</apll> (because of the inexactness of the Midpoint <apll>1.2</apll>) which displays as <apll>1.2<_pom/>2.2<_E/>¯16</apll>. Although the input notation supports it, the need to override the default Radius should be '''quite rare'''.

==New Functions and Operator==

===Functions To Convert To And From Ball Arithmetic Numbers===

A Ball Arithmetic constant number may be created through notation (as seen above). Otherwise, a Ball Arithmetic number also may be created from its constituent parts under program control with a new primitive function '''Contract''', and separated into its constituent parts with a new primitive function '''Distract'''.

====Contract====

'''Contract''' is a monadic function (<apll>≤R</apll>) whose right argument is a Real array. It may be used to create a Ball Arithmetic number from the one- or two-elements of its Midpoint and (optional) Radius. The Axis Operator may specify any length <apll>1</apll> or <apll>2</apll> coordinate along which the conversion takes place. It converts the first element along the specified dimension to the Midpoint of the result and (if present) the second element to its Radius. If the second element is not present (i.e., length <apll>1</apll> coordinate), the Radius used is calculated from the precision of the Midpoint.

<apll><pre>
≤23
23
≤1.2
1.2<_pom/>1.1<_E/>¯16
≤1.3 1<_E/>¯20
1.3<_pom/>1<_E/>¯20
≤1r2 1<_E/>¯20
0.5<_pom/>1<_E/>¯20
≤○1
3.141592653589793<_pom/>1.1<_E/>¯16
≤(○1) 1<_E/>¯20
3.141592653589793<_pom/>1<_E/>¯20
</pre></apll>

====Distract====

'''Distract''' is a monadic function (<apll>≥R</apll>) whose right argument is a Real array. It may be used to separate out a Ball Arithmetic's Midpoint and Radius into two distinct values. The Axis Operator may specify the place where the new length <apll>2</apll> coordinate is to be inserted. The pair of values along this new coordinate are the Midpoint and Radius of the corresponding number from the right argument. If the number can be represented exactly, the Radius is zero.

<apll><pre>
⎕PP←10
≥1.2
1.2 1.1102230246251565<_E/>¯16
≥1.2<_x/>
6<_r/>5 0
≥1.2<_v/>
1.2 1.110223024625157<_E/>¯16
≥1.2<_pom/>
1.2 2.220446049250313<_E/>¯16
≥1.5
1.5 0
≥1.5<_x/>
3<_r/>2 0
≥1.5<_v/>
1.5 1.110223024625157<_E/>¯16
≥1.5<_pom/>
1.5 0
a←(○1<_pom/>)×⍳4
a,'|',≥a
3.141592654<_pom/>4.4<_E/>¯16 | 3.141592654 4.440892148<_E/>¯16
6.283185307<_pom/>8.9<_E/>¯16 | 6.283185307 8.881784296<_E/>¯16
9.424777961<_pom/>1.3<_E/>¯15 | 9.424777961 1.332267643<_E/>¯15
12.56637061<_pom/>1.8<_E/>¯15 | 12.56637061 1.776356859<_E/>¯15
</pre></apll>

These numbers can be represented exactly:

<apll><pre>
≥⍳4<_pom/>
1 0
2 0
3 0
4 0
≥[2]⍳4<_pom/>
1 0
2 0
3 0
4 0
≥[1]⍳4<_pom/>
1 2 3 4
0 0 0 0
</pre></apll>

====Identities====

The Contract of the Distract of a number both along the same axis, returns the same number. The Distract of the Contract of a number both along the same length two axis, returns the same number.

===Operator To Provide Arbitrary Precision Specific Derived Functions===

The DoubleTilde monadic operator (<apll>≈</apll>, Alt-'Q' or Alt-shift-'q') takes a primitive function as its left operand and produces an arbitrary-precision-specific derived function.

This operator should be considered preliminary.

====Conversion Between Balls and Intervals====

At the moment, the monadic derived function is defined for two primitive functions.

Because Ball Arithmetic has a close relationship with Interval Arithmetic, it may be desirable to convert between the two forms: Endpoint to Midpoint-Radius and its inverse. The resulting derived function is not scalar.

For example,

<apll><pre>
x←0<_pom/>∞ 1<_pom/>2
≥≈x
¯∞ ∞
¯1 3
≥[1]≈x
¯∞ ¯1
∞ 3
≥[2]≈x
¯∞ ∞
¯1 3
≤≈≥≈x
0<_pom/>∞ 1<_pom/>2
≤≈≥[2]≈x
0<_pom/>∞ 1<_pom/>2
≤[2]≈≥≈x
0<_pom/>∞ 1<_pom/>2
≤[2]≈≥[2]≈x
0<_pom/>∞ 1<_pom/>2
≤[1]≈≥[1]≈x
0<_pom/>∞ 1<_pom/>2
</pre></apll>

The domain of these functions is limited to Real numbers. Hypercomplex numbers signal a <apll>DOMAIN ERROR</apll>. To see the Endpoints of a Hypercomplex number, first '''Dilate''' it so as to separate out the coefficients and then apply <apll>≥≈</apll> to extract the Endpoints as in <apll>≥≈>R</apll>.

For example,

<apll><pre> ≥≈ 2±0i0±1
DOMAIN ERROR
≥≈ 2±0i0±1
∧
≥≈>2±0i0±1
2 2
¯1 1
<≤≈≥≈>2±0i0±1
2±0i0±1
</pre></apll>

====Union and Intersection of Balls====

The left operand of the DoubleTilde operator is defined for the '''Union''' and '''Intersection''' functions (<apll>∪∩</apll>) as follows:

Because Ball arithmetic numbers represent a range of floating point numbers, at times it is desirable to produce their Union or Intersection.

For example,

<apll><pre>
x←0.5<pom/>0.5 ⍝ The interval [0,1]
x∪≈x+2 ⍝ The interval [0,1] ∪ [2,3] ←→ [0,3]
1.5<pom/>1.5
x∪≈x+1 ⍝ The interval [0,1] ∪ [1,2] ←→ [0,2]
1<pom/>1
x∩≈x+1 ⍝ The interval [0,1] ∩ [1,2] ←→ [1,1]
1
</pre></apll>

If two Balls have no points in common, the result is a NaN (<apll>∅</apll>):

<apll><pre>
x∩≈x+1.1 ⍝ The interval [0,1] ∩ [1.1,2.1] ←→ ∅
∅
</pre></apll>

The identity element for <apll>∩≈</apll> is <apll>0±∞</apll>.

The derived functions are scalar.

====Comparison Functions====

The left operand of the DoubleTilde operator is defined for the six comparison functions (<apll><≤=≥>≠</apll>) as follows:

For <apll>f</apll> as one of the above six functions, <apll>Lf≈R</apll> is TRUE iff for all <apll>t∊L</apll> and <apll>u∊R</apll>, then <apll>  t f u  </apll> is TRUE. These comparisons are all done with <apll>⎕CT←0</apll> because, in a sense, Ball comparisons carry their own Fuzz.

For example,

<apll><pre> (L R)←4<_pom/>1 6<_pom/>1
L<≈R
0
L≤≈R
1
</pre></apll>

The derived functions are scalar.

====Containing Functions====

The left operand of the DoubleTilde operator is defined for two containing functions (<apll>⊂⊃</apll>) as follows:

The expression <apll>L⊂≈R</apll> is TRUE iff for all <apll>t∊L</apll> it is TRUE that <apll>t∊R</apll> and that the containment is strict — that is, none of the endpoints match.

The expression <apll>L⊃≈R</apll> is TRUE iff <apll>R⊂≈L</apll> is TRUE.

For example,

<apll><pre>
(L R)←6±4 6±2
L⊃≈R
1
R⊂≈L
1
</pre></apll>

The derived functions are scalar.

====Positional Relationships====

For <apll>L</apll> and <apll>R</apll> as finite intervals with non-zero Radii (i.e., multi-valued), here are the possible '''Relative Positions''' between the two Intervals, where <apll>L0</apll>, <apll>L1</apll> and <apll>R0</apll>, <apll>R1</apll> are the endpoints of <apll>L</apll> and <apll>R</apll>, respectively — <apll>(L0 L1)←≥≈L</apll> and <apll>(R0 R1)←≥≈R</apll>:

<table border="1" cellspacing="0" rowspacing="5" style="font-family:APL385, SImPL;">
<th>Relative Positions</th>
<th>Endpoint Relationships</th>
<th colspan="2">Defining Expressions</th>
<th>Sample Data L   R</th>
<th>Distinct Ball {⍺∘≈⍵}</th>

<tr>
<td> ├──L──┤  ├──R──┤ </td>
<td> L0 L1≡R0 R1 </td>
<td> L=≈R </td>
<td> </td>
<td> 6<_pom/>4 6<_pom/>4 </td>
<td>  0 </td>
</tr>

<tr>
<td> ├──L──┤  ├──R─┤ </td>
<td> (L0=R0)∧L1>R1 </td>
<td> </td>
<td> </td>
<td> 6<_pom/>4 5<_pom/>3 </td>
<td>  1 </td>
</tr>

<tr>
<td> ├──L─┤  ├──R──┤ </td>
<td> (L0=R0)∧L1<R1 </td>
<td> </td>
<td> </td>
<td> 6<_pom/>4 7<_pom/>5 </td>
<td> ¯1 </td>
</tr>

<tr>
<td> ├──L──┤   ├──R─┤ </td>
<td> (L0<R0)∧L1=R1 </td>
<td> </td>
<td> </td>
<td> 6<_pom/>4 7<_pom/>3 </td>
<td>  2 </td>
</tr>

<tr>
<td>  ├──L─┤  ├──R──┤ </td>
<td> (L0>R0)∧L1=R1 </td>
<td> </td>
<td> </td>
<td> 6<_pom/>4 5<_pom/>5 </td>
<td> ¯2 </td>
</tr>

<tr>
<td> ├──L──┤   ├─R─┤ </td>
<td> (L0<R0)∧L1>R1 </td>
<td> L⊃≈R </td>
<td> </td>
<td> 6<_pom/>4 6<_pom/>2 </td>
<td>  3 </td>
</tr>

<tr>
<td>  ├─L─┤  ├──R──┤ </td>
<td> (L0>R0)∧L1<R1 </td>
<td> L⊂≈R </td>
<td> </td>
<td> 6<_pom/>4 6<_pom/>6 </td>
<td> ¯3 </td>
</tr>

<tr>
<td> ├──L──┤   ├──R──┤ </td>
<td> (L0<R0)∧(L1<R1)∧L1>R0 </td>
<td> </td>
<td> </td>
<td> 6<_pom/>4 8<_pom/>4 </td>
<td>  4 </td>
</tr>

<tr>
<td>  ├──L──┤  ├──R──┤ </td>
<td> (L0>R0)∧(L1>R1)∧L0<R1 </td>
<td> </td>
<td> </td>
<td> 6<_pom/>4 4<_pom/>4 </td>
<td> ¯4 </td>
</tr>

<tr>
<td> ├──L──┼──R──┤ </td>
<td> L1=R0 </td>
<td> (L≤≈R)>L<≈R </td>
<td> L=∀R </td>
<td> 6<_pom/>4 14<_pom/>4 </td>
<td>  5 </td>
</tr>

<tr>
<td> ├──R──┼──L──┤ </td>
<td> L0=R1 </td>
<td> (L≥≈R)>L>≈R </td>
<td> L=∃R </td>
<td> 6<_pom/>4 0<_pom/>2 </td>
<td> ¯5 </td>
</tr>

<tr>
<td> ├──L──┤ ├──R──┤ </td>
<td> L1<R0 </td>
<td> L<≈R </td>
<td> L<∀R </td>
<td> 6<_pom/>4 17<_pom/>5 </td>
<td>  6 </td>
</tr>

<tr>
<td> ├──R──┤ ├──L──┤ </td>
<td> L0>R1 </td>
<td> L>≈R </td>
<td> L>∃R </td>
<td> 6<_pom/>4 0<_pom/>1 </td>
<td> ¯6 </td>
</tr>
</table>

The above table can be interpreted as follows:

<ul>
<li>The first column lists all possible '''Relative Positions''' two finite non-zero Radii (i.e., multi-valued) Balls may assume.</li>
<li>The second column lists the '''Endpoint Relationship''' that defines the corresponding '''Relative Position'''.</li>
<li>The third column lists a '''Defining Expression''' using the '''DoubleTilde''' operator that is TRUE for the corresponding '''Relative Position''' and FALSE for all of the other '''Relative Positions''' in different rows.</li>
<li>The fourth column lists a '''Defining Expression''' using the '''Always''' and '''Possibly''' operators [http://forum.nars2000.org/ball-arithmetic-t614.html#p1803 proposed] by David Rabenhorst.</li>
<li>The fifth column contains '''Sample Data''' that defines the corresponding '''Relative Position'''.</li>
<li>The sixth column gives the value of the '''Distinct Ball''' function <apll>{⍺∘≈⍵}</apll> when executed on the '''Sample Data'''.</li>
</ul>

Note the table is organized such that apart from the first row (whose '''Endpoint Relationship''' is Symmetric in <apll>L</apll> and <apll>R</apll>) the rows occur in six pairs where the '''Endpoint Relationship''' for each row in the pair can be obtained from the other by swapping <apll>L</apll> and <apll>R</apll>. The pairs of rows are ordered roughly by the amount of overlap between the two Balls.

As you can see, we are missing '''Defining Expressions''' for several of the '''Relative Positions'''. I was hoping that this table would be filled in with more detail so it would serve as a measure of whether or not we had a complete set of functions/operators which could be used to tell apart all of the '''Relative Positions''' without having to resort to the '''Distinct Ball''' function <apll>{⍺∘≈⍵}</apll> constructed explicitly for this purpose. Please add to this table if you find a missing expression, or wish to propose whole new functions and/or operators.

If you would like to try a '''Defining Expression''' on the '''Sample Data''', collect the data into a <apll>13</apll> by <apll>2</apll> matrix, put your '''Defining Expression''' into (say) an Anonymous Function and reduce the same data matrix along the last coordinate — if you have a successful function, you should see a <apll>13</apll>-element Boolean vector with a single one in it. For example,

<apll><pre>
sd←6<_pom/>4,⍪6<_pom/>4 5<_pom/>3 7<_pom/>5 7<_pom/>3 5<_pom/>5 6<_pom/>2 6<_pom/>6 8<_pom/>4 4<_pom/>4 14<_pom/>4 0<_pom/>2 17<_pom/>5 0<_pom/>1
{(⍺≤≈⍵)>⍺<≈⍵}/sd
0 0 0 0 0 0 0 0 0 1 0 0 0
{(⍺≥≈⍵)>⍺>≈⍵}/sd
0 0 0 0 0 0 0 0 0 0 1 0 0
{⍺<≈⍵}/sd
0 0 0 0 0 0 0 0 0 0 0 1 0
{⍺>≈⍵}/sd
0 0 0 0 0 0 0 0 0 0 0 0 1
</pre></apll>

Finally, the following example shows that the '''Sample Data''' are all distinct:

<apll><pre>
{⍺∘≈⍵}/sd
0 1 ¯1 2 ¯2 3 ¯3 4 ¯4 5 ¯5 6 ¯6
</pre></apll>

==Hypercomplex Ball Arithmetic==

Ball Arithmetic also extends to all Hypercomplex numbers (Complex, Quaternion, and Octonion) such as <apll>1.2<_pom/><_i/>2.2</apll> which represents a Complex Ball Arithmetic number whose '''Real''' part is <apll>1.2<_pom/></apll> and whose '''Imaginary''' part is <apll>2.2<_pom/></apll>. Note that the display of that number (at <apll>53</apll> bits of precision)

<apll><pre>
1.2<_pom/><_i/>2.2
1.2<_pom/>2.2<_E/>¯16<_i/>2.2<_pom/>4.4<_E/>¯16
</pre></apll>

demonstrates the inexactness of both floating point coefficients.

==Small v. Large Radii==

Ball Arithmetic works best when given small Radii, that is, small relative to the Midpoint. Calculations on large Radii numbers such as <apll>0<_i/>0.5<_pom/>0.5</apll> don't always produce useful results. For example,

<apll><pre>
⎕PP←8 ⍝ To reduce fractional clutter
x←0<_i/>0.5<_pom/>0.5
√x
0.4267767<_pom/>0.57<_i/>0<_pom/>1.4
</pre></apll>

Unfortunately in this case, because the Radius is so large it overwhelms the Midpoint, as can be seen when we Square the Square Root:

<apll><pre>
(√x)*2
0.18213835<_pom/>2.8<_i/>0<_pom/>2.8
x=(√x)*2
0
</pre></apll>

These calculations work better using smaller Radii, such as

<apll><pre>
y←0<_i/>0.5<_pom/>0.5<_E/>¯30
√y
0.5<_pom/>2.5<_E/>¯31<_i/>0.5<_pom/>7.5<_E/>¯31
(√y)*2
0<_pom/>1<_E/>¯30<_i/>0.5<_pom/>1<_E/>¯30
y=(√y)*2
1
</pre></apll>

where Square Root now returns a meaningful result, and its Square returns the original value.

==Mixed Notation==

The notation for Ball Arithmetic may be combined with other Point Notations in many ways, such as

<apll><pre>
1<_p/>1<_pom/>
3.141592653589793<_pom/>4.4<_E/>¯16
1<_x/>1<_pom/>
2.718281828459045<_pom/>4.5<_E/>¯16
</pre></apll>

For more details about combining different forms of Point Notation, see [[Point_Notation|Point Notation]] and in particular the section on [[Point_Notation#Mixed_Notation|Mixed Notation]].

==Implementation==

===Datataypes===

This feature adds four new datatypes to the language: one each for the four Division Algebras. Note that marking any one of the Hypercomplex coefficients as a Ball automatically marks them all as Balls.
{|
|-
|Real||    <apll>23<_pom/></apll> or <apll>1.2<_pom/></apll>
|-
|Complex||    <apll>23<_i/>1.2<_pom/></apll>
|-
|Quaternion||    <apll>23<_i/>1.2<_j/>0.5<_k/>¯3.14<_pom/></apll>
|-
|Octonion||    <apll>23<_i/>1.2<_j/>0.5<_k/>¯3.14<l/>¯3<_ij/>19<_jk/>¯3<_kl/>6<_pom/></apll>
|}

===Primitive Functions===

As described above, this feature adds two new monadic primitive functions '''Contract''' (<apll>≤R</apll>) and its inverse '''Distract''' (<apll>≥R</apll>), along with one new monadic primitive operator '''DoubleTilde''' <apll>f≈ R</apll> and <apll>L f≈ R</apll>.

All primitive and system functions and operators are extended to work with Ball Arithmetic, so you can test out this feature on all algorithms, with the exception of trying to find the Factorial of a multiple-precision/Ball Arithmetic Hypercomplex number, in which case you are temporarily out of luck.

==Workspace Version Number==

For the moment, if you save a workspace which contains a Ball Arithmetic number in a variable, the workspace version number is <apll>0.03</apll>. Otherwise, the version number is <apll>0.02</apll>. Workspaces with Ball Arithmetic Numbers may not be loaded by an earlier version of the interpreter.

==Pending Design Decisions==

<ul>
<li>The direct comparison functions (<apll><≤=≥>≠≡≢</apll>) and indirect ones (<apll>⍳⍸∊⍷⍒⍋</apll>) are currently implemented to compare the Midpoints of the arguments ignoring the Radius — this may change.</li>
<li>The Floor function and its related functions (Ceiling, Residue, GCD, LCM, and Encode) are currently implemented to execute on the Midpoints of the arguments ignoring the Radius — this may change.</li>
<li>We need new primitives to describe the different ways in which Balls can intersect/overlap, possibly through a new monadic operator (perhaps ≈) on the existing direct and indirect comparison and Floor-related functions.</li>
<li>Certain primitives are defined on limited ranges, and sometimes multiple ranges where it might use different algorithms for different ranges (e.g. <apll>L!R</apll>). There are open questions about how to handle a Ball that spans two or more of these ranges as to how to calculate the result.</li>
<li>Infinity and NaN have not been given any special treatment as yet. The results are whatever ARB has chosen which we may need to override.</li>
</ul>

==Download==

This feature is available as an Alpha release — it may be downloaded from the [http://www.nars2000.org/download/binaries/alpha/ Alpha] release directory.

==Acknowledgements==

<ul>
<li>The Arb team has contributed their Ball Arithmetic<ref name="BallArithmetic" /> code which made this feature possible.</li>
<li>David Rabenhorst has contributed many very useful ideas and many hours debugging the implementation.</li>
</ul>

==References==

<references />

Ball Arithmetic

2018-10-26T10:16:06Z

Drab:

==Introduction==

Ball<ref name="BallArithmetic">NARS2000 Wiki, [http://arblib.org Ball Arithmetic]</ref> and Interval Arithmetic<ref name="IntervalArithmetic">Wikipedia, [https://en.wikipedia.org/wiki/Interval_arithmetic Interval Arithmetic]</ref> both refer to the same idea of calculating with arbitrary precision inexact numbers — that is, with a range of floating point values instead of a single representative. This allows such numbers to reflect their inexactness explicitly.

"The most common use is to keep track of and handle rounding errors directly during the calculation and of uncertainties in the knowledge of the exact values of physical and technical parameters. The latter often arise from measurement errors and tolerances for components or '''due to limits on computational accuracy'''. Interval arithmetic also helps find reliable and guaranteed solutions to equations and optimization problems."<ref name="IntervalArithmetic" />

The difference between the two is in how each approach represents its values: Interval Arithmetic uses the two ends of the interval <apll>[a, b]</apll> whereas Ball Arithmetic uses a Ball with a Midpoint and Radius. Each has its advantages. Interval Arithmetic is implemented by several libraries such as MPFI<ref name="MPFI">[http://perso.ens-lyon.fr/nathalie.revol/software.html MPFI]</ref> and Ball Arithmetic is implemented by ARB<ref name="BallArithmetic" />.

This feature uses the ARB<ref name="BallArithmetic" /> library to implement Ball Arithmetic.

==Environment==

In the examples below, unless otherwise specified, the precision of multi-precision floating-point numbers is set to <apll>53</apll> bits via <apll>⎕FPC←53</apll> which is the same precision as standard double-precision numbers. The Printing Precision is set to <apll>17</apll> digits via <apll>⎕PP←17</apll> so as to display the maximum number of significant digits at <apll>53</apll> bits of precision. The Comparison Tolerance is set via <apll>⎕CT←1<_E/>¯10</apll>.

==Calculating With Ball Arithmetic==

If you have ever entered an expression such as

<apll><pre>
÷3
0.3333333333333333
</pre></apll>

undoubtedly, you are aware that the result is inexact. Moreover, as you incorporate into your code the many other functions that produce inexact numbers, you should be wondering "How does the accumulation of inexactness affect my results?".

The answer to that question lies in the province of [https://en.wikipedia.org/wiki/Numerical_analysis Numerical Analysis] where it is described as "the study of algorithms that use numerical approximation". In general, applying the techniques of Numerical Analysis often involves time-consuming analysis. Fortunately, Ball Arithmetic is an easy to use alternative to that. That is, instead of analyzing the algorithm, run the algorithm with all calculations done with Ball Arithmetic.

Here is the same expression in Ball Arithmetic (use Alt-‘q’ or Ctrl-‘q’, depending on your keyboard layout, to enter the <apll><_pom/></apll> symbol):

<apll><pre>
÷3<_pom/>
0.3333333333333333<_pom/>5.6<_E/>¯17
</pre></apll>

Essentially, each Ball Arithmetic number carries with it an explicit Radius that describes the interval which contains the true value. Because Ball Arithmetic data propagates throughout a calculation and is never demoted to a lower datatype, its advantage is that after each operation the resulting Ball reliably includes the true answer. Across an entire algorithm, the final result is a Ball that indicates explicitly how precise it is, thus eliminating the guesswork as to the effect on the algorithm of the inexactness of the input and any other inexact constants used.

For example,

<apll><pre>
sqrt←{⍺←= ⋄ {0.5×⍵+⍺÷⍵}⍣⍺⍨⍵} ⍝ Calculate square root using Newton's method
√2<_pom/>
1.414213562373095<_pom/>2.2<_E/>¯16
sqrt 2<_pom/>
1.414213562373095<_pom/>1.1<_E/>¯15
10 sqrt 2<_pom/>
1.414213562373095<_pom/>2.8<_E/>¯15
100 sqrt 2<_pom/>
1.414213562373095<_pom/>3.3<_E/>¯14
1000 sqrt 2<_pom/>
1.414213562373095<_pom/>3.3<_E/>¯13
10000 sqrt 2<_pom/>
1.414213562373095<_pom/>3.5<_E/>¯12
</pre></apll>

As you can see, after a few iterations, the accuracy of the result doesn't improve. However, because it is an iterative algorithm, the precision worsens as the number of iterations increases. Only by using Ball Arithmetic can you see just how the precision of the result is affected.

To demonstrate with an extreme example, consider the following function:

<apll><pre>
f←{⍪⍵-3√(⍵*3)-⍵*2}
f 10*3×⍳5
0.33344450621439137
0.3333334452472627
0.33333444595336914
0.3348388671875
2.25
</pre></apll>

Without knowing much about this ([https://en.wikipedia.org/wiki/Condition_number ill-conditioned]) algorithm, you have no reason to suspect these results. As it turns out, there's a mathematical proof that as the argument increases to infinity, the algorithm converges to <apll>÷3</apll>. Only when you calculate with Ball Arithmetic do you see that there is a problem:

<apll><pre>
f 10*3×⍳5<_pom/>
0.3334445062153009<_pom/>2.2<_E/>¯12
0.3333334461785853<_pom/>4.3<_E/>¯9
0.333337664604187<_pom/>0.0000062
0.3385009765625<_pom/>0.0085
1.25<_pom/>5
</pre></apll>

As you can see, the precision worsens quickly to the point where the result is meaningless. The last line even says that the answer lies somewhere between <apll>¯3.75</apll> and <apll>6.25</apll>. In other words, that answer is so '''imprecise''' as to be '''inaccurate'''! At higher powers of ten, the answers only worsen.

The good news with Ball Arithmetic is that if you need greater precision, just increase the value of <apll>⎕FPC</apll>. For example,

<apll><pre>
⎕FPC←128 ⋄ ⎕PP←40
sqrt 2<_pom/>
1.41421356237309504880168872420969807857<_pom/>5.9<_E/>¯39
</pre></apll>

In the case of an ill-conditioned algorithm, we just push the problem off a bit, but not very far — the problem is that the algorithm requires infinite precision:

<apll><pre>
f 10*3×⍳10<_pom/>
0.333444506214021971364848444464298690637<_pom/>3.4<_E/>¯35
0.333333444444506172880658466392407138179<_pom/>1.1<_E/>¯31
0.333333333444444444506172839632668381652<_pom/>1.6<_E/>¯28
0.333333333333444444444444547458328123363<_pom/>2.2<_E/>¯25
0.333333333333333444444619474517508476812<_pom/>3.2<_E/>¯22
0.333333333333333333552432522356445687706<_pom/>1.3<_E/>¯19
0.333333333333333491771410805881714622956<_pom/>4.5<_E/>¯16
0.333333333333417414223731611855328083038<_pom/>4.5<_E/>¯13
0.33333333363771089352667331695556640625<_pom/>4.5<_E/>¯10
0.3333334587514400482177734375<_pom/>2.6<_E/>¯7
</pre></apll>

While your floating point code might not be so ill-conditioned, the only way to be sure of its results is one of

<ul>
<li>Hire a Numerical Analyst</li>
<li>Become a Numerical Analyst</li>
<li>Use Ball Arithmetic</li>
</ul>

==Ball Arithmetic Numbers==

A Ball Arithmetic number consist of two parts: a Midpoint and Radius. The Midpoint value represents the middle of the Ball and the Radius represents the distance around the Ball of possible values. There is no way to tell which value in the Ball is the correct answer — you should assume that all values are potential candidates.

The precision at which the Midpoint and Radius are stored depends upon the current value of <apll>⎕FPC</apll>. For example, the following display

<apll><pre>
⎕FPC←53
÷3<_pom/>
0.3333333333333333<_pom/>5.6<_E/>¯17
</pre></apll>

indicates that the number <apll>÷3<_pom/></apll> has a Midpoint of <apll>0.3333333333333333</apll> at <apll>53</apll>-bit precision and a Radius of <apll>5.6<_E/>¯17</apll>. The Midpoint is stored as a variable precision number, and the Radius is stored as a low-precision number.

Alternatively,

<apll><pre>
⎕FPC←128 ⋄ ⎕PP←40
÷3<_pom/>
0.333333333333333333333333333333333333332<_pom/>1.5<_E/>¯39
</pre></apll>

shows the corresponding value at a higher precision and correspondingly narrower Radius.


The approximate number of digits of precision for a given value of <apll>⎕FPC</apll> is about <apll>10⍟2*⎕FPC</apll>. The value of <apll>⎕FPC</apll> for a given required precision is about <apll>⌈2⍟10*digits</apll>.



The current <apll>⎕FPC</apll> precision
should be sufficient to cover both the midpoint and the radius
of all ball numbers being used.
Otherwise,
some results might be misleading.



For example,

<apll>
2v*⎕fpc←128
3.402823669E38
⍟ ⎕← * 1±1e¯1000
2.718281828±1.2E¯38
1±7.3E¯39
</apll>


What happened? -- <apll>⎕FPC</apll> was too small!


<apll>
2v*⎕fpc←4096
1.044388881E1233
⍟ ⎕← * 1±1e¯1000
2.718281828±2.7E¯1000
1±1E¯1000
</apll>

 

A Ball describes a range of '''variable-precision floating point numbers''', and so the expression <apll>23<_pom/>1</apll> represents the set of Real numbers between <apll>22</apll> and <apll>24</apll>, inclusive, '''not''' the three integers <apll>22</apll>, <apll>23</apll>, <apll>24</apll>.

In general, suffix notation on Ball Arithmetic constants should be sufficient for '''almost all purposes'''. That is, the Radius is already correctly calculated for constants such as <apll>23<_pom/></apll> whose Radius is zero, and <apll>1.2<_pom/></apll> whose Radius at <apll>53</apll> bits of precision is on the order of <apll>2*¯53</apll> (because of the inexactness of the Midpoint <apll>1.2</apll>) which displays as <apll>1.2<_pom/>2.2<_E/>¯16</apll>. Although the input notation supports it, the need to override the default Radius should be '''quite rare'''.

==New Functions and Operator==

===Functions To Convert To And From Ball Arithmetic Numbers===

A Ball Arithmetic constant number may be created through notation (as seen above). Otherwise, a Ball Arithmetic number also may be created from its constituent parts under program control with a new primitive function '''Contract''', and separated into its constituent parts with a new primitive function '''Distract'''.

====Contract====

'''Contract''' is a monadic function (<apll>≤R</apll>) whose right argument is a Real array. It may be used to create a Ball Arithmetic number from the one- or two-elements of its Midpoint and (optional) Radius. The Axis Operator may specify any length <apll>1</apll> or <apll>2</apll> coordinate along which the conversion takes place. It converts the first element along the specified dimension to the Midpoint of the result and (if present) the second element to its Radius. If the second element is not present (i.e., length <apll>1</apll> coordinate), the Radius used is calculated from the precision of the Midpoint.

<apll><pre>
≤23
23
≤1.2
1.2<_pom/>1.1<_E/>¯16
≤1.3 1<_E/>¯20
1.3<_pom/>1<_E/>¯20
≤1r2 1<_E/>¯20
0.5<_pom/>1<_E/>¯20
≤○1
3.141592653589793<_pom/>1.1<_E/>¯16
≤(○1) 1<_E/>¯20
3.141592653589793<_pom/>1<_E/>¯20
</pre></apll>

====Distract====

'''Distract''' is a monadic function (<apll>≥R</apll>) whose right argument is a Real array. It may be used to separate out a Ball Arithmetic's Midpoint and Radius into two distinct values. The Axis Operator may specify the place where the new length <apll>2</apll> coordinate is to be inserted. The pair of values along this new coordinate are the Midpoint and Radius of the corresponding number from the right argument. If the number can be represented exactly, the Radius is zero.

<apll><pre>
⎕PP←10
≥1.2
1.2 1.1102230246251565<_E/>¯16
≥1.2<_x/>
6<_r/>5 0
≥1.2<_v/>
1.2 1.110223024625157<_E/>¯16
≥1.2<_pom/>
1.2 2.220446049250313<_E/>¯16
≥1.5
1.5 0
≥1.5<_x/>
3<_r/>2 0
≥1.5<_v/>
1.5 1.110223024625157<_E/>¯16
≥1.5<_pom/>
1.5 0
a←(○1<_pom/>)×⍳4
a,'|',≥a
3.141592654<_pom/>4.4<_E/>¯16 | 3.141592654 4.440892148<_E/>¯16
6.283185307<_pom/>8.9<_E/>¯16 | 6.283185307 8.881784296<_E/>¯16
9.424777961<_pom/>1.3<_E/>¯15 | 9.424777961 1.332267643<_E/>¯15
12.56637061<_pom/>1.8<_E/>¯15 | 12.56637061 1.776356859<_E/>¯15
</pre></apll>

These numbers can be represented exactly:

<apll><pre>
≥⍳4<_pom/>
1 0
2 0
3 0
4 0
≥[2]⍳4<_pom/>
1 0
2 0
3 0
4 0
≥[1]⍳4<_pom/>
1 2 3 4
0 0 0 0
</pre></apll>

====Identities====

The Contract of the Distract of a number both along the same axis, returns the same number. The Distract of the Contract of a number both along the same length two axis, returns the same number.

===Operator To Provide Arbitrary Precision Specific Derived Functions===

The DoubleTilde monadic operator (<apll>≈</apll>, Alt-'Q' or Alt-shift-'q') takes a primitive function as its left operand and produces an arbitrary-precision-specific derived function.

This operator should be considered preliminary.

====Conversion Between Balls and Intervals====

At the moment, the monadic derived function is defined for two primitive functions.

Because Ball Arithmetic has a close relationship with Interval Arithmetic, it may be desirable to convert between the two forms: Endpoint to Midpoint-Radius and its inverse. The resulting derived function is not scalar.

For example,

<apll><pre>
x←0<_pom/>∞ 1<_pom/>2
≥≈x
¯∞ ∞
¯1 3
≥[1]≈x
¯∞ ¯1
∞ 3
≥[2]≈x
¯∞ ∞
¯1 3
≤≈≥≈x
0<_pom/>∞ 1<_pom/>2
≤≈≥[2]≈x
0<_pom/>∞ 1<_pom/>2
≤[2]≈≥≈x
0<_pom/>∞ 1<_pom/>2
≤[2]≈≥[2]≈x
0<_pom/>∞ 1<_pom/>2
≤[1]≈≥[1]≈x
0<_pom/>∞ 1<_pom/>2
</pre></apll>

The domain of these functions is limited to Real numbers. Hypercomplex numbers signal a <apll>DOMAIN ERROR</apll>. To see the Endpoints of a Hypercomplex number, first '''Dilate''' it so as to separate out the coefficients and then apply <apll>≥≈</apll> to extract the Endpoints as in <apll>≥≈>R</apll>.

For example,

<apll><pre> ≥≈ 2±0i0±1
DOMAIN ERROR
≥≈ 2±0i0±1
∧
≥≈>2±0i0±1
2 2
¯1 1
<≤≈≥≈>2±0i0±1
2±0i0±1
</pre></apll>

====Union and Intersection of Balls====

The left operand of the DoubleTilde operator is defined for the '''Union''' and '''Intersection''' functions (<apll>∪∩</apll>) as follows:

Because Ball arithmetic numbers represent a range of floating point numbers, at times it is desirable to produce their Union or Intersection.

For example,

<apll><pre>
x←0.5<pom/>0.5 ⍝ The interval [0,1]
x∪≈x+2 ⍝ The interval [0,1] ∪ [2,3] ←→ [0,3]
1.5<pom/>1.5
x∪≈x+1 ⍝ The interval [0,1] ∪ [1,2] ←→ [0,2]
1<pom/>1
x∩≈x+1 ⍝ The interval [0,1] ∩ [1,2] ←→ [1,1]
1
</pre></apll>

If two Balls have no points in common, the result is a NaN (<apll>∅</apll>):

<apll><pre>
x∩≈x+1.1 ⍝ The interval [0,1] ∩ [1.1,2.1] ←→ ∅
∅
</pre></apll>

The identity element for <apll>∩≈</apll> is <apll>0±∞</apll>.

The derived functions are scalar.

====Comparison Functions====

The left operand of the DoubleTilde operator is defined for the six comparison functions (<apll><≤=≥>≠</apll>) as follows:

For <apll>f</apll> as one of the above six functions, <apll>Lf≈R</apll> is TRUE iff for all <apll>t∊L</apll> and <apll>u∊R</apll>, then <apll>  t f u  </apll> is TRUE. These comparisons are all done with <apll>⎕CT←0</apll> because, in a sense, Ball comparisons carry their own Fuzz.

For example,

<apll><pre> (L R)←4<_pom/>1 6<_pom/>1
L<≈R
0
L≤≈R
1
</pre></apll>

The derived functions are scalar.

====Containing Functions====

The left operand of the DoubleTilde operator is defined for two containing functions (<apll>⊂⊃</apll>) as follows:

The expression <apll>L⊂≈R</apll> is TRUE iff for all <apll>t∊L</apll> it is TRUE that <apll>t∊R</apll> and that the containment is strict — that is, none of the endpoints match.

The expression <apll>L⊃≈R</apll> is TRUE iff <apll>R⊂≈L</apll> is TRUE.

For example,

<apll><pre>
(L R)←6±4 6±2
L⊃≈R
1
R⊂≈L
1
</pre></apll>

The derived functions are scalar.

====Positional Relationships====

For <apll>L</apll> and <apll>R</apll> as finite intervals with non-zero Radii (i.e., multi-valued), here are the possible '''Relative Positions''' between the two Intervals, where <apll>L0</apll>, <apll>L1</apll> and <apll>R0</apll>, <apll>R1</apll> are the endpoints of <apll>L</apll> and <apll>R</apll>, respectively — <apll>(L0 L1)←≥≈L</apll> and <apll>(R0 R1)←≥≈R</apll>:

<table border="1" cellspacing="0" rowspacing="5" style="font-family:APL385, SImPL;">
<th>Relative Positions</th>
<th>Endpoint Relationships</th>
<th colspan="2">Defining Expressions</th>
<th>Sample Data L   R</th>
<th>Distinct Ball {⍺∘≈⍵}</th>

<tr>
<td> ├──L──┤  ├──R──┤ </td>
<td> L0 L1≡R0 R1 </td>
<td> L=≈R </td>
<td> </td>
<td> 6<_pom/>4 6<_pom/>4 </td>
<td>  0 </td>
</tr>

<tr>
<td> ├──L──┤  ├──R─┤ </td>
<td> (L0=R0)∧L1>R1 </td>
<td> </td>
<td> </td>
<td> 6<_pom/>4 5<_pom/>3 </td>
<td>  1 </td>
</tr>

<tr>
<td> ├──L─┤  ├──R──┤ </td>
<td> (L0=R0)∧L1<R1 </td>
<td> </td>
<td> </td>
<td> 6<_pom/>4 7<_pom/>5 </td>
<td> ¯1 </td>
</tr>

<tr>
<td> ├──L──┤   ├──R─┤ </td>
<td> (L0<R0)∧L1=R1 </td>
<td> </td>
<td> </td>
<td> 6<_pom/>4 7<_pom/>3 </td>
<td>  2 </td>
</tr>

<tr>
<td>  ├──L─┤  ├──R──┤ </td>
<td> (L0>R0)∧L1=R1 </td>
<td> </td>
<td> </td>
<td> 6<_pom/>4 5<_pom/>5 </td>
<td> ¯2 </td>
</tr>

<tr>
<td> ├──L──┤   ├─R─┤ </td>
<td> (L0<R0)∧L1>R1 </td>
<td> L⊃≈R </td>
<td> </td>
<td> 6<_pom/>4 6<_pom/>2 </td>
<td>  3 </td>
</tr>

<tr>
<td>  ├─L─┤  ├──R──┤ </td>
<td> (L0>R0)∧L1<R1 </td>
<td> L⊂≈R </td>
<td> </td>
<td> 6<_pom/>4 6<_pom/>6 </td>
<td> ¯3 </td>
</tr>

<tr>
<td> ├──L──┤   ├──R──┤ </td>
<td> (L0<R0)∧(L1<R1)∧L1>R0 </td>
<td> </td>
<td> </td>
<td> 6<_pom/>4 8<_pom/>4 </td>
<td>  4 </td>
</tr>

<tr>
<td>  ├──L──┤  ├──R──┤ </td>
<td> (L0>R0)∧(L1>R1)∧L0<R1 </td>
<td> </td>
<td> </td>
<td> 6<_pom/>4 4<_pom/>4 </td>
<td> ¯4 </td>
</tr>

<tr>
<td> ├──L──┼──R──┤ </td>
<td> L1=R0 </td>
<td> (L≤≈R)>L<≈R </td>
<td> L=∀R </td>
<td> 6<_pom/>4 14<_pom/>4 </td>
<td>  5 </td>
</tr>

<tr>
<td> ├──R──┼──L──┤ </td>
<td> L0=R1 </td>
<td> (L≥≈R)>L>≈R </td>
<td> L=∃R </td>
<td> 6<_pom/>4 0<_pom/>2 </td>
<td> ¯5 </td>
</tr>

<tr>
<td> ├──L──┤ ├──R──┤ </td>
<td> L1<R0 </td>
<td> L<≈R </td>
<td> L<∀R </td>
<td> 6<_pom/>4 17<_pom/>5 </td>
<td>  6 </td>
</tr>

<tr>
<td> ├──R──┤ ├──L──┤ </td>
<td> L0>R1 </td>
<td> L>≈R </td>
<td> L>∃R </td>
<td> 6<_pom/>4 0<_pom/>1 </td>
<td> ¯6 </td>
</tr>
</table>

The above table can be interpreted as follows:

<ul>
<li>The first column lists all possible '''Relative Positions''' two finite non-zero Radii (i.e., multi-valued) Balls may assume.</li>
<li>The second column lists the '''Endpoint Relationship''' that defines the corresponding '''Relative Position'''.</li>
<li>The third column lists a '''Defining Expression''' using the '''DoubleTilde''' operator that is TRUE for the corresponding '''Relative Position''' and FALSE for all of the other '''Relative Positions''' in different rows.</li>
<li>The fourth column lists a '''Defining Expression''' using the '''Always''' and '''Possibly''' operators [http://forum.nars2000.org/ball-arithmetic-t614.html#p1803 proposed] by David Rabenhorst.</li>
<li>The fifth column contains '''Sample Data''' that defines the corresponding '''Relative Position'''.</li>
<li>The sixth column gives the value of the '''Distinct Ball''' function <apll>{⍺∘≈⍵}</apll> when executed on the '''Sample Data'''.</li>
</ul>

Note the table is organized such that apart from the first row (whose '''Endpoint Relationship''' is Symmetric in <apll>L</apll> and <apll>R</apll>) the rows occur in six pairs where the '''Endpoint Relationship''' for each row in the pair can be obtained from the other by swapping <apll>L</apll> and <apll>R</apll>. The pairs of rows are ordered roughly by the amount of overlap between the two Balls.

As you can see, we are missing '''Defining Expressions''' for several of the '''Relative Positions'''. I was hoping that this table would be filled in with more detail so it would serve as a measure of whether or not we had a complete set of functions/operators which could be used to tell apart all of the '''Relative Positions''' without having to resort to the '''Distinct Ball''' function <apll>{⍺∘≈⍵}</apll> constructed explicitly for this purpose. Please add to this table if you find a missing expression, or wish to propose whole new functions and/or operators.

If you would like to try a '''Defining Expression''' on the '''Sample Data''', collect the data into a <apll>13</apll> by <apll>2</apll> matrix, put your '''Defining Expression''' into (say) an Anonymous Function and reduce the same data matrix along the last coordinate — if you have a successful function, you should see a <apll>13</apll>-element Boolean vector with a single one in it. For example,

<apll><pre>
sd←6<_pom/>4,⍪6<_pom/>4 5<_pom/>3 7<_pom/>5 7<_pom/>3 5<_pom/>5 6<_pom/>2 6<_pom/>6 8<_pom/>4 4<_pom/>4 14<_pom/>4 0<_pom/>2 17<_pom/>5 0<_pom/>1
{(⍺≤≈⍵)>⍺<≈⍵}/sd
0 0 0 0 0 0 0 0 0 1 0 0 0
{(⍺≥≈⍵)>⍺>≈⍵}/sd
0 0 0 0 0 0 0 0 0 0 1 0 0
{⍺<≈⍵}/sd
0 0 0 0 0 0 0 0 0 0 0 1 0
{⍺>≈⍵}/sd
0 0 0 0 0 0 0 0 0 0 0 0 1
</pre></apll>

Finally, the following example shows that the '''Sample Data''' are all distinct:

<apll><pre>
{⍺∘≈⍵}/sd
0 1 ¯1 2 ¯2 3 ¯3 4 ¯4 5 ¯5 6 ¯6
</pre></apll>

==Hypercomplex Ball Arithmetic==

Ball Arithmetic also extends to all Hypercomplex numbers (Complex, Quaternion, and Octonion) such as <apll>1.2<_pom/><_i/>2.2</apll> which represents a Complex Ball Arithmetic number whose '''Real''' part is <apll>1.2<_pom/></apll> and whose '''Imaginary''' part is <apll>2.2<_pom/></apll>. Note that the display of that number (at <apll>53</apll> bits of precision)

<apll><pre>
1.2<_pom/><_i/>2.2
1.2<_pom/>2.2<_E/>¯16<_i/>2.2<_pom/>4.4<_E/>¯16
</pre></apll>

demonstrates the inexactness of both floating point coefficients.

==Small v. Large Radii==

Ball Arithmetic works best when given small Radii, that is, small relative to the Midpoint. Calculations on large Radii numbers such as <apll>0<_i/>0.5<_pom/>0.5</apll> don't always produce useful results. For example,

<apll><pre>
⎕PP←8 ⍝ To reduce fractional clutter
x←0<_i/>0.5<_pom/>0.5
√x
0.4267767<_pom/>0.57<_i/>0<_pom/>1.4
</pre></apll>

Unfortunately in this case, because the Radius is so large it overwhelms the Midpoint, as can be seen when we Square the Square Root:

<apll><pre>
(√x)*2
0.18213835<_pom/>2.8<_i/>0<_pom/>2.8
x=(√x)*2
0
</pre></apll>

These calculations work better using smaller Radii, such as

<apll><pre>
y←0<_i/>0.5<_pom/>0.5<_E/>¯30
√y
0.5<_pom/>2.5<_E/>¯31<_i/>0.5<_pom/>7.5<_E/>¯31
(√y)*2
0<_pom/>1<_E/>¯30<_i/>0.5<_pom/>1<_E/>¯30
y=(√y)*2
1
</pre></apll>

where Square Root now returns a meaningful result, and its Square returns the original value.

==Mixed Notation==

The notation for Ball Arithmetic may be combined with other Point Notations in many ways, such as

<apll><pre>
1<_p/>1<_pom/>
3.141592653589793<_pom/>4.4<_E/>¯16
1<_x/>1<_pom/>
2.718281828459045<_pom/>4.5<_E/>¯16
</pre></apll>

For more details about combining different forms of Point Notation, see [[Point_Notation|Point Notation]] and in particular the section on [[Point_Notation#Mixed_Notation|Mixed Notation]].

==Implementation==

===Datataypes===

This feature adds four new datatypes to the language: one each for the four Division Algebras. Note that marking any one of the Hypercomplex coefficients as a Ball automatically marks them all as Balls.
{|
|-
|Real||    <apll>23<_pom/></apll> or <apll>1.2<_pom/></apll>
|-
|Complex||    <apll>23<_i/>1.2<_pom/></apll>
|-
|Quaternion||    <apll>23<_i/>1.2<_j/>0.5<_k/>¯3.14<_pom/></apll>
|-
|Octonion||    <apll>23<_i/>1.2<_j/>0.5<_k/>¯3.14<l/>¯3<_ij/>19<_jk/>¯3<_kl/>6<_pom/></apll>
|}

===Primitive Functions===

As described above, this feature adds two new monadic primitive functions '''Contract''' (<apll>≤R</apll>) and its inverse '''Distract''' (<apll>≥R</apll>), along with one new monadic primitive operator '''DoubleTilde''' <apll>f≈ R</apll> and <apll>L f≈ R</apll>.

All primitive and system functions and operators are extended to work with Ball Arithmetic, so you can test out this feature on all algorithms, with the exception of trying to find the Factorial of a multiple-precision/Ball Arithmetic Hypercomplex number, in which case you are temporarily out of luck.

==Workspace Version Number==

For the moment, if you save a workspace which contains a Ball Arithmetic number in a variable, the workspace version number is <apll>0.03</apll>. Otherwise, the version number is <apll>0.02</apll>. Workspaces with Ball Arithmetic Numbers may not be loaded by an earlier version of the interpreter.

==Pending Design Decisions==

<ul>
<li>The direct comparison functions (<apll><≤=≥>≠≡≢</apll>) and indirect ones (<apll>⍳⍸∊⍷⍒⍋</apll>) are currently implemented to compare the Midpoints of the arguments ignoring the Radius — this may change.</li>
<li>The Floor function and its related functions (Ceiling, Residue, GCD, LCM, and Encode) are currently implemented to execute on the Midpoints of the arguments ignoring the Radius — this may change.</li>
<li>We need new primitives to describe the different ways in which Balls can intersect/overlap, possibly through a new monadic operator (perhaps ≈) on the existing direct and indirect comparison and Floor-related functions.</li>
<li>Certain primitives are defined on limited ranges, and sometimes multiple ranges where it might use different algorithms for different ranges (e.g. <apll>L!R</apll>). There are open questions about how to handle a Ball that spans two or more of these ranges as to how to calculate the result.</li>
<li>Infinity and NaN have not been given any special treatment as yet. The results are whatever ARB has chosen which we may need to override.</li>
</ul>

==Download==

This feature is available as an Alpha release — it may be downloaded from the [http://www.nars2000.org/download/binaries/alpha/ Alpha] release directory.

==Acknowledgements==

<ul>
<li>The Arb team has contributed their Ball Arithmetic<ref name="BallArithmetic" /> code which made this feature possible.</li>
<li>David Rabenhorst has contributed many very useful ideas and many hours debugging the implementation.</li>
</ul>

==References==

<references />

Ball Arithmetic

2018-10-26T10:14:23Z

Drab:

Ball Arithmetic

2018-10-26T10:12:44Z

Drab:

==Introduction==

Ball<ref name="BallArithmetic">NARS2000 Wiki, [http://arblib.org Ball Arithmetic]</ref> and Interval Arithmetic<ref name="IntervalArithmetic">Wikipedia, [https://en.wikipedia.org/wiki/Interval_arithmetic Interval Arithmetic]</ref> both refer to the same idea of calculating with arbitrary precision inexact numbers — that is, with a range of floating point values instead of a single representative. This allows such numbers to reflect their inexactness explicitly.

"The most common use is to keep track of and handle rounding errors directly during the calculation and of uncertainties in the knowledge of the exact values of physical and technical parameters. The latter often arise from measurement errors and tolerances for components or '''due to limits on computational accuracy'''. Interval arithmetic also helps find reliable and guaranteed solutions to equations and optimization problems."<ref name="IntervalArithmetic" />

The difference between the two is in how each approach represents its values: Interval Arithmetic uses the two ends of the interval <apll>[a, b]</apll> whereas Ball Arithmetic uses a Ball with a Midpoint and Radius. Each has its advantages. Interval Arithmetic is implemented by several libraries such as MPFI<ref name="MPFI">[http://perso.ens-lyon.fr/nathalie.revol/software.html MPFI]</ref> and Ball Arithmetic is implemented by ARB<ref name="BallArithmetic" />.

This feature uses the ARB<ref name="BallArithmetic" /> library to implement Ball Arithmetic.

==Environment==

In the examples below, unless otherwise specified, the precision of multi-precision floating-point numbers is set to <apll>53</apll> bits via <apll>⎕FPC←53</apll> which is the same precision as standard double-precision numbers. The Printing Precision is set to <apll>17</apll> digits via <apll>⎕PP←17</apll> so as to display the maximum number of significant digits at <apll>53</apll> bits of precision. The Comparison Tolerance is set via <apll>⎕CT←1<_E/>¯10</apll>.

==Calculating With Ball Arithmetic==

If you have ever entered an expression such as

<apll><pre>
÷3
0.3333333333333333
</pre></apll>

undoubtedly, you are aware that the result is inexact. Moreover, as you incorporate into your code the many other functions that produce inexact numbers, you should be wondering "How does the accumulation of inexactness affect my results?".

The answer to that question lies in the province of [https://en.wikipedia.org/wiki/Numerical_analysis Numerical Analysis] where it is described as "the study of algorithms that use numerical approximation". In general, applying the techniques of Numerical Analysis often involves time-consuming analysis. Fortunately, Ball Arithmetic is an easy to use alternative to that. That is, instead of analyzing the algorithm, run the algorithm with all calculations done with Ball Arithmetic.

Here is the same expression in Ball Arithmetic (use Alt-‘q’ or Ctrl-‘q’, depending on your keyboard layout, to enter the <apll><_pom/></apll> symbol):

<apll><pre>
÷3<_pom/>
0.3333333333333333<_pom/>5.6<_E/>¯17
</pre></apll>

Essentially, each Ball Arithmetic number carries with it an explicit Radius that describes the interval which contains the true value. Because Ball Arithmetic data propagates throughout a calculation and is never demoted to a lower datatype, its advantage is that after each operation the resulting Ball reliably includes the true answer. Across an entire algorithm, the final result is a Ball that indicates explicitly how precise it is, thus eliminating the guesswork as to the effect on the algorithm of the inexactness of the input and any other inexact constants used.

For example,

<apll><pre>
sqrt←{⍺←= ⋄ {0.5×⍵+⍺÷⍵}⍣⍺⍨⍵} ⍝ Calculate square root using Newton's method
√2<_pom/>
1.414213562373095<_pom/>2.2<_E/>¯16
sqrt 2<_pom/>
1.414213562373095<_pom/>1.1<_E/>¯15
10 sqrt 2<_pom/>
1.414213562373095<_pom/>2.8<_E/>¯15
100 sqrt 2<_pom/>
1.414213562373095<_pom/>3.3<_E/>¯14
1000 sqrt 2<_pom/>
1.414213562373095<_pom/>3.3<_E/>¯13
10000 sqrt 2<_pom/>
1.414213562373095<_pom/>3.5<_E/>¯12
</pre></apll>

As you can see, after a few iterations, the accuracy of the result doesn't improve. However, because it is an iterative algorithm, the precision worsens as the number of iterations increases. Only by using Ball Arithmetic can you see just how the precision of the result is affected.

To demonstrate with an extreme example, consider the following function:

<apll><pre>
f←{⍪⍵-3√(⍵*3)-⍵*2}
f 10*3×⍳5
0.33344450621439137
0.3333334452472627
0.33333444595336914
0.3348388671875
2.25
</pre></apll>

Without knowing much about this ([https://en.wikipedia.org/wiki/Condition_number ill-conditioned]) algorithm, you have no reason to suspect these results. As it turns out, there's a mathematical proof that as the argument increases to infinity, the algorithm converges to <apll>÷3</apll>. Only when you calculate with Ball Arithmetic do you see that there is a problem:

<apll><pre>
f 10*3×⍳5<_pom/>
0.3334445062153009<_pom/>2.2<_E/>¯12
0.3333334461785853<_pom/>4.3<_E/>¯9
0.333337664604187<_pom/>0.0000062
0.3385009765625<_pom/>0.0085
1.25<_pom/>5
</pre></apll>

As you can see, the precision worsens quickly to the point where the result is meaningless. The last line even says that the answer lies somewhere between <apll>¯3.75</apll> and <apll>6.25</apll>. In other words, that answer is so '''imprecise''' as to be '''inaccurate'''! At higher powers of ten, the answers only worsen.

The good news with Ball Arithmetic is that if you need greater precision, just increase the value of <apll>⎕FPC</apll>. For example,

<apll><pre>
⎕FPC←128 ⋄ ⎕PP←40
sqrt 2<_pom/>
1.41421356237309504880168872420969807857<_pom/>5.9<_E/>¯39
</pre></apll>

In the case of an ill-conditioned algorithm, we just push the problem off a bit, but not very far — the problem is that the algorithm requires infinite precision:

<apll><pre>
f 10*3×⍳10<_pom/>
0.333444506214021971364848444464298690637<_pom/>3.4<_E/>¯35
0.333333444444506172880658466392407138179<_pom/>1.1<_E/>¯31
0.333333333444444444506172839632668381652<_pom/>1.6<_E/>¯28
0.333333333333444444444444547458328123363<_pom/>2.2<_E/>¯25
0.333333333333333444444619474517508476812<_pom/>3.2<_E/>¯22
0.333333333333333333552432522356445687706<_pom/>1.3<_E/>¯19
0.333333333333333491771410805881714622956<_pom/>4.5<_E/>¯16
0.333333333333417414223731611855328083038<_pom/>4.5<_E/>¯13
0.33333333363771089352667331695556640625<_pom/>4.5<_E/>¯10
0.3333334587514400482177734375<_pom/>2.6<_E/>¯7
</pre></apll>

While your floating point code might not be so ill-conditioned, the only way to be sure of its results is one of

<ul>
<li>Hire a Numerical Analyst</li>
<li>Become a Numerical Analyst</li>
<li>Use Ball Arithmetic</li>
</ul>

==Ball Arithmetic Numbers==

A Ball Arithmetic number consist of two parts: a Midpoint and Radius. The Midpoint value represents the middle of the Ball and the Radius represents the distance around the Ball of possible values. There is no way to tell which value in the Ball is the correct answer — you should assume that all values are potential candidates.

The precision at which the Midpoint and Radius are stored depends upon the current value of <apll>⎕FPC</apll>. For example, the following display

<apll><pre>
⎕FPC←53
÷3<_pom/>
0.3333333333333333<_pom/>5.6<_E/>¯17
</pre></apll>

indicates that the number <apll>÷3<_pom/></apll> has a Midpoint of <apll>0.3333333333333333</apll> at <apll>53</apll>-bit precision and a Radius of <apll>5.6<_E/>¯17</apll>. The Midpoint is stored as a variable precision number, and the Radius is stored as a low-precision number.

Alternatively,

<apll><pre>
⎕FPC←128 ⋄ ⎕PP←40
÷3<_pom/>
0.333333333333333333333333333333333333332<_pom/>1.5<_E/>¯39
</pre></apll>

shows the corresponding value at a higher precision and correspondingly narrower Radius.


The approximate number of digits of precision for a given value of <apll>⎕FPC</apll> is about <apll>10⍟2*⎕FPC</apll>. The value of <apll>⎕FPC</apll> for a given required precision is about <apll>⌈2⍟10*digits</apll>.



The current <apll>⎕FPC</apll> precision
should be sufficient to cover both the midpoint and the radius
of all ball numbers being used.
Otherwise,
some results might be misleading.



For example,
<apll>
2v*⎕fpc←128
3.402823669E38
⍟ ⎕← * 1±1e¯1000
2.718281828±1.2E¯38
1±7.3E¯39
</apll>



What happened? -- <apll>⎕FPC</apll> was too small!


<apll>
2v*⎕fpc←4096
1.044388881E1233
⍟ ⎕← * 1±1e¯1000
2.718281828±2.7E¯1000
1±1E¯1000
</apll>

 

A Ball describes a range of '''variable-precision floating point numbers''', and so the expression <apll>23<_pom/>1</apll> represents the set of Real numbers between <apll>22</apll> and <apll>24</apll>, inclusive, '''not''' the three integers <apll>22</apll>, <apll>23</apll>, <apll>24</apll>.

In general, suffix notation on Ball Arithmetic constants should be sufficient for '''almost all purposes'''. That is, the Radius is already correctly calculated for constants such as <apll>23<_pom/></apll> whose Radius is zero, and <apll>1.2<_pom/></apll> whose Radius at <apll>53</apll> bits of precision is on the order of <apll>2*¯53</apll> (because of the inexactness of the Midpoint <apll>1.2</apll>) which displays as <apll>1.2<_pom/>2.2<_E/>¯16</apll>. Although the input notation supports it, the need to override the default Radius should be '''quite rare'''.

==New Functions and Operator==

===Functions To Convert To And From Ball Arithmetic Numbers===

A Ball Arithmetic constant number may be created through notation (as seen above). Otherwise, a Ball Arithmetic number also may be created from its constituent parts under program control with a new primitive function '''Contract''', and separated into its constituent parts with a new primitive function '''Distract'''.

====Contract====

'''Contract''' is a monadic function (<apll>≤R</apll>) whose right argument is a Real array. It may be used to create a Ball Arithmetic number from the one- or two-elements of its Midpoint and (optional) Radius. The Axis Operator may specify any length <apll>1</apll> or <apll>2</apll> coordinate along which the conversion takes place. It converts the first element along the specified dimension to the Midpoint of the result and (if present) the second element to its Radius. If the second element is not present (i.e., length <apll>1</apll> coordinate), the Radius used is calculated from the precision of the Midpoint.

<apll><pre>
≤23
23
≤1.2
1.2<_pom/>1.1<_E/>¯16
≤1.3 1<_E/>¯20
1.3<_pom/>1<_E/>¯20
≤1r2 1<_E/>¯20
0.5<_pom/>1<_E/>¯20
≤○1
3.141592653589793<_pom/>1.1<_E/>¯16
≤(○1) 1<_E/>¯20
3.141592653589793<_pom/>1<_E/>¯20
</pre></apll>

====Distract====

'''Distract''' is a monadic function (<apll>≥R</apll>) whose right argument is a Real array. It may be used to separate out a Ball Arithmetic's Midpoint and Radius into two distinct values. The Axis Operator may specify the place where the new length <apll>2</apll> coordinate is to be inserted. The pair of values along this new coordinate are the Midpoint and Radius of the corresponding number from the right argument. If the number can be represented exactly, the Radius is zero.

<apll><pre>
⎕PP←10
≥1.2
1.2 1.1102230246251565<_E/>¯16
≥1.2<_x/>
6<_r/>5 0
≥1.2<_v/>
1.2 1.110223024625157<_E/>¯16
≥1.2<_pom/>
1.2 2.220446049250313<_E/>¯16
≥1.5
1.5 0
≥1.5<_x/>
3<_r/>2 0
≥1.5<_v/>
1.5 1.110223024625157<_E/>¯16
≥1.5<_pom/>
1.5 0
a←(○1<_pom/>)×⍳4
a,'|',≥a
3.141592654<_pom/>4.4<_E/>¯16 | 3.141592654 4.440892148<_E/>¯16
6.283185307<_pom/>8.9<_E/>¯16 | 6.283185307 8.881784296<_E/>¯16
9.424777961<_pom/>1.3<_E/>¯15 | 9.424777961 1.332267643<_E/>¯15
12.56637061<_pom/>1.8<_E/>¯15 | 12.56637061 1.776356859<_E/>¯15
</pre></apll>

These numbers can be represented exactly:

<apll><pre>
≥⍳4<_pom/>
1 0
2 0
3 0
4 0
≥[2]⍳4<_pom/>
1 0
2 0
3 0
4 0
≥[1]⍳4<_pom/>
1 2 3 4
0 0 0 0
</pre></apll>

====Identities====

The Contract of the Distract of a number both along the same axis, returns the same number. The Distract of the Contract of a number both along the same length two axis, returns the same number.

===Operator To Provide Arbitrary Precision Specific Derived Functions===

The DoubleTilde monadic operator (<apll>≈</apll>, Alt-'Q' or Alt-shift-'q') takes a primitive function as its left operand and produces an arbitrary-precision-specific derived function.

This operator should be considered preliminary.

====Conversion Between Balls and Intervals====

At the moment, the monadic derived function is defined for two primitive functions.

Because Ball Arithmetic has a close relationship with Interval Arithmetic, it may be desirable to convert between the two forms: Endpoint to Midpoint-Radius and its inverse. The resulting derived function is not scalar.

For example,

<apll><pre>
x←0<_pom/>∞ 1<_pom/>2
≥≈x
¯∞ ∞
¯1 3
≥[1]≈x
¯∞ ¯1
∞ 3
≥[2]≈x
¯∞ ∞
¯1 3
≤≈≥≈x
0<_pom/>∞ 1<_pom/>2
≤≈≥[2]≈x
0<_pom/>∞ 1<_pom/>2
≤[2]≈≥≈x
0<_pom/>∞ 1<_pom/>2
≤[2]≈≥[2]≈x
0<_pom/>∞ 1<_pom/>2
≤[1]≈≥[1]≈x
0<_pom/>∞ 1<_pom/>2
</pre></apll>

The domain of these functions is limited to Real numbers. Hypercomplex numbers signal a <apll>DOMAIN ERROR</apll>. To see the Endpoints of a Hypercomplex number, first '''Dilate''' it so as to separate out the coefficients and then apply <apll>≥≈</apll> to extract the Endpoints as in <apll>≥≈>R</apll>.

For example,

<apll><pre> ≥≈ 2±0i0±1
DOMAIN ERROR
≥≈ 2±0i0±1
∧
≥≈>2±0i0±1
2 2
¯1 1
<≤≈≥≈>2±0i0±1
2±0i0±1
</pre></apll>

====Union and Intersection of Balls====

The left operand of the DoubleTilde operator is defined for the '''Union''' and '''Intersection''' functions (<apll>∪∩</apll>) as follows:

Because Ball arithmetic numbers represent a range of floating point numbers, at times it is desirable to produce their Union or Intersection.

For example,

<apll><pre>
x←0.5<pom/>0.5 ⍝ The interval [0,1]
x∪≈x+2 ⍝ The interval [0,1] ∪ [2,3] ←→ [0,3]
1.5<pom/>1.5
x∪≈x+1 ⍝ The interval [0,1] ∪ [1,2] ←→ [0,2]
1<pom/>1
x∩≈x+1 ⍝ The interval [0,1] ∩ [1,2] ←→ [1,1]
1
</pre></apll>

If two Balls have no points in common, the result is a NaN (<apll>∅</apll>):

<apll><pre>
x∩≈x+1.1 ⍝ The interval [0,1] ∩ [1.1,2.1] ←→ ∅
∅
</pre></apll>

The identity element for <apll>∩≈</apll> is <apll>0±∞</apll>.

The derived functions are scalar.

====Comparison Functions====

The left operand of the DoubleTilde operator is defined for the six comparison functions (<apll><≤=≥>≠</apll>) as follows:

For <apll>f</apll> as one of the above six functions, <apll>Lf≈R</apll> is TRUE iff for all <apll>t∊L</apll> and <apll>u∊R</apll>, then <apll>  t f u  </apll> is TRUE. These comparisons are all done with <apll>⎕CT←0</apll> because, in a sense, Ball comparisons carry their own Fuzz.

For example,

<apll><pre> (L R)←4<_pom/>1 6<_pom/>1
L<≈R
0
L≤≈R
1
</pre></apll>

The derived functions are scalar.

====Containing Functions====

The left operand of the DoubleTilde operator is defined for two containing functions (<apll>⊂⊃</apll>) as follows:

The expression <apll>L⊂≈R</apll> is TRUE iff for all <apll>t∊L</apll> it is TRUE that <apll>t∊R</apll> and that the containment is strict — that is, none of the endpoints match.

The expression <apll>L⊃≈R</apll> is TRUE iff <apll>R⊂≈L</apll> is TRUE.

For example,

<apll><pre>
(L R)←6±4 6±2
L⊃≈R
1
R⊂≈L
1
</pre></apll>

The derived functions are scalar.

====Positional Relationships====

For <apll>L</apll> and <apll>R</apll> as finite intervals with non-zero Radii (i.e., multi-valued), here are the possible '''Relative Positions''' between the two Intervals, where <apll>L0</apll>, <apll>L1</apll> and <apll>R0</apll>, <apll>R1</apll> are the endpoints of <apll>L</apll> and <apll>R</apll>, respectively — <apll>(L0 L1)←≥≈L</apll> and <apll>(R0 R1)←≥≈R</apll>:

<table border="1" cellspacing="0" rowspacing="5" style="font-family:APL385, SImPL;">
<th>Relative Positions</th>
<th>Endpoint Relationships</th>
<th colspan="2">Defining Expressions</th>
<th>Sample Data L   R</th>
<th>Distinct Ball {⍺∘≈⍵}</th>

<tr>
<td> ├──L──┤  ├──R──┤ </td>
<td> L0 L1≡R0 R1 </td>
<td> L=≈R </td>
<td> </td>
<td> 6<_pom/>4 6<_pom/>4 </td>
<td>  0 </td>
</tr>

<tr>
<td> ├──L──┤  ├──R─┤ </td>
<td> (L0=R0)∧L1>R1 </td>
<td> </td>
<td> </td>
<td> 6<_pom/>4 5<_pom/>3 </td>
<td>  1 </td>
</tr>

<tr>
<td> ├──L─┤  ├──R──┤ </td>
<td> (L0=R0)∧L1<R1 </td>
<td> </td>
<td> </td>
<td> 6<_pom/>4 7<_pom/>5 </td>
<td> ¯1 </td>
</tr>

<tr>
<td> ├──L──┤   ├──R─┤ </td>
<td> (L0<R0)∧L1=R1 </td>
<td> </td>
<td> </td>
<td> 6<_pom/>4 7<_pom/>3 </td>
<td>  2 </td>
</tr>

<tr>
<td>  ├──L─┤  ├──R──┤ </td>
<td> (L0>R0)∧L1=R1 </td>
<td> </td>
<td> </td>
<td> 6<_pom/>4 5<_pom/>5 </td>
<td> ¯2 </td>
</tr>

<tr>
<td> ├──L──┤   ├─R─┤ </td>
<td> (L0<R0)∧L1>R1 </td>
<td> L⊃≈R </td>
<td> </td>
<td> 6<_pom/>4 6<_pom/>2 </td>
<td>  3 </td>
</tr>

<tr>
<td>  ├─L─┤  ├──R──┤ </td>
<td> (L0>R0)∧L1<R1 </td>
<td> L⊂≈R </td>
<td> </td>
<td> 6<_pom/>4 6<_pom/>6 </td>
<td> ¯3 </td>
</tr>

<tr>
<td> ├──L──┤   ├──R──┤ </td>
<td> (L0<R0)∧(L1<R1)∧L1>R0 </td>
<td> </td>
<td> </td>
<td> 6<_pom/>4 8<_pom/>4 </td>
<td>  4 </td>
</tr>

<tr>
<td>  ├──L──┤  ├──R──┤ </td>
<td> (L0>R0)∧(L1>R1)∧L0<R1 </td>
<td> </td>
<td> </td>
<td> 6<_pom/>4 4<_pom/>4 </td>
<td> ¯4 </td>
</tr>

<tr>
<td> ├──L──┼──R──┤ </td>
<td> L1=R0 </td>
<td> (L≤≈R)>L<≈R </td>
<td> L=∀R </td>
<td> 6<_pom/>4 14<_pom/>4 </td>
<td>  5 </td>
</tr>

<tr>
<td> ├──R──┼──L──┤ </td>
<td> L0=R1 </td>
<td> (L≥≈R)>L>≈R </td>
<td> L=∃R </td>
<td> 6<_pom/>4 0<_pom/>2 </td>
<td> ¯5 </td>
</tr>

<tr>
<td> ├──L──┤ ├──R──┤ </td>
<td> L1<R0 </td>
<td> L<≈R </td>
<td> L<∀R </td>
<td> 6<_pom/>4 17<_pom/>5 </td>
<td>  6 </td>
</tr>

<tr>
<td> ├──R──┤ ├──L──┤ </td>
<td> L0>R1 </td>
<td> L>≈R </td>
<td> L>∃R </td>
<td> 6<_pom/>4 0<_pom/>1 </td>
<td> ¯6 </td>
</tr>
</table>

The above table can be interpreted as follows:

<ul>
<li>The first column lists all possible '''Relative Positions''' two finite non-zero Radii (i.e., multi-valued) Balls may assume.</li>
<li>The second column lists the '''Endpoint Relationship''' that defines the corresponding '''Relative Position'''.</li>
<li>The third column lists a '''Defining Expression''' using the '''DoubleTilde''' operator that is TRUE for the corresponding '''Relative Position''' and FALSE for all of the other '''Relative Positions''' in different rows.</li>
<li>The fourth column lists a '''Defining Expression''' using the '''Always''' and '''Possibly''' operators [http://forum.nars2000.org/ball-arithmetic-t614.html#p1803 proposed] by David Rabenhorst.</li>
<li>The fifth column contains '''Sample Data''' that defines the corresponding '''Relative Position'''.</li>
<li>The sixth column gives the value of the '''Distinct Ball''' function <apll>{⍺∘≈⍵}</apll> when executed on the '''Sample Data'''.</li>
</ul>

Note the table is organized such that apart from the first row (whose '''Endpoint Relationship''' is Symmetric in <apll>L</apll> and <apll>R</apll>) the rows occur in six pairs where the '''Endpoint Relationship''' for each row in the pair can be obtained from the other by swapping <apll>L</apll> and <apll>R</apll>. The pairs of rows are ordered roughly by the amount of overlap between the two Balls.

As you can see, we are missing '''Defining Expressions''' for several of the '''Relative Positions'''. I was hoping that this table would be filled in with more detail so it would serve as a measure of whether or not we had a complete set of functions/operators which could be used to tell apart all of the '''Relative Positions''' without having to resort to the '''Distinct Ball''' function <apll>{⍺∘≈⍵}</apll> constructed explicitly for this purpose. Please add to this table if you find a missing expression, or wish to propose whole new functions and/or operators.

If you would like to try a '''Defining Expression''' on the '''Sample Data''', collect the data into a <apll>13</apll> by <apll>2</apll> matrix, put your '''Defining Expression''' into (say) an Anonymous Function and reduce the same data matrix along the last coordinate — if you have a successful function, you should see a <apll>13</apll>-element Boolean vector with a single one in it. For example,

<apll><pre>
sd←6<_pom/>4,⍪6<_pom/>4 5<_pom/>3 7<_pom/>5 7<_pom/>3 5<_pom/>5 6<_pom/>2 6<_pom/>6 8<_pom/>4 4<_pom/>4 14<_pom/>4 0<_pom/>2 17<_pom/>5 0<_pom/>1
{(⍺≤≈⍵)>⍺<≈⍵}/sd
0 0 0 0 0 0 0 0 0 1 0 0 0
{(⍺≥≈⍵)>⍺>≈⍵}/sd
0 0 0 0 0 0 0 0 0 0 1 0 0
{⍺<≈⍵}/sd
0 0 0 0 0 0 0 0 0 0 0 1 0
{⍺>≈⍵}/sd
0 0 0 0 0 0 0 0 0 0 0 0 1
</pre></apll>

Finally, the following example shows that the '''Sample Data''' are all distinct:

<apll><pre>
{⍺∘≈⍵}/sd
0 1 ¯1 2 ¯2 3 ¯3 4 ¯4 5 ¯5 6 ¯6
</pre></apll>

==Hypercomplex Ball Arithmetic==

Ball Arithmetic also extends to all Hypercomplex numbers (Complex, Quaternion, and Octonion) such as <apll>1.2<_pom/><_i/>2.2</apll> which represents a Complex Ball Arithmetic number whose '''Real''' part is <apll>1.2<_pom/></apll> and whose '''Imaginary''' part is <apll>2.2<_pom/></apll>. Note that the display of that number (at <apll>53</apll> bits of precision)

<apll><pre>
1.2<_pom/><_i/>2.2
1.2<_pom/>2.2<_E/>¯16<_i/>2.2<_pom/>4.4<_E/>¯16
</pre></apll>

demonstrates the inexactness of both floating point coefficients.

==Small v. Large Radii==

Ball Arithmetic works best when given small Radii, that is, small relative to the Midpoint. Calculations on large Radii numbers such as <apll>0<_i/>0.5<_pom/>0.5</apll> don't always produce useful results. For example,

<apll><pre>
⎕PP←8 ⍝ To reduce fractional clutter
x←0<_i/>0.5<_pom/>0.5
√x
0.4267767<_pom/>0.57<_i/>0<_pom/>1.4
</pre></apll>

Unfortunately in this case, because the Radius is so large it overwhelms the Midpoint, as can be seen when we Square the Square Root:

<apll><pre>
(√x)*2
0.18213835<_pom/>2.8<_i/>0<_pom/>2.8
x=(√x)*2
0
</pre></apll>

These calculations work better using smaller Radii, such as

<apll><pre>
y←0<_i/>0.5<_pom/>0.5<_E/>¯30
√y
0.5<_pom/>2.5<_E/>¯31<_i/>0.5<_pom/>7.5<_E/>¯31
(√y)*2
0<_pom/>1<_E/>¯30<_i/>0.5<_pom/>1<_E/>¯30
y=(√y)*2
1
</pre></apll>

where Square Root now returns a meaningful result, and its Square returns the original value.

==Mixed Notation==

The notation for Ball Arithmetic may be combined with other Point Notations in many ways, such as

<apll><pre>
1<_p/>1<_pom/>
3.141592653589793<_pom/>4.4<_E/>¯16
1<_x/>1<_pom/>
2.718281828459045<_pom/>4.5<_E/>¯16
</pre></apll>

For more details about combining different forms of Point Notation, see [[Point_Notation|Point Notation]] and in particular the section on [[Point_Notation#Mixed_Notation|Mixed Notation]].

==Implementation==

===Datataypes===

This feature adds four new datatypes to the language: one each for the four Division Algebras. Note that marking any one of the Hypercomplex coefficients as a Ball automatically marks them all as Balls.
{|
|-
|Real||    <apll>23<_pom/></apll> or <apll>1.2<_pom/></apll>
|-
|Complex||    <apll>23<_i/>1.2<_pom/></apll>
|-
|Quaternion||    <apll>23<_i/>1.2<_j/>0.5<_k/>¯3.14<_pom/></apll>
|-
|Octonion||    <apll>23<_i/>1.2<_j/>0.5<_k/>¯3.14<l/>¯3<_ij/>19<_jk/>¯3<_kl/>6<_pom/></apll>
|}

===Primitive Functions===

As described above, this feature adds two new monadic primitive functions '''Contract''' (<apll>≤R</apll>) and its inverse '''Distract''' (<apll>≥R</apll>), along with one new monadic primitive operator '''DoubleTilde''' <apll>f≈ R</apll> and <apll>L f≈ R</apll>.

All primitive and system functions and operators are extended to work with Ball Arithmetic, so you can test out this feature on all algorithms, with the exception of trying to find the Factorial of a multiple-precision/Ball Arithmetic Hypercomplex number, in which case you are temporarily out of luck.

==Workspace Version Number==

For the moment, if you save a workspace which contains a Ball Arithmetic number in a variable, the workspace version number is <apll>0.03</apll>. Otherwise, the version number is <apll>0.02</apll>. Workspaces with Ball Arithmetic Numbers may not be loaded by an earlier version of the interpreter.

==Pending Design Decisions==

<ul>
<li>The direct comparison functions (<apll><≤=≥>≠≡≢</apll>) and indirect ones (<apll>⍳⍸∊⍷⍒⍋</apll>) are currently implemented to compare the Midpoints of the arguments ignoring the Radius — this may change.</li>
<li>The Floor function and its related functions (Ceiling, Residue, GCD, LCM, and Encode) are currently implemented to execute on the Midpoints of the arguments ignoring the Radius — this may change.</li>
<li>We need new primitives to describe the different ways in which Balls can intersect/overlap, possibly through a new monadic operator (perhaps ≈) on the existing direct and indirect comparison and Floor-related functions.</li>
<li>Certain primitives are defined on limited ranges, and sometimes multiple ranges where it might use different algorithms for different ranges (e.g. <apll>L!R</apll>). There are open questions about how to handle a Ball that spans two or more of these ranges as to how to calculate the result.</li>
<li>Infinity and NaN have not been given any special treatment as yet. The results are whatever ARB has chosen which we may need to override.</li>
</ul>

==Download==

This feature is available as an Alpha release — it may be downloaded from the [http://www.nars2000.org/download/binaries/alpha/ Alpha] release directory.

==Acknowledgements==

<ul>
<li>The Arb team has contributed their Ball Arithmetic<ref name="BallArithmetic" /> code which made this feature possible.</li>
<li>David Rabenhorst has contributed many very useful ideas and many hours debugging the implementation.</li>
</ul>

==References==

<references />

System Variable FPC

2018-10-26T09:58:16Z

Drab:

<h1><aplxx>⎕FPC</aplxx> - Floating Point Control</h1>
{{BoxStart|<apll>Z←⎕FPC</apll>
|returns the current value of the Floating Point Control precision.}}
{{BoxEnd|<apll>Z</apll> is an integer scalar that represents the current value of the number of bits of precision used for [[Rational and VFP Numbers|Multiple-Precision Floating Point numbers]].}}
 

{| border="1" cellpadding="5" cellspacing="0" rules="none" summary=""
|-
||<apll>⎕FPC←R</apll>
||sets the value of the Floating Point Control precision to <apll>R</apll>.
|-
||<apll>⎕FPC←⍬</apll>
||sets the value of the Floating Point Control precision to its system default state.
|}
 


The approximate number of digits of precision for a given value of <apll>⎕FPC</apll> is about <apll>10⍟2*⎕FPC</apll>. The value of <apll>⎕FPC</apll> for a given required precision is about <apll>⌈2⍟10*digits</apll>.


The minimum value for <apll>⎕FPC</apll> is <apll>1</apll>; the maximum value is the largest positive 32-bit integer <apll>2147483647</apll>.

The system default value is <apll>128</apll>.

For example, (note: these numbers are normally one long string of digits; they are broken here, with the break shown by '''(continued)''' for your convenience in reading)

<apll><pre>
⎕PP←240
⎕FPC←64 ⋄ *1x
2.7182818284590452354
⎕FPC←128 ⋄ *1x
2.71828182845904523536028747135266249776
⎕FPC←256 ⋄ *1x
2.71828182845904523536028747135266249775724709369995957496696762 (Continued)
772407663035355
⎕FPC←512 ⋄ *1x
2.71828182845904523536028747135266249775724709369995957496696762 (Continued)
7724076630353547594571382178525166427427466391932003059921817413 (Continued)
5966290435729003342952605956
⎕FPC←1024 ⋄ *1x
2.71828182845904523536028747135266249775724709369995957496696762 (Continued)
7724076630353547594571382178525166427427466391932003059921817413 (Continued)
5966290435729003342952605956307381323286279434907632338298807531 (Continued)
9525101901157383418793070215408914993488416750924</pre></apll>

Whereas, according to Wolfram Alpha, the exact value to 230 digits is

<apll><pre>2.71828182845904523536028747135266249775724709369995957496696762 (Continued)
7724076630353547594571382178525166427427466391932003059921817413 (Continued)
5966290435729003342952605956307381323286279434907632338298807531 (Continued)
9525101901157383418793070215408914993488...</pre></apll>

'''See Also: [[System_Variable_CT|Comparison Tolerance ⎕CT]]'''

{{System Variables}}
{{Article footer}}

System Variable FPC

2018-10-26T09:56:39Z

Drab:

<h1><aplxx>⎕FPC</aplxx> - Floating Point Control</h1>
{{BoxStart|<apll>Z←⎕FPC</apll>
|returns the current value of the Floating Point Control precision.}}
{{BoxEnd|<apll>Z</apll> is an integer scalar that represents the current value of the number of bits of precision used for [[Rational and VFP Numbers|Multiple-Precision Floating Point numbers]].}}
 

{| border="1" cellpadding="5" cellspacing="0" rules="none" summary=""
|-
||<apll>⎕FPC←R</apll>
||sets the value of the Floating Point Control precision to <apll>R</apll>.
|-
||<apll>⎕FPC←⍬</apll>
||sets the value of the Floating Point Control precision to its system default state.
|}
 


The approximate number of digits of precision for a given value of ⎕FPC is about 10⍟2*⎕FPC. The value of ⎕FPC for a given required precision is about ⌈2⍟10*digits.


The minimum value for <apll>⎕FPC</apll> is <apll>1</apll>; the maximum value is the largest positive 32-bit integer <apll>2147483647</apll>.

The system default value is <apll>128</apll>.

For example, (note: these numbers are normally one long string of digits; they are broken here, with the break shown by '''(continued)''' for your convenience in reading)

<apll><pre>
⎕PP←240
⎕FPC←64 ⋄ *1x
2.7182818284590452354
⎕FPC←128 ⋄ *1x
2.71828182845904523536028747135266249776
⎕FPC←256 ⋄ *1x
2.71828182845904523536028747135266249775724709369995957496696762 (Continued)
772407663035355
⎕FPC←512 ⋄ *1x
2.71828182845904523536028747135266249775724709369995957496696762 (Continued)
7724076630353547594571382178525166427427466391932003059921817413 (Continued)
5966290435729003342952605956
⎕FPC←1024 ⋄ *1x
2.71828182845904523536028747135266249775724709369995957496696762 (Continued)
7724076630353547594571382178525166427427466391932003059921817413 (Continued)
5966290435729003342952605956307381323286279434907632338298807531 (Continued)
9525101901157383418793070215408914993488416750924</pre></apll>

Whereas, according to Wolfram Alpha, the exact value to 230 digits is

<apll><pre>2.71828182845904523536028747135266249775724709369995957496696762 (Continued)
7724076630353547594571382178525166427427466391932003059921817413 (Continued)
5966290435729003342952605956307381323286279434907632338298807531 (Continued)
9525101901157383418793070215408914993488...</pre></apll>

'''See Also: [[System_Variable_CT|Comparison Tolerance ⎕CT]]'''

{{System Variables}}
{{Article footer}}

Find

2016-12-04T15:51:16Z

Drab:

<apll>
Z←L⍷R
</apll>

 Find looks for patterns <apll>L</apll> in <apll>R</apll>.
 It returns a simple boolean array with the same shape as <apll>R</apll>.
 
<apll>L</apll> and <apll>R</apll> are arbitrary arrays.
But to be useful,
the rank of <apll>L</apll>
should not be more than the rank of <apll>R</apll>,
and both arguments should have the same type
(either both character or both numeric).
 
If <apll>L</apll> is a scalar,
then <apll>Z←L⍷R ←→ R∊L</apll>.
 
If <apll>L</apll> is not a scalar,
then it represents a contiguous pattern,
which may or may not be found in one or more places in <apll>R</apll>.
 
If the arguments are numeric,
then the result is sensitive to <apll>⎕CT</apll>.

 
Examples:
 
<apll>

        R←3 8⍴'thirteenfourteenfifteen '
        R
 thirteen
 fourteen
 fifteen

 
        't'⍷R
 1 0 0 0 1 0 0 0
 0 0 0 0 1 0 0 0
 0 0 0 1 0 0 0 0

 
        'teen'⍷R
 0 0 0 0 1 0 0 0
 0 0 0 0 1 0 0 0
 0 0 0 1 0 0 0 0

</apll>

Find

2016-12-03T13:12:00Z

Drab: Created page with "<apll> Z←L⍷R </apll> Find looks for patterns <apll>R</apll> in <apll>L</apll>. It returns a simple boolean array with the same shape as <apll>L</apll>..."

<apll>
Z←L⍷R
</apll>

 Find looks for patterns <apll>R</apll> in <apll>L</apll>.
 It returns a simple boolean array with the same shape as <apll>L</apll>.
 
<apll>L</apll> and <apll>R</apll> are arbitrary arrays.
But to be useful,
the rank of <apll>R</apll>
should not be more than the rank of <apll>L</apll>,
and both arguments should have the same type
(either both character or both numeric).
 
If <apll>R</apll> is a scalar,
then <apll>Z←L⍷R ←→ L∊R</apll>.
 
If <apll>R</apll> is not a scalar,
then it represents a contiguous pattern,
which may or may not be found in one or more places in <apll>L</apll>.
 
If the arguments are numeric,
then the result is sensitive to <apll>⎕CT</apll>.

 
Examples:
 
<apll>

        L←3 8⍴'thirteenfourteenfifteen '
        L
 thirteen
 fourteen
 fifteen

 
        L⍷'t'
 1 0 0 0 1 0 0 0
 0 0 0 0 1 0 0 0
 0 0 0 1 0 0 0 0

 
        L⍷'teen'
 0 0 0 0 1 0 0 0
 0 0 0 0 1 0 0 0
 0 0 0 1 0 0 0 0

</apll>

Compose

2016-12-01T13:56:58Z

Drab:

In the following descriptions, <apll>f</apll> and <apll>g</apll> represent functions and <apll>a</apll> and <apll>b</apll> represent variables.

* The form <apll>f∘g</apll> may be used both monadically and dyadically.
:
:Monadic: <apll>Z←f∘gR</apll> is identical to <apll>Z←fgR</apll>.
:
:Dyadic: <apll>Z←Lf∘gR</apll> is identical to <apll>Z←LfgR</apll>.

* The form <apll>f∘b</apll> may be used monadically only.
:
:Monadic: <apll>Z←(f∘b) R</apll> is identical to <apll>Z←R f b</apll>.
:Note that parentheses are required around the function to avoid interpreting <apll>b R</apll> as a strand.

* The form <apll>a∘g</apll> may be used monadically only.
:
:Monadic: <apll>Z←a∘gR</apll> is identical to <apll>Z←a g R</apll>.

* The derived function from the form <apll>a∘b</apll> always signals a <apll>SYNTAX ERROR</apll>.

Compose can be useful for function assignment (but enclosing parentheses are necessary).
 
For example 
<apll>
      p1←(1∘+¯2∘π) ⍝ 1 more than the Nth prime number 
      p1 ⍳9 
3 4 6 8 12 14 18 20 24
</apll>

Rational and VFP Numbers

2016-12-01T12:20:30Z

Drab: /* Conversions */

===Introduction===

Rational numbers complement the existing 64-bit integer datatype to provide infinite precision (or <apll>WS FULL</apll>) at the cost of some performance. Similarly, Variable-precision Floating Point (VFP) numbers complement the existing 64-bit floating point (IEEE-754) datatype to provide more precision (as much as the user cares to specify) again with some cost in performance.

There is no separate datatype for infinite precision Integers. Instead, they are represented as a special case of rational numbers. In the discussion below, the phrase rational integer means a rational number whose denominator is one.

Throughout this discussion the similarity between integer and rational numbers as well as floating point and VFP numbers will become apparent.

===Rationale===

Simply put: precision. If the 53 bits of precision in the floating point result of, say, <apll>2÷3</apll> is not enough, you now have two more choices: one as an exact number, and one as a binary floating point number with as much precision as you care to specify.

For example,

<apll>      ⎕PP←60 ⋄ ⎕FPC←512 
      2÷3 
0.6666666666666666 
      2÷3x 
2r3 
      2÷3v 
0.666666666666666666666666666666666666666666666666666666666667 
      2*200 
1.6069380442589903E60 
      2*200x 
1606938044258990275541962092341162602522202993782792835301376 
      ○1 
3.141592653589793 
      ○1v 
3.14159265358979323846264338327950288419716939937510582097494</apll>

===Constants===

Rational constants may be entered by suffixing an integer constant with an <apll>x</apll> (for a rational integer) or separating the numerator and denominator with an <apll>r</apll> (for a rational number) as follows:

* <apll>123x</apll> for the constant <apll>123</apll> (the suffix <apll>x</apll> is but a shorthand for <apll>r1</apll>)
* <apll>123r4567</apll> for the constant <apll>123÷4567</apll>

VFP constants may be entered by suffixing the integer or floating point constant with a <apll>v</apll> as follows:

* <apll>123v</apll> for the constant <apll>123</apll>
* <apll>123.4567v</apll> for the constant <apll>123.4567</apll>
* <apll>123.4567e3v</apll> for the constant <apll>123.4567e3</apll>

The above formats for constants (except for the suffix <apll>x</apll>) may be used in other constants such as <apll>1r4p2</apll> to generate a shorter and more accurate value for <apll>π2/4</apll> than, say, <apll>((○1)*2)÷4</apll>.

===Precision===

Rational numbers have infinite precision. They are stored with a separate numerator and denominator, both of which are exact numbers in the sense that their size (and hence precision) grows limited only by the available workspace.

VFP numbers have user-controlled variable precision, and each number may have a different precision. The default precision at startup is controlled by the value of the system variable <apll>⎕FPC</apll>. The system default of this value is <apll>128</apll> in units of bits of precision of the mantissa (the digits) of the number, not counting the exponent (which is of fixed size). The current precision may be changed as needed by assigning a new value to the system variable. All newly created VFP numbers will have the new precision – the precision of VFP numbers already present in the workspace does not change.

Generally, precision is set once for a particular application and unchanged thereafter. Although not recommended, it is possible to mix VFP numbers of different precisions in a single array – presumably you really know what you are doing. The system function <apll>0 ⎕DR</apll> may be used to display an array's precision(s).

===Datatype Propagation===

Generally, the datatype of constants propagates through a calculation. That is, if you start with a rational number and don't calculate with irrational or transcendental functions, you'll end up with a rational result, and if you start with a VFP number, you'll end up with a VFP result.

An example from the programming problems site [http://projecteuler.net/ ProjectEuler.net] illustrates this point. [http://projecteuler.net/problem=48 Problem #48] asks what are the low-order ten digits of the sum of the first thousand instances of NN?

The obvious expression <apll>¯10↑⍕+/*⍨⍳1000</apll> at first sight seems to solve the problem until you realize that it quickly runs afoul of the limited precision of 64-bit integer and floating point numbers. Clearly, this is a problem for the infinite precision of rational integers.

As <apll>⍳1000</apll> generates the first thousand integers as an integer datatype (actually an Arithmetic Progression Array), <apll>⍳1000x</apll> generates the same values as rational integers. Next, <apll>*⍨⍳1000x</apll> generates the first thousand instances of NN as exact rational integers, and unlike its integer counterpart, there is no overflow to floating point, just an increase in precision (as well as space used in the workspace). Then, <apll>+/*⍨⍳1000x</apll> sums them into a single 3001-digit rational integer, and finally <apll>¯10↑⍕+/*⍨⍳1000x</apll> converts the large integer to characters and extracts the low-order ten digits — <apll>9110846700</apll> — all in a small number of milliseconds.

Note how we started with an obvious expression that failed because of its limited precision, and made a single change to suffix the constant <apll>1000</apll> with an <apll>x</apll> to convert it to a rational integer which then propagates through the calculation with infinite precision to yield the correct result.

===Display===

Rational integers are displayed as an integer with no special adornment; rational non-integers are displayed as a numerator and denominator separated by an <apll>r</apll> as in <apll>34r9</apll>. As with the integer datatype, the numerator and denominator of a rational number are displayed exactly, unaffected by the current setting for Printing Precision (<apll>⎕PP</apll>).
 
<apll>      !40 
8.159152832478977E47 
      !40x 
815915283247897734345611269596115894272000000000 
      +\÷⍳10x 
1 3r2 11r6 25r12 137r60 49r20 363r140 761r280 7129r2520 7381r2520 
</apll>

VFP numbers are displayed as decimal numbers to the precision inherent in the number or <apll>⎕PP</apll>, whichever is smaller, just as floating point numbers are displayed. For example,

<apll>      ⎕PP←100 
      ⎕FPC←64 
      ○1 
3.141592653589793 
      ○1x 
3.14159265358979323851 
      ⎕FPC←128 
      ○1x 
3.141592653589793238462643383279502884195</apll> 

where both of the above displays were limited by the precision of the number, not <apll>⎕PP</apll>.

However, the first of the following displays is limited by <apll>⎕PP</apll>:

<apll>      ⎕FPC←128 
      ⎕PP←20 
      !40v 
81591528324789773435____________________________ 
      ⎕PP←80 
      !40v 
8159152832478977343456112695961158942720________</apll> 

In the last display, the current setting of Printing Precision is large enough, but the current setting of the Floating Point Control (<apll>⎕FPC</apll>) whose value is in bits is too small, so the display is truncated.

===Formatted Display===

The system function <apll>⎕FMT</apll> has been enhanced to allow formatting of rational numbers via the (new) <apll>R</apll>-format specifier. For example,

<apll>      'R4.2' ⎕FMT ∘.÷⍨⍳6x 
1   1r2 1r3 1r4 1r5 1r6 
2   1   2r3 1r2 2r5 1r3 
3   3r2 1   3r4 3r5 1r2 
4   2   4r3 1   4r5 2r3 
5   5r2 5r3 5r4 1   5r6 
6   3   2   3r2 6r5 1 </apll>

Moreover, the Symbol Substitution (<apll>S<…></apll>) feature of <apll>⎕FMT</apll> allows you to substitute a different symbol for the default <apll>r</apll> used to separate the numerator and denominator of a rational number, as in

<apll>      'S<r/>R4.2' ⎕FMT ∘.÷⍨⍳6x 
1   1/2 1/3 1/4 1/5 1/6 
2   1   2/3 1/2 2/5 1/3 
3   3/2 1   3/4 3/5 1/2 
4   2   4/3 1   4/5 2/3 
5   5/2 5/3 5/4 1   5/6 
6   3   2   3/2 6/5 1 </apll>

===Datatype Promotion===

For the most part, rational numbers beget rational numbers and VFP numbers beget VFP numbers. However, when irrational, transcendental, and certain other functions are used, rational numbers beget VFP numbers. For example,

<apll>      *1 
2.718281828459045 
      *1x 
2.718281828459045235360287471352662497759</apll> 

where the datatype of the two results are floating point and VFP, respectively. That is, in a manner similar to how some primitive functions with integer arguments may return floating point results, when a rational number is used as an argument to a primitive function that can't return a result with infinite precision, it returns a VFP number.

The reason irrational, transcendental, and certain other functions on rational numbers do not return rational numbers is that, by definition, the result of such a function is, in general, not representable as a rational number; instead, VFP numbers are better suited to represent irrational results where the end user may control exactly how much precision is desired in an obviously inexact number.

Two special functions are the prime decomposition (<apll>πR</apll>)/number theoretic (<apll>LπR</apll>) functions. In these cases, fractional or VFP right arguments are converted to integers or rational integers, respectively, which is the datatype of the result except for <apll>0πR</apll> (Primality Test) which always returns a Boolean result regardless of the type of <apll>R</apll>.

Ignoring purely structural functions, the list of functions that produce VFP numbers given rational numbers is as follows:

* Power: <apll>*R</apll> and <apll>L*R</apll> (except when <apll>R</apll> is a 32-bit integer, in which case the result is a rational number)
* Logarithm: <apll>⍟R</apll> and <apll>L⍟R</apll>
* Pi Times and Circle functions: <apll>○R</apll> and <apll>L○R</apll>
* Root: <apll>√R</apll> and <apll>L√R</apll>
* Factorial and Binomial: <apll>!R</apll> and <apll>L!R</apll> (except when the arguments are rational integers, in which case the result is a rational integer)

Beyond the ones mentioned above, the list of functions that don't produce a rational or VFP result given those argument(s) is as follows:

* Depth: <apll>≡</apll> (Integer)
* Dyadic Comparison: <apll>= ≠ < ≤ ≥ > ≡ ≢</apll> (Boolean)
* Nand and Nor: <apll>⍲ ⍱</apll> (Boolean)
* Grade Up/Down: <apll>⍋ ⍒</apll> (Integer)
* Index Of: <apll>⍳</apll> (Integer)
* Member Of: <apll>∊</apll> (Boolean)
* Find: <apll>⍷</apll> (Boolean)
* Subset and Superset: <apll>⊆ ⊇</apll> (Boolean)
* Format: <apll>⍕</apll> (Character)

Otherwise, rational argument(s) produce rational result(s) and VFP argument(s) produce VFP result(s).

===Datatype Demotion===

It is common in APL implementations to demote datatypes where appropriate. For example, the constant <apll>1.0</apll> might actually be represented as an integer or even Boolean datatype. The idea is there is no loss of precision and the storage is typically smaller which might lead to a more efficient algorithm when next used, so why not?

With rational and VFP numbers those reasons no longer apply. While the constant <apll>1x</apll> might have the same precision as the constant <apll>1.0</apll>, the difference in latent precision between the two is vast. In fact, in order for datatype propagation of rational and VFP numbers to work at all, we must be careful not to demote them automatically to a smaller datatype. Otherwise, it would require an intolerable degree of analysis on the part of the programmer to ensure that the desired datatype (rational or VFP) remains in effect throughout a calculation.

===Conversions===

To convert manually from one datatype to another,

* Integer to rational, use <apll>0x+</apll>
* Integer, floating point, or rational number to VFP, use <apll>0v+</apll>
* Rational integer to a 64-bit integer (possibly losing some precision), use <apll>⍎⍕</apll>
* Rational non-integer to a floating point number (possibly losing some precision), use <apll>0.0+</apll>
* VFP number to a floating point number (possibly losing some precision), use <apll>⍎⍕</apll>
* VFP number to rational — at the moment, no APL method exists although the underlying libraries do have such a function

===Comparisons===

<ul><li>Comparisons between two rational numbers or a rational number and any other integer is exact — just as they are between integers.</li>

<li>Comparisons between a rational number and a floating point number convert both arguments to VFP numbers and compare the two as below.</li>

<li>Comparisons between a VFP number and any other number is sensitive to the current setting of Comparison Tolerance (<apll>⎕CT</apll>) — just as they are between floating point numbers.</li>

</ul>

That is, comparisons continue the analogy between integers and rationals as well as floats and VFPs.

===Integer Tolerance===

Both rational and VFP numbers may be used where the system ordinarily requires an integer (such axis coordinates, indexing, left argument to structural primitives, etc.) just as the system tolerates floating point numbers in those contexts if they are sufficiently near an integer. In all cases, the system attempts to convert the non-integer to an integer using the fixed system comparison tolerance (at the moment, <apll>3E¯15</apll>).

===Infinities===

Support for <apll>±∞</apll> has been extended to rational and VFP numbers in the same manner as it applies to 64-bit integers and 64-bit floats. That is, the same cases covered by the system variable <apll>⎕IC</apll> (Indeterminate Control) also apply to infinite rational and VFP numbers. Moreover, infinite numeric constants may be entered, for example, as

* <apll>∞x</apll>
* <apll>∞r1</apll>
* <apll>∞v</apll>
* <apll>∞v0</apll>

Also constants such as <apll>2r∞</apll> resolve to <apll>0x</apll>.

===New And/Or Different Behavior===

<ul>
<li>Both roll (<apll>?R</apll>) and deal (<apll>L?R</apll>) on rational integers use a built-in random number generator so as to use the entire range of rational integers – this algorithm uses its own internal seeds that are much more complicated than the simple integer seed that is <apll>⎕RL</apll> (Random Link). Thus <apll>⎕RL</apll> is unchanged by these functions on rationals.

For example, if you need really large random numbers

<apll>      ?10*60x 
370857192605742854709703007683731949504799559659692534573173</apll>
</li>
<li>Matrix inverse (<apll>⌹R</apll>) and matrix division (<apll>L⌹R</apll>) on rational or VFP arguments each have two limitations above and beyond that of normal conformability:
<ul>
<li>for a square right argument that it be non-singular, and</li>
<li>for an overdetermined (<apll>>/⍴R</apll>) right argument that the symmetric matrix <apll>(⍉R)+.×R</apll> be non-singular.</li>
</ul>
These limitations are due to the algorithm ([http://en.wikipedia.org/wiki/Gauss%E2%80%93Jordan_elimination Gauss-Jordan Elimination]) used to implement Matrix Inverse/Divide on rational and VFP numbers.

Integer and floating point arguments are not subject to these limitations because they use a more general algorithm ([http://en.wikipedia.org/wiki/Singular_value_decomposition Singular Value Decomposition]) that produces a unique result even for singular arguments (e.g., <apll>⌹5 3⍴0</apll>).</li>
</ul>

===Conclusions===

The new datatypes offer several benefits:

* They extend the precision of existing integer and floating point datatypes to a much greater level.
* As integer blows up to floating point, rational blows up to VFP, providing a natural parallel progression for irrational and transcendental primitive functions.
* There is a close similarity between integer and rational numbers as well as floating point and VFP numbers.
* Datatype propagation without demotion allows one to code an algorithm in either of the new types easily and without the need for detailed analysis of the datatype in intermediate results.
* All primitives extend naturally to encompass the new types as numbers.
* The notation for constants builds on existing point notation formats.

===Acknowledgments===

The designers of J are thanked for having the foresight to include [http://www.jsoftware.com/help/dictionary/dictg.htm rational numbers] as a separate datatype.

The following LGPL libraries have been used to provide support for these datatypes:

* MPIR (Multiple Precision Integers and Rationals) at [http://www.mpir.org mpir.org].
* MPFR (Multiple Precision Floating-Point Reliable Library ) at [http://www.mpfr.org mpfr.org].

===References===

For a PDF version of this page, view it [http://www.sudleyplace.com/APL/Rational%20&%20Variable-Precision%20FP.pdf here].

Commute-Duplicate

2016-12-01T11:15:03Z

Drab:

<table border="1" cellpadding="5" cellspacing="0" rules="none" summary="">
<tr>
<td>
<table border="0" cellpadding="5" cellspacing="0" summary="">
<tr>
<td valign="top">Duplicate: <apll>Z←f⍨ R</apll></td>
<td></td>
<td>    </td>
<td>returns <apll>R f R</apll>.</td>
</tr>
</table>
</td>
</tr>

<tr>
<td>
<table border="0" cellpadding="5" cellspacing="0" summary="">
<tr>
<td valign="top">Commute: <apll>Z←L f⍨ R</apll></td>
<td></td>
<td></td>
<td>returns <apll>R f L</apll>.</td>
</tr>
</table>
</td>
</tr>

<tr>
<td><apll>L</apll> and <apll>R</apll> are arrays.</td>
</tr>

<tr>
<td><apll>f</apll> is a function.</td>
</tr>
</table>

For example:

<apll>
      ×⍨⍳5 ⍝ squares 
1 4 9 16 25 
      ∘.=⍨⍳4 ⍝ Identity matrix 
1 0 0 0 
0 1 0 0 
0 0 1 0 
0 0 0 1 
</apll>

<apll>
      2*⍨⍳5 ⍝ squares 
1 4 9 16 25 
      √+/2*⍨6 2⍴3 4, 5 12, 8 15, 7 24, 20 21, 12 35 ⍝ Pythagorean hypotenuses 
5 13 17 25 29 37 
</apll>