diff --git a/Package.swift b/Package.swift index 911a1da8..cf8c7cd6 100644 --- a/Package.swift +++ b/Package.swift @@ -1,9 +1,9 @@ -// swift-tools-version:5.5 +// swift-tools-version:5.9 //===--- Package.swift ----------------------------------------*- swift -*-===// // // This source file is part of the Swift Numerics open source project // -// Copyright (c) 2019-2021 Apple Inc. and the Swift Numerics project authors +// Copyright (c) 2019-2025 Apple Inc. and the Swift Numerics project authors // Licensed under Apache License v2.0 with Runtime Library Exception // // See https://swift.org/LICENSE.txt for license information @@ -81,19 +81,6 @@ let package = Package( name: "RealTests", dependencies: ["_TestSupport"], exclude: ["CMakeLists.txt"] - ), - - // MARK: - Test executables - .executableTarget( - name: "ComplexLog", - dependencies: ["Numerics", "_TestSupport"], - path: "Tests/Executable/ComplexLog" - ), - - .executableTarget( - name: "ComplexLog1p", - dependencies: ["Numerics", "_TestSupport"], - path: "Tests/Executable/ComplexLog1p" ) ] ) diff --git a/Sources/IntegerUtilities/DivideWithRounding.swift b/Sources/IntegerUtilities/DivideWithRounding.swift index f1d8fa08..caf0b177 100644 --- a/Sources/IntegerUtilities/DivideWithRounding.swift +++ b/Sources/IntegerUtilities/DivideWithRounding.swift @@ -116,22 +116,6 @@ extension BinaryInteger { // If q is already odd, we have the correct result. if q._lowWord & 1 == 1 { return q } - case .stochastically: - let bmag = other.magnitude - let rmag = r.magnitude - var bhi: UInt64 - var rhi: UInt64 - if other.magnitude <= UInt64.max { - bhi = UInt64(bmag) - rhi = UInt64(rmag) - } else { - let shift = bmag._msb - 63 - bhi = UInt64(truncatingIfNeeded: bmag >> shift) - rhi = UInt64(truncatingIfNeeded: rmag >> shift) - } - let (sum, car) = rhi.addingReportingOverflow(.random(in: 0 ..< bhi)) - if sum < bhi && !car { return q } - case .requireExact: preconditionFailure("Division was not exact.") } @@ -299,22 +283,6 @@ extension SignedInteger { // If q is already odd, we have the correct result. if q._lowWord & 1 == 1 { return (q, r) } - case .stochastically: - let bmag = other.magnitude - let rmag = r.magnitude - var bhi: UInt64 - var rhi: UInt64 - if other.magnitude <= UInt64.max { - bhi = UInt64(bmag) - rhi = UInt64(rmag) - } else { - let shift = bmag._msb - 63 - bhi = UInt64(truncatingIfNeeded: bmag >> shift) - rhi = UInt64(truncatingIfNeeded: rmag >> shift) - } - let (sum, car) = rhi.addingReportingOverflow(.random(in: 0 ..< bhi)) - if sum < bhi && !car { return (q, r) } - case .requireExact: preconditionFailure("Division was not exact.") } @@ -347,7 +315,7 @@ extension SignedInteger { /// is not representable. /// /// - Returns: `(quotient, remainder)`, with `0 <= remainder < b.magnitude`. -func euclideanDivision(_ a: T, _ b: T) -> (quotient: T, remainder: T) +public func euclideanDivision(_ a: T, _ b: T) -> (quotient: T, remainder: T) where T: SignedInteger { a.divided(by: b, rounding: a >= 0 ? .towardZero : .awayFromZero) diff --git a/Sources/IntegerUtilities/RoundingRule.swift b/Sources/IntegerUtilities/RoundingRule.swift index 317b8b17..60f0f15d 100644 --- a/Sources/IntegerUtilities/RoundingRule.swift +++ b/Sources/IntegerUtilities/RoundingRule.swift @@ -2,21 +2,13 @@ // // This source file is part of the Swift Numerics open source project // -// Copyright (c) 2021-2024 Apple Inc. and the Swift Numerics project authors +// Copyright (c) 2021-2025 Apple Inc. and the Swift Numerics project authors // Licensed under Apache License v2.0 with Runtime Library Exception // // See https://swift.org/LICENSE.txt for license information // //===----------------------------------------------------------------------===// -// TODO: it's unfortunate that we can't specify a custom random source -// for the stochastic rounding rule, but I don't see a nice way to have -// that share the API with the other rounding rules, because we'd then -// have to take either the rule in-out or have an additional RNG/state -// parameter. The same problem applies to rounding with dithering and -// any other stateful rounding method. We should consider adding a -// stateful rounding API down the road to support those use cases. - /// A rule that defines how to select one of the two representable results /// closest to a given value. /// @@ -59,22 +51,22 @@ /// 2.0 | 2 | 2 | 2 | 2 | 2 | /// -------+----------+----------+----------+----------+----------+ /// -/// Specialized rounding rules +/// Specialized rounding rules /// -/// value | toOdd | stochastically | requireExact | -/// =======+==============+=======================+================+ -/// -1.5 | -1 | 50% -2, 50% -1 | trap | -/// -------+--------------+-----------------------+----------------+ -/// -0.5 | -1 | 50% -1, 50% 0 | trap | -/// -------+--------------+-----------------------+----------------+ -/// 0.5 | 1 | 50% 0, 50% 1 | trap | -/// -------+--------------+-----------------------+----------------+ -/// 0.7 | 1 | 30% 0, 70% 1 | trap | -/// -------+--------------+-----------------------+----------------+ -/// 1.2 | 1 | 80% 1, 20% 2 | trap | -/// -------+--------------+-----------------------+----------------+ -/// 2.0 | 2 | 2 | 2 | -/// -------+--------------+-----------------------+----------------+ +/// value | toOdd | requireExact | +/// =======+==============+================+ +/// -1.5 | -1 | trap | +/// -------+--------------+----------------+ +/// -0.5 | -1 | trap | +/// -------+--------------+----------------+ +/// 0.5 | 1 | trap | +/// -------+--------------+----------------+ +/// 0.7 | 1 | trap | +/// -------+--------------+----------------+ +/// 1.2 | 1 | trap | +/// -------+--------------+----------------+ +/// 2.0 | 2 | 2 | +/// -------+--------------+----------------+ /// ``` public enum RoundingRule { /// Produces the closest representable value that is less than or equal @@ -200,45 +192,6 @@ public enum RoundingRule { /// because 5/2 = 2.5 is equally close to 2 and 3, and 2 is even. case toNearestOrEven - /// Adds a uniform random value from [0, d) to the value being rounded, - /// where d is the distance between the two closest representable values, - /// then rounds the sum downwards. - /// - /// Unlike all the other rounding modes, this mode is _not deterministic_; - /// repeated calls to rounding operations with this mode will generally - /// produce different results. There is a tradeoff implicit in using this - /// mode: you can sacrifice _reproducible_ results to get _more accurate_ - /// results in aggregate. For a contrived but illustrative example, consider - /// the following: - /// ``` - /// let data = Array(repeating: 1, count: 100) - /// let result = data.reduce(0) { - /// $0 + $1.divided(by: 3, rounding: rule) - /// } - /// ``` - /// because 1/3 is always the same value between 0 and 1, any - /// deterministic rounding rule must produce either 0 or 100 for - /// this computation. But rounding `stochastically` will - /// produce a value close to 33. The _error_ of the computation - /// is smaller, but the result will now change between runs of the - /// program. - /// - /// For this simple case a better solution would be to add the - /// values first, and then divide. This gives a result that is both - /// reproducible _and_ accurate: - /// ``` - /// let result = data.reduce(0, +)/3 - /// ``` - /// but this isn't always possible in more sophisticated scenarios, - /// and in those cases this rounding rule may be useful. - /// - /// Examples: - /// - `(-4).divided(by: 3, rounding: .stochastically)` - /// will be –1 with probability 2/3 and –2 with probability 1/3. - /// - `5.shifted(rightBy: 1, rounding: .stochastically)` - /// will be 2 with probability 1/2 and 3 with probability 1/2. - case stochastically - /// If the value being rounded is representable, that value is returned. /// Otherwise, a precondition failure occurs. /// @@ -304,15 +257,6 @@ extension FloatingPoint { // which way that rounds, then select the other value. let even = (trunc + one/2).rounded(.toNearestOrEven) return trunc == even ? trunc + one : trunc - case .stochastically: - let trunc = rounded(.towardZero) - if trunc == self { return trunc } - // We have eliminated all large values at this point; add dither in - // ±[0,1) and then truncate. - let bits = Swift.min(-Self.ulpOfOne.exponent, 32) - let random = Self(UInt32.random(in: 0 ... (1 << bits &- 1))) - let dither = Self(sign: sign, exponent: -bits, significand: random) - return (self + dither).rounded(.towardZero) case .requireExact: let trunc = rounded(.towardZero) precondition(isInfinite || trunc == self, "\(self) is not an exact integer.") diff --git a/Sources/IntegerUtilities/ShiftWithRounding.swift b/Sources/IntegerUtilities/ShiftWithRounding.swift index 0493d049..b83f1e64 100644 --- a/Sources/IntegerUtilities/ShiftWithRounding.swift +++ b/Sources/IntegerUtilities/ShiftWithRounding.swift @@ -130,27 +130,6 @@ extension BinaryInteger { return floor + Self((round + lost) >> count) case .toOdd: return floor | (lost == 0 ? 0 : 1) - case .stochastically: - // In theory, u01 should be Self.random(in: 0 ..< onesBit), but the - // random(in:) method does not exist on BinaryInteger. This is - // (arguably) good, though, because there's actually no reason to - // generate large amounts of randomness just to implement stochastic - // rounding; 32b suffices for almost all purposes, and 64b is more - // than enough. - var g = SystemRandomNumberGenerator() - let u01 = g.next() - if count < 64 { - // count is small, so mask and lost are representable as both - // UInt64 and Self, regardless of what type Self actually is. - return floor + Self(((u01 & UInt64(mask)) + UInt64(lost)) >> count) - } else { - // count is large, so lost may not be representable as UInt64; pre- - // shift by count-64 to isolate the high 64b of the fraction, then - // add u01 and carry-out to round. - let highWord = UInt64(truncatingIfNeeded: lost >> (Int(count) - 64)) - let (_, carry) = highWord.addingReportingOverflow(u01) - return floor + (carry ? 1 : 0) - } case .requireExact: precondition(lost == 0, "shift was not exact.") return floor diff --git a/Tests/Executable/ComplexLog/main.swift b/Tests/Executable/ComplexLog/main.swift deleted file mode 100644 index bf39da7a..00000000 --- a/Tests/Executable/ComplexLog/main.swift +++ /dev/null @@ -1,108 +0,0 @@ -//===--- main.swift -------------------------------------------*- swift -*-===// -// -// This source file is part of the Swift Numerics open source project -// -// Copyright (c) 2020 Apple Inc. and the Swift Numerics project authors -// Licensed under Apache License v2.0 with Runtime Library Exception -// -// See https://swift.org/LICENSE.txt for license information -// See https://swift.org/CONTRIBUTORS.txt for the list of Swift project authors -// -//===----------------------------------------------------------------------===// - -import Numerics -import _TestSupport - -#if DEBUG -fatalError("Run this test in Release configuration") -#else - -var componentError = Double(Float.ulpOfOne/2) -var complexError = Double(Float.ulpOfOne/2) -var componentMaxInput = Complex.zero -var complexMaxInput = Complex.zero - -func test(_ z: Complex) { - let tst = Complex.log(z) - // TODO: we _should_ be able to say Complex(z), but that goes through - // the slow, generic path for BinaryFloatingPoint conversion; once we get - // that resolved in the standard library, we can replace the conversion on - // this and the next line. Currently it dominate testing time if we do that. - let ref = Complex.log(Complex(Double(z.real), Double(z.imaginary))) - if tst == Complex(Float(ref.real), Float(ref.imaginary)) { return } - let thisError = relativeError(tst, ref) - if thisError > complexError { - complexMaxInput = z - complexError = thisError - } - let thisComponentError = componentwiseError(tst, ref) - if thisComponentError > componentError { - componentMaxInput = z - componentError = thisComponentError - } -} - -func testWithSymmetries(_ x: Float, _ y: Float) { - test(Complex( x, y)) - test(Complex( y, x)) - test(Complex(-y, x)) - test(Complex(-x, y)) - test(Complex(-x,-y)) - test(Complex(-y,-x)) - test(Complex( y,-x)) - test(Complex( x,-y)) -} - -// The hardest to evaluate cases for log are those close to the unit circle, -// where log(z) is nearly zero. We want to have plausibly dense test coverage -// close to the circle. -let radii: [Float] = [0.9, - 0.95, - 0.99, - 0.999, - 0.9999, - 1 - 10 * .ulpOfOne, - 1, - 1 + 10 * .ulpOfOne, - 1.0001, - 1.001, - 1.01, - 1.05, - 1.1] - -for r in radii { - for x in Interval(from: 1/Float.sqrt(2), to: r) { - // Generate the two y values that put us closest to the circle of radius r - let base = Double.sqrt(Double(r).addingProduct(Double(-x), Double(x))) - let a, b: Float - if Double(Float(base)) < base { a = Float(base); b = a.nextUp } - else { b = Float(base); a = b.nextDown } - testWithSymmetries(x, a) - testWithSymmetries(x, b) - } -} - -// Away from the unit circle is "easy" but we still want to get some coverage. -// Generate 1 million uniform random points inside the circle, and then test -// both z and 1/z. -var g = SystemRandomNumberGenerator() -var count = 0 -while count < 1_000_000 { - let z = Complex(.random(in: -1 ... 1, using: &g), .random(in: -1 ... 1, using: &g)) - if z.length > 1 { continue } - count += 1 - test(z) - test(1/z) -} - -print("Worst complex norm error seen for log was \(complexError)") -print("For input \(complexMaxInput).") -print("Reference result: \(Complex.log(Complex(complexMaxInput)))") -print(" Observed result: \(Complex.log(complexMaxInput))") - -print("Worst componentwise error seen for log was \(componentError)") -print("For input \(componentMaxInput).") -print("Reference result: \(Complex.log(Complex(componentMaxInput)))") -print(" Observed result: \(Complex.log(componentMaxInput))") - -#endif diff --git a/Tests/Executable/ComplexLog1p/main.swift b/Tests/Executable/ComplexLog1p/main.swift deleted file mode 100644 index 51a43975..00000000 --- a/Tests/Executable/ComplexLog1p/main.swift +++ /dev/null @@ -1,108 +0,0 @@ -//===--- main.swift -------------------------------------------*- swift -*-===// -// -// This source file is part of the Swift Numerics open source project -// -// Copyright (c) 2020 Apple Inc. and the Swift Numerics project authors -// Licensed under Apache License v2.0 with Runtime Library Exception -// -// See https://swift.org/LICENSE.txt for license information -// See https://swift.org/CONTRIBUTORS.txt for the list of Swift project authors -// -//===----------------------------------------------------------------------===// - -import Numerics -import _TestSupport - -#if DEBUG -fatalError("Run this test in Release configuration") -#else - -var componentError = Double(Float.ulpOfOne/2) -var complexError = Double(Float.ulpOfOne/2) -var componentMaxInput = Complex.zero -var complexMaxInput = Complex.zero - -func test(_ z: Complex) { - let tst = Complex.log(onePlus: z) - // TODO: we _should_ be able to say Complex(z), but that goes through - // the slow, generic path for BinaryFloatingPoint conversion; once we get - // that resolved in the standard library, we can replace the conversion on - // this and the next line. Currently it dominate testing time if we do that. - let ref = Complex.log(onePlus: Complex(Double(z.real), Double(z.imaginary))) - if tst == Complex(Float(ref.real), Float(ref.imaginary)) { return } - let thisError = relativeError(tst, ref) - if thisError > complexError { - complexMaxInput = z - complexError = thisError - } - let thisComponentError = componentwiseError(tst, ref) - if thisComponentError > componentError { - componentMaxInput = z - componentError = thisComponentError - } -} - -func testWithSymmetry(_ x: Float, _ y: Float) { - test(Complex( x-1, y)) - test(Complex( y-1, x)) - test(Complex(-y-1, x)) - test(Complex(-x-1, y)) - test(Complex(-x-1,-y)) - test(Complex(-y-1,-x)) - test(Complex( y-1,-x)) - test(Complex( x-1,-y)) -} - -// The hardest to evaluate cases for log are those close to the unit circle -// centered at -1, where log(1+z) is nearly zero. We want to have plausibly -// dense test coverage close to the circle. -let radii: [Float] = [0.9, - 0.95, - 0.99, - 0.999, - 0.9999, - 1 - 10 * .ulpOfOne, - 1, - 1 + 10 * .ulpOfOne, - 1.0001, - 1.001, - 1.01, - 1.05, - 1.1] - -for r in radii { - for x in Interval(from: 1/Float.sqrt(2), to: r) { - // Generate the two y values that put us closest to the circle of radius r - let base = Double.sqrt(Double(r).addingProduct(Double(-x), Double(x))) - let a, b: Float - if Double(Float(base)) < base { a = Float(base); b = a.nextUp } - else { b = Float(base); a = b.nextDown } - testWithSymmetry(x, a) - testWithSymmetry(x, b) - } -} - -// Away from the unit circle is "easy" but we still want to get some coverage. -// Generate 1 million uniform random points inside the circle, and then test -// both z and 1/z. -var g = SystemRandomNumberGenerator() -var count = 0 -while count < 1_000_000 { - let z = Complex(.random(in: -1 ... 1, using: &g), .random(in: -1 ... 1, using: &g)) - if z.length > 1 { continue } - count += 1 - test(z - 1) - test(1/z - 1) -} - -print("Worst complex norm error seen for log(onePlus:) was \(complexError)") -print("For input \(complexMaxInput).") -print("Reference result: \(Complex.log(onePlus: Complex(complexMaxInput)))") -print(" Observed result: \(Complex.log(onePlus: complexMaxInput))") - -print("Worst componentwise error seen for log(onePlus:) was \(componentError)") -print("For input \(componentMaxInput).") -print("Reference result: \(Complex.log(onePlus: Complex(componentMaxInput)))") -print(" Observed result: \(Complex.log(onePlus: componentMaxInput))") - -#endif diff --git a/Tests/IntegerUtilitiesTests/DivideTests.swift b/Tests/IntegerUtilitiesTests/DivideTests.swift index 5150393a..70292fc0 100644 --- a/Tests/IntegerUtilitiesTests/DivideTests.swift +++ b/Tests/IntegerUtilitiesTests/DivideTests.swift @@ -374,17 +374,6 @@ final class IntegerUtilitiesDivideTests: XCTestCase { } } - func testDivideStochastic(_ values: [T]) { - for a in values { - for b in values where b != 0 { - // Skip any SignedInt.min / -1 cases, because those will trap. - if a == .min && b == -1 { continue } - let (q, r) = a.divided(by: b, rounding: .stochastically) - let _ = divisionRuleHolds(a, b, q, r) - } - } - } - func testDivideExact(_ values: [T]) { for a in values { for b in values where b != 0 { @@ -412,7 +401,6 @@ final class IntegerUtilitiesDivideTests: XCTestCase { testDivideToNearestOrAway(values) testDivideToNearestOrEven(values) testDivideToOdd(values) - testDivideStochastic(values) testDivideExact(values) } @@ -433,7 +421,6 @@ final class IntegerUtilitiesDivideTests: XCTestCase { testDivideToNearestOrAway(values) testDivideToNearestOrEven(values) testDivideToOdd(values) - testDivideStochastic(values) testDivideExact(values) } @@ -454,7 +441,6 @@ final class IntegerUtilitiesDivideTests: XCTestCase { testDivideToNearestOrAway(values) testDivideToNearestOrEven(values) testDivideToOdd(values) - testDivideStochastic(values) testDivideExact(values) } @@ -489,74 +475,10 @@ final class IntegerUtilitiesDivideTests: XCTestCase { } } - // stochastically rounding doesn't have a deterministic "expected" answer, - // but we know that the result must be either the floor or the ceiling. - // The above tests ensure that, but that's not a very strong guarantee; - // an implementation could just implement it as self / other and pass - // that test. - // - // Here we round the _same_ value many times, compute the average, and - // check that it is acceptably close to the exact expected value; simple - // use of any deterministic rounding rule will not achieve this. - func testStochasticDivide(_ a: T, _ b: T) -> Bool { - let trunc = a/b - let sum = (0..<1024).reduce(into: 0.0) { sum, _ in - let rounding = a.divided(by: b, rounding: .stochastically) - trunc - sum += Double(rounding) - } - let expected = 1024*Double(a%b)/Double(b) - let difference = abs(sum - expected) - // Waving my hands slightly instead of giving a precise explanation - // here, the expectation is that difference should be about - // 1/2 sqrt(1024). If we're more than a couple standard deviations - // off, we should flag that. This isn't _necessarily_ a problem, - // but if you see a repeated failure, that's almost surely a real bug. - // - // TODO: precise justification of thresholds - XCTAssertLessThanOrEqual(difference, 64, - "Accumulated error (\(difference)) was unexpectedly large in \(a).divided(by: \(b))" - ) - return difference > 16 - } - - func testDivideStochasticInt8() { - var values = [Int8](repeating: 0, count: 32) - for i in 0 ..< values.count { - while values[i] == 0 { - values[i] = .random(in: .min ... .max) - } - } - var fails = 0 - for a in values { - for b in values { - if a == .min && b == -1 { continue } - fails += testStochasticDivide(a, b) ? 1 : 0 - } - } - XCTAssertLessThanOrEqual(fails, 32*16) - } - - func testDivideStochasticUInt32() { - var values = [UInt32](repeating: 0, count: 32) - for i in 0 ..< values.count { - while values[i] == 0 { - values[i] = .random(in: .min ... .max) - } - } - var fails = 0 - for a in values { - for b in values { - fails += testStochasticDivide(a, b) ? 1 : 0 - } - } - XCTAssertLessThanOrEqual(fails, 32*16) - } - func testRemainderByMinusOne() { // These would trap if implemented as a - bq or similar, even though // the remainder is well-defined. XCTAssertEqual(0, Int.min.remainder(dividingBy: -1)) XCTAssertEqual(0, Int.min.remainder(dividingBy: -1, rounding: .up)) - XCTAssertEqual(0, Int.min.remainder(dividingBy: -1, rounding: .stochastically)) } } diff --git a/Tests/IntegerUtilitiesTests/ShiftTests.swift b/Tests/IntegerUtilitiesTests/ShiftTests.swift index 47f62f54..2291f5dd 100644 --- a/Tests/IntegerUtilitiesTests/ShiftTests.swift +++ b/Tests/IntegerUtilitiesTests/ShiftTests.swift @@ -77,18 +77,6 @@ final class IntegerUtilitiesShiftTests: XCTestCase { case 0b11: expected = ceiling default: preconditionFailure() } - case .stochastically: - // Just test that it's floor if exact, otherwise either floor - // or ceiling. - let observed = value.shifted(rightBy: count, rounding: rule) - if observed != floor && observed != ceiling { - print("Error found in \(T.self).shifted(rightBy: \(count), rounding: \(rule)).") - print(" Value: \(String(value, radix: 2))") - print("Expected: \(String(floor, radix: 2)) or \(String(ceiling, radix: 2))") - print("Observed: \(String(observed, radix: 2))") - XCTFail() - } - return case .requireExact: preconditionFailure() } @@ -130,7 +118,6 @@ final class IntegerUtilitiesShiftTests: XCTestCase { testRoundingShift(Int8.self, rounding: .toNearestOrAway) testRoundingShift(Int8.self, rounding: .toNearestOrEven) testRoundingShift(Int8.self, rounding: .toOdd) - testRoundingShift(Int8.self, rounding: .stochastically) testRoundingShift(UInt8.self, rounding: .down) testRoundingShift(UInt8.self, rounding: .up) @@ -142,7 +129,6 @@ final class IntegerUtilitiesShiftTests: XCTestCase { testRoundingShift(UInt8.self, rounding: .toNearestOrAway) testRoundingShift(UInt8.self, rounding: .toNearestOrEven) testRoundingShift(UInt8.self, rounding: .toOdd) - testRoundingShift(UInt8.self, rounding: .stochastically) testRoundingShift(Int.self, rounding: .down) testRoundingShift(Int.self, rounding: .up) @@ -154,7 +140,6 @@ final class IntegerUtilitiesShiftTests: XCTestCase { testRoundingShift(Int.self, rounding: .toNearestOrAway) testRoundingShift(Int.self, rounding: .toNearestOrEven) testRoundingShift(Int.self, rounding: .toOdd) - testRoundingShift(Int.self, rounding: .stochastically) testRoundingShift(UInt.self, rounding: .down) testRoundingShift(UInt.self, rounding: .up) @@ -166,62 +151,5 @@ final class IntegerUtilitiesShiftTests: XCTestCase { testRoundingShift(UInt.self, rounding: .toNearestOrAway) testRoundingShift(UInt.self, rounding: .toNearestOrEven) testRoundingShift(UInt.self, rounding: .toOdd) - testRoundingShift(UInt.self, rounding: .stochastically) - } - - // Stochastic rounding doesn't have a deterministic "expected" answer, - // but we know that the result must be either the floor or the ceiling. - // The above tests ensure that, but that's not a very strong guarantee; - // an implementation could just implement it as self >> count and pass - // that test. - // - // Here we round the _same_ value many times, compute the average, and - // check that it is acceptably close to the exact expected value; simple - // use of any deterministic rounding rule will not achieve this. - func testStochasticAverage(_ value: T) { - var fails = 0 - for count in 1 ... T.bitWidth { - let sum = (0..<256).reduce(into: DoubleWidth.zero) { sum, _ in - let rounded = value.shifted(rightBy: count, rounding: .stochastically) - sum += DoubleWidth(rounded) - } - let expected = DoubleWidth(value) << (8 - count) - let difference = sum >= expected ? sum - expected : expected - sum - // Waving my hands slightly instead of giving a precise explanation - // here, the expectation is that difference should be about - // 1/2 sqrt(256). If it's repeatedly bigger than that, we _may_ - // have a problem, but it's OK for this to fail occasionally. - // - // TODO: precise justification of thresholds - if difference > 8 { fails += 1 } - // On the other hand, if we're more than a couple standard deviations - // off, we should flag that. This still isn't _necessarily_ a problem, - // but if you see a repeated failure for a given shift count, that's - // almost surely a real bug. - XCTAssertLessThanOrEqual(difference, 32, - "Accumulated error (\(difference)) was unexpectedly large in \(value).shifted(rightBy: \(count))" - ) - } - // Threshold chosen so that this is expected to _usually_ pass, but - // it will fail sporadically even with a correct implementation. This is - // not a great fit for CI workflows, sorry. Basically ignore one-off - // failures, but a repeated failure here is an indication that a bug - // exists. - XCTAssertLessThanOrEqual(fails, T.bitWidth/2, - "Accumulated error was large more often than expected for \(value).shifted(rightBy:)" - ) - } - - func testStochasticShifts() { - testStochasticAverage(Int8.random(in: .min ... .max)) - testStochasticAverage(Int16.random(in: .min ... .max)) - testStochasticAverage(Int32.random(in: .min ... .max)) - testStochasticAverage(UInt8.random(in: .min ... .max)) - testStochasticAverage(UInt16.random(in: .min ... .max)) - testStochasticAverage(UInt32.random(in: .min ... .max)) - testStochasticAverage(Int64.random(in: .min ... .max)) - testStochasticAverage(UInt64.random(in: .min ... .max)) - testStochasticAverage(DoubleWidth.random(in: .min ... .max)) - testStochasticAverage(DoubleWidth.random(in: .min ... .max)) } } diff --git a/Tests/RealTests/IntegerExponentTests.swift b/Tests/RealTests/IntegerExponentTests.swift index f618fb6f..d41113e3 100644 --- a/Tests/RealTests/IntegerExponentTests.swift +++ b/Tests/RealTests/IntegerExponentTests.swift @@ -157,10 +157,10 @@ extension Double { #endif assertClose( 1.7976931348623151738531864721534215e308, Double.pow(-u, 3196577161300664268), allowedError: tol) assertClose(-1.7976931348623155730212483790972209e308, Double.pow(-u, 3196577161300664269), allowedError: tol) - assertClose( 1.7976931348623159721893102860411089e308, Double.pow(-u, 3196577161300664270), allowedError: tol) + assertClose( 1.7976931348623159721893102860411089e308, Double.pow(-u, 3196577161300664270), allowedError: tol) // warning expected on non-x86 assertClose( 1.7976931348623157075547244136070910e308, Double.pow(-d, -6393154322601327474), allowedError: tol) - assertClose(-1.7976931348623159071387553670790721e308, Double.pow(-d, -6393154322601327475), allowedError: tol) - assertClose( 1.7976931348623161067227863205510754e308, Double.pow(-d, -6393154322601327476), allowedError: tol) + assertClose(-1.7976931348623159071387553670790721e308, Double.pow(-d, -6393154322601327475), allowedError: tol) // warning expected on non-x86 + assertClose( 1.7976931348623161067227863205510754e308, Double.pow(-d, -6393154322601327476), allowedError: tol) // warning expected on non-x86 // Exponents close to underflow boundary. assertClose( 2.4703282292062334560337346683707907e-324, Double.pow(-u, -3355781687888880946)) assertClose(-2.4703282292062329075106789791206172e-324, Double.pow(-u, -3355781687888880947))