Probability Distribution and Statistical Functions -- Uniform Distribution Module #272
+1,234 −3
chg. complex numbers with kinds
- Loading branch information
commit b58681c58ffcaaf8699cc01a0d5ad743a18cd5ff
There are no files selected for viewing
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters. Learn more about bidirectional Unicode characters
Add this suggestion to a batch that can be applied as a single commit. This suggestion is invalid because no changes were made to the code. Suggestions cannot be applied while the pull request is closed. Suggestions cannot be applied while viewing a subset of changes. Only one suggestion per line can be applied in a batch. Add this suggestion to a batch that can be applied as a single commit. Applying suggestions on deleted lines is not supported. You must change the existing code in this line in order to create a valid suggestion. Outdated suggestions cannot be applied. This suggestion has been applied or marked resolved. Suggestions cannot be applied from pending reviews. Suggestions cannot be applied on multi-line comments. Suggestions cannot be applied while the pull request is queued to merge. Suggestion cannot be applied right now. Please check back later.
Uh oh!
There was an error while loading. Please reload this page.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Note that Goualard mentions dividing integers by a Mersenne number to make a floating-point number to make a point that doing so can lead to non-uniformity issues, and not necessarily to show that libraries should generate random floating-point numbers this way. He makes the same point with dividing integers by a power of 2.
Also, there are ways to implement uniform floating-point number generation that covers all numbers the
realtype covers in [0, 1], but in general they are far from trivial which is why I don't demand it here. Examples include the Rademacher Floating-Point Library by @christophe-conrads as well as an algorithm I give.Uh oh!
There was an error while loading. Please reload this page.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Not quite understand your point. Yes, dividing Mersenne integer leads to non-uniformity. Multiplying the inverse of the same number reduces the issue significantly.
Uh oh!
There was an error while loading. Please reload this page.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
My point is indeed to reduce non-uniformity, especially since this code is intended to form part of the standard library and thus be used by thousands of applications with different requirements, and while some of them may not care about the statistical quality of pseudorandom numbers, others may. In this sense, of the methods that involve multiplying or dividing an integer by a constant to create a binary64, the "best" one is to multiply a binary64 integer in [1, 2^53) (which is representable in binary64) by the binary64 constant 1/2^53. This is "better" than multiplying an integer by a Mersenne number or its inverse, because it "better" ensures the binary64 values are evenly spaced in [0, 1] as the result of the multiplication is representable in binary64. This difference is illustrated if we multiply 5201 by either 1/2^53 or 1/2^53-1 (in exact rational arithmetic), then convert the result to binary64 and back. (Notice that I put "best" and "better" in quotes here, in part because, in my opinion, floating-point random numbers should be avoided whenever random integers or random rational numbers suffice.)
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Based on Goualard paper, the best one is a binary64 integer in [1, 2^53) multiplying a float of 1/(2^53-1) as shown in Fig 7. in the paper. This is exactly what I implemented in the module. Dividing binary64 integer by 2^53 leads to more non-uniformity as GNU Fortran 9.2.0 does.