-
Notifications
You must be signed in to change notification settings - Fork 85
Implement exp, log, power and ** in ruby #347
New issue
Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.
By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.
Already on GitHub? Sign in to your account
Conversation
21e90c5 to
5c4d747
Compare
674e7f1 to
b3212f4
Compare
| # Detect overflow/underflow before consuming infinite memory | ||
| if (xn.exponent.abs - 1) * int_part / n > 1000000000000000000 | ||
| return ((xn.exponent > 0) ^ neg ? BigDecimal::INFINITY : BigDecimal(0)) * (int_part.even? || x > 0 ? 1 : -1) | ||
| end |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
This part is tested as a spec in test_bigdecimal.rb
assert_positive_infinite(BigDecimal(3) ** (2**100))
|
|
||
| # exp(x * 10**cnt) = exp(x)**(10**cnt) | ||
| cnt = x > 1 ? x.exponent : 0 | ||
| prec2 = prec + BigDecimal.double_fig + cnt |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
y**10 makes relative error 10 times larger. (y*(1+e))**10 = y**10 * (1+10*e+...)
In the last step of this method cnt.times{y = y**10}, error will be 10**cnt times larger.
Calculation of exp(x) needs cnt more digits.
|
|
||
| # Taylor series for log(x) around 1 | ||
| # log(x) = -log((1 + X) / (1 - X)) where X = (x - 1) / (x + 1) | ||
| # log(x) = 2 * (X + X**3 / 3 + X**5 / 5 + X**7 / 7 + ...) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
This is what the original log implementation in C was doing.
Converts x: 0.1..10 → (x-1)/(x+1): -0.818..0.818 to make Taylor series converge. (but converge of 0.818**n was very slow. This is the slowness reason of log(10))
| prec += BigDecimal.double_fig | ||
|
|
||
| # log(x) = log(sqrt(sqrt(sqrt(sqrt(x))))) * 2**sqrt_steps | ||
| sqrt_steps = [2 * Integer.sqrt(prec) + 3 * x_minus_one_exponent, 0].max |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
To minimize the computation cost, we need to minimize this:
sqrt_cost * sqrt_steps + mult_cost * taylor_steps
Taylor steps can be approximated as:
taylor_steps = const1 * prec / (-log10(|x-1| / 2**sqrt_steps))
Assuming sqrt_cost and mult_cost are ideally the same order, (note that currently it is not)
We need to minimize something like this:
sqrt_steps + const1 * prec / (log10(|x-1|) + log10(2) * sqrt_step)
sqrt_steps that minimizes the cost should be
sqrt_steps = const2 * sqrt(prec) + log10(|x-1|) / log10(2)
≒ const2 * sqrt(prec) + 3 * x_minus_one_exponent
5f90e04 to
94deded
Compare
94deded to
dc66026
Compare
dc66026 to
29c15a5
Compare
Co-authored-by: Kenta Murata <[email protected]>
29c15a5 to
00aee7c
Compare
Writing in Ruby is easier than C.
This will bring both maintainability and performance (because implementing faster algorithm is easy).
Eliminates problem like #340 in the long run.(already fixed)Use case of BigMath.log BigMath.exp are high-precision calculation (if not, using Float is enough).
Cost of multiplications/divisions are dominant. There is not much need to implement it in C.
Fixes #83, #345, #352
Fixes power,
log and exp part of #340(already fixed)Benchmark
Makes
BigMath.log(10, prec)which is the bottleneck ofBigDecimal#**fasterMakes BigMath.exp(large_x) faster (#82)