1# cython: language_level=3 2 3""" 4Mathematical Constants 5 6Numeric, Arithmetic, or Symbolic constants like Pi, E, or Infinity. 7""" 8 9import math 10import mpmath 11import numpy 12import sympy 13 14from mathics.version import __version__ # noqa used in loading to check consistency. 15 16from mathics.builtin.base import Predefined, SympyObject 17from mathics.core.expression import ( 18 MachineReal, 19 PrecisionReal, 20 Symbol, 21 strip_context, 22) 23from mathics.core.numbers import get_precision, PrecisionValueError, machine_precision 24 25 26def mp_constant(fn: str, d=None) -> mpmath.mpf: 27 """ 28 Return the mpmath constant _fn_ with integer precision _d_. 29 """ 30 if d is None: 31 return getattr(mpmath, fn)() 32 else: 33 # TODO: In some functions like Pi, you can 34 # ask for a certain number of digits, but the 35 # accuracy will be less than that. Figure out 36 # what's up and compensate somehow. 37 mpmath.mp.dps = int_d = int(d * 3.321928) 38 return getattr(mpmath, fn)(prec=int_d) 39 40 41def mp_convert_constant(obj, **kwargs): 42 if isinstance(obj, mpmath.ctx_mp_python._constant): 43 prec = kwargs.get("prec", None) 44 if prec is not None: 45 return sympy.Float(obj(prec=prec)) 46 return sympy.Float(obj) 47 return obj 48 49 50def numpy_constant(name: str, d=None) -> float: 51 if d: 52 # by mmatera: Here I have a question: 53 # 0.0123`2 should be rounded to 54 # 0.01 or to 0.0123? 55 # (absolute versus relative accuracy) 56 val = getattr(numpy, name) 57 val = numpy.round(val, d) 58 return val 59 else: 60 return getattr(numpy, name) 61 62 63def sympy_constant(fn, d=None): 64 return getattr(sympy, fn).evalf(n=d) 65 66 67class _Constant_Common(Predefined): 68 69 attributes = ("Constant", "Protected", "ReadProtected") 70 nargs = 0 71 options = {"Method": "Automatic"} 72 73 def apply_N(self, precision, evaluation, options={}): 74 "N[%(name)s, precision_?NumericQ, OptionsPattern[%(name)s]]" 75 76 preference = self.get_option(options, "Method", evaluation).get_string_value() 77 if preference == "Automatic": 78 return self.get_constant(precision, evaluation) 79 else: 80 return self.get_constant(precision, evaluation, preference) 81 82 def apply_N2(self, evaluation, options={}): 83 "N[%(name)s, OptionsPattern[%(name)s]]" 84 return self.apply_N(None, evaluation, options) 85 86 def is_constant(self) -> bool: 87 return True 88 89 def get_constant(self, precision, evaluation, preference=None): 90 # first, determine the precision 91 machine_d = int(0.30103 * machine_precision) 92 d = None 93 if precision: 94 try: 95 d = get_precision(precision, evaluation) 96 except PrecisionValueError: 97 pass 98 99 if d is None: 100 d = machine_d 101 102 # If preference not especified, determine it 103 # from the precision. 104 if preference is None: 105 if d <= machine_d: 106 preference = "numpy" 107 else: 108 preference = "mpmath" 109 # If preference is not valid, send a message and return. 110 if not (preference in ("sympy", "numpy", "mpmath")): 111 evaluation.message(f'{preference} not in ("sympy", "numpy", "mpmath")') 112 return 113 # Try to determine the numeric value 114 value = None 115 if preference == "mpmath" and not hasattr(self, "mpmath_name"): 116 preference = "numpy" 117 elif preference == "sympy" and not hasattr(self, "sympy_name"): 118 preference = "numpy" 119 120 if preference == "numpy" and not hasattr(self, "numpy_name"): 121 if hasattr(self, "sympy_name"): 122 preference = "sympy" 123 elif hasattr(self, "mpmath_name"): 124 preference = "mpmath" 125 else: 126 preference = "" 127 if preference == "numpy": 128 value = numpy_constant(self.numpy_name) 129 if d == machine_d: 130 return MachineReal(value) 131 if preference == "sympy": 132 value = sympy_constant(self.sympy_name, d + 2) 133 if preference == "mpmath": 134 value = mp_constant(self.mpmath_name, d * 2) 135 if value: 136 return PrecisionReal(sympy.Float(str(value), d)) 137 # If the value is not available, return none 138 # and keep it unevaluated. 139 return 140 141 142class MPMathConstant(_Constant_Common): 143 """Representation of a constant in mpmath, e.g. Pi, E, I, etc.""" 144 145 # Subclasses should define this. 146 mpmath_name = None 147 148 mathics_to_mpmath = {} 149 150 def __init__(self, *args, **kwargs): 151 super().__init__(*args, **kwargs) 152 if self.mpmath_name is None: 153 self.mpmath_name = strip_context(self.get_name()).lower() 154 self.mathics_to_mpmath[self.__class__.__name__] = self.mpmath_name 155 156 def to_mpmath(self, args): 157 if self.mpmath_name is None or len(args) != 0: 158 return None 159 return getattr(mpmath, self.mpmath_name) 160 161 162class NumpyConstant(_Constant_Common): 163 """Representation of a constant in numpy, e.g. Pi, E, etc.""" 164 165 # Subclasses should define this. 166 numpy_name = None 167 168 mathics_to_numpy = {} 169 170 def __init__(self, *args, **kwargs): 171 super().__init__(*args, **kwargs) 172 if self.numpy_name is None: 173 self.numpy_name = strip_context(self.get_name()).lower() 174 self.mathics_to_numpy[self.__class__.__name__] = self.numpy_name 175 176 def to_numpy(self, args): 177 if self.numpy_name is None or len(args) != 0: 178 return None 179 return self.get_constant() 180 181 182class SympyConstant(_Constant_Common, SympyObject): 183 """Representation of a constant in Sympy, e.g. Pi, E, I, Catalan, etc.""" 184 185 # Subclasses should define this. 186 sympy_name = None 187 188 def to_sympy(self, expr=None, **kwargs): 189 if expr is None or expr.is_atom(): 190 result = getattr(sympy, self.sympy_name) 191 if kwargs.get("evaluate", False): 192 result = mp_convert_constant(result, **kwargs) 193 return result 194 else: 195 # there is no "native" SymPy expression for e.g. E[x] 196 return None 197 198 199class Catalan(MPMathConstant, SympyConstant): 200 """ 201 <dl> 202 <dt>'Catalan' 203 <dd>is Catalan's constant with numerical value \u2243 0.915966. 204 </dl> 205 206 >> Catalan // N 207 = 0.915965594177219 208 209 >> N[Catalan, 20] 210 = 0.91596559417721901505 211 """ 212 213 mpmath_name = "catalan" 214 # numpy_name = "catalan" ## This is not defined in numpy 215 sympy_name = "Catalan" 216 217 218class ComplexInfinity(SympyConstant): 219 """ 220 <dl> 221 <dt>'ComplexInfinity' 222 <dd>represents an infinite complex quantity of undetermined direction. 223 </dl> 224 225 >> 1 / ComplexInfinity 226 = 0 227 >> ComplexInfinity * Infinity 228 = ComplexInfinity 229 >> FullForm[ComplexInfinity] 230 = DirectedInfinity[] 231 232 ## Issue689 233 #> ComplexInfinity + ComplexInfinity 234 : Indeterminate expression ComplexInfinity + ComplexInfinity encountered. 235 = Indeterminate 236 #> ComplexInfinity + Infinity 237 : Indeterminate expression ComplexInfinity + Infinity encountered. 238 = Indeterminate 239 """ 240 241 sympy_name = "zoo" 242 243 rules = { 244 "ComplexInfinity": "DirectedInfinity[]", 245 } 246 247 248class Degree(MPMathConstant, NumpyConstant, SympyConstant): 249 """ 250 <dl> 251 <dt>'Degree' 252 <dd>is the number of radians in one degree. It has a numerical value of \u03c0 / 180. 253 </dl> 254 >> Cos[60 Degree] 255 = 1 / 2 256 257 Degree has the value of Pi / 180 258 >> Degree == Pi / 180 259 = True 260 261 #> Cos[Degree[x]] 262 = Cos[Degree[x]] 263 264 ## Issue 274 265 #> \\[Degree] == ° == Degree 266 = True 267 268 #> N[Degree] 269 = 0.0174533 270 #> N[Degree, 30] 271 = 0.0174532925199432957692369076849 272 """ 273 274 mpmath_name = "degree" 275 276 def to_sympy(self, expr=None, **kwargs): 277 if expr == Symbol("System`Degree"): 278 # return mpmath.degree 279 return sympy.pi / 180 280 281 def to_numpy(self, expr=None, **kwargs): 282 if expr == Symbol("System`Degree"): 283 # return mpmath.degree 284 return numpy.pi / 180 285 286 def apply_N(self, precision, evaluation, options={}): 287 "N[Degree, precision_, OptionsPattern[%(name)s]]" 288 try: 289 if precision: 290 d = get_precision(precision, evaluation) 291 else: 292 d = get_precision(Symbol("System`MachinePrecision"), evaluation) 293 except PrecisionValueError: 294 return 295 296 # FIXME: There are all sorts of interactions between in the trig functions, 297 # that are expected to work out right. Until we have convertion between 298 # mpmath and sympy worked out so that values can be made the to the same 299 # precision and compared. we have to not use mpmath right now. 300 # return self.get_constant(precision, evaluation, preference="mpmath") 301 302 if d is None: 303 return MachineReal(math.pi / 180) 304 else: 305 return PrecisionReal((sympy.pi / 180).n(d)) 306 307 308class E(MPMathConstant, NumpyConstant, SympyConstant): 309 """ 310 <dl> 311 <dt>'E'</dt> 312 <dd>is the constant \u2147 with numerical value \u2243 2.71828. 313 </dl> 314 315 >> N[E] 316 = 2.71828 317 >> N[E, 50] 318 = 2.7182818284590452353602874713526624977572470937000 319 320 #> 5. E 321 = 13.5914 322 """ 323 324 mpmath_name = "e" 325 numpy_name = "e" 326 sympy_name = "E" 327 328 def apply_N(self, precision, evaluation, options={}): 329 "N[E, precision_, OptionsPattern[%(name)s]]" 330 return self.get_constant(precision, evaluation) 331 332 333class EulerGamma(MPMathConstant, NumpyConstant, SympyConstant): 334 """ 335 <dl> 336 <dt>'EulerGamma'</dt> 337 <dd>is Euler's constant \u03b3 with numerial value \u2243 0.577216. 338 </dl> 339 340 >> EulerGamma // N 341 = 0.577216 342 343 >> N[EulerGamma, 40] 344 = 0.5772156649015328606065120900824024310422 345 """ 346 347 mpmath_name = "euler" 348 numpy_name = "euler_gamma" 349 sympy_name = "EulerGamma" 350 351 352class Glaisher(MPMathConstant): 353 """ 354 <dl> 355 <dt>'Glaisher'</dt> 356 <dd>is Glaisher's constant, with numerical value \u2243 1.28243. 357 </dl> 358 359 >> N[Glaisher] 360 = 1.28242712910062 361 >> N[Glaisher, 50] 362 = 1.2824271291006226368753425688697917277676889273250 363 # 1.2824271291006219541941391071304678916931152343750 364 """ 365 366 mpmath_name = "glaisher" 367 368 369class GoldenRatio(MPMathConstant, SympyConstant): 370 """ 371 <dl> 372 <dt>'GoldenRatio' 373 <dd>is the golden ratio, \u03D5 = (1+Sqrt[5])/2. 374 </dl> 375 376 >> GoldenRatio // N 377 = 1.61803398874989 378 >> N[GoldenRatio, 40] 379 = 1.618033988749894848204586834365638117720 380 """ 381 382 sympy_name = "GoldenRatio" 383 mpmath_name = "phi" 384 385 386class Indeterminate(SympyConstant): 387 """ 388 <dl> 389 <dt>'Indeterminate'</dt> 390 <dd>represents an indeterminate result. 391 </dl> 392 393 >> 0^0 394 : Indeterminate expression 0 ^ 0 encountered. 395 = Indeterminate 396 397 >> Tan[Indeterminate] 398 = Indeterminate 399 """ 400 401 sympy_name = "nan" 402 403 404class Infinity(SympyConstant): 405 """ 406 <dl> 407 <dt>'Infinity' 408 <dd>represents an infinite real quantity. 409 </dl> 410 411 >> 1 / Infinity 412 = 0 413 >> Infinity + 100 414 = Infinity 415 416 Use 'Infinity' in sum and limit calculations: 417 >> Sum[1/x^2, {x, 1, Infinity}] 418 = Pi ^ 2 / 6 419 420 #> FullForm[Infinity] 421 = DirectedInfinity[1] 422 #> (2 + 3.5*I) / Infinity 423 = 0. 424 #> Infinity + Infinity 425 = Infinity 426 #> Infinity / Infinity 427 : Indeterminate expression 0 Infinity encountered. 428 = Indeterminate 429 """ 430 431 sympy_name = "oo" 432 numpy_name = "Inf" 433 mpmath_name = "inf" 434 python_equivalent = math.inf 435 436 rules = { 437 "Infinity": "DirectedInfinity[1]", 438 "MakeBoxes[Infinity, f:StandardForm|TraditionalForm]": ('"\\[Infinity]"'), 439 } 440 441 442class Khinchin(MPMathConstant): 443 """ 444 <dl> 445 <dt>'Khinchin'</dt> 446 <dd>is Khinchin's constant, with numerical value \u2243 2.68545. 447 </dl> 448 449 >> N[Khinchin] 450 = 2.68545200106531 451 >> N[Khinchin, 50] 452 = 2.6854520010653064453097148354817956938203822939945 453 # = 2.6854520010653075701156922150403261184692382812500 454 """ 455 456 mpmath_name = "khinchin" 457 458 459class Pi(MPMathConstant, SympyConstant): 460 """ 461 <dl> 462 <dt>'Pi'</dt> 463 <dd>is the constant \u03c0. 464 </dl> 465 466 >> N[Pi] 467 = 3.14159 468 469 Force using the value given from numpy to compute Pi. 470 >> N[Pi, Method->"numpy"] 471 = 3.14159 472 473 Force using the value given from sympy to compute Pi to 3 places, 474 two places after the decimal point. 475 476 Note that sympy is the default method. 477 >> N[Pi, 3, Method->"sympy"] 478 = 3.14 479 480 >> N[Pi, 50] 481 = 3.1415926535897932384626433832795028841971693993751 482 >> Attributes[Pi] 483 = {Constant, Protected, ReadProtected} 484 """ 485 486 sympy_name = "pi" 487 mpmath_name = "pi" 488 numpy_name = "pi" 489