0

Does 1 / x / (1 / x) always give 1 for floating point calculation if I use highest optimization of compilers (with simd) for x < 0 <= 1.0 and 1/x doesn't overflow and becomes inf or -inf.

Or, I have to ensure that vector x and vector y in the following program have the same memory alignment?

Program:

#include <vector>
#include <iostream>
#include <cstdlib>
#include <time.h>
int main() {
    int a = 1;
    int length = 65537;
    srand(time(NULL));

    double c = 0;
    while (c == 0) c = (double)(rand()) / (RAND_MAX); 
    std::cout << c << "\n";
    std::vector<double> x(length, a); 
    std::vector<double> y(length, a); 
    for (auto& a: x) a /= c;
    for (auto& a: y) a /= c;
    for (int i = 0; i < length; ++i) x[i] /= y[i];
    for (int i = 0; i < length; ++i) {
        if (x[i] != 1) std::cout << "X";
    }   
    std::cout<<"\n";
}

Result (no X):

0.22659
hamster on wheels
  • 2,771
  • 17
  • 50
  • Floating point arithmetic is not exact in C++ (or Java, Python, etc.). – Tim Biegeleisen Aug 14 '17 at 01:49
  • Incorrect optimizations and the `x < 0 <= 1.0` typo notwithstanding, no. `1 / 5e-324 * 5e-324` is infinity. http://coliru.stacked-crooked.com/a/fb14f74f791b5ac6 – Ry- Aug 14 '17 at 01:51
  • If you visit the link, you’ll see that it’s not. – Ry- Aug 14 '17 at 01:51
  • ok. `1/(5e-324)` probably overflows, then `inf`/`inf` is `nan`. http://coliru.stacked-crooked.com/a/3dd48b27fa00c209 – hamster on wheels Aug 14 '17 at 01:53
  • @Ryan what if we are talking about normal numbers and no overflow, underflow, +-inf, or nan. no crazies. – hamster on wheels Aug 14 '17 at 01:55
  • Okay, you changed the question from `1 / x * x` to `(1 / x) / (1 / x)`. Yes, that will always be `1.0` when `1 / x` has a finite value, because it will be equal to itself. 5 × 10^−324 is a normal number, though, not a “crazy” one. – Ry- Aug 14 '17 at 01:56
  • thank god. my sanity is back. – hamster on wheels Aug 14 '17 at 01:56
  • Possible duplicate of [Is floating point math broken?](https://stackoverflow.com/questions/588004/is-floating-point-math-broken) – Henri Menke Aug 14 '17 at 02:08
  • @TimBiegeleisen: Floating point math is not _guaranteed_ to be exact. It often is. IEEE754 requires it for the 5 basic operations `+-*/ sqrt` – MSalters Aug 14 '17 at 08:00
  • @MSalters How do you define "exact"? – Patricia Shanahan Aug 14 '17 at 11:48
  • @PatriciaShanahan: as if the math is done with infinite precision, and then correctly rounded to the number of bits in the result. – MSalters Aug 14 '17 at 11:55
  • @MSalters I thought that might be what you meant by it. I've more usually seen it used to mean the real number arithmetic result of the operation, without rounding. – Patricia Shanahan Aug 14 '17 at 12:00
  • @PatriciaShanahan: That's obviously not going to fit in a fixed-size `float` or `double`. In fact, `sqrt(2)` isn't going to fit in _anything_. – MSalters Aug 14 '17 at 12:18
  • 1
    @MSalters What term would you use for the real number result before rounding? – Patricia Shanahan Aug 14 '17 at 12:23
  • an alternative is to use symbolic computation. like using a matrix to represent an algebraic number. not quite sure if that way can represent all real numbers though... – hamster on wheels Aug 14 '17 at 14:34

0 Answers0