search

msg-byu / symlib

Spacegroup finder. Includes symmetry-related routines for cluster expansion and other codes that rely on symmetries of lattices and crystals.
MIT License
12 stars 15 forks source link

Numerical precision issues #4

Closed wsmorgan closed 9 years ago

wsmorgan commented 9 years ago

As I've been running unit tests there have been some precision issues that have come up in the code. For example in the tests for the cross_product subroutine in vector_matrix_utilities about 25% of the tests differ, from those independently generated, in the last one or two digits of precision that we are using in fortpy (so around the 11th and 12th decimal places). These differences are small and can be accommodated in the unit tests by having a tolerance of 2E-11 for the comparisons. My question though is if this is what we want or expect?

These kinds of errors could just be machine noise but they do propagate, the tolerance for the volume subroutines output, in the same module, has to be set to 2E-10 in order to pass 100% of the time because of it's dependence on the cross_product who's errors get compounded when multiplied by other values.

glwhart commented 9 years ago

The cross product is essentially a determinant. I suspect that when two vectors are nearly parallel that the precision suffers. Is there something in common with those that fail the tests? Can you poke around and see if there are algorithms for doing this better? Maybe we should just up the tolerance, but I don't want to if there is a more stable way of doing the computation.

glwhart commented 9 years ago

I'm surprised at how large the errors can be. It's much bigger than the normal finite precision errors I'm used to.

wsmorgan commented 9 years ago

I'll look around for a new approach and try to see which tests are failing without a tolerance and if there is a seemingly good reason for why.

Wiley Morgan http://hagrid.byu.edu/Wiley#Bio

On Fri, May 29, 2015 at 9:53 AM, Gus Hart notifications@github.com wrote:

I'm surprised at how large the errors can be. It's much bigger than the normal finite precision errors I'm used to.

— Reply to this email directly or view it on GitHub https://github.com/msg-byu/symlib/issues/4#issuecomment-106855024.

wsmorgan commented 9 years ago

I can't find a new way to find the cross product but I think I've nailed down the issue. When fortran generated the random numbers used to make the tests it created numbers that had 15 decimal places to work with and so it used those 15 decimals to compute the cross product, then fortpy only saved out 12 of those decimals, for the input and the output, and so the differences we're seeing could well be due to trailing decimal places that fortpy cutoff. They contributed to the original answer but not to the unit tests answer.

This would explain the numerical issues popping up in other places as well. I'll see what I can do to verify this. Probably have mathematica read in the inputs and save new outputs with the correct number of trailing decimals.

Wiley Morgan http://hagrid.byu.edu/Wiley#Bio

On Fri, May 29, 2015 at 9:53 AM, Gus Hart notifications@github.com wrote:

I'm surprised at how large the errors can be. It's much bigger than the normal finite precision errors I'm used to.

— Reply to this email directly or view it on GitHub https://github.com/msg-byu/symlib/issues/4#issuecomment-106855024.

glwhart commented 9 years ago

Good catch Wiley. A much better answer than my original stab in the dark. You probably don't need to go to mathematica to correct it. Just a short python script, reprinting only to 12 places will do the trick.

wsmorgan commented 9 years ago

I managed to get the output files to agree with the input files now. So I'm closing this issue.

glwhart commented 9 years ago

Fabulous. Love to see the process working. We're getting better at software!