wrightomega.py 955 B

1234567891011121314151617181920212223242526272829303132333435363738394041
  1. import numpy as np
  2. try:
  3. import mpmath
  4. except ImportError:
  5. pass
  6. def mpmath_wrightomega(x):
  7. return mpmath.lambertw(mpmath.exp(x), mpmath.mpf('-0.5'))
  8. def wrightomega_series_error(x):
  9. series = x
  10. desired = mpmath_wrightomega(x)
  11. return abs(series - desired) / desired
  12. def wrightomega_exp_error(x):
  13. exponential_approx = mpmath.exp(x)
  14. desired = mpmath_wrightomega(x)
  15. return abs(exponential_approx - desired) / desired
  16. def main():
  17. desired_error = 2 * np.finfo(float).eps
  18. print('Series Error')
  19. for x in [1e5, 1e10, 1e15, 1e20]:
  20. with mpmath.workdps(100):
  21. error = wrightomega_series_error(x)
  22. print(x, error, error < desired_error)
  23. print('Exp error')
  24. for x in [-10, -25, -50, -100, -200, -400, -700, -740]:
  25. with mpmath.workdps(100):
  26. error = wrightomega_exp_error(x)
  27. print(x, error, error < desired_error)
  28. if __name__ == '__main__':
  29. main()