Fermat's factorisation in C++

Question

For fun, I've been implementing some maths stuff in C++, and I've been attempting to implement Fermats Factorisation Method, however, I don't know that I understand what it's supposed to return. This implementation I have, returns `105` for the example number 5959 given in the Wikipedia article. The pseudocode in Wikipedia looks like this: One tries various values of a, hoping that is a square. ``` FermatFactor(N): // N should be odd a → ceil(sqrt(N)) b2 → a*a - N while b2 isn't a square: a → a + 1 // equivalently: b2 → b2 + 2*a + 1 b2 → a*a - N // a → a + 1 endwhile return a - sqrt(b2) // or a + sqrt(b2) ``` My C++ implementation, look like this: ``` int FermatFactor(int oddNumber) { double a = ceil(sqrt(static_cast<double>(oddNumber))); double b2 = a*a - oddNumber; std::cout << "B2: " << b2 << "a: " << a << std::endl; double tmp = sqrt(b2); tmp = round(tmp,1); while (compare_doubles(tmp*tmp, b2)) //does this line look correct? { a = a + 1; b2 = a*a - oddNumber; std::cout << "B2: " <…

Accepted Answer

Do note that you should be doing all those calculations on integer types, not on floating point types. It would be much, much simpler (and possibly more correct).

Your `compare_doubles` function is wrong. `diff` should be a `double`.

And once you fix that, you'll need to fix your test line. `compare_doubles` will return true if its inputs are "nearly equal". You need to loop while they are "not nearly equal".

So:

```
bool compare_doubles(double a, double b)
{
    double diff = std::fabs(a - b);
    return diff < std::numeric_limits<double>::epsilon();
}
```

And:

```
while (!compare_doubles(tmp*tmp, b2))  // now it is
{
```

And you will get you the correct result (`101`) for this input.

You'll also need to call your `round` function with `0` as "places" as vhallac points out - you shouldn't be rounding to one digit after the decimal point.

The Wikipedia article you link has the equation that allows you to identify `b` from `N` and `a-b`.

Fermat's factorisation in C++

Solution

Problem