1
  2
  3
  4
  5
  6
  7
  8
  9
 10
 11
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
// Copyright 2014-2018 Optimal Computing (NZ) Ltd.
// Licensed under the MIT license.  See LICENSE for details.

use std::cmp::Ordering;
use std::num::{FpCategory};
use super::{Ulps, ApproxEqUlps};

/// ApproxOrdUlps is for sorting floating point values where approximate equality
/// is considered equal.
pub trait ApproxOrdUlps: ApproxEqUlps {
    /// This method returns an ordering between `self` and `other` values
    /// if one exists, where Equal is returned if they are approximately
    /// equal within `ulps` floating point representations.  See module
    /// documentation for an understanding of `ulps`
    fn approx_cmp_ulps(&self, other: &Self, ulps: <<Self as ApproxEqUlps>::Flt as Ulps>::U)
                       -> Ordering;

    /// This method tests less than (for `self` < `other`), where values
    /// within `ulps` of each other are not less than.  See module
    /// documentation for an understanding of `ulps`.
    #[inline]
    fn approx_lt_ulps(&self, other: &Self, ulps: <<Self as ApproxEqUlps>::Flt as Ulps>::U)
                      -> bool
    {
        match self.approx_cmp_ulps(other, ulps) {
            Ordering::Less => true,
            _ => false,
        }
    }

    /// This method tests less than or equal to (for `self` <= `other`)
    /// where values within `ulps` are equal.  See module documentation
    /// for an understanding of `ulps`.
    #[inline]
    fn approx_le_ulps(&self, other: &Self, ulps: <<Self as ApproxEqUlps>::Flt as Ulps>::U)
                      -> bool
    {
        match self.approx_cmp_ulps(other, ulps) {
            Ordering::Less | Ordering::Equal => true,
            _ => false,
        }
    }

    /// This method tests greater than (for `self` > `other`)
    /// where values within `ulps` are not greater than.  See module
    /// documentation for an understanding of `ulps`
    #[inline]
    fn approx_gt_ulps(&self, other: &Self, ulps: <<Self as ApproxEqUlps>::Flt as Ulps>::U)
                      -> bool
    {
        match self.approx_cmp_ulps(other, ulps) {
            Ordering::Greater => true,
            _ => false,
        }
    }

    /// This method tests greater than or equal to (for `self` > `other`)
    /// where values within `ulps` are equal.  See module documentation
    /// for an understanding of `ulps`.
    #[inline]
    fn approx_ge_ulps(&self, other: &Self, ulps: <<Self as ApproxEqUlps>::Flt as Ulps>::U)
                      -> bool
    {
        match self.approx_cmp_ulps(other, ulps) {
            Ordering::Greater | Ordering::Equal => true,
            _ => false,
        }
    }
}

impl ApproxOrdUlps for f32 {
    fn approx_cmp_ulps(&self, other: &f32, ulps: <<Self as ApproxEqUlps>::Flt as Ulps>::U)
                       -> Ordering
    {
        let selfclass = self.classify();
        let otherclass = other.classify();

        // -0 and +0 are drastically far in ulps terms, so
        // we need a special case for that.
        if selfclass==FpCategory::Zero && otherclass==FpCategory::Zero {
            return Ordering::Equal;
        }

        // Handle differing signs as a special case, even if they are very
        // close, most people consider them unequal.
        let self_pos = self.is_sign_positive();
        let other_pos = other.is_sign_positive();

        let udiff: i32 = match (self_pos, other_pos) {
            (true, false) => return Ordering::Greater,
            (false, true) => return Ordering::Less,
            (true, true) => self.ulps(other),
            (false, false) => other.ulps(self), // invert ulps for negatives
        };

        match udiff {
            x if x < -ulps => Ordering::Less,
            x if x >= -ulps && x <= ulps => Ordering::Equal,
            x if x > ulps => Ordering::Greater,
            _ => unreachable!()
        }
    }
}

#[test]
fn f32_approx_cmp_test1() {
    let f: f32 = 0.1_f32;
    let mut sum: f32 = 0.0_f32;
    for _ in 0_isize..10_isize { sum += f; }
    let product: f32 = f * 10.0_f32;
    assert!(sum != product); // Should not be directly equal:
    println!("Ulps Difference: {}",sum.ulps(&product));
    assert!(sum.approx_cmp_ulps(&product,1) == Ordering::Equal); // But should be close
    assert!(sum.approx_cmp_ulps(&product,0) != Ordering::Equal);
    assert!(product.approx_cmp_ulps(&sum,0) != Ordering::Equal);
}
#[test]
fn f32_approx_cmp_test2() {
    let x: f32 = 1000000_f32;
    let y: f32 = 1000000.1_f32;
    assert!(x != y); // Should not be directly equal
    println!("Ulps Difference: {}",x.ulps(&y));
    assert!(x.approx_cmp_ulps(&y,2) == Ordering::Equal);
    assert!(x.approx_cmp_ulps(&y,1) == Ordering::Less);
    assert!(y.approx_cmp_ulps(&x,1) == Ordering::Greater);
}
#[test]
fn f32_approx_cmp_negatives() {
    let x: f32 = -1.0;
    let y: f32 = -2.0;
    assert!(x.approx_cmp_ulps(&y, 2) == Ordering::Greater);
}

// In all cases, approx_cmp() should be the same as cmp() if ulps=0
#[test]
fn f32_approx_cmp_vs_partial_cmp() {
    use std::mem;

    let testcases: [u32; 20] = [
        0,          // +0
        0x80000000, // -0
        0x00000101, // + denormal
        0x80000101, // - denormal
        0x00001101, // + denormal
        0x80001101, // - denormal
        0x01000101, // + normal
        0x81000101, // - normal
        0x01001101, // + normal
        0x81001101, // - normal
        0x7F800000, // +infinity
        0xFF800000, // -infinity
        0x7F801101, // SNaN
        0xFF801101, // SNaN
        0x7FC01101, // QNaN
        0xFFC01101, // QNaN
        0x7F801110, // SNaN
        0xFF801110, // SNaN
        0x7FC01110, // QNaN
        0xFFC01110, // QNaN

    ];

    let mut xf: f32;
    let mut yf: f32;
    for xbits in testcases.iter() {
        for ybits in testcases.iter() {
            xf = unsafe { mem::transmute::<u32,f32>(*xbits) };
            yf = unsafe { mem::transmute::<u32,f32>(*ybits) };
            if let Some(ordering) = xf.partial_cmp(&yf) {
                if ordering != xf.approx_cmp_ulps(&yf, 0) {
                    panic!("{} ({:x}) vs {} ({:x}): partial_cmp gives {:?} \
                            approx_cmp_ulps gives {:?}",
                           xf, xbits, yf, ybits, ordering , xf.approx_cmp_ulps(&yf, 0));
                }
            }
        }
        print!(".");
    }
}

impl ApproxOrdUlps for f64 {
    fn approx_cmp_ulps(&self, other: &f64, ulps: <<Self as ApproxEqUlps>::Flt as Ulps>::U)
                       -> Ordering
    {
        let selfclass = self.classify();
        let otherclass = other.classify();

        // -0 and +0 are drastically far in ulps terms, so
        // we need a special case for that.
        if selfclass==FpCategory::Zero && otherclass==FpCategory::Zero {
            return Ordering::Equal;
        }

        // Handle differing signs as a special case, even if they are very
        // close, most people consider them unequal.
        let self_pos = self.is_sign_positive();
        let other_pos = other.is_sign_positive();

        let udiff: i64 = match (self_pos, other_pos) {
            (true, false) => return Ordering::Greater,
            (false, true) => return Ordering::Less,
            (true, true) => self.ulps(other),
            (false, false) => other.ulps(self), // invert ulps for negatives
        };

        match udiff {
            x if x < -ulps => Ordering::Less,
            x if x >= -ulps && x <= ulps => Ordering::Equal,
            x if x > ulps => Ordering::Greater,
            _ => unreachable!()
        }
    }
}

#[test]
fn f64_approx_cmp_ulps_test1() {
    let f: f64 = 0.000000001_f64;
    let mut sum: f64 = 0.0_f64;
    for _ in 0_isize..10_isize { sum += f; }
    let product: f64 = f * 10.0_f64;
    assert!(sum != product); // Should not be directly equal:
    println!("Ulps Difference: {}",sum.ulps(&product));
    assert!(sum.approx_cmp_ulps(&product,1) == Ordering::Equal); // But should be close
    assert!(sum.approx_cmp_ulps(&product,0) != Ordering::Equal);
    assert!(product.approx_cmp_ulps(&sum,0) != Ordering::Equal);
}
#[test]
fn f64_approx_cmp_ulps_test2() {
    let x: f64 = 1000000_f64;
    let y: f64 = 1000000.0000000003_f64;
    assert!(x != y); // Should not be directly equal
    println!("Ulps Difference: {}",x.ulps(&y));
    assert!(x.approx_cmp_ulps(&y,3) == Ordering::Equal);
    assert!(x.approx_cmp_ulps(&y,2) == Ordering::Less);
    assert!(y.approx_cmp_ulps(&x,2) == Ordering::Greater);
}
#[test]
fn f64_approx_cmp_ulps_negatives() {
    let x: f64 = -1.0;
    let y: f64 = -2.0;
    assert!(x.approx_cmp_ulps(&y, 2) == Ordering::Greater);
}

// In all cases, approx_cmp_ulps() should be the same as cmp() if ulps=0
#[test]
fn f64_approx_cmp_ulps_vs_partial_cmp() {
    use std::mem;

    let testcases: [u64; 20] = [
        0,                   // +0
        0x80000000_00000000, // -0
        0x00000010_10000000, // + denormal
        0x80000010_10000000, // - denormal
        0x00000110_10000000, // + denormal
        0x80000110_10000000, // - denormal
        0x01000101_00100100, // + normal
        0x81000101_00100100, // - normal
        0x01001101_00100100, // + normal
        0x81001101_00100100, // - normal
        0x7FF00000_00000000, // +infinity
        0xFFF00000_00000000, // -infinity
        0x7FF01101_00100100, // SNaN
        0xFFF01101_00100100, // SNaN
        0x7FF81101_00100100, // QNaN
        0xFFF81101_00100100, // QNaN
        0x7FF01110_00100100, // SNaN
        0xFFF01110_00100100, // SNaN
        0x7FF81110_00100100, // QNaN
        0xFFF81110_00100100, // QNaN
    ];

    let mut xf: f64;
    let mut yf: f64;
    for xbits in testcases.iter() {
        for ybits in testcases.iter() {
            xf = unsafe { mem::transmute::<u64,f64>(*xbits) };
            yf = unsafe { mem::transmute::<u64,f64>(*ybits) };
            if let Some(ordering) = xf.partial_cmp(&yf) {
                if ordering != xf.approx_cmp_ulps(&yf, 0) {
                    panic!("{} ({:x}) vs {} ({:x}): partial_cmp gives {:?} \
                            approx_cmp_ulps gives {:?}",
                           xf, xbits, yf, ybits, ordering , xf.approx_cmp_ulps(&yf, 0));
                }
            }
        }
        print!(".");
    }
}