{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Computational Errors notebook" ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "import numpy as np #initialization\n", "import matplotlib.pyplot as plt\n", "import sys # Note the new module!\n", "from scipy.integrate import odeint\n", " \n", "%matplotlib inline " ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Implementation information for floating point and integer values" ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Help on built-in module sys:\n", "\n", "NAME\n", " sys\n", "\n", "MODULE REFERENCE\n", " https://docs.python.org/3.7/library/sys\n", " \n", " The following documentation is automatically generated from the Python\n", " source files. It may be incomplete, incorrect or include features that\n", " are considered implementation detail and may vary between Python\n", " implementations. When in doubt, consult the module reference at the\n", " location listed above.\n", "\n", "DESCRIPTION\n", " This module provides access to some objects used or maintained by the\n", " interpreter and to functions that interact strongly with the interpreter.\n", " \n", " Dynamic objects:\n", " \n", " argv -- command line arguments; argv[0] is the script pathname if known\n", " path -- module search path; path[0] is the script directory, else ''\n", " modules -- dictionary of loaded modules\n", " \n", " displayhook -- called to show results in an interactive session\n", " excepthook -- called to handle any uncaught exception other than SystemExit\n", " To customize printing in an interactive session or to install a custom\n", " top-level exception handler, assign other functions to replace these.\n", " \n", " stdin -- standard input file object; used by input()\n", " stdout -- standard output file object; used by print()\n", " stderr -- standard error object; used for error messages\n", " By assigning other file objects (or objects that behave like files)\n", " to these, it is possible to redirect all of the interpreter's I/O.\n", " \n", " last_type -- type of last uncaught exception\n", " last_value -- value of last uncaught exception\n", " last_traceback -- traceback of last uncaught exception\n", " These three are only available in an interactive session after a\n", " traceback has been printed.\n", " \n", " Static objects:\n", " \n", " builtin_module_names -- tuple of module names built into this interpreter\n", " copyright -- copyright notice pertaining to this interpreter\n", " exec_prefix -- prefix used to find the machine-specific Python library\n", " executable -- absolute path of the executable binary of the Python interpreter\n", " float_info -- a struct sequence with information about the float implementation.\n", " float_repr_style -- string indicating the style of repr() output for floats\n", " hash_info -- a struct sequence with information about the hash algorithm.\n", " hexversion -- version information encoded as a single integer\n", " implementation -- Python implementation information.\n", " int_info -- a struct sequence with information about the int implementation.\n", " maxsize -- the largest supported length of containers.\n", " maxunicode -- the value of the largest Unicode code point\n", " platform -- platform identifier\n", " prefix -- prefix used to find the Python library\n", " thread_info -- a struct sequence with information about the thread implementation.\n", " version -- the version of this interpreter as a string\n", " version_info -- version information as a named tuple\n", " __stdin__ -- the original stdin; don't touch!\n", " __stdout__ -- the original stdout; don't touch!\n", " __stderr__ -- the original stderr; don't touch!\n", " __displayhook__ -- the original displayhook; don't touch!\n", " __excepthook__ -- the original excepthook; don't touch!\n", " \n", " Functions:\n", " \n", " displayhook() -- print an object to the screen, and save it in builtins._\n", " excepthook() -- print an exception and its traceback to sys.stderr\n", " exc_info() -- return thread-safe information about the current exception\n", " exit() -- exit the interpreter by raising SystemExit\n", " getdlopenflags() -- returns flags to be used for dlopen() calls\n", " getprofile() -- get the global profiling function\n", " getrefcount() -- return the reference count for an object (plus one :-)\n", " getrecursionlimit() -- return the max recursion depth for the interpreter\n", " getsizeof() -- return the size of an object in bytes\n", " gettrace() -- get the global debug tracing function\n", " setcheckinterval() -- control how often the interpreter checks for events\n", " setdlopenflags() -- set the flags to be used for dlopen() calls\n", " setprofile() -- set the global profiling function\n", " setrecursionlimit() -- set the max recursion depth for the interpreter\n", " settrace() -- set the global debug tracing function\n", "\n", "FUNCTIONS\n", " __breakpointhook__ = breakpointhook(...)\n", " breakpointhook(*args, **kws)\n", " \n", " This hook function is called by built-in breakpoint().\n", " \n", " __displayhook__ = displayhook(...)\n", " displayhook(object) -> None\n", " \n", " Print an object to sys.stdout and also save it in builtins._\n", " \n", " __excepthook__ = excepthook(...)\n", " excepthook(exctype, value, traceback) -> None\n", " \n", " Handle an exception by displaying it with a traceback on sys.stderr.\n", " \n", " breakpointhook(...)\n", " breakpointhook(*args, **kws)\n", " \n", " This hook function is called by built-in breakpoint().\n", " \n", " call_tracing(...)\n", " call_tracing(func, args) -> object\n", " \n", " Call func(*args), while tracing is enabled. The tracing state is\n", " saved, and restored afterwards. This is intended to be called from\n", " a debugger from a checkpoint, to recursively debug some other code.\n", " \n", " callstats(...)\n", " callstats() -> tuple of integers\n", " \n", " Return a tuple of function call statistics, if CALL_PROFILE was defined\n", " when Python was built. Otherwise, return None.\n", " \n", " When enabled, this function returns detailed, implementation-specific\n", " details about the number of function calls executed. The return value is\n", " a 11-tuple where the entries in the tuple are counts of:\n", " 0. all function calls\n", " 1. calls to PyFunction_Type objects\n", " 2. PyFunction calls that do not create an argument tuple\n", " 3. PyFunction calls that do not create an argument tuple\n", " and bypass PyEval_EvalCodeEx()\n", " 4. PyMethod calls\n", " 5. PyMethod calls on bound methods\n", " 6. PyType calls\n", " 7. PyCFunction calls\n", " 8. generator calls\n", " 9. All other calls\n", " 10. Number of stack pops performed by call_function()\n", " \n", " exc_info(...)\n", " exc_info() -> (type, value, traceback)\n", " \n", " Return information about the most recent exception caught by an except\n", " clause in the current stack frame or in an older stack frame.\n", " \n", " exit(...)\n", " exit([status])\n", " \n", " Exit the interpreter by raising SystemExit(status).\n", " If the status is omitted or None, it defaults to zero (i.e., success).\n", " If the status is an integer, it will be used as the system exit status.\n", " If it is another kind of object, it will be printed and the system\n", " exit status will be one (i.e., failure).\n", " \n", " get_asyncgen_hooks(...)\n", " get_asyncgen_hooks()\n", " \n", " Return a namedtuple of installed asynchronous generators hooks (firstiter, finalizer).\n", " \n", " get_coroutine_origin_tracking_depth()\n", " Check status of origin tracking for coroutine objects in this thread.\n", " \n", " get_coroutine_wrapper(...)\n", " get_coroutine_wrapper()\n", " \n", " Return the wrapper for coroutine objects set by sys.set_coroutine_wrapper.\n", " \n", " getallocatedblocks(...)\n", " getallocatedblocks() -> integer\n", " \n", " Return the number of memory blocks currently allocated, regardless of their\n", " size.\n", " \n", " getcheckinterval(...)\n", " getcheckinterval() -> current check interval; see setcheckinterval().\n", " \n", " getdefaultencoding(...)\n", " getdefaultencoding() -> string\n", " \n", " Return the current default string encoding used by the Unicode \n", " implementation.\n", " \n", " getdlopenflags(...)\n", " getdlopenflags() -> int\n", " \n", " Return the current value of the flags that are used for dlopen calls.\n", " The flag constants are defined in the os module.\n", " \n", " getfilesystemencodeerrors(...)\n", " getfilesystemencodeerrors() -> string\n", " \n", " Return the error mode used to convert Unicode filenames in\n", " operating system filenames.\n", " \n", " getfilesystemencoding(...)\n", " getfilesystemencoding() -> string\n", " \n", " Return the encoding used to convert Unicode filenames in\n", " operating system filenames.\n", " \n", " getprofile(...)\n", " getprofile()\n", " \n", " Return the profiling function set with sys.setprofile.\n", " See the profiler chapter in the library manual.\n", " \n", " getrecursionlimit(...)\n", " getrecursionlimit()\n", " \n", " Return the current value of the recursion limit, the maximum depth\n", " of the Python interpreter stack. This limit prevents infinite\n", " recursion from causing an overflow of the C stack and crashing Python.\n", " \n", " getrefcount(...)\n", " getrefcount(object) -> integer\n", " \n", " Return the reference count of object. The count returned is generally\n", " one higher than you might expect, because it includes the (temporary)\n", " reference as an argument to getrefcount().\n", " \n", " getsizeof(...)\n", " getsizeof(object, default) -> int\n", " \n", " Return the size of object in bytes.\n", " \n", " getswitchinterval(...)\n", " getswitchinterval() -> current thread switch interval; see setswitchinterval().\n", " \n", " gettrace(...)\n", " gettrace()\n", " \n", " Return the global debug tracing function set with sys.settrace.\n", " See the debugger chapter in the library manual.\n", " \n", " intern(...)\n", " intern(string) -> string\n", " \n", " ``Intern'' the given string. This enters the string in the (global)\n", " table of interned strings whose purpose is to speed up dictionary lookups.\n", " Return the string itself or the previously interned string object with the\n", " same value.\n", " \n", " is_finalizing(...)\n", " is_finalizing()\n", " Return True if Python is exiting.\n", " \n", " set_asyncgen_hooks(...)\n", " set_asyncgen_hooks(*, firstiter=None, finalizer=None)\n", " \n", " Set a finalizer for async generators objects.\n", " \n", " set_coroutine_origin_tracking_depth(depth)\n", " Enable or disable origin tracking for coroutine objects in this thread.\n", " \n", " Coroutine objects will track 'depth' frames of traceback information about\n", " where they came from, available in their cr_origin attribute. Set depth of 0\n", " to disable.\n", " \n", " set_coroutine_wrapper(...)\n", " set_coroutine_wrapper(wrapper)\n", " \n", " Set a wrapper for coroutine objects.\n", " \n", " setcheckinterval(...)\n", " setcheckinterval(n)\n", " \n", " Tell the Python interpreter to check for asynchronous events every\n", " n instructions. This also affects how often thread switches occur.\n", " \n", " setdlopenflags(...)\n", " setdlopenflags(n) -> None\n", " \n", " Set the flags used by the interpreter for dlopen calls, such as when the\n", " interpreter loads extension modules. Among other things, this will enable\n", " a lazy resolving of symbols when importing a module, if called as\n", " sys.setdlopenflags(0). To share symbols across extension modules, call as\n", " sys.setdlopenflags(os.RTLD_GLOBAL). Symbolic names for the flag modules\n", " can be found in the os module (RTLD_xxx constants, e.g. os.RTLD_LAZY).\n", " \n", " setprofile(...)\n", " setprofile(function)\n", " \n", " Set the profiling function. It will be called on each function call\n", " and return. See the profiler chapter in the library manual.\n", " \n", " setrecursionlimit(...)\n", " setrecursionlimit(n)\n", " \n", " Set the maximum depth of the Python interpreter stack to n. This\n", " limit prevents infinite recursion from causing an overflow of the C\n", " stack and crashing Python. The highest possible limit is platform-\n", " dependent.\n", " \n", " setswitchinterval(...)\n", " setswitchinterval(n)\n", " \n", " Set the ideal thread switching delay inside the Python interpreter\n", " The actual frequency of switching threads can be lower if the\n", " interpreter executes long sequences of uninterruptible code\n", " (this is implementation-specific and workload-dependent).\n", " \n", " The parameter must represent the desired switching delay in seconds\n", " A typical value is 0.005 (5 milliseconds).\n", " \n", " settrace(...)\n", " settrace(function)\n", " \n", " Set the global debug tracing function. It will be called on each\n", " function call. See the debugger chapter in the library manual.\n", "\n", "DATA\n", " __stderr__ = <_io.TextIOWrapper name='' mode='w' encoding='UTF...\n", " __stdin__ = <_io.TextIOWrapper name='' mode='r' encoding='UTF-8...\n", " __stdout__ = <_io.TextIOWrapper name='' mode='w' encoding='UTF...\n", " abiflags = 'm'\n", " api_version = 1013\n", " argv = ['/Users/nemenman/anaconda3/lib/python3.7/site-packages/ipykern...\n", " base_exec_prefix = '/Users/nemenman/anaconda3'\n", " base_prefix = '/Users/nemenman/anaconda3'\n", " builtin_module_names = ('_abc', '_ast', '_codecs', '_collections', '_f...\n", " byteorder = 'little'\n", " copyright = 'Copyright (c) 2001-2018 Python Software Foundati...ematis...\n", " displayhook = \n", " dont_write_bytecode = False\n", " exec_prefix = '/Users/nemenman/anaconda3'\n", " executable = '/Users/nemenman/anaconda3/bin/python'\n", " flags = sys.flags(debug=0, inspect=0, interactive=0, opt...ation=1, is...\n", " float_info = sys.float_info(max=1.7976931348623157e+308, max_...epsilo...\n", " float_repr_style = 'short'\n", " hash_info = sys.hash_info(width=64, modulus=2305843009213693...iphash2...\n", " hexversion = 50790896\n", " implementation = namespace(_multiarch='darwin', cache_tag='cpytho...in...\n", " int_info = sys.int_info(bits_per_digit=30, sizeof_digit=4)\n", " maxsize = 9223372036854775807\n", " maxunicode = 1114111\n", " meta_path = [, , \n", " stdin = <_io.TextIOWrapper name='' mode='r' encoding='UTF-8'>\n", " stdout = \n", " thread_info = sys.thread_info(name='pthread', lock='mutex+cond', versi...\n", " version = '3.7.1 (default, Dec 14 2018, 13:28:58) \\n[Clang 4.0.1 (tags...\n", " version_info = sys.version_info(major=3, minor=7, micro=1, releaseleve...\n", " warnoptions = []\n", "\n", "FILE\n", " (built-in)\n", "\n", "\n", "namespace(_multiarch='darwin', cache_tag='cpython-37', hexversion=50790896, name='cpython', version=sys.version_info(major=3, minor=7, micro=1, releaselevel='final', serial=0))\n", "sys.int_info(bits_per_digit=30, sizeof_digit=4)\n", "sys.float_info(max=1.7976931348623157e+308, max_exp=1024, max_10_exp=308, min=2.2250738585072014e-308, min_exp=-1021, min_10_exp=-307, dig=15, mant_dig=53, epsilon=2.220446049250313e-16, radix=2, rounds=1)\n" ] } ], "source": [ "dir(sys)\n", "help(sys)\n", "print(sys.implementation)\n", "print(sys.int_info)\n", "print(sys.float_info)\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Pitfalls of doing floating point aruthmetics on a computer" ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "1.0000000000000002 1.0\n", "10.0\n", "1000.0000000000002 1000.0\n", "2.2250738585072014e-308 2.225073858507e-311 2.223e-321 0.0\n", "1.7976931348623157e+308 inf inf inf\n" ] } ], "source": [ "# epislon is the smallest possible number you can add to 1 an still see a difference\n", "print(1+sys.float_info.epsilon, 1+sys.float_info.epsilon/2)\n", "print(10+sys.float_info.epsilon)\n", "# epsilon is a measure of the relative error -- you can see a difference by adding x*epsilon to any x \n", "print(1e3+1e3*sys.float_info.epsilon, 1e3+1e3*sys.float_info.epsilon/4)\n", "\n", "# note that the handling of the minimum and the maximum numbers is different\n", "# This is because the minimum number is 1.11111 (binary)... 2^(min_exp); this number can be further decreased, \n", "# but the accuracy will be lost. In contrast, the maximum number is 1.1111(binary) 2^{max_exp}, \n", "# and this cannot be increased\n", "print(sys.float_info.min, sys.float_info.min/1e3, sys.float_info.min/1e13, sys.float_info.min/1e16)\n", "print(sys.float_info.max, sys.float_info.max*(1+ sys.float_info.epsilon), sys.float_info.max*1e3, sys.float_info.max*1e16)\n" ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "1e+20\n", "1.00000000000001e+20\n" ] } ], "source": [ "# it's important to not add small numbers to large numbers.\n", "x = np.ones(int(1e6))\n", "x[0] = 1e20\n", "print(np.sum(x)) # this results in completely neglecting the million entries after the first one\n", "# because each entry is smaller than epsilon compared to the first one \n", "# in contrast, in the following line, one first sums up the small numbersm making them bigger, and then adds \n", "# them to the big number\n", "print(np.sum(np.sort(x)))" ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "6.000000000000002\n", "False\n" ] } ], "source": [ "# Comparing floating point numbers will lead to errors: because of truncation, numbers\n", "# that you want to be the same will probably not be the same\n", "x=0.002*np.ones(3000)\n", "y= np.sum(x)\n", "print(y)\n", "print(y==6.0)" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [] } ], "metadata": { "kernelspec": { "display_name": "Python 3", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.7.1" } }, "nbformat": 4, "nbformat_minor": 2 }