- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content
There is a nice set of Daubechies wavelets coded in old Fortran by John Burkardt from Pittsburg. I have a building that is shaking itself to pieces slowly, but surely. I was thinking of looking at the various acceleration signals from the building using these wavelets. A blasted long shot I am sure.
I enclose a zip of the code changed to F90. I started with a 1600 long random set with D2. The issue is firstly the input, as shown in the figure generated using the code random method.
The mean is 0.496, the stdev is 0.289, 0.5/0.289 is 1.73 so it is a truncated Gaussian, the skew is 0.00176 and the kurt is -1.17
For instance, the uniform distribution (i.e. one that is uniformly finite over some bound and zero elsewhere) is platykurtic. As this one is.
And if I look at the first 50 numbers as a graph it is clearly a simple Fourier series.
The code for the RNG, it is clearly not really random.
!*********************************************************************72
!
! R8VEC_UNIFORM_01 returns a unit pseudorandom R8VEC.
!
! Discussion:
!
! An R8VEC is a vector of R8's.
!
!
! Reference:
!
! Paul Bratley, Bennett Fox, Linus Schrage,
! A Guide to Simulation,
! Springer Verlag, pages 201-202, 1983.
!
! Bennett Fox,
! Algorithm 647:
! Implementation and Relative Efficiency of Quasirandom
! Sequence Generators,
! ACM Transactions on Mathematical Software,
! Volume 12, Number 4, pages 362-376, 1986.
!
! Peter Lewis, Allen Goodman, James Miller,
! A Pseudo-Random Number Generator for the System/360,
! IBM Systems Journal,
! Volume 8, pages 136-143, 1969.
!
! Parameters:
!
! Input, integer N, the number of entries in the vector.
!
! Input/output, integer SEED, the "seed" value, which should NOT be 0.
! On output, SEED has been updated.
!
! Output, double precision R(N), the vector of pseudorandom values.
!
implicit none
integer n
integer i
integer k
integer seed
double precision r(n)
do i = 1, n
k = seed / 127773
seed = 16807 * ( seed - k * 127773 ) - k * 2836
if ( seed .lt. 0 ) then
seed = seed + 2147483647
end if
r(i) = dble ( seed ) * 4.656612875D-10
end do
return
end
A Fourier series is not a random set, it cannot be, even if you pass the random number test.
@mecej4 or @jimdempseythecove or one of the other gurus, please tell me why I am wrong?
Below is the D2 results.
I can change to the MKL randoms, but I cannot see that the result over a bounded width willnot be the same, I will always have a platykurtic. distribution as we do from the acceleration data and FFT's for some components.
Link Copied
- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content
You see, one goes along for years on this forum, and the same people provide the guiding light for a poor programmer. And then they just go quiet. Even though you have never met them, it leaves a hole in your soul.
--------------------------------------------------------------------------------------------------------------------------------------------------------
The problem with using Fortran as a RNG is that the language has strict structure, it is designed to solve engineering problems, not generate RNGs. It repeats in RNG's, there is always a repeat, this has plagued Fortran for the last 80 years. Basica always had a better built in random number generator, one only has to look at the most basic Genetic Algorithm from the SciAmer in the 1980's to see that.
A real RNG from nature is thermal or similar, in the picture everything up to 2.5 Hz is purely random and Gaussian, you cannot get any information on the underlying structure, which is what you want in an RNG, but above that we can determine structure and hence it is not a RNG generator. The thermal line in the FFT data is always the same shape up to a point of structural interaction.
I realize this is not about the Intel Fortran, but an underlying issue.
- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content
John,
If you are going to roll your own RNG, try something along the lines of this:
program Prime
implicit none
integer(8) :: seed = 0, i8temp
integer(8), parameter :: largeCircularPrime = 999331_8, mask=Z'FFFFFFFF'
integer, parameter :: n=10000
real(8) r(n), rTemp, rDivisor
integer :: i
rDivisor = mask + 1
do i=1,n
if(seed == 0) seed = 1
seed = seed * largeCircularPrime
i8temp = iand(ishft(seed, 16), mask) ! random number between 2^32-1 and 0
rTemp = real(i8temp,8)
r(i) = rTemp / rDivisor ! r(i) >= 0, r(i) < 1.0
end do
end program Prime
Note, the above written for clarity of understanding not performance. You can choose a different largeCircular prime and/or change the mask width (with appropriate ishft bits). IOW the above produces pseudo random numbers of 32 bits of precision.
To get more (say 48 bits), widen the mask by 16 bits and reduce the shift by 8 bits.
Jim Dempsey
- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content
Your method has an underlying structure that repeats, the lines should all be along the bottom as a solid block, the popups are the repeats, this is in 4096 numbers, that may appear random but are not. If you use this series looking for a Monte Carlo Analysis you will run into problems with the structure. First 26 points show the structure.
If I measure a vibrating string and I get a non mode frequency < 1 or 2 depending on your luck, you have a random Gaussian, if the vibration is from the sun. If you have an underlying math structure I can find the pattern. It looks random, it may even pass a so called random test, but it is not truly random.
- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content
If you run Gardners genetic algorithm with Fortran and any RNG it always fails on the ones I have tested.
BASICA however does not fail, I just do not know the form of their RNG.
- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content
>> the popups are the repeats, this is in 4096 numbers
Thanks. This would warrant investigation.
What did you use to chart the pop-ups?
This may be due to bad choice of multiplier 999331, can you try a non-circular prime number.
Here is a list of different types of prime numbers:
https://en.wikipedia.org/wiki/List_of_prime_numbers
Try a balanced prime first. (e.g. 5387)
Jim
- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content
The plotting is on EXCEL, for simple stuff it is easier.
My problem is I look for patterns in data. Random is never random.
So your system gives a usable set of a Monte Carlo analysis, but it is a long way from truly random.
If I take a bridge I have, it has a theoretical lowest frequency of 3 Hz, we find that nicely if we listen for an hour. Non Gaussian in a Gaussian stream.
After three years of listening, I found it had a frequency of 1 Hz, that occurred 10 times in 3 million records. If I pull out those ten records the 1 Hz data is Gaussian Random, but once a month it is non Gaussian. 10 in 3 million is a low probability but for my stuff it matters. We now know we have low frequencies that are hard to find.
A prime number generator is a set of circular steps, anyone who says Prime numbers are random, really does not understand the relationship between the primes and the integer number line.
Think of this number system
2 1
3 10
5 100
7 1000
11 10000
So an composite number can come from the second list 1, 10, 100 etc. it is just a weird multiplication, one could say the next prime is 100000 and it is simple to find.
so remove all factors of 2, then 3 pops up, remove all factors of 3 and 6 etc. this goes all the way. There should be a representation system that is not the integer number line that makes it easy. So 49 is just 2000, 9 is 20 etc. , 6 is just 11 and 25 is 200
- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content
This simple set proves the primes are infinite, but they get further apart in a non random way.

- Subscribe to RSS Feed
- Mark Topic as New
- Mark Topic as Read
- Float this Topic for Current User
- Bookmark
- Subscribe
- Printer Friendly Page