diff --git a/.formatter.exs b/.formatter.exs index cc79f69..56f192e 100644 --- a/.formatter.exs +++ b/.formatter.exs @@ -1,3 +1,3 @@ [ - inputs: ["{bench,mix,.formatter}.exs", "{config,lib,test}/**/*.{ex,exs}"] + inputs: ["*.exs", "{config,lib,test}/**/*.{ex,exs}"] ] diff --git a/CHANGELOG.md b/CHANGELOG.md index 2317c2d..ff903db 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -2,6 +2,20 @@ ## Unreleased +### Enhancements + +* Reduce the cost of every operation that goes through the context: results + that signal nothing no longer copy the context or write it back to the + process dictionary, the coefficient's digit count is computed once per + operation instead of up to four times, and digit counting no longer walks a + ladder of bignum comparisons before falling back to its estimate. Division + is ~1.9x faster, `add`/`sub`/`mult`/`round`/`normalize` ~1.3x, comparison of + same-scale values ~1.4x, parsing ~1.2x, all with 20-30% less allocation. + +* Make `Decimal.to_float/1` scale the operand with a computed shift instead of + one bit at a time: ~1.6x faster for typical values and ~11x for exponents + near the ends of the double range, with 96% fewer collections. + ### Bug fixes * Fix `Decimal.div/2` rounding the wrong way on inexact results. The long diff --git a/bench.exs b/bench.exs index addf6ac..1584b13 100644 --- a/bench.exs +++ b/bench.exs @@ -1,11 +1,22 @@ -Mix.install([ - {:decimal, path: ".", override: true}, - {:benchee, "~> 1.0"}, - {:benchee_html, "~> 1.0"} -]) - -{head, 0} = System.cmd("git", ["symbolic-ref", "--short", "HEAD"]) -{hash, 0} = System.cmd("git", ["rev-parse", "--short", "HEAD"]) +Code.require_file("bench_helper.exs", __DIR__) + +# Measure production-compiled code: +# +# MIX_ENV=prod elixir bench.exs +# +# To compare two checkouts, measure one and load its results into the run of +# the other, which makes Benchee print both per job: +# +# DECIMAL_PATH=../decimal-main MIX_ENV=prod elixir bench.exs +# LOAD=benchmarks/HEAD-0c0f72c.benchee MIX_ENV=prod elixir bench.exs +# +decimal_path = BenchHelper.install!([{:benchee, "~> 1.0"}, {:benchee_html, "~> 1.0"}]) + +# The tag names the code under test, which is not necessarily this checkout. +# `--abbrev-ref` reports `HEAD` rather than failing when that checkout has no +# branch, as a worktree at a bare commit does. +{head, 0} = System.cmd("git", ["-C", decimal_path, "rev-parse", "--abbrev-ref", "HEAD"]) +{hash, 0} = System.cmd("git", ["-C", decimal_path, "rev-parse", "--short", "HEAD"]) tag = "#{String.trim(head)}-#{String.trim(hash)}" @@ -73,6 +84,20 @@ each = fn decimals, fun -> fn -> Enum.each(decimals, fun) end end +money_strings = BenchHelper.money_strings() +money = Enum.map(money_strings, &Decimal.new/1) +money_pairs = Enum.zip(money, Enum.reverse(money)) + +floats = for i <- 1..200, do: i * 1.37 + +# Conversions to float scale the operand into the significand of a double, so +# operands near the ends of the exponent range do the most work. Exponents stay +# inside the double range, which `to_float/1` requires. +float_range_decimals = + for coef <- [1, 15, 1_234_567_890_123_456], exp <- [-290, -30, -1, 0, 1, 30, 290] do + struct(Decimal, %{sign: 1, coef: coef, exp: exp}) + end + jobs = %{ "compare" => each_pair.(decimal_pairs, &Decimal.compare/2), "compare same scale" => each_pair.(same_scale_pairs, &Decimal.compare/2), @@ -87,12 +112,31 @@ jobs = %{ "normalize" => each.(decimals, &Decimal.normalize/1), "sqrt" => each.(positive_decimals, &Decimal.sqrt/1), "to_string scientific" => each.(decimals, &Decimal.to_string(&1, :scientific)), - "to_string normal" => each.(decimals, &Decimal.to_string(&1, :normal)) + "to_string normal" => each.(decimals, &Decimal.to_string(&1, :normal)), + "new from string" => each.(money_strings, &Decimal.new/1), + "from_float" => each.(floats, &Decimal.from_float/1), + "to_float" => each.(float_range_decimals, &Decimal.to_float/1), + "inspect" => each.(money, &inspect/1), + "money add" => each_pair.(money_pairs, &Decimal.add/2), + "money mult" => each_pair.(money_pairs, &Decimal.mult/2), + "money div" => each_pair.(money_pairs, &Decimal.div/2), + "money compare" => each_pair.(money_pairs, &Decimal.compare/2), + "money round" => each.(money, &Decimal.round(&1, 2)) } +load = System.get_env("LOAD", "") |> String.split(",", trim: true) + +# Benchee resolves load paths with Path.wildcard, which turns a mistyped path +# into an empty list and a comparison run into a plain one, with no error. +for path <- load, Path.wildcard(path) == [] do + IO.puts(:stderr, "LOAD path #{path} matches no saved benchmark") + System.halt(1) +end + Benchee.run(jobs, time: 10, memory_time: 2, save: [path: "benchmarks/#{tag}.benchee", tag: tag], + load: load, formatters: [Benchee.Formatters.Console] ) diff --git a/bench_helper.exs b/bench_helper.exs new file mode 100644 index 0000000..6b17830 --- /dev/null +++ b/bench_helper.exs @@ -0,0 +1,27 @@ +defmodule BenchHelper do + @moduledoc false + + # Both bench scripts must measure the same checkout the same way, so the + # install preamble lives here. Returns the path of the code under test. + def install!(extra_deps \\ []) do + decimal_path = System.get_env("DECIMAL_PATH", ".") + + Mix.install([{:decimal, path: decimal_path, override: true} | extra_deps]) + + if Mix.env() != :prod do + IO.puts(:stderr, "refusing to benchmark a #{Mix.env()} build; rerun with MIX_ENV=prod") + System.halt(1) + end + + decimal_path + end + + # Values at a fixed scale with small coefficients: the shape of monetary + # amounts. The cents are zero-padded so every value has exponent -2 and the + # same-exponent paths of `add/2` and `compare/2` are what run. + def money_strings do + for i <- 1..200 do + "#{i * 37}.#{String.pad_leading("#{Integer.mod(i * 13, 100)}", 2, "0")}" + end + end +end diff --git a/bench_resources.exs b/bench_resources.exs new file mode 100644 index 0000000..3ec0d77 --- /dev/null +++ b/bench_resources.exs @@ -0,0 +1,374 @@ +Code.require_file("bench_helper.exs", __DIR__) + +# Benchmarks must measure production-compiled code: +# +# MIX_ENV=prod elixir bench_resources.exs +# +# To compare two builds, save one run and compare the next against it. The +# saved file holds the measurements themselves, so nothing has to parse the +# printed table back: +# +# DECIMAL_PATH=../decimal-main SAVE=/tmp/main.bench MIX_ENV=prod \ +# elixir bench_resources.exs +# COMPARE=/tmp/main.bench MIX_ENV=prod elixir bench_resources.exs +# +BenchHelper.install!() + +defmodule Resources do + @moduledoc """ + Resource accounting for a workload, measured in a dedicated process. + + Wall time alone is not a verdict: an implementation that is faster but burns + more CPU or allocates more is usually a bad trade. So every workload reports + + * `wall/op` - wall clock per operation + * `cpu/op` - VM CPU time per operation, summed over *all* schedulers + (a wall-time win paid for with extra CPU shows up here) + * `reds/op` - reductions per operation: deterministic work units, immune + to machine noise + * `alloc/op` - heap words allocated per operation (allocation churn) + * `copy/op` - heap words copied by the collector per operation. This is + the GC's actual cost, and it is the metric that punishes + allocating many small objects when they are *not* garbage + by the time the collector runs + * `gc/1k` - collections per 1k operations, minor + major + * `peak KB` - peak process footprint during the run (the memory spike) + + Every workload is measured twice: + + * `discard` - results are dropped immediately. Garbage dies young, minor + collections are cheap, and allocation looks almost free. This is what a + throwaway benchmark process measures. + * `retain` - a live set is held across the whole run and every result is + kept in a list. Now each collection has to trace and copy live data, and + produced terms get promoted to the old heap, so the true cost of + producing many objects is visible. + """ + + def measure(fun, reps, mode, live_set) do + parent = self() + + pid = + spawn(fn -> + # `live_set` is referenced after the loop so it stays live for the + # whole run and has to be traced by every collection. + receive do: (:go -> :ok) + acc = loop(fun, reps, mode, []) + send(parent, {:done, self()}) + # The hash keeps `acc` and the live set referenced past the loop. It + # runs only after the parent has taken every measurement, so its cost + # lands in none of them. + receive do: (:measured -> :ok) + send(parent, {:hash, self(), :erlang.phash2({acc, live_set})}) + receive do: (:stop -> :ok) + end) + + :erlang.trace(pid, true, [:garbage_collection]) + info0 = Process.info(pid, [:reductions, :garbage_collection_info]) + sched0 = scheduler_active() + t0 = System.monotonic_time(:nanosecond) + send(pid, :go) + + receive do + {:done, ^pid} -> :ok + end + + t1 = System.monotonic_time(:nanosecond) + sched1 = scheduler_active() + :erlang.trace(pid, false, [:garbage_collection]) + info = Process.info(pid, [:reductions, :garbage_collection_info]) + send(pid, :measured) + + receive do + {:hash, ^pid, _hash} -> :ok + end + + send(pid, :stop) + + events = drain([]) + gc = summarize(events, info, used_words(info0[:garbage_collection_info] || [])) + + %{ + wall_ns: (t1 - t0) / reps, + cpu_ns: (sched1 - sched0) / reps, + reds: (info[:reductions] - info0[:reductions]) / reps, + alloc_w: gc.allocated / reps, + copy_w: gc.copied / reps, + minor: gc.minor * 1000 / reps, + major: gc.major * 1000 / reps, + peak_kb: gc.peak_words * 8 / 1024 + } + end + + defp loop(_fun, 0, _mode, acc), do: acc + + defp loop(fun, n, :discard, acc) do + fun.() + loop(fun, n - 1, :discard, acc) + end + + defp loop(fun, n, :retain, acc) do + loop(fun, n - 1, :retain, [fun.() | acc]) + end + + # Sum of active time over all schedulers, normal and dirty: the VM's total + # CPU consumption, not just the wall time of one process. + defp scheduler_active do + :erlang.statistics(:scheduler_wall_time) + |> Enum.reduce(0, fn {_id, active, _total}, sum -> sum + active end) + end + + defp drain(acc) do + receive do + {:trace, _pid, event, info} -> drain([{event, info} | acc]) + after + 0 -> :lists.reverse(acc) + end + end + + # Allocation is the growth of the used heap between the end of one collection + # and the start of the next, on top of the heap the process was spawned with. + # Copying is what a collection moves: a minor collection copies the surviving + # young heap plus what it promotes and leaves the rest of the old heap where + # it is, while a major collection rebuilds the whole live set. + defp summarize(events, info, spawned_words) do + init = %{ + allocated: 0, + copied: 0, + last_end: spawned_words, + old_at_start: 0, + peak_words: 0, + minor: 0, + major: 0 + } + + acc = + Enum.reduce(events, init, fn {event, gc}, acc -> + used = used_words(gc) + + case event do + :gc_minor_start -> + %{ + acc + | allocated: acc.allocated + max(used - acc.last_end, 0), + old_at_start: Keyword.get(gc, :old_heap_size, 0), + peak_words: max(acc.peak_words, held(gc)), + minor: acc.minor + 1 + } + + :gc_major_start -> + %{ + acc + | allocated: acc.allocated + max(used - acc.last_end, 0), + old_at_start: Keyword.get(gc, :old_heap_size, 0), + peak_words: max(acc.peak_words, held(gc)), + major: acc.major + 1 + } + + :gc_minor_end -> + copied = + Keyword.get(gc, :heap_size, 0) + + (Keyword.get(gc, :old_heap_size, 0) - acc.old_at_start) + + %{acc | copied: acc.copied + max(copied, 0), last_end: used} + + :gc_major_end -> + %{acc | copied: acc.copied + used, last_end: used} + + _ -> + acc + end + end) + + final = info[:garbage_collection_info] || [] + + %{ + acc + | allocated: acc.allocated + max(used_words(final) - acc.last_end, 0), + peak_words: max(acc.peak_words, held(final)) + } + end + + defp used_words(gc) do + Keyword.get(gc, :heap_size, 0) + Keyword.get(gc, :old_heap_size, 0) + end + + defp held(gc) do + Keyword.get(gc, :heap_block_size, 0) + Keyword.get(gc, :old_heap_block_size, 0) + + Keyword.get(gc, :mbuf_size, 0) + Keyword.get(gc, :stack_size, 0) + end + + # Reported in this order, each as `{label, key, precision}`. + @metrics [ + {"wall/op", :wall_ns, 1}, + {"cpu/op", :cpu_ns, 1}, + {"reds/op", :reds, 1}, + {"alloc/op", :alloc_w, 1}, + {"copy/op", :copy_w, 1}, + {"minor/1k", :minor, 2}, + {"major/1k", :major, 2}, + {"peak KB", :peak_kb, 1} + ] + + @name_width 18 + @mode_width 8 + @value_width 10 + @delta_width 6 + + def header(comparing?) do + width = if comparing?, do: @value_width + @delta_width, else: @value_width + + @metrics + |> Enum.map(fn {label, _key, _precision} -> String.pad_leading(label, width) end) + |> emit_row("operation", "mode") + end + + # `baseline` is the matching measurement from the run being compared against, + # or nil when this workload was not in it. + def row(name, mode, measurement, baseline, comparing?) do + @metrics + |> Enum.map(&cell(&1, measurement, baseline, comparing?)) + |> emit_row(name, Atom.to_string(mode)) + end + + defp emit_row(cells, name, mode) do + name = String.pad_trailing(name, @name_width) + mode = String.pad_trailing(mode, @mode_width) + IO.puts(Enum.join([name, mode | cells], " ")) + end + + defp cell({_label, key, precision}, measurement, baseline, comparing?) do + value = Map.fetch!(measurement, key) + formatted = value |> :erlang.float_to_binary(decimals: precision) |> pad(@value_width) + + cond do + not comparing? -> formatted + baseline == nil -> formatted <> pad("new", @delta_width) + true -> formatted <> delta(value, Map.fetch!(baseline, key)) + end + end + + # A percentage against a zero baseline says nothing. + defp delta(_value, baseline) when baseline == 0, do: pad("-", @delta_width) + + defp delta(value, baseline) do + percent = (value - baseline) / baseline * 100 + sign = if percent < 0, do: "-", else: "+" + pad("#{sign}#{:erlang.float_to_binary(abs(percent), decimals: 0)}%", @delta_width) + end + + defp pad(string, width), do: String.pad_leading(string, width) + + # A workload runs `count` operations per repetition, so everything except the + # peak footprint - which belongs to the process, not to one operation - is + # divided down to a single operation. + @per_operation [:wall_ns, :cpu_ns, :reds, :alloc_w, :copy_w, :minor, :major] + + def per_operation(measurement, count) do + Enum.reduce(@per_operation, measurement, fn key, acc -> + Map.update!(acc, key, &(&1 / count)) + end) + end + + def load(nil), do: %{} + + def load(path) do + measurements = path |> File.read!() |> :erlang.binary_to_term() + keys = Enum.map(@metrics, fn {_label, key, _precision} -> key end) + + valid? = + is_map(measurements) and + Enum.all?(Map.values(measurements), fn measurement -> + is_map(measurement) and Enum.all?(keys, &Map.has_key?(measurement, &1)) + end) + + if valid? do + measurements + else + raise "#{path} does not hold measurements for #{inspect(keys)}; re-save it with SAVE=" + end + end + + def save(nil, _rows), do: :ok + + def save(path, rows) do + File.write!(path, :erlang.term_to_binary(rows)) + IO.puts("\nsaved #{map_size(rows)} measurements to #{path}") + end +end + +## Workloads ################################################################## + +:erlang.system_flag(:scheduler_wall_time, true) + +money_strings = BenchHelper.money_strings() +money = Enum.map(money_strings, &Decimal.new/1) +money_pairs = Enum.zip(money, Enum.reverse(money)) + +wide = + Enum.map(1..50, fn i -> + Decimal.new(1, 1_234_567_890_123_456_789_012_345_678_901_234 + i, -6 + rem(i, 5)) + end) + +wide_pairs = Enum.zip(wide, Enum.reverse(wide)) +tiny = Enum.map(1..50, &Decimal.new("#{&1}e-300")) +floats = Enum.map(1..200, &(&1 * 1.37)) + +# A live set the workload process holds for its entire lifetime, so the +# collector has real work to do on every pass, the way it does in a long-lived +# process that owns state. +live_set = Enum.map(1..20_000, &Decimal.new(1, &1 * 7919, -2)) + +pairs = fn list, fun -> fn -> Enum.map(list, fn {a, b} -> fun.(a, b) end) end end +over = fn list, fun -> fn -> Enum.map(list, fun) end end + +workloads = [ + {"add money", pairs.(money_pairs, &Decimal.add/2), 150, 200}, + {"sub money", pairs.(money_pairs, &Decimal.sub/2), 150, 200}, + {"mult money", pairs.(money_pairs, &Decimal.mult/2), 150, 200}, + {"div money", pairs.(money_pairs, &Decimal.div/2), 60, 200}, + {"compare money", pairs.(money_pairs, &Decimal.compare/2), 300, 200}, + {"round money 2", over.(money, &Decimal.round(&1, 2)), 150, 200}, + {"normalize money", over.(money, &Decimal.normalize/1), 300, 200}, + {"new money string", over.(money_strings, &Decimal.new/1), 150, 200}, + {"to_string money", over.(money, &Decimal.to_string/1), 150, 200}, + {"to_float money", over.(money, &Decimal.to_float/1), 60, 200}, + {"from_float", over.(floats, &Decimal.from_float/1), 60, 200}, + {"add wide", pairs.(wide_pairs, &Decimal.add/2), 300, 50}, + {"mult wide", pairs.(wide_pairs, &Decimal.mult/2), 300, 50}, + {"div wide", pairs.(wide_pairs, &Decimal.div/2), 150, 50}, + {"sqrt wide", over.(wide, &Decimal.sqrt/1), 60, 50}, + {"to_float tiny", over.(tiny, &Decimal.to_float/1), 15, 50} +] + +compare_to = System.get_env("COMPARE") +baseline = Resources.load(compare_to) +save_to = System.get_env("SAVE") +rounds = String.to_integer(System.get_env("ROUNDS", "3")) + +IO.puts("Decimal beam: #{:code.which(Decimal)}") +IO.puts("live set: #{length(live_set)} decimals held across every run") +if compare_to, do: IO.puts("comparing against #{compare_to}") +IO.puts("") + +Resources.header(compare_to != nil) + +rows = + for {name, fun, reps, per_rep} <- workloads, mode <- [:discard, :retain], into: %{} do + # one untimed pass so the JIT has compiled everything + Resources.measure(fun, max(div(reps, 10), 1), mode, live_set) + + # Wall and CPU time are noisy on a shared machine while reductions and + # allocation are deterministic, so take the fastest of several runs: the + # minimum is the estimate least polluted by unrelated system activity. + measurement = + Enum.map(1..rounds, fn _ -> Resources.measure(fun, reps, mode, live_set) end) + |> Enum.min_by(& &1.wall_ns) + |> Resources.per_operation(per_rep) + + Resources.row(name, mode, measurement, baseline[{name, mode}], compare_to != nil) + + {{name, mode}, measurement} + end + +Resources.save(save_to, rows) diff --git a/lib/decimal.ex b/lib/decimal.ex index ab55948..6070c39 100644 --- a/lib/decimal.ex +++ b/lib/decimal.ex @@ -168,6 +168,10 @@ defmodule Decimal do @default_max_exponent 6_144 @default_to_string_max_digits 6_178 + # A literal map is a compile-time constant, so the limits `parse/1`, `new/1` + # and `cast/1` always use cost no allocation. + @default_parse_limits %{max_digits: @default_max_digits, max_exponent: @default_max_exponent} + # Below 10^2000 the BIF `:erlang.integer_to_binary/1` is fast enough; for # larger integers `integer_to_decimal_iodata/3` recursively splits on a # power of 10 (down to chunks of `@decimal_conversion_leaf_digits` digits) @@ -376,15 +380,21 @@ defmodule Decimal do coef2 == 0 -> add_zero(num2, num1, ctx) + # The bounded path below guards against coefficient amplification from an + # exponent gap (CVE-2026-32686). Equal exponents need no alignment, so + # there is nothing to amplify: the bounded path would settle on + # `base_exp == exp1`, scale both coefficients by 10^0 and compute this + # very sum. Reaching that conclusion costs counting both coefficients' + # digits, which is all `add_bounded?/3` does. + exp1 == exp2 -> + add_coefs(sign1, coef1, sign2, coef2, exp1, ctx) + add_bounded?(num1, num2, ctx) -> add_bounded(num1, num2, ctx) true -> {coef1, coef2} = add_align(coef1, exp1, coef2, exp2) - coef = sign1 * coef1 + sign2 * coef2 - exp = Kernel.min(exp1, exp2) - sign = add_sign(sign1, sign2, coef, ctx) - context(%Decimal{sign: sign, coef: Kernel.abs(coef), exp: exp}, [], false, ctx) + add_coefs(sign1, coef1, sign2, coef2, Kernel.min(exp1, exp2), ctx) end end @@ -504,6 +514,17 @@ defmodule Decimal do def compare(%Decimal{sign: 1}, %Decimal{sign: -1}), do: :gt def compare(%Decimal{sign: -1}, %Decimal{sign: 1}), do: :lt + # With equal exponents the adjusted exponents differ exactly as the + # coefficient lengths do, so the coefficients decide it on their own and no + # digits need counting. + def compare(%Decimal{sign: sign, coef: coef1, exp: exp}, %Decimal{coef: coef2, exp: exp}) do + cond do + coef1 == coef2 -> :eq + coef1 < coef2 -> if sign == 1, do: :lt, else: :gt + true -> if sign == 1, do: :gt, else: :lt + end + end + def compare(%Decimal{} = num1, %Decimal{} = num2) do adjusted_exp1 = adjust_exp(num1) adjusted_exp2 = adjust_exp(num2) @@ -537,29 +558,47 @@ defmodule Decimal do compare(decimal(num1), decimal(num2)) end + @compile {:inline, adjust_exp: 1} + defp adjust_exp(%Decimal{coef: coef, exp: exp}) do coef_adjustment = coef_length(coef) exp + coef_adjustment - 1 end - defp coef_length(0), do: 1 - defp coef_length(coef) when coef < 10, do: 1 - defp coef_length(coef) when coef < 100, do: 2 - defp coef_length(coef) when coef < 1_000, do: 3 - defp coef_length(coef) when coef < 10_000, do: 4 - defp coef_length(coef) when coef < 100_000, do: 5 - defp coef_length(coef) when coef < 1_000_000, do: 6 - defp coef_length(coef) when coef < 10_000_000, do: 7 - defp coef_length(coef) when coef < 100_000_000, do: 8 - defp coef_length(coef) when coef < 1_000_000_000, do: 9 - defp coef_length(coef) when coef < 10_000_000_000, do: 10 - defp coef_length(coef) when coef < 100_000_000_000, do: 11 - defp coef_length(coef) when coef < 1_000_000_000_000, do: 12 - defp coef_length(coef) when coef < 10_000_000_000_000, do: 13 - defp coef_length(coef) when coef < 100_000_000_000_000, do: 14 - defp coef_length(coef) when coef < 1_000_000_000_000_000, do: 15 - defp coef_length(coef) when coef < 10_000_000_000_000_000, do: 16 - defp coef_length(coef) when coef < 100_000_000_000_000_000, do: 17 + # The ladder compares only against literals that fit a machine word + # (2^59 - 1 is the largest), so each rung is a register compare. Bignums + # leave it after two rungs, because comparing them costs an order of + # magnitude more than the bit-length estimate they fall through to. + defp coef_length(coef) when coef < 1_000_000_000 do + cond do + coef < 10 -> 1 + coef < 100 -> 2 + coef < 1_000 -> 3 + coef < 10_000 -> 4 + coef < 100_000 -> 5 + coef < 1_000_000 -> 6 + coef < 10_000_000 -> 7 + coef < 100_000_000 -> 8 + true -> 9 + end + end + + defp coef_length(coef) when coef <= 576_460_752_303_423_487 do + cond do + coef < 10_000_000_000 -> 10 + coef < 100_000_000_000 -> 11 + coef < 1_000_000_000_000 -> 12 + coef < 10_000_000_000_000 -> 13 + coef < 100_000_000_000_000 -> 14 + coef < 1_000_000_000_000_000 -> 15 + coef < 10_000_000_000_000_000 -> 16 + coef < 100_000_000_000_000_000 -> 17 + true -> 18 + end + end + + # One more rung for the bignums that are still 18 digits: the estimate costs + # about ten comparisons, so this one pays for itself. defp coef_length(coef) when coef < 1_000_000_000_000_000_000, do: 18 defp coef_length(coef), do: integer_decimal_digit_count(coef) @@ -796,7 +835,7 @@ defmodule Decimal do else ctx = Context.get() {coef1, coef2, adjust} = div_adjust(coef1, coef2) - {coef, adjust, rem, signals} = div_calc(coef1, coef2, adjust, ctx.precision) + {coef, adjust, rem, signals, digits} = div_calc(coef1, coef2, adjust, ctx.precision) # `rem` is the leftover of the division below the digits we kept. # It must be carried into rounding as the sticky bit: a nonzero `rem` @@ -804,7 +843,13 @@ defmodule Decimal do # so a guard digit of 5 is not an exact tie (`:half_even`/`:half_down`) # and a guard digit of 0 is still nonzero for `:ceiling`/`:floor`/`:up`. # Without it, ~5% of inexact divisions round the wrong way. - context(%Decimal{sign: sign, coef: coef, exp: exp1 - exp2 - adjust}, signals, rem != 0, ctx) + context( + %Decimal{sign: sign, coef: coef, exp: exp1 - exp2 - adjust}, + signals, + rem != 0, + ctx, + digits + ) end end @@ -1278,11 +1323,18 @@ defmodule Decimal do %{num | exp: 0} end - def normalize(%Decimal{sign: sign, coef: coef, exp: exp}) do - if coef == 0 do - %Decimal{sign: sign, coef: 0, exp: 0} - else - %{do_normalize(coef, exp) | sign: sign} |> context + def normalize(%Decimal{sign: sign, coef: coef, exp: exp} = num) do + cond do + coef == 0 -> + %Decimal{sign: sign, coef: 0, exp: 0} + + # A coefficient that does not end in zero is already minimal, so the + # input struct is the normalized value - no need to rebuild it. + Kernel.rem(coef, 10) != 0 -> + context(num) + + true -> + %{do_normalize(coef, exp) | sign: sign} |> context end end @@ -1556,7 +1608,7 @@ defmodule Decimal do """ @spec cast(term) :: {:ok, t} | :error - def cast(term), do: cast_with_limits(term, default_parse_limits()) + def cast(term), do: cast_with_limits(term, @default_parse_limits) @doc """ Creates a new decimal number from an integer, string, float, or existing decimal @@ -1621,7 +1673,7 @@ defmodule Decimal do """ @spec parse(binary()) :: {t(), binary()} | :error def parse(binary) when is_binary(binary) do - parse_with_limits(binary, default_parse_limits()) + parse_with_limits(binary, @default_parse_limits) end @doc """ @@ -1818,8 +1870,7 @@ defmodule Decimal do end defp integer_decimal_digit_count(int) do - bits = int |> :binary.encode_unsigned() |> bit_length() - digits = Kernel.div((bits - 1) * @log10_2_num, @log10_2_den) + 1 + digits = Kernel.div((bit_length(int) - 1) * @log10_2_num, @log10_2_den) + 1 integer_decimal_digit_count(int, digits) end @@ -1836,7 +1887,13 @@ defmodule Decimal do end end - defp bit_length(<>) do + # Index of the most significant set bit of a positive integer. The VM has no + # BIF for it, so go through the integer's shortest big-endian encoding. + defp bit_length(int) do + int |> :binary.encode_unsigned() |> binary_bit_length() + end + + defp binary_bit_length(<>) do byte_size(rest) * 8 + byte_bit_length(byte) end @@ -1903,7 +1960,7 @@ defmodule Decimal do defp to_string_digit_count(%Decimal{coef: coef}, _type) when coef in [:NaN, :inf], do: 0 defp to_string_digit_count(%Decimal{coef: coef, exp: exp}, :normal), - do: normal_digit_count(coef, exp) + do: normal_digit_count(coef_length(coef), exp) defp to_string_digit_count(%Decimal{coef: coef, exp: exp}, :xsd), do: xsd_digit_count(coef, exp) @@ -1919,14 +1976,12 @@ defmodule Decimal do cond do exp == 0 -> digits - exp < 0 and adjusted >= -6 -> normal_digit_count(coef, exp) + exp < 0 and adjusted >= -6 -> normal_digit_count(digits, exp) true -> digits + integer_digit_count(adjusted) end end - defp normal_digit_count(coef, exp) do - digits = coef_length(coef) - + defp normal_digit_count(digits, exp) do if exp >= 0 do digits + exp else @@ -1944,11 +1999,12 @@ defmodule Decimal do defp xsd_digit_count(coef, exp) do %Decimal{coef: coef, exp: exp} = do_normalize(coef, exp) + digits = coef_length(coef) if exp >= 0 do - coef_length(coef) + exp + 1 + digits + exp + 1 else - normal_digit_count(coef, exp) + normal_digit_count(digits, exp) end end @@ -2072,16 +2128,35 @@ defmodule Decimal do @spec scale(t) :: non_neg_integer() def scale(%Decimal{exp: exp}), do: Kernel.max(0, -exp) + # Shifts `num` up until it reaches `den`, so that dividing the two yields the + # 53 bits of a double's significand. The distance is the difference of the + # bit lengths, which lands within one bit, so a single comparison finishes + # it. defp scale_up(num, den, exp) when num >= den, do: {num, exp} - defp scale_up(num, den, exp), do: scale_up(num <<< 1, den, exp - 1) + defp scale_up(num, den, exp) do + shift = bit_length(den) - bit_length(num) + scaled = num <<< shift + + if scaled >= den do + {scaled, exp - shift} + else + {scaled <<< 1, exp - shift - 1} + end + end + + # The other direction: `den` grows until `num < 2 * den`. The caller passes it + # pre-multiplied by 2^52, which the result gives back. Landing one bit past + # the target means undoing that bit, hence the extra shift in the first + # branch. defp scale_down(num, den, exp) do - new_den = den <<< 1 + shift = Kernel.max(bit_length(num) - bit_length(den), 1) + scaled = den <<< shift - if num < new_den do - {den >>> 52, exp} + if scaled > num do + {scaled >>> 53, exp + shift - 1} else - scale_down(num, new_den, exp + 1) + {scaled >>> 52, exp + shift} end end @@ -2143,6 +2218,12 @@ defmodule Decimal do ## ARITHMETIC ## + defp add_coefs(sign1, coef1, sign2, coef2, exp, ctx) do + coef = sign1 * coef1 + sign2 * coef2 + sign = add_sign(sign1, sign2, coef, ctx) + context(%Decimal{sign: sign, coef: Kernel.abs(coef), exp: exp}, [], false, ctx) + end + defp add_align(coef1, exp1, coef2, exp2) when exp1 == exp2, do: {coef1, coef2} defp add_align(coef1, exp1, coef2, exp2) when exp1 > exp2, @@ -2293,6 +2374,11 @@ defmodule Decimal do # matching the exit conditions of the digit-at-a-time loop this replaces # (including the loop's inexact-shaped signals for exact quotients whose # adjust stays negative). + # + # The same invariant pins the quotient's digit count to `precision + 1`, so + # it is returned for rounding to use instead of counting the digits of a + # number just produced. Stripping trailing zeros changes that length, so the + # exact branch reports `nil` and its digits are counted as before. defp div_calc(coef1, coef2, adjust, precision) do scaled = coef1 * pow10(precision) coef = Kernel.div(scaled, coef2) @@ -2300,16 +2386,18 @@ defmodule Decimal do cond do rem != 0 -> - signals = if base10?(rem), do: [:rounded], else: [:inexact, :rounded] - {coef, adjust + precision, rem, signals} + # A nonzero remainder is the sticky bit, and precision/5 signals + # :inexact whenever it is set, so the result is always inexact here + # regardless of the remainder's shape. + {coef, adjust + precision, rem, [:inexact, :rounded], precision + 1} adjust + precision < 0 -> - {coef, adjust + precision, 0, [:inexact, :rounded]} + {coef, adjust + precision, 0, [:inexact, :rounded], precision + 1} true -> {stripped, zeros} = strip_trailing_zeros(coef, 0) strip = Kernel.min(zeros, adjust + precision) - {stripped * pow10(zeros - strip), adjust + precision - strip, 0, []} + {stripped * pow10(zeros - strip), adjust + precision - strip, 0, [], nil} end end @@ -2386,12 +2474,10 @@ defmodule Decimal do defp ratio(coef, exp) when exp >= 0, do: {coef * pow10(exp), 1} defp ratio(coef, exp) when exp < 0, do: {coef, pow10(-exp)} - pow10_max = - Enum.reduce(0..104, 1, fn int, acc -> - defp pow10(unquote(int)), do: unquote(acc) - defp base10?(unquote(acc)), do: true - acc * 10 - end) + Enum.reduce(0..104, 1, fn int, acc -> + defp pow10(unquote(int)), do: unquote(acc) + acc * 10 + end) # Binary powering (square-and-multiply): O(log n) multiplications instead # of a linear chain that repeatedly multiplies an ever-growing bignum. @@ -2406,16 +2492,6 @@ defmodule Decimal do end end - defp base10?(num) when num >= unquote(pow10_max) do - if Kernel.rem(num, unquote(pow10_max)) == 0 do - base10?(Kernel.div(num, unquote(pow10_max))) - else - false - end - end - - defp base10?(_num), do: false - ## ROUNDING ## defp do_round(sign, coef, exp, target_exp, rounding) do @@ -2438,16 +2514,24 @@ defmodule Decimal do end end - defp precision(%Decimal{coef: :NaN} = num, _precision, _rounding, _sticky?) do - {num, []} + # The result's digit count comes back too, since the caller needs it for the + # exponent limits and it is known here either way. + defp precision(%Decimal{coef: :NaN} = num, _digits, _precision, _rounding, _sticky?) do + {num, [], 0} end - defp precision(%Decimal{coef: :inf} = num, _precision, _rounding, _sticky?) do - {num, []} + defp precision(%Decimal{coef: :inf} = num, _digits, _precision, _rounding, _sticky?) do + {num, [], 0} end - defp precision(%Decimal{sign: sign, coef: coef, exp: exp} = num, precision, rounding, sticky?) do - num_digits = coef_length(coef) + defp precision( + %Decimal{sign: sign, coef: coef, exp: exp} = num, + digits, + precision, + rounding, + sticky? + ) do + num_digits = digits || coef_length(coef) cond do num_digits > precision -> @@ -2457,7 +2541,7 @@ defmodule Decimal do do_precision(sign, coef, num_digits, exp, num_digits, rounding, sticky?) true -> - {num, []} + {num, [], num_digits} end end @@ -2483,7 +2567,10 @@ defmodule Decimal do exp = exp + drop + carry dec = %Decimal{sign: sign, coef: signif, exp: exp} - {dec, signals} + # Dropping `drop` of `num_digits` digits leaves exactly `precision` of + # them, and the carry above restores that length when the increment + # lengthened it. + {dec, signals, precision} end # Splits `coef` into the leading digits that survive dropping the `drop` @@ -2493,6 +2580,12 @@ defmodule Decimal do # then a leading zero and all of `coef` lands in the rest. defp split_digits(coef, 0, sticky?), do: {coef, 0, sticky?} + # Dropping one digit, what every division does with its guard digit, needs + # no powers of ten at all. + defp split_digits(coef, 1, sticky?) do + {Kernel.div(coef, 10), Kernel.rem(coef, 10), sticky?} + end + defp split_digits(coef, drop, sticky?) do guard_pow = pow10(drop - 1) divisor = guard_pow * 10 @@ -2534,17 +2627,38 @@ defmodule Decimal do defp context(num, signals, sticky?), do: context(num, signals, sticky?, Context.get()) defp context(num, signals, sticky?, %Context{} = context) do - {result, prec_signals} = precision(num, context.precision, context.rounding, sticky?) - {result, exp_signals} = exponent_limits(result, context) - signals = signals |> put_uniq(prec_signals) |> put_uniq(exp_signals) - error(signals, nil, result, context) + context(num, signals, sticky?, context, nil) end - defp exponent_limits(%Decimal{coef: coef} = num, _context) when coef in [:NaN, :inf, 0], - do: {num, []} + # `digits` is the coefficient's digit count when the caller already knows it, + # `nil` when it has to be counted. + defp context(num, signals, sticky?, %Context{} = context, digits) do + {result, prec_signals, digits} = + precision(num, digits, context.precision, context.rounding, sticky?) + + {result, exp_signals} = exponent_limits(result, digits, context) + error(merge_signals(signals, prec_signals, exp_signals), nil, result, context) + end - defp exponent_limits(%Decimal{} = num, %Context{} = context) do - adjusted_exp = adjust_exp(num) + # Flags are recorded in the order the signals are merged, so the order is + # kept exactly. The first clauses only skip the membership scanning for the + # two shapes that cover virtually every operation: the caller's signals + # already cover the rounding ones, as in an inexact division, and rounding + # being the only thing that signalled, as in any rounded result. + defp merge_signals(signals, [], []), do: signals + defp merge_signals(signals, prec_signals, []) when prec_signals == signals, do: signals + defp merge_signals([], prec_signals, []), do: :lists.reverse(prec_signals) + + defp merge_signals(signals, prec_signals, exp_signals) do + signals |> put_uniq(prec_signals) |> put_uniq(exp_signals) + end + + defp exponent_limits(%Decimal{coef: coef} = num, _digits, _context) + when coef in [:NaN, :inf, 0], + do: {num, []} + + defp exponent_limits(%Decimal{exp: exp} = num, digits, %Context{} = context) do + adjusted_exp = exp + digits - 1 cond do above_emax?(adjusted_exp, context.emax) -> @@ -2581,35 +2695,31 @@ defmodule Decimal do defp overflow_to_infinity?(:ceiling, sign), do: sign == 1 defp overflow_to_infinity?(_rounding, _sign), do: true - defp put_uniq(list, elems) when is_list(elems) do - Enum.reduce(elems, list, &put_uniq(&2, &1)) + defp put_uniq(list, []), do: list + + defp put_uniq(list, [elem | elems]) do + list |> put_uniq(elem) |> put_uniq(elems) end defp put_uniq(list, elem) do - if elem in list, do: list, else: [elem | list] + if :lists.member(elem, list), do: list, else: [elem | list] end ## PARSING ## + defp parse_limits!([]), do: @default_parse_limits + defp parse_limits!(opts) do - Enum.reduce( - opts, - %{max_digits: @default_max_digits, max_exponent: @default_max_exponent}, - fn - {:max_digits, value}, acc -> - %{acc | max_digits: limit!(:max_digits, value)} - - {:max_exponent, value}, acc -> - %{acc | max_exponent: limit!(:max_exponent, value)} - - {key, _value}, _acc -> - raise ArgumentError, "unknown option #{inspect(key)}" - end - ) - end + Enum.reduce(opts, @default_parse_limits, fn + {:max_digits, value}, acc -> + %{acc | max_digits: limit!(:max_digits, value)} - defp default_parse_limits do - %{max_digits: @default_max_digits, max_exponent: @default_max_exponent} + {:max_exponent, value}, acc -> + %{acc | max_exponent: limit!(:max_exponent, value)} + + {key, _value}, _acc -> + raise ArgumentError, "unknown option #{inspect(key)}" + end) end defp limit!(_key, :infinity), do: :infinity @@ -2621,36 +2731,34 @@ defmodule Decimal do "#{inspect(key)} must be a non-negative integer or :infinity, got: #{inspect(value)}" end - defp parse_digits_count(<>, acc, count, leading_zeros) - when count == leading_zeros do - parse_digits_count(rest, acc, count + 1, leading_zeros + 1) - end - - defp parse_digits_count(<>, acc, count, leading_zeros) - when digit in ?0..?9 do - parse_digits_count(rest, [digit | acc], count + 1, leading_zeros) - end + # One pass counts the digits, which is what the limits are checked against, + # and accumulates their value. Leading zeros carry no value, so only the + # significant digits count against the accumulator: up to `@accum_digits` of + # them the accumulator stays inside a machine word and the scan yields the + # coefficient itself, with no intermediate list or binary. Beyond that it + # would become a bignum and grow quadratically, so longer runs are counted + # only and converted in one step - after the limit check, so that an + # over-long input is still rejected without its coefficient ever being + # materialized. + @accum_digits 17 - defp parse_digits_count(rest, acc, count, leading_zeros) do - {acc, count, leading_zeros, rest} + defp parse_digits_count(<>, count, leading_zeros, acc) + when count == leading_zeros do + parse_digits_count(rest, count + 1, leading_zeros + 1, acc) end - defp digits_acc_to_integer([], _size), do: 0 - defp digits_acc_to_integer(acc, _size), do: :erlang.list_to_integer(:lists.reverse(acc)) - - defp parse_exp(<>) - when e in [?e, ?E] and sign in [?+, ?-] and digit in ?0..?9 do - {digits, rest} = parse_digits(rest) - {[sign, digit | digits], rest} + defp parse_digits_count(<>, count, leading_zeros, acc) + when digit in ?0..?9 and count - leading_zeros < @accum_digits do + parse_digits_count(rest, count + 1, leading_zeros, acc * 10 + (digit - ?0)) end - defp parse_exp(<>) when e in [?e, ?E] and digit in ?0..?9 do - {digits, rest} = parse_digits(rest) - {[digit | digits], rest} + defp parse_digits_count(<>, count, leading_zeros, acc) + when digit in ?0..?9 do + parse_digits_count(rest, count + 1, leading_zeros, acc) end - defp parse_exp(bin) do - {[], bin} + defp parse_digits_count(rest, count, leading_zeros, acc) do + {count, leading_zeros, acc, rest} end defp parse_unsign(<>, _limits) @@ -2681,15 +2789,12 @@ defmodule Decimal do end defp parse_unsign(bin, limits) do - {int_rev, int_size, leading_zeros, after_int} = parse_digits_count(bin, [], 0, 0) + {int_size, leading_zeros, acc, after_int} = parse_digits_count(bin, 0, 0, 0) - {coef_rev, total_size, leading_zeros, after_float} = + {total_size, leading_zeros, acc, after_float} = case after_int do - <> -> - parse_digits_count(after_dot, int_rev, int_size, leading_zeros) - - _ -> - {int_rev, int_size, leading_zeros, after_int} + <> -> parse_digits_count(after_dot, int_size, leading_zeros, acc) + _ -> {int_size, leading_zeros, acc, after_int} end cond do @@ -2700,14 +2805,12 @@ defmodule Decimal do :error true -> - {exp, rest} = parse_exp(after_float) - exp_chars = if exp == [], do: ~c"0", else: exp float_size = total_size - int_size - case bounded_exponent(exp_chars, float_size, limits.max_exponent) do - {:ok, exp_int} -> - coef = digits_acc_to_integer(coef_rev, total_size) - {%Decimal{coef: coef, exp: exp_int}, rest} + case parse_exp(after_float, float_size, limits.max_exponent) do + {:ok, exp, rest} -> + coef = parse_coef(bin, int_size, float_size, total_size - leading_zeros, acc) + {%Decimal{coef: coef, exp: exp}, rest} :error -> :error @@ -2715,76 +2818,101 @@ defmodule Decimal do end end - defp decimal_within_limits?(%Decimal{coef: coef, exp: exp}, limits) do - not exceeds_limit?(decimal_digit_count(coef), limits.max_digits) and - within_exponent_limit?(exp, limits.max_exponent) + # Short coefficients came out of the scan already. Longer ones convert each + # digit run separately and shift the integer part up, rather than copying the + # two runs into one binary to convert together. A fraction can only exist + # when a dot was consumed, so its digits sit at `int_size + 1` in the input. + defp parse_coef(_bin, _int_size, _float_size, significant, acc) + when significant <= @accum_digits, + do: acc + + defp parse_coef(bin, int_size, 0, _significant, _acc) do + :erlang.binary_to_integer(binary_part(bin, 0, int_size)) end - defp decimal_digit_count(coef) when coef in [:NaN, :inf], do: 0 - defp decimal_digit_count(coef), do: coef_length(coef) + defp parse_coef(bin, 0, float_size, _significant, _acc) do + :erlang.binary_to_integer(binary_part(bin, 1, float_size)) + end - defp exceeds_limit?(_value, :infinity), do: false - defp exceeds_limit?(value, limit), do: value > limit + defp parse_coef(bin, int_size, float_size, _significant, _acc) do + int = :erlang.binary_to_integer(binary_part(bin, 0, int_size)) + frac = :erlang.binary_to_integer(binary_part(bin, int_size + 1, float_size)) + int * pow10(float_size) + frac + end - defp within_exponent_limit?(_exp, :infinity), do: true - defp within_exponent_limit?(exp, max_exponent), do: Kernel.abs(exp) <= max_exponent + # `e` notation. The digits are checked against the limit before being turned + # into an integer, so `1e` is rejected without ever being + # materialized. With no exponent, or no digits after the marker - in which + # case the marker is not part of the number - only the fraction digits count. + defp parse_exp(<> = bin, float_size, max_exponent) when e in [?e, ?E] do + {negative?, digits} = + case rest do + <> -> {true, tail} + <> -> {false, tail} + _ -> {false, rest} + end - defp bounded_exponent(chars, float_digits, :infinity) do - {:ok, List.to_integer(chars) - float_digits} - end + case parse_digits_count(digits, 0, 0, 0) do + {0, _leading_zeros, _acc, _rest} -> + no_exp(bin, float_size, max_exponent) - defp bounded_exponent(chars, float_digits, max_exponent) do - with {:ok, exp} <- bounded_integer(chars, max_exponent + float_digits) do - exp = exp - float_digits - if within_exponent_limit?(exp, max_exponent), do: {:ok, exp}, else: :error + {size, leading_zeros, acc, rest} -> + exp_value(digits, size, leading_zeros, acc, negative?, rest, float_size, max_exponent) end end - defp bounded_integer([?- | digits], bound) do - with {:ok, int} <- bounded_non_neg_integer(digits, bound), do: {:ok, -int} - end - - defp bounded_integer([?+ | digits], bound), do: bounded_non_neg_integer(digits, bound) - defp bounded_integer(digits, bound), do: bounded_non_neg_integer(digits, bound) + defp parse_exp(bin, float_size, max_exponent), do: no_exp(bin, float_size, max_exponent) - defp bounded_non_neg_integer(digits, bound) do - digits = trim_leading_zeroes(digits) - bound_digits = integer_to_charlist(bound) - digits_length = length(digits) - bound_length = length(bound_digits) - - cond do - digits == [] -> - {:ok, 0} + defp no_exp(rest, float_size, max_exponent) do + exp = -float_size + if within_exponent_limit?(exp, max_exponent), do: {:ok, exp, rest}, else: :error + end - digits_length > bound_length -> - :error + defp exp_value(digits, size, leading_zeros, acc, negative?, rest, float_size, :infinity) do + value = exp_digits_value(digits, size, leading_zeros, acc) + {:ok, exponent(value, negative?, float_size), rest} + end - digits_length == bound_length and digits_gt?(digits, bound_digits) -> - :error + defp exp_value(digits, size, leading_zeros, acc, negative?, rest, float_size, max_exponent) do + # More digits than the largest allowed exponent has cannot be within the + # limit, and rejecting on the count alone keeps a hostile exponent from + # being turned into an integer at all. What survives is at most one digit + # longer than the limit, so the limit check below settles it. + if size - leading_zeros > coef_length(max_exponent + float_size) do + :error + else + value = exp_digits_value(digits, size, leading_zeros, acc) + exp = exponent(value, negative?, float_size) - true -> - {:ok, List.to_integer(digits)} + if within_exponent_limit?(exp, max_exponent), do: {:ok, exp, rest}, else: :error end end - defp trim_leading_zeroes([?0 | rest]), do: trim_leading_zeroes(rest) - defp trim_leading_zeroes(digits), do: digits - - defp digits_gt?([digit | rest1], [digit | rest2]), do: digits_gt?(rest1, rest2) - defp digits_gt?([digit1 | _], [digit2 | _]), do: digit1 > digit2 - defp digits_gt?([], []), do: false + defp exponent(value, true, float_size), do: -value - float_size + defp exponent(value, false, float_size), do: value - float_size - defp parse_digits(bin), do: parse_digits(bin, []) + defp exp_digits_value(_digits, size, leading_zeros, acc) + when size - leading_zeros <= @accum_digits, + do: acc - defp parse_digits(<>, acc) when digit in ?0..?9 do - parse_digits(rest, [digit | acc]) + defp exp_digits_value(digits, size, leading_zeros, _acc) do + :erlang.binary_to_integer(binary_part(digits, leading_zeros, size - leading_zeros)) end - defp parse_digits(rest, acc) do - {:lists.reverse(acc), rest} + defp decimal_within_limits?(%Decimal{coef: coef, exp: exp}, limits) do + not exceeds_limit?(decimal_digit_count(coef), limits.max_digits) and + within_exponent_limit?(exp, limits.max_exponent) end + defp decimal_digit_count(coef) when coef in [:NaN, :inf], do: 0 + defp decimal_digit_count(coef), do: coef_length(coef) + + defp exceeds_limit?(_value, :infinity), do: false + defp exceeds_limit?(value, limit), do: value > limit + + defp within_exponent_limit?(_exp, :infinity), do: true + defp within_exponent_limit?(exp, max_exponent), do: Kernel.abs(exp) <= max_exponent + # Util defp decimal(%Decimal{} = num), do: num @@ -2796,37 +2924,68 @@ defmodule Decimal do "implicit conversion of #{inspect(other)} to Decimal is not allowed. Use Decimal.from_float/1" end - defp handle_error(signals, reason, result, context) do + # Most operations signal nothing, and with no signals there are no flags to + # add and no trap to fire: the context is unchanged, so the copy and the + # process dictionary write are skipped entirely. + defp handle_error([], _reason, result, _context), do: {:ok, result} + + defp handle_error(signals, reason, result, context) when is_list(signals) do + do_handle_error(signals, reason, result, context) + end + + defp handle_error(signal, reason, result, context) do + do_handle_error([signal], reason, result, context) + end + + defp do_handle_error(signals, reason, result, context) do context = context || Context.get() - signals = List.wrap(signals) - flags = Enum.reduce(signals, context.flags, &put_uniq(&2, &1)) - Context.set(%{context | flags: flags}) - error_signal = Enum.find(signals, &(&1 in context.traps)) + flags = put_uniq(context.flags, signals) - if error_signal do - error = [signal: error_signal, reason: reason] - {:error, error} - else - {:ok, result} + # Flags are sticky, so an already recorded signal leaves the context + # untouched: in a loop that keeps signalling the same thing, only the first + # operation writes. + if flags !== context.flags do + Context.set(%{context | flags: flags}) end - end - defp fix_float_exp(digits) do - fix_float_exp(digits, []) + case find_trap(signals, context.traps) do + nil -> {:ok, result} + error_signal -> {:error, [signal: error_signal, reason: reason]} + end end - defp fix_float_exp([?e | rest], [?0 | [?. | result]]) do - fix_float_exp(rest, [?e | result]) - end + defp find_trap([], _traps), do: nil - defp fix_float_exp([digit | rest], result) do - fix_float_exp(rest, [digit | result]) + defp find_trap([signal | signals], traps) do + if :lists.member(signal, traps) do + signal + else + find_trap(signals, traps) + end end - defp fix_float_exp([], result), do: :lists.reverse(result) + # `:io_lib_format.fwrite_g/1` renders exponent notation with a redundant + # fraction: `1.0e5`. Dropping it keeps `from_float/1` from reading that as a + # coefficient of 10 with the exponent one lower. + defp fix_float_exp([?., ?0, ?e | rest]), do: [?e | fix_float_exp(rest)] + defp fix_float_exp([char | rest]), do: [char | fix_float_exp(rest)] + defp fix_float_exp([]), do: [] + + # An adjusted exponent strictly inside ±308 puts the value between + # 10^adjusted and 10^(adjusted+1), clear of both DBL_MAX and DBL_MIN. Only + # the extremes need the exact comparisons below. + defp check_dbl_min_max(%Decimal{coef: 0} = num), do: num - defp check_dbl_min_max(%Decimal{sign: 1} = num) do + defp check_dbl_min_max(%Decimal{} = num) do + if Kernel.abs(adjust_exp(num)) < 308 do + num + else + check_dbl_range(num) + end + end + + defp check_dbl_range(%Decimal{sign: 1} = num) do cond do Decimal.gt?(num, dbl_max(1)) -> raise Error, reason: "number bigger than DBL_MAX: #{inspect(num)}" @@ -2839,7 +2998,7 @@ defmodule Decimal do end end - defp check_dbl_min_max(num) do + defp check_dbl_range(num) do cond do Decimal.lt?(num, dbl_max(-1)) -> raise Error, reason: "negative number smaller than DBL_MAX: #{inspect(num)}" @@ -2855,12 +3014,6 @@ defmodule Decimal do defp dbl_min(sign), do: %Decimal{sign: sign, coef: 22_250_738_585_072_014, exp: -324} defp zero(sign), do: %Decimal{sign: sign, coef: 0, exp: 0} defp dbl_max(sign), do: %Decimal{sign: sign, coef: 17_976_931_348_623_158, exp: 292} - - if Version.compare(System.version(), "1.3.0") == :lt do - defp integer_to_charlist(string), do: Integer.to_char_list(string) - else - defp integer_to_charlist(string), do: Integer.to_charlist(string) - end end defimpl Inspect, for: Decimal do diff --git a/test/decimal/property_test.exs b/test/decimal/property_test.exs index d64a9bd..2ed4d71 100644 --- a/test/decimal/property_test.exs +++ b/test/decimal/property_test.exs @@ -293,6 +293,81 @@ defmodule Decimal.PropertyTest do end end + describe "against independent oracles" do + property "compare/2 agrees with comparing the operands as scaled integers" do + check all(a <- decimal(), b <- decimal(), max_runs: 200) do + assert Decimal.compare(a, b) == scaled_integer_compare(a, b) + end + end + + property "compare/2 agrees with scaled integers at equal exponents" do + # Equal exponents skip the adjusted-exponent comparison entirely, so + # pin that path against the oracle as well. + check all(a <- decimal(), b <- decimal(), max_runs: 200) do + b = %{b | exp: a.exp} + assert Decimal.compare(a, b) == scaled_integer_compare(a, b) + end + end + + property "add/2 agrees with exact integer addition when the sum fits the precision" do + # Coefficients up to 16 digits at equal exponents sum without rounding, + # so the result must be the exact integer sum. + check all( + a <- decimal(coef_max: 9_999_999_999_999_999, exp_min: -20, exp_max: 20), + b <- decimal(coef_max: 9_999_999_999_999_999, exp_min: -20, exp_max: 20), + max_runs: 200 + ) do + b = %{b | exp: a.exp} + sum = a.sign * a.coef + b.sign * b.coef + result = Decimal.add(a, b) + + assert result.sign * result.coef == sum + assert result.exp == a.exp + end + end + + property "to_float/1 agrees with the runtime's own decimal to float conversion" do + # `String.to_float/1` is correctly rounded, and the exponent range here + # stays inside the double range that `to_float/1` accepts. + check all( + a <- + positive_decimal( + coef_max: 9_999_999_999_999_999_999_999_999_999_999_999, + exp_min: -60, + exp_max: 60 + ), + max_runs: 200 + ) do + expected = String.to_float("#{a.coef}.0e#{a.exp}") + + assert Decimal.to_float(a) == expected + assert Decimal.to_float(%{a | sign: -1}) == -expected + end + end + + property "parse/1 agrees with building the coefficient from the digits" do + check all( + int_digits <- StreamData.string(?0..?9, min_length: 1, max_length: 17), + frac_digits <- StreamData.string(?0..?9, max_length: 17), + exponent <- StreamData.integer(-40..40), + max_runs: 200 + ) do + string = "#{int_digits}.#{frac_digits}e#{exponent}" + {parsed, ""} = Decimal.parse(string, max_digits: :infinity) + + assert parsed.coef == String.to_integer(int_digits <> frac_digits) + assert parsed.exp == exponent - String.length(frac_digits) + end + end + end + + defp scaled_integer_compare(a, b) do + scale = min(a.exp, b.exp) + left = a.sign * a.coef * Integer.pow(10, a.exp - scale) + right = b.sign * b.coef * Integer.pow(10, b.exp - scale) + term_compare(left, right) + end + defp to_dec(float) when is_float(float), do: Decimal.from_float(float) defp to_dec(other), do: Decimal.new(other) diff --git a/test/decimal_test.exs b/test/decimal_test.exs index 0197510..fb7802f 100644 --- a/test/decimal_test.exs +++ b/test/decimal_test.exs @@ -1178,6 +1178,29 @@ defmodule DecimalTest do end) end + @tag timeout: @bounded_smoke_timeout + test "add/2 with a very large coefficient at the same exponent stays bounded" do + # Equal exponents skip the bounded path, so pin that this cannot amplify: + # with nothing to align, the coefficient is whatever the caller already + # built, and the result is rounded to the context precision. + # kept under `emax` so the result is a number rather than an overflow + digits = 5_000 + + big = %Decimal{ + sign: 1, + coef: :erlang.binary_to_integer(String.duplicate("9", digits)), + exp: 0 + } + + small = %Decimal{sign: 1, coef: 1, exp: 0} + + assert_runs_quickly("add large coef at equal exp", fn -> + # 5_000 nines plus one is 10^5_000, which has 5_001 digits and rounds to + # the 34 the default context allows + assert Decimal.add(big, small) == d(1, Integer.pow(10, 33), digits + 1 - 34) + end) + end + @tag timeout: @bounded_smoke_timeout test "add/2 with very large coefficient and small addend" do big = %Decimal{ @@ -1526,11 +1549,21 @@ defmodule DecimalTest do c19 = 9_999_999_999_999_999_999 pow10 = fn n -> String.to_integer("1" <> String.duplicate("0", n)) end - # 18 -> 19 digits crosses from the guard-chain clauses into the - # bit-length estimate; neither should be touched at default precision + c17 = 99_999_999_999_999_999 + + # 17 -> 18 digits crosses from the comparison ladder into the bit-length + # estimate, 18 -> 19 crosses the machine word; none should be touched at + # default precision + assert Decimal.apply_context(Decimal.new(1, c17, 0)) == d(1, c17, 0) + assert Decimal.apply_context(Decimal.new(1, c17 + 1, 0)) == d(1, c17 + 1, 0) assert Decimal.apply_context(Decimal.new(1, c18, 0)) == d(1, c18, 0) assert Decimal.apply_context(Decimal.new(1, c19, 0)) == d(1, c19, 0) + # the same boundary through scaling and comparison + assert Decimal.compare(Decimal.new(1, c17, 0), Decimal.new(1, c17 + 1, 0)) == :lt + assert Decimal.add(Decimal.new(1, c17, 0), Decimal.new(1, 1, 0)) == d(1, c17 + 1, 0) + assert Decimal.to_string(Decimal.new(1, c17 + 1, 0)) == "100000000000000000" + # 34 digits fits the default precision exactly; 35 digits rounds assert Decimal.apply_context(Decimal.new(1, pow10.(33), 0)) == d(1, pow10.(33), 0) assert Decimal.apply_context(Decimal.new(1, pow10.(34), 0)) == d(1, pow10.(33), 1) @@ -1560,6 +1593,192 @@ defmodule DecimalTest do "1" <> String.duplicate("0", 120) <> ".0" end + test "compare/2 with equal exponents" do + assert Decimal.compare(d(1, 12_345, -2), d(1, 12_346, -2)) == :lt + assert Decimal.compare(d(1, 12_346, -2), d(1, 12_345, -2)) == :gt + assert Decimal.compare(d(1, 12_345, -2), d(1, 12_345, -2)) == :eq + + # a negative sign inverts the coefficient order + assert Decimal.compare(d(-1, 12_345, -2), d(-1, 12_346, -2)) == :gt + assert Decimal.compare(d(-1, 12_346, -2), d(-1, 12_345, -2)) == :lt + + # coefficients of different lengths at the same exponent + assert Decimal.compare(d(1, 9, 0), d(1, 1_000, 0)) == :lt + assert Decimal.compare(d(-1, 9, 0), d(-1, 1_000, 0)) == :gt + + # 34-digit coefficients at the same scale + coef = 1_234_567_890_123_456_789_012_345_678_901_234 + assert Decimal.compare(d(1, coef, -5), d(1, coef + 1, -5)) == :lt + assert Decimal.compare(d(-1, coef, -5), d(-1, coef + 1, -5)) == :gt + assert Decimal.compare(d(1, coef, -5), d(1, coef, -5)) == :eq + + # zero and NaN are still decided before the coefficients are compared + assert Decimal.compare(d(1, 0, -2), d(1, 0, -2)) == :eq + assert Decimal.compare(d(1, 0, -2), d(1, 5, -2)) == :lt + assert Decimal.compare(d(-1, 0, -2), d(-1, 5, -2)) == :gt + assert Decimal.compare(d(1, 5, -2), d(-1, 5, -2)) == :gt + + assert_raise Error, fn -> Decimal.compare(d(1, :NaN, 0), d(1, 5, 0)) end + end + + test "add/2 with equal exponents and a coefficient length gap over the precision" do + # An exponent gap is what makes alignment expensive; with equal exponents + # there is nothing to align, so the exact sum is computed and rounded. + wide = String.to_integer("1" <> String.duplicate("0", 40)) + rounded = String.to_integer("1" <> String.duplicate("0", 33)) + + Context.with(%Context{traps: []}, fn -> + assert Decimal.add(d(1, wide, 0), d(1, 1, 0)) == d(1, rounded, 7) + assert :inexact in Context.get().flags + assert :rounded in Context.get().flags + end) + + # subtractive cancellation at equal exponents + assert Decimal.add(d(1, wide, 0), d(-1, wide, 0)) == d(1, 0, 0) + assert Decimal.add(d(1, 12_345, -2), d(-1, 12_346, -2)) == d(-1, 1, -2) + + # in-precision sums are exact + assert Decimal.add(d(1, 12_345, -2), d(1, 1, -2)) == d(1, 12_346, -2) + assert Decimal.sub(d(1, 12_345, -2), d(1, 12_345, -2)) == d(1, 0, -2) + end + + test "to_float/1 across the double range" do + # subnormals are below DBL_MIN and rejected, as documented + assert_raise Error, fn -> Decimal.to_float(~d"5e-324") end + + assert Decimal.to_float(~d"2.2250738585072014e-308") == 2.2250738585072014e-308 + assert Decimal.to_float(~d"1.7976931348623157e308") == 1.7976931348623157e308 + assert Decimal.to_float(~d"-1.7976931348623157e308") == -1.7976931348623157e308 + + # `String.to_float/1` is correctly rounded, so it is an independent oracle + # for the scaling `to_float/1` does with integer arithmetic + for coef <- [1, 3, 7, 9, 15, 123, 999_999_999_999_999, 1_234_567_890_123_456], + exp <- [-30, -17, -7, -3, -1, 0, 1, 3, 7, 17, 30] do + decimal = Decimal.new(1, coef, exp) + expected = String.to_float("#{coef}.0e#{exp}") + + assert Decimal.to_float(decimal) == expected + assert Decimal.to_float(Decimal.new(-1, coef, exp)) == -expected + end + + # every power of ten in range survives the round trip + for exp <- -300..300 do + float = :math.pow(10.0, exp) + assert float |> Decimal.from_float() |> Decimal.to_float() == float + end + end + + test "from_float/1 drops the redundant fraction from exponent notation" do + assert Decimal.from_float(1.0e5) == d(1, 1, 5) + assert Decimal.from_float(100_000.0) == d(1, 1, 5) + assert Decimal.from_float(1.0e-5) == d(1, 1, -5) + assert Decimal.from_float(-1.0e300) == d(-1, 1, 300) + assert Decimal.from_float(5.0e-324) == d(1, 5, -324) + + # values rendered without an exponent are untouched + assert Decimal.from_float(1.5) == d(1, 15, -1) + assert Decimal.from_float(0.1) == d(1, 1, -1) + assert Decimal.from_float(-3.14) == d(-1, 314, -2) + assert Decimal.from_float(0.0) == d(1, 0, -1) + end + + test "parse/2 exponent digits" do + assert Decimal.parse("1e0000") == {d(1, 1, 0), ""} + assert Decimal.parse("1e00000000000000000000005") == {d(1, 1, 5), ""} + assert Decimal.parse("1.5e-2x") == {d(1, 15, -3), "x"} + + # the marker is only part of the number when digits follow it + assert Decimal.parse("1e") == {d(1, 1, 0), "e"} + assert Decimal.parse("1e+") == {d(1, 1, 0), "e+"} + assert Decimal.parse("1e-") == {d(1, 1, 0), "e-"} + assert Decimal.parse("1E") == {d(1, 1, 0), "E"} + + # an exponent past the limit is rejected without being materialized + assert Decimal.parse("1e" <> String.duplicate("9", 10_000)) == :error + assert Decimal.parse("1e-" <> String.duplicate("9", 10_000)) == :error + assert Decimal.parse("1e6144") == {d(1, 1, 6144), ""} + assert Decimal.parse("1e6145") == :error + + # without a limit it is parsed in full + digits = String.duplicate("9", 40) + + assert Decimal.parse("1e" <> digits, max_exponent: :infinity) == + {d(1, 1, String.to_integer(digits)), ""} + end + + test "parse/1 coefficients around the digit accumulator boundary" do + for length <- 15..21 do + digits = String.duplicate("9", length) + assert Decimal.parse(digits) == {d(1, String.to_integer(digits), 0), ""} + + with_point = String.duplicate("9", length) <> "." <> String.duplicate("7", 5) + + assert Decimal.parse(with_point) == + {d(1, String.to_integer(String.duplicate("9", length) <> "77777"), -5), ""} + end + + # leading zeros are not significant digits and do not stop the accumulator + assert Decimal.parse("0000000000000000000000000000000000000001.5") == {d(1, 15, -1), ""} + assert Decimal.parse("0.0000000000000000000001") == {d(1, 1, -22), ""} + assert Decimal.parse(String.duplicate("0", 40)) == {d(1, 0, 0), ""} + + # past the boundary the digits convert from the input, leading zeros and all + nines = String.duplicate("9", 18) + assert Decimal.parse("000" <> nines) == {d(1, String.to_integer(nines), 0), ""} + + assert Decimal.parse("000" <> nines <> ".55") == + {d(1, String.to_integer(nines <> "55"), -2), ""} + + assert Decimal.parse("0." <> String.duplicate("0", 21) <> nines) == + {d(1, String.to_integer(nines), -39), ""} + + # exponent digits go through the same scan + assert Decimal.parse("1e" <> String.duplicate("0", 20) <> "5") == {d(1, 1, 5), ""} + end + + test "flags accumulate across operations" do + Context.with(%Context{traps: []}, fn -> + assert Context.get().flags == [] + + Decimal.add(~d"1", ~d"2") + assert Context.get().flags == [] + + # The order is asserted, not just the membership: operations merge their + # own signals with the ones rounding raises, and the merge must not + # reshuffle what is already recorded. + Decimal.div(~d"1", ~d"3") + assert Context.get().flags == [:rounded, :inexact] + + # re-signalling the same condition leaves the flags as they are + Decimal.div(~d"2", ~d"7") + Decimal.add(~d"1", ~d"2") + assert Context.get().flags == [:rounded, :inexact] + + # a new signal is still recorded on top + Decimal.mult(~d"1e6000", ~d"1e6000") + flags = Context.get().flags + assert :overflow in flags + assert :inexact in flags + assert :rounded in flags + end) + + # flags are per context, so a fresh one starts clean + Context.with(%Context{}, fn -> assert Context.get().flags == [] end) + end + + test "round/3 dropping exactly one digit" do + assert Decimal.round(~d"1.25", 1, :half_even) == d(1, 12, -1) + assert Decimal.round(~d"1.35", 1, :half_even) == d(1, 14, -1) + assert Decimal.round(~d"1.25", 1, :half_up) == d(1, 13, -1) + assert Decimal.round(~d"-1.25", 1, :half_up) == d(-1, 13, -1) + assert Decimal.round(~d"1.25", 1, :half_down) == d(1, 12, -1) + assert Decimal.round(~d"1.29", 1, :down) == d(1, 12, -1) + assert Decimal.round(~d"1.21", 1, :up) == d(1, 13, -1) + assert Decimal.round(~d"1.21", 1, :ceiling) == d(1, 13, -1) + assert Decimal.round(~d"1.29", 1, :floor) == d(1, 12, -1) + assert Decimal.round(~d"-1.21", 1, :floor) == d(-1, 13, -1) + end + defp assert_runs_quickly(name, fun) do {elapsed_us, _result} = :timer.tc(fun) diff --git a/verify_diff.exs b/verify_diff.exs new file mode 100644 index 0000000..0c458cf --- /dev/null +++ b/verify_diff.exs @@ -0,0 +1,414 @@ +decimal_path = System.get_env("DECIMAL_PATH", ".") +out = System.get_env("OUT") || raise "set OUT=" + +Mix.install([{:decimal, path: decimal_path, override: true}]) + +# Dumps the observable behaviour of the whole public surface - results, struct +# representations, raised errors and context flags - over a deterministic case +# matrix. Run it against two builds and diff the output: any behaviour change +# shows up as a diff line. +# +# OUT=/tmp/a DECIMAL_PATH=/tmp/decimal-baseline MIX_ENV=prod elixir verify_diff.exs +# OUT=/tmp/b MIX_ENV=prod elixir verify_diff.exs +# diff /tmp/a /tmp/b + +# `Decimal` is only available at runtime under Mix.install, so build structs +# through constructors instead of struct literals. +ctx = fn fields -> struct(Decimal.Context, fields) end +default_ctx = ctx.([]) + +coefs = [ + 0, + 1, + 2, + 5, + 9, + 10, + 11, + 50, + 99, + 100, + 999, + 1_000, + 12_345, + 999_999_999, + 1_000_000_000, + 99_999_999_999_999_999, + 100_000_000_000_000_000, + 999_999_999_999_999_999, + 1_000_000_000_000_000_000, + 1_000_000_000_000_000_001, + 1_234_567_890_123_456_789_012_345_678_901_234, + 9_999_999_999_999_999_999_999_999_999_999_999, + 10_000_000_000_000_000_000_000_000_000_000_000, + 12_345_678_901_234_567_890_123_456_789_012_345_678_901_234, + :NaN, + :inf +] + +exps = [-6200, -320, -300, -40, -35, -34, -20, -7, -6, -2, -1, 0, 1, 2, 6, 7, 20, 34, 300, 6200] + +decimals = + for sign <- [1, -1], coef <- coefs, exp <- exps do + Decimal.new(sign, coef, exp) + end + +# every 23rd decimal keeps the pair matrix at a few thousand cases +pair_sample = Enum.take_every(decimals, 23) +pairs = for a <- pair_sample, b <- pair_sample, do: {a, b} + +raw = fn d -> if is_struct(d, Decimal), do: {d.sign, d.coef, d.exp}, else: d end + +# Each case runs in a pristine context so the flags reported belong to it. +run = fn fun -> + Decimal.Context.with(default_ctx, fn -> + result = + try do + {:ok, fun.()} + rescue + e -> {:raised, e.__struct__, Exception.message(e)} + end + + {result, Decimal.Context.get().flags} + end) +end + +format = fn + {{:ok, {a, b}}, flags} when is_struct(a, Decimal) and is_struct(b, Decimal) -> + "#{inspect({raw.(a), raw.(b)})} #{inspect(flags)}" + + {{:ok, value}, flags} -> + "#{inspect(raw.(value))} #{inspect(flags)}" + + {{:raised, mod, msg}, flags} -> + "raise #{inspect(mod)} #{inspect(msg)} #{inspect(flags)}" +end + +file = File.open!(out, [:write, :raw, :delayed_write]) +emit = fn label, fun -> IO.binwrite(file, [label, " => ", format.(run.(fun)), "\n"]) end + +for {a, b} <- pairs do + key = "#{inspect(raw.(a))} #{inspect(raw.(b))}" + emit.("add #{key}", fn -> Decimal.add(a, b) end) + emit.("sub #{key}", fn -> Decimal.sub(a, b) end) + emit.("mult #{key}", fn -> Decimal.mult(a, b) end) + emit.("div #{key}", fn -> Decimal.div(a, b) end) + emit.("div_int #{key}", fn -> Decimal.div_int(a, b) end) + emit.("rem #{key}", fn -> Decimal.rem(a, b) end) + emit.("div_rem #{key}", fn -> Decimal.div_rem(a, b) end) + emit.("compare #{key}", fn -> Decimal.compare(a, b) end) + emit.("equal? #{key}", fn -> Decimal.equal?(a, b) end) + emit.("max #{key}", fn -> Decimal.max(a, b) end) + emit.("min #{key}", fn -> Decimal.min(a, b) end) +end + +for d <- decimals do + key = inspect(raw.(d)) + emit.("normalize #{key}", fn -> Decimal.normalize(d) end) + emit.("abs #{key}", fn -> Decimal.abs(d) end) + emit.("negate #{key}", fn -> Decimal.negate(d) end) + emit.("apply_context #{key}", fn -> Decimal.apply_context(d) end) + emit.("sqrt #{key}", fn -> Decimal.sqrt(d) end) + emit.("integer? #{key}", fn -> Decimal.integer?(d) end) + emit.("scale #{key}", fn -> Decimal.scale(d) end) + emit.("positive? #{key}", fn -> Decimal.positive?(d) end) + emit.("negative? #{key}", fn -> Decimal.negative?(d) end) + emit.("to_integer #{key}", fn -> Decimal.to_integer(d) end) + emit.("to_float #{key}", fn -> Decimal.to_float(d) end) + emit.("inspect #{key}", fn -> inspect(d) end) + + for type <- [:scientific, :normal, :xsd, :raw] do + emit.("to_string #{type} #{key}", fn -> Decimal.to_string(d, type) end) + end + + for places <- [-40, -3, 0, 1, 2, 7, 40], + mode <- [:down, :half_up, :half_even, :ceiling, :floor, :half_down, :up] do + emit.("round #{places} #{mode} #{key}", fn -> Decimal.round(d, places, mode) end) + end +end + +# rounding and precision behaviour under non-default contexts +for precision <- [1, 2, 5, 7, 34], + rounding <- [:down, :half_up, :half_even, :ceiling, :floor, :half_down, :up], + {a, b} <- Enum.take_every(pairs, 97) do + case_ctx = ctx.(precision: precision, rounding: rounding) + key = "#{precision} #{rounding} #{inspect(raw.(a))} #{inspect(raw.(b))}" + + emit.("ctx add #{key}", fn -> Decimal.Context.with(case_ctx, fn -> Decimal.add(a, b) end) end) + emit.("ctx mult #{key}", fn -> Decimal.Context.with(case_ctx, fn -> Decimal.mult(a, b) end) end) + emit.("ctx div #{key}", fn -> Decimal.Context.with(case_ctx, fn -> Decimal.div(a, b) end) end) + emit.("ctx sqrt #{key}", fn -> Decimal.Context.with(case_ctx, fn -> Decimal.sqrt(a) end) end) +end + +for emax <- [:infinity, 6144, 20, 2], emin <- [:infinity, -6143, -20, -2] do + case_ctx = ctx.(emax: emax, emin: emin, traps: []) + + for {a, b} <- Enum.take_every(pairs, 211) do + key = "#{inspect(emax)} #{inspect(emin)} #{inspect(raw.(a))} #{inspect(raw.(b))}" + emit.("lim add #{key}", fn -> Decimal.Context.with(case_ctx, fn -> Decimal.add(a, b) end) end) + + emit.("lim mult #{key}", fn -> + Decimal.Context.with(case_ctx, fn -> Decimal.mult(a, b) end) + end) + + emit.("lim div #{key}", fn -> Decimal.Context.with(case_ctx, fn -> Decimal.div(a, b) end) end) + end +end + +strings = [ + "0", + "-0", + "0.0", + "3.14", + "-3.14", + "+3.14", + ".5", + "5.", + "1e10", + "1E-10", + "1e+10", + "-1.1e3", + "0.0000000001", + "1234567890123456789012345678901234", + "12345678901234567890123456789012345", + "0.000000000000000000000000000000000000001", + "1e6144", + "1e6145", + "1e-6143", + "1e-6200", + "00000000000000000000000000000000000000001.5", + "inf", + "Infinity", + "-inf", + "nan", + "NaN", + "-NaN", + "bad", + "1.2.3", + "1e", + "1e+", + "", + "1_000" +] + +for s <- strings do + emit.("parse #{inspect(s)}", fn -> Decimal.parse(s) end) + emit.("new #{inspect(s)}", fn -> Decimal.new(s) end) + emit.("cast #{inspect(s)}", fn -> Decimal.cast(s) end) + + for opts <- [ + [], + [max_digits: 5], + [max_digits: :infinity], + [max_exponent: 10], + [max_exponent: :infinity] + ] do + emit.("parse2 #{inspect(s)} #{inspect(opts)}", fn -> Decimal.parse(s, opts) end) + emit.("cast2 #{inspect(s)} #{inspect(opts)}", fn -> Decimal.cast(s, opts) end) + end +end + +# Flags are sticky and accumulate across operations, so run sequences of +# operations inside one context and record the flags after every step. A +# single-operation-per-context check would never exercise re-signalling. +for {label, ops} <- [ + {"inexact chain", + [ + {:div, "1", "3"}, + {:div, "1", "3"}, + {:add, "1", "0.0000000000000000000000000000000000001"}, + {:mult, "1.5", "2"}, + {:div, "10", "2"}, + {:div, "2", "7"}, + {:sqrt, "2", nil}, + {:add, "1", "2"} + ]}, + {"overflow chain", + [ + {:mult, "1e6000", "1e6000"}, + {:add, "1", "2"}, + {:mult, "1e-6000", "1e-6000"}, + {:div, "1", "3"}, + {:mult, "1e6000", "1e6000"} + ]}, + {"clean chain", [{:add, "1", "2"}, {:mult, "3", "4"}, {:sub, "5", "6"}, {:div, "8", "2"}]} + ] do + emit.("flag sequence #{label}", fn -> + Decimal.Context.with(ctx.(traps: []), fn -> + Enum.map(ops, fn + {:sqrt, a, _} -> + result = Decimal.sqrt(a) + {Decimal.to_string(result, :raw), Decimal.Context.get().flags} + + {op, a, b} -> + result = apply(Decimal, op, [a, b]) + {Decimal.to_string(result, :raw), Decimal.Context.get().flags} + end) + end) + end) +end + +# The same sequences with traps enabled, so a trapped signal raises mid-chain +for signal <- [:inexact, :rounded, :overflow, :underflow] do + emit.("trap #{signal}", fn -> + Decimal.Context.with(ctx.(traps: [signal]), fn -> + Enum.map([{:div, "1", "3"}, {:mult, "1e6000", "1e6000"}, {:add, "1", "2"}], fn {op, a, b} -> + try do + {Decimal.to_string(apply(Decimal, op, [a, b]), :raw), Decimal.Context.get().flags} + rescue + e -> {:raised, Exception.message(e), Decimal.Context.get().flags} + end + end) + end) + end) +end + +# Deterministic fuzz over the parser: digits, dots, signs, exponents with +# leading zeros and huge magnitudes, garbage prefixes and suffixes. +:rand.seed(:exsss, {17, 42, 4711}) + +digits = fn n -> for _ <- 1..n, into: "", do: <> end + +fuzz_strings = + for _ <- 1..4000 do + int_len = Enum.random(0..40) + frac_len = Enum.random(0..40) + zeros = String.duplicate("0", Enum.random(0..8)) + + int = if int_len == 0, do: "", else: zeros <> digits.(int_len) + frac = if frac_len == 0, do: "", else: digits.(frac_len) + + # a trailing point with no fraction digits is valid input too + dot = + cond do + frac != "" -> "." + Enum.random(1..4) == 1 -> "." + true -> "" + end + + exp = + case Enum.random(1..8) do + 1 -> "" + 2 -> "e" <> digits.(Enum.random(1..3)) + 3 -> "E-" <> digits.(Enum.random(1..3)) + 4 -> "e+" <> digits.(Enum.random(1..5)) + 5 -> "e" <> String.duplicate("0", Enum.random(1..10)) <> digits.(Enum.random(1..4)) + 6 -> "e-" <> digits.(Enum.random(6..25)) + 7 -> "e" + 8 -> "E+" + end + + sign = Enum.random(["", "-", "+"]) + tail = Enum.random(["", "x", ".5", "e3", " ", "-"]) + + sign <> int <> dot <> frac <> exp <> tail + end + +fuzz_limits = [ + [], + [max_digits: 0], + [max_digits: 1], + [max_digits: 5], + [max_digits: 34], + [max_digits: :infinity], + [max_exponent: 0], + [max_exponent: 1], + [max_exponent: 10], + [max_exponent: 6144], + [max_exponent: :infinity], + [max_digits: :infinity, max_exponent: :infinity] +] + +for s <- fuzz_strings do + emit.("fuzz parse #{inspect(s)}", fn -> Decimal.parse(s) end) + + opts = Enum.random(fuzz_limits) + emit.("fuzz parse2 #{inspect(s)} #{inspect(opts)}", fn -> Decimal.parse(s, opts) end) + emit.("fuzz cast #{inspect(s)} #{inspect(opts)}", fn -> Decimal.cast(s, opts) end) +end + +for s <- Enum.take(fuzz_strings, 400), opts <- fuzz_limits do + emit.("fuzz grid #{inspect(s)} #{inspect(opts)}", fn -> Decimal.parse(s, opts) end) +end + +floats = [ + 0.0, + -0.0, + 1.0, + 1.5, + -1.5, + 0.1, + 3.14, + 100_000.0, + 1.0e-5, + 1.0e-300, + 1.0e300, + 2.2250738585072014e-308, + 1.7976931348623157e308, + 123_456_789.123456, + 1.0e16, + 1.0e17 +] + +for f <- floats do + emit.("from_float #{inspect(f)}", fn -> Decimal.from_float(f) end) + emit.("cast_float #{inspect(f)}", fn -> Decimal.cast(f) end) + emit.("roundtrip #{inspect(f)}", fn -> f |> Decimal.from_float() |> Decimal.to_float() end) +end + +# `from_float/1` reformats the rendered float, so sweep the shapes that +# rendering can take: every power of ten, integral values, subnormals, and a +# large random sample across the exponent range. +# products that overflow the double range raise, so keep only what survives +finite = fn thunk -> + try do + [thunk.()] + rescue + _ -> [] + end +end + +float_fuzz = + Enum.flat_map(-323..308, fn e -> + base = :math.pow(10.0, e) + + Enum.flat_map( + [ + fn -> base end, + fn -> base * 1.5 end, + fn -> base * 9.87654321 end, + fn -> -base end, + fn -> base * 2.0 end, + fn -> base * 1.0000000000000002 end + ], + finite + ) + end) ++ + Enum.map(1..40, &(&1 * 1.0)) ++ + Enum.map(1..40, &(1.0 / &1)) ++ + [5.0e-324, 1.0e-323, 2.2250738585072011e-308, 1.7976931348623157e308] ++ + Enum.flat_map(1..6000, fn i -> + finite.(fn -> :rand.uniform() * :math.pow(10.0, rem(i * 7, 600) - 300) end) + end) + +for f <- float_fuzz do + emit.("float fuzz from #{inspect(f)}", fn -> Decimal.from_float(f) end) + emit.("float fuzz raw #{inspect(f)}", fn -> Decimal.to_string(Decimal.from_float(f), :raw) end) + + emit.("float fuzz trip #{inspect(f)}", fn -> f |> Decimal.from_float() |> Decimal.to_float() end) +end + +for i <- [0, 1, -1, 42, -42, 10_000_000_000_000_000_000, -10_000_000_000_000_000_000] do + emit.("new_int #{i}", fn -> Decimal.new(i) end) + emit.("cast_int #{i}", fn -> Decimal.cast(i) end) + emit.("to_float_int #{i}", fn -> i |> Decimal.new() |> Decimal.to_float() end) +end + +# every to_float that a double can represent, across the whole exponent range +for exp <- -330..310, coef <- [1, 15, 1_234_567_890_123_456, 9_999_999_999_999_999] do + emit.("to_float_sweep #{coef} #{exp}", fn -> Decimal.to_float(Decimal.new(1, coef, exp)) end) + emit.("to_float_sweep- #{coef} #{exp}", fn -> Decimal.to_float(Decimal.new(-1, coef, exp)) end) +end + +File.close(file) +IO.puts("wrote #{out}")