1# -*- coding: utf-8 -*- 2# cython: language_level=3 3 4import sympy 5import mpmath 6 7from math import log, ceil 8import string 9 10import typing 11 12C = log(10, 2) # ~ 3.3219280948873626 13 14 15# Number of bits of machine precision 16machine_precision = 53 17 18 19machine_epsilon = 2 ** (1 - machine_precision) 20 21 22def reconstruct_digits(bits) -> int: 23 """ 24 Number of digits needed to reconstruct a number with given bits of precision. 25 26 >>> reconstruct_digits(53) 27 17 28 """ 29 return int(ceil(bits / C) + 1) 30 31 32class PrecisionValueError(Exception): 33 pass 34 35 36class SpecialValueError(Exception): 37 def __init__(self, name) -> None: 38 self.name = name 39 40 41def _get_float_inf(value, evaluation) -> typing.Optional[float]: 42 value = value.evaluate(evaluation) 43 if value.has_form("DirectedInfinity", 1): 44 if value.leaves[0].get_int_value() == 1: 45 return float("inf") 46 elif value.leaves[0].get_int_value() == -1: 47 return float("-inf") 48 else: 49 return None 50 return value.round_to_float(evaluation) 51 52 53def get_precision(value, evaluation) -> typing.Optional[int]: 54 if value.get_name() == "System`MachinePrecision": 55 return None 56 else: 57 from mathics.core.expression import Symbol, MachineReal 58 59 dmin = _get_float_inf(Symbol("$MinPrecision"), evaluation) 60 dmax = _get_float_inf(Symbol("$MaxPrecision"), evaluation) 61 d = value.round_to_float(evaluation) 62 assert dmin is not None and dmax is not None 63 if d is None: 64 evaluation.message("N", "precbd", value) 65 elif d < dmin: 66 dmin = int(dmin) 67 evaluation.message("N", "precsm", value, MachineReal(dmin)) 68 return dmin 69 elif d > dmax: 70 dmax = int(dmax) 71 evaluation.message("N", "preclg", value, MachineReal(dmax)) 72 return dmax 73 else: 74 return d 75 raise PrecisionValueError() 76 77 78def get_type(value) -> typing.Optional[str]: 79 if isinstance(value, sympy.Integer): 80 return "z" 81 elif isinstance(value, sympy.Rational): 82 return "q" 83 elif isinstance(value, sympy.Float) or isinstance(value, mpmath.mpf): 84 return "f" 85 elif ( 86 isinstance(value, sympy.Expr) and value.is_number and not value.is_real 87 ) or isinstance(value, mpmath.mpc): 88 return "c" 89 else: 90 return None 91 92 93def sameQ(v1, v2) -> bool: 94 """Mathics SameQ""" 95 return get_type(v1) == get_type(v2) and v1 == v2 96 97 98def dps(prec) -> int: 99 return max(1, int(round(int(prec) / C - 1))) 100 101 102def prec(dps) -> int: 103 return max(1, int(round((int(dps) + 1) * C))) 104 105 106def min_prec(*args): 107 result = None 108 for arg in args: 109 prec = arg.get_precision() 110 if result is None or (prec is not None and prec < result): 111 result = prec 112 return result 113 114 115def pickle_mp(value): 116 return (get_type(value), str(value)) 117 118 119def unpickle_mp(value): 120 type, value = value 121 if type == "z": 122 return sympy.Integer(value) 123 elif type == "q": 124 return sympy.Rational(value) 125 elif type == "f": 126 return sympy.Float(value) 127 else: 128 return value 129 130 131# algorithm based on 132# http://stackoverflow.com/questions/5110177/how-to-convert-floating-point-number-to-base-3-in-python # nopep8 133 134 135def convert_base(x, base, precision=10) -> str: 136 sign = -1 if x < 0 else 1 137 x *= sign 138 139 length_of_int = 0 if x == 0 else int(log(x, base)) 140 iexps = list(range(length_of_int, -1, -1)) 141 digits = string.digits + string.ascii_lowercase 142 143 if base > len(digits): 144 raise ValueError 145 146 def convert(x, base, exponents): 147 out = [] 148 for e in exponents: 149 d = int(x / (base ** e)) 150 x -= d * (base ** e) 151 out.append(digits[d]) 152 if x == 0 and e < 0: 153 break 154 return out 155 156 int_part = convert(int(x), base, iexps) 157 if sign == -1: 158 int_part.insert(0, "-") 159 160 if isinstance(x, (float, sympy.Float)): 161 fexps = list(range(-1, -int(precision + 1), -1)) 162 real_part = convert(x - int(x), base, fexps) 163 164 return "%s.%s" % ("".join(int_part), "".join(real_part)) 165 elif isinstance(x, int): 166 return "".join(int_part) 167 else: 168 raise TypeError(x) 169 170 171def convert_int_to_digit_list(x, base) -> typing.List[int]: 172 if x == 0: 173 return [0] 174 175 x = abs(x) 176 177 length_of_int = int(log(x, base)) + 1 178 iexps = list(range(length_of_int, -1, -1)) 179 180 def convert(x, base, exponents): 181 out = [] 182 for e in exponents: 183 d = int(x // (base ** e)) 184 x -= d * (base ** e) 185 if out or d != 0: # drop any leading zeroes 186 out.append(d) 187 if x == 0 and e < 0: 188 break 189 return out 190 191 return convert(x, base, iexps) 192