I am solving a numerical optimization problem with my own L-BFGS (implemented in c++). The problem has $\approx 10^6$ optimization parameters.

To find the objective function gradient, I am taking a discrete Fourier transform by using FFTW3 library in double precision, all my other calculations/variables are also in double precision.

Now, the problem is that with the same initial condition and all parameters, the algorithm sometimes converges with $\approx 2000$ iterations, sometimes with $\approx 3000$.

  • Can this non-deterministic behavior be caused by accumulating floating point error (or FFTW error)?
  • OR does it mean that my code has a bug and I am probably somewhere accessing unallocated memory (or some other bug which causes non-deterministic behavior)?

Non-reproducible behaviors in computing amidst different runs can involve several mechanisms, sometimes mixed. They can be especially sensitive when one iterates calculations on large sets of data, like inverse 3D tomography.

  • soft errors, caused by defects, interaction with high-energy particles in spacecrafts (let's rule them out for now)
  • number representations, different systems and OSes and not sufficiently specified rules for floats like rounding, overflow,
  • code and compiler options: sometimes x=a+b differs from x=b+a (non-commutativity),
  • hardware management, especially related to forms of parallelization on heterogeneous collections of processors (CPU, GPU): sometimes, x=a+b+c+d does not give the same results when computed (on two units) either as x=(a+b)+(c+d) or x=(a+c)+(b+d) (non-associativity).

The non-deterministic aspect might arise, mostly, from the fourth type, as from one run to the other, a different set, or a different grouping can be performed.

[EDIT] Specifically to FFTW3, the documentation says:

Question 3.8. FFTW gives different results between runs

If you use FFTW_MEASURE or FFTW_PATIENT mode, then the algorithm FFTW employs is not deterministic: it depends on runtime performance measurements. This will cause the results to vary slightly from run to run. However, the differences should be slight, on the order of the floating-point precision, and therefore should have no practical impact on most applications. If you use saved plans (wisdom) or FFTW_ESTIMATE mode, however, then the algorithm is deterministic and the results should be identical between runs.

It is for instance discussed in Not So Fast - The Hacker Factor Blog

Some of the shortcuts are very fast but are non-deterministic

This could be checked in you case.

A few quick sources:

  • 1
    $\begingroup$ Yes, I was using a non-deterministic FFTW mode. Switched to FFTW_ESTIMATE and every value is exactly the same, even if I print e.g. 30 digits of a double precision variable. Therefore, I am assuming that point 1 in your answer might not be correct. $\endgroup$ – eimrek Aug 24 '16 at 14:05
  • 1
    $\begingroup$ Points 1, 2 and 3 only lead to 'non-deterministic' behavior if you change hardware or compiler (flags). $\endgroup$ – Jannis Teunissen Aug 24 '16 at 14:27
  • $\begingroup$ @Jannis Teunissen Indeed, I have rewritten this part, hoping it is less incorrect now. My knowledge is too shallow to write bold statements $\endgroup$ – Laurent Duval Aug 24 '16 at 14:29

FFTW includes some adaptivity to select the fastest algorithms on your hardware. I'm no expert on the internals of FFTW, but it seems unlikely that this library causes the large differences you observe, in particular because it is widely used. Therefore I would suggest:

  • Try to compile with as much run-time error checking as possible. For example, with gfortran the flag -fcheck=all can be very helpful when debugging. Unfortunately g++ does not seem to offer such checks, although -Wall -Wextra could sometimes help. If available, try different compilers.

  • Run your program using a tool like valgrind, perhaps with a reduced problem size, to find out about unallocated, uninitialized or unfreed memory.

You don't mention whether your code runs in parallel, if that is the case, then the parallel load balancing can be a source of non-deterministic behavior.


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.