The study is dedicated to considering some examples of computer calculations that are completely correct from the point of view of the standards of computer work with numbers and arithmetic operations, the result is incorrect to a human performer. It's demonstrated that mathematically correct calculations are often impossible to reproduce correctly within the IEEE-754 standard due to the fundamental limitations of the memory allocated to represent a specific number. The authors focused on the study of the specified problem using the example of the numpy module, as one of the main software tools for modern scientific calculations. The work has methodological value in the professional training of specialists in Computer Science.
The limitation on the amount of memory allocated by computer technology for number processing imposes fundamental limitations on the possibilities and logic of the organization of computer calculations. It’s also requires the involvement of specific mathematical research methods.
The standard and most popular tool for computer calculations, especially in machine
learning tasks, is the Python programming language, in particular with its powerful extension
NumPy module. It allows you to efficiently work with large multidimensional arrays and
matrices, along with a large library of mathematical functions. The NumPy module provides
much wider possibilities for working with mathematical functions and objects compared to the
math module. The NumPy is an effective tool for providing complex mathematical calculations
not only in the tasks of machine learning and data analysis, but is also an indispensable
component in the implementation of a large volume of computer calculations in a wide class of
modern engineering and scientific tasks [
The investigation considers the main aspects related to the problems of computer
calculations within the IEEE-754 [
The history of the development of computing and the basics of the theory of computing with
floating-point numbers are described in [
The presented work is a continuation of the authors' research [20-21]. Previously, the
authors investigated rather complex cases of violation of the correctness of computer
calculations, the example of Ramp [22] and the sequence of Müller [
Let's consider the numeric data types of the NumPy library, one of the most popular Python libraries for scientific computing. Figure 1 shows the hierarchical diagram of the numpy.generic data class, the base class for numpy scalar types. We not consider components of numpy.generic like nympy.bool_, nympy.object_, nympy.datetime64 or nympy.flexible.
Floating point numbers according to the IEEE 754 standard are presented in the form = (−1) × 1. × 2 where S is number sign (signum), M is mantissa, E is number order (exponent). Figure 2 shows the distribution of bits of a single-precision number in the binary32 format (the most significant bit is the sign, the order of the number is preserved using 8 bits, 23 bits are allocated to the mantissa – a total of 32 bits).
The mantissa is normalized (the integer part is equal to one), as we can see from the formal representation of the number n. To calculate the order, a constant shift exponent is used, which for single-precision numbers is 127 (or 01111111 in binary code).
Let's consider for the clarity, what will happen to the representation of the rational number 1⁄3 according to the described principles of IEEE 754. Mathematically, in the decimal and binary number systems, we have
1⁄3 = 0. (3)10 = 0. (01)2.
For a single-precision number, we can write
1⁄3 = (−1)0 ∙ 1.01010101010101010101011∙ 2−10, i.e S = 0 (number is positive), M = 01010101010101010101011 (since the first discarded bit was a unit, according to the rules of the standard, we add 1 to the mantissa), E = -10 (the order is -2 in the decimal number system).
We have float32-approximation to the 1⁄3 in the mathematical form (0.01010101010101010101011)2 = 0 ∙ 2−1 + 1 ∙ 2−2 + 0 ∙ 2−3 + 1 ∙ 2−4 + 1 ∙ 2−5 … = 12 =1
1 1 = ∑ 22 + 225 = (0.3333333432674407958984375)10
The obtained approximate value of the number is greater than the real value 1⁄3. Everything is still correct in the mathematical sense
0. (01)2 = ∑∞=1 212 = 1 −1⁄14⁄4 = 1⁄3
Let's consider an important concept of the theory of computer calculations – machine epsilon [23]. This is numerical value below which it is impossible to specify the relative accuracy for any algorithm that produces real numbers.
A close concept, but not identical, is the concept of machine zero. That is a numerical value with such a negative order that the computer perceives it as zero.
In fact, to get the machine epsilon value in Python you can use the function numpy.finfo() for specific floating point number formats. But it seems appropriate to give a software implementation (Listing #1) of calculating the values of machine epsilon and machine zero, in order to understand the essence of these concepts in a constructive way.
Listing #1 1 from numpy import * 2 def machineEpsilon(float_format = float): where the exponent and the mantissa determine the number of bits for the corresponding data type.
Let's consider one interesting example of working with integers when overflowing a bit grid.
This problem itself is well known, but when training specialists, it is important to demonstrate that the results obtained are incorrect to a human programmer, but are nevertheless absolutely correct from the point of view of executing the program within the framework of existing standards. A simple example is the calculation of the value of the factorial of some natural number n. With implementing directly in Python (the factorial function from the math module), there are no problems with accuracy. Python implements the BigInt [24] byte short intc int_ ubyte ushort uintc uint n
We obtained a paradoxical result, that the product of natural numbers is equal to a negative number.
In integers in the int and long long formats, everything is calculated correctly. 21! = (51090942171709440000)10, 21! = (101100010100000111011111010011011010111000110001000000000000000000)2
An entry in binary code is 66 bits long, while the int type in C is 32 bits long (the most significant bit is the sign), and the long long type is 64 bits long. This means that in binary code, discarding the higher bits, we will actually get int factorial(21) = (10111000110001000000000000000000)2 long long (21) = (1100010100000111011111010011011010111000110001000000000000000000)2 Taking into account the fact that we are working with sign types and the presented results are written in additive code, the reverse conversion to the decimal numbering system will give the numbers – 1195114496 and –4249290049419214848, respectively, that is, the calculation is correct not from a mathematical point of view, but according to standards and algorithms work with integer types.
In order to increase the accuracy of calculations and obtain correct calculation results, it is often recommended to use floating-point number types instead of integer data types. Sometimes it can give a certain effect, but it is also significantly limited, and the correctness of the obtained results also requires mandatory verification.
We will use the floating point data types presented in Table 1 instead of integer types in the described example of calculating the factorial of a number.
The float32 type ceases to guarantee an exact result already at n = 14:
float32(14!) = (87178289152)10 while the correct value 14! = 87178291200. In terms of the IEEE-754 standard, Python calculates again without errors.
After converting to binary code, we have (87178291200)10 = (1010001001100001110110010100000000000)2 float32(87178291200) = 0 sign 10100011 01000100110000111011001
exponent mantissa that is, we discard "0100000000000" in the lower bits, which is exactly 2048 in the decimal number system. The number 2048 is the error between the actual calculation value of 14! and the result of the calculation is float32(14!).
Starting with n = 23, the type of floating point numbers float64 starts to return an incorrect result. The value of the factorial for n > 25 is incorrect even for the longdouble type, while we will get a value of longdouble(26!) greater than the real value of 26!, which also follows from the peculiarities of the application of the IEEE-754 standard:
longdouble (403291461126605635584000000) = In the decimal code, 403291461126605635592388608.
0 100000001010111 sign exponent 0100110110011000010010011110101000110111111011101010110010010010 mantissa the specified number is converted as 3.2. Homer Simpson’s numbers Let's consider another example of manipulating integer formats. The images in Figure 3 are stills from the episode "The Wizard of Evergreen Terrace" [26] and the episode "Treehouse of Horror VI. Homer3 (Homer Cubed)" [27] of the popular animated series The Simpsons.
Mathematical statements relating to Fermat's Great Theorem will be true when checked with a calculator [28]. Let's consider the first example when working with numeric data types from the numpy module.
Of course, the actual numerical values do not match
Note that an attempt to convert any number to int32, int64 (uint32, uint64) types generates an OverflowError error.
Direct conversion to float32 also generates an OverflowError, so let's try to consider a mathematically identical expression
398.712 + 436.512 ? 447.212
Now the result of the operation in Python (398.7**12 + 436.5**12 == 447.2**12) will be True. Indeed, it is fashionable to write within the IEEE-754 standard
Arithmetic expression: Correct value: Float32 format: Binary format: Arithmetic expression: Correct value: Float32 format: Binary format:
That is, the reason lies in discarding the lower bits when converting a number to the float32 format. Let's consider the same problem, but staying within the framework of integer formats.
For example, for n = 5 for the int32 format, more than 1400 triples of natural numbers (A, B, C) can be found, so that when calculating int32(A**5) + int32(B**5) == int32(C**5 ) would return a boolean value of True.
The maximum value of the parameter C: C = 6208, and as can be seen from the Table 4, several pairs (A, B) can correspond to one value of C. The explanation for this is simple: different large integers (more than 4 bytes) due to the discarding of the high-order bits in the int32 format can have the same value. For triples from the Table 4 (32, B, 6208) we have the obtained result is incorrect from a mathematical point of view. But you can choose an example when an untrained user can mistakenly perceive the result as correct, for example, implementing C++ programs for the standard type int (4 bytes). From our point of view, this is the example of A = 1504, B = 2960, C = 4496, n = 5. At the same time, it is interesting that the values of the expressions A5 and B5 in the int32 format are generally negative numbers: int32(1504**5) = -570425344, int32(2960**5) = -1584398336,
In the additional code in int32 format, we can write (-570425344)10 = (11011110000000000000000000000000)2 (-1584398336)10 = (10100001100100000000000000000000)2 In the binary code, the sum of these numbers will have a length of 33 bits and is equal to (-570425344)10 + (-1584398336)10 = (101111111100100000000000000000000)2 Discarding the most significant bit (33 bits), we get the result
(01111111100100000000000000000000)2 = (2140143616)10.
Let's move on to consider some typical examples of incorrect behavior of computational algorithms about working with floating point numbers.
Consider the following problem: it is necessary to find such a minimum value of n that 1 1 + =
1 2020 , ∈ ∈ , >
This problem is an example of Diophantine equation. In a general sense, the Diophantine equation is an equation in which only integer solutions are allowed.
The program for solving this equation is presented in listing #2.
Listing #2 1 2 3 4 5 6 7 8 9 10 import numpy as np def mathExpr (num_format): value = num_format(1/2020) n = 1 while True: n1 = num_format(1/n) for m in range(1, n):
return n, m n += 1
if n1 + num_format(1/m) == value: 11 float_formats = (np.float16, np.float32, np.float64, np.float128) for form in float_formats:
print(f"{form}:", mathExpr (form))
As a result, we got four different answers depending on the format of the floating point numbers. At the same time, the pairs (4545, 3636) and (6060, 3030) are the correct solution to the equation.
It is easy to show that n = 4545, if we look for the minimum value of n. That is, when using the float16 and float32 formats, we will get incorrect answers, and in the case of the float128 format, we will get a correct answer from a mathematical point of view, but we will skip over the minimum possible value of n (n = 4545), which satisfies the problem. We can engage the specialized apparatus of computer computation to improve the accuracy of computation. Let's try to increase the accuracy of calculations in Python using a special decimal module [30] by adding a decimal.Decimal element to the float_formats tuple. However, when trying to execute the program, we will get into an infinite loop, and the task of finding an exact match of values can be considered as failed.
We can correct this situation, for example, using the isclose method from the math module. To do this, in listing #2, line 8 must be replaced with the command if math.isclose(n1+num_format(1/m), value, rel_tol=1e-9), and in this case, the float128 and decimal.Decimal format will already have n = 4545, m = 3636. Let's change line 3 in listing #1 to the command value = num_format(1/10), and after starting the program, we will get the result: <class 'numpy.float16'>: n = 30, m = 15 <class 'numpy.float32'>: n = 35, m = 14 <class 'numpy.float64'>: n = 30, m = 15 <class 'numpy.longdouble'>: Infinite loop <class decimal.Decimal'>: Infinite loop
Obviously, the answer to the problem is the pair n = 30, m = 15 (the pair n = 35, m = 14 is also a correct result, but the minimum value of the parameter n is still equal to 30). With the double-precision format, the program falls into an infinite loop, that is, the result cannot be calculated.
If we calculate the specified example using the isclose function, then for the case of rel_tol=1e-7 for all formats of floating point numbers the result will be n = 30, m = 15, and for the value of the rel_tol parameter of greater precision (for example , already rel_tol=1e-8) – for float16, float64, float128, decimal.Decimal formats the result will be n = 30, m = 15, and for float32 format the result will still be n = 35, m = 14.
We can make an interesting intermediate conclusion: in some examples of computer calculations with floating point numbers, a given lower precision of calculations can give a more correct result than the procedure of increasing the precision This effect is precisely what we observe in the described example: the correct answer is provided by the calculation with the float16 data type (as opposed to float32), the correct result is achieved with the value of the parameter rel_tol=1e-7, and not with the value 0< rel_tol 1e-8.
Finally, we will consider the comparison of calculating the values of the square root of integers using two operations from the math module: sqrt and pow. Despite the simplicity of the problem statement, calculating the square root of a number, the whole root of a number, and especially the inverse square root of a natural number [29] is a relevant task both to correctness of calculations and the time efficiency of algorithms. For example, consider the range of numbers from 1 to 10000, the data type of floating point numbers is standard float (numpy.float64). We will iterate through the values until we find the number value such that math.sqrt(number) != math.pow(number). Table 5 shows the following numbers (for example, less than 10000) depending on the version of Python and the compiler. The table shows only typical examples of the results of computer calculations of the proposed problem. For some versions (for example, for Python 3.6.5 [GCC 7.3.1], there are no such values under the specified restrictions).
We see from analysis of the binary representation of calculation results in float64 format that the difference is in the least significant bit of the mantissa. The reason is that the first discarded bit in the representation of the number is 1, and with the calculating math.sqrt(2921) this fact is taken into account. And 1 is added to the least significant bit of the mantissa, but not for calculating math.pow(2921, 0.5). For example, both values are calculated with the increment of the least significant bit taken into account in Python 3.6.5 [GCC 7.3.1].
The work presents the results of the authors' research in the field of methodological foundations of teaching academic disciplines for students of the Computer Science specialty. The authors focused their attention on the topic of accuracy of computer calculations of various problems. This problem is relevant not only when performing engineering and scientific calculations, but also, as shown in the article, when solving the simple training problems. The authors focused on the study of the specified problem using the example of the numpy module, as one of the main software tools for modern scientific calculations.
The authors of the presented study focused on considering fairly simple examples of computer calculations that are completely correct from the point of view of the standards of computer work with numbers and arithmetic operations, the results are incorrect to a human performer. It has been demonstrated that mathematically correct calculations are often impossible to reproduce correctly within the IEEE-754 standard due to the fundamental limitations of the memory allocated to represent a specific number.
Using specific examples, it has been shown how the use of special Python computation features, such as the decimal module, can nevertheless lead to incorrect results. The use of external modules such as mpmath requires further research.
The work also gives paradoxical examples, when less accurate calculations allow you to get a more correct result. Such a statement contradicts the classical practices of calculation methods, and requires a more thorough study of the specifics of computer calculations.
It’s important, for the training future specialists in the field of Computer Science during the teaching of subjects related to the accuracy of computer calculations and the fundamental limitations of computer calculations, to present all facts and principles from a single point of view so that students do not there was a false opinion about the arbitrariness of the results of computer calculations and the randomness of the corresponding errors of computer calculations. For a Computer Science professional, it is important not only to understand the basic principles of working with the floating-point number format from the standpoint of the generally accepted IEEE standard, but also to know, even at a basic level, possible modern alternatives, such as, for example, the concept of Posit.
All computer calculations and implementation of computer-mathematical models were supported by the computational capacities of the Laboratory of Machine Learning and Intellectual Data Analysis of the Vasyl’ Stus Donetsk National University. The laboratory was founded in 2020 by Artem Baiev to conduct academic research and educational work. [17] D.M. Russinoff, Formal Verification of Floating-Point Hardware Design, Springer
Cham, 2019. doi:10.1007/978-3-319-95513-1 [18] S. Boldo and G. Melquiond, Computer Arithmetic and Rormal Proofs, ISTE Press
Elsevier, 2017. doi.org:10.1016/C2015-0-01301-6 [19] M.T. Nakao, M. Plum, and Y. Watanabe, "Numerical Verification Methods and Computer-Assisted Proofs for Partial Differential Equations", Springer Singapore, 2019. doi:10.1007/978-981-13-7669-6 [20] O. Vietrov and R. Bilous, "Special methods of increasing the accuracy of computer calculations", 2022 IEEE 3rd KhPI Week on Advanced Technology (KhPIWeek) (2022), 1–5. doi:10.1109/KhPIWeek57572.2022.9916383 [21] O. Vietrov and R. Bilous, "Study of the Convergence of Muller's Sequence Computer Calculations", 2021 IEEE 3rd Ukraine Conference on Electrical and Computer Engineering (UKRCON) (2021), 547–551. doi:10.1109/UKRCON53503.2021.9575546 [22] S.M. Rump, "Algorithms for Verified Inclusions: Theory and Practice", Reliability in Computing: the Role of Interval Methods in Scientific Computing (1988), 109–126. doi:10.15480/882.316 [23] L.G. de la Fraga, "Differential Evolution under FixedPoint Arithmetic and FP16
Numbers", Math. Comput. Appl., 26 (1), 13 (2021). doi:10.3390/mca26010013 [24] Integer Objects, 2024. URL:https://docs.python.org/3/c-api/long.html [25] Data types, 2023. URL:https://numpy.org/doc/stable/user/basics.types.html [26] R. Gray, Did Homer Simpson discover the HIGGS BOSON?, 2015. URL: https://www.dailymail.co.uk/sciencetech/article-2975606/Did-Homer-Simpsondiscover-HIGGS-BOSON-Maths-1998-episode-predicts-particle-s-mass-14-yearsCERN.html [27] D. Snierson, The Simpsons: A somewhat complete history of 'Homer?' from 'Treehouse of Horror VI', 2018. URL:https://ew.com/tv/2018/10/19/the-simpsonshomer-cubed-treehouse-of-horror-vi/ [28] S. Singh, The Simpsons and Their Mathematical Secrets, Bloomsbury USA, 2013. [29] C.J. Walczyk, L. Moroz, and J. Cieslinski, "A Modification of the Fast Inverse Square Root Algorithm", Computation, 7(3):41 (2019). doi:10.3390/computation7030041 [30] Decimal fixed point and floating point arithmetic, 2024.
URL:https://docs.python.org/3/library/decimal.html