Intel® oneAPI Math Kernel Library
Ask questions and share information with other developers who use Intel® Math Kernel Library.

MKL random number stream equivalent to Matlab default RandStream

zuch
Beginner
832 Views

Is it possible to generate random numbers with MKL that are equivalent to Matlab random numbers? 

I use the following Matlab/ Fortran codes but the results are different

Matlab
stream=RandStream('mt19937ar','Seed',0);
RandStream.setGlobalStream(stream);
reset(stream,0)

fid = fopen('~/Desktop/iseed.bin', 'w');
fwrite(fid, stream.State,'int32');
fclose(fid);

rand(10,1)

% 0.8147
% 0.9058
% 0.1270
% 0.9134
% 0.6324
% 0.0975
% 0.2785
% 0.5469
% 0.9575
% 0.9649

Fortran

include 'mkl_vsl.f90'
program main
USE MKL_VSL_TYPE
USE MKL_VSL
implicit none
integer :: params(625)

real(kind=8) r(10) ! buffer for random numbers
TYPE (VSL_STREAM_STATE) :: stream
integer(kind=4) :: errcode
integer(kind=4) :: i,j
integer :: brng,method,seed,n

open(1, file='~/Desktop/iseed.bin', form='binary')
read(1) params
close(1)

n = 10
brng=VSL_BRNG_MT19937
method=VSL_RNG_METHOD_UNIFORM_STD
seed=0
errcode=vslnewstreamex(stream, brng, 625, params)
! alternative option without matlab stream state
! errcode=vslnewstream( stream, brng, seed )
errcode=vdrnguniform(method, stream, n, r, 0.0d0, 1.0d0)
write(*,'(f)')  r(1:10)
errcode=vsldeletestream( stream )
end

! result
0.1327075199224055
0.3464016632642597
0.7798899088520557
0.4143710811622441
0.4759427784010768
0.4244252541102469
0.0815817557740957
0.9338225021492690
0.5113811327610165
0.5184877812862396
0 Kudos
9 Replies
Ying_H_Intel
Employee
832 Views

Hi Zuch, 

Could you please also attach the file ~/Desktop/iseed.bin'.  ? 

Thanks

Ying 

0 Kudos
zuch
Beginner
832 Views

Hi Ying, I'am attaching seed state in binary form. 

Thanks

0 Kudos
Andrey_N_Intel
Employee
832 Views

Hi Zuch,

What value of Seed do you provide into Matlab version of MT19937? Did you provide the same seed into MKL version of the generator?

Initialization of the generator relies on the algorithm which converts the seed (a scalar or array) into its internal state. Thus, if you provide the array of 625 numbers into the generator, the initialization does not take the array as its state, but transforms the array into another state.

Additional details about initialization of MT19937 in MKL are available in VSL Notes, https://software.intel.com/en-us/mkl_11.1_vslnotes, the section "Testing of Basic Random Number Generators\ Basic Random Generator Properties and Testing Results \ MT19937".

Please, let us know, if you have more questions and/or results of your further experiments with Matlab and MKL versions of MT19937 generator.

Thanks,

Andrey

 

 

0 Kudos
zuch
Beginner
832 Views

Thanks Andrey, I will read the details in VSL notes.  

Re: seed value, I use 0 for both MKL and Matlab

 

 

0 Kudos
Andrey_N_Intel
Employee
832 Views

Hi Zuch,

If I correctly understand, you see different outputs between Intel MKL and Matlab versions of MT19937 also in the case when you provide the same seed (=0) into vslNewStream routine. If this is the case, the possible reason might be in the initialization of the generator. MKL version of MT19937 relies on the scheme which corresponds to the initialization function init_by_array() in the original code http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/MT2002/CODES/mt19937ar.c. That is, when you provide a scalar seed, it would be treated as array of size 1, and the state of the generator is initialized properly. While I'm not aware about the initialization of the Matlab version of the generator for case of a scalar seed, it might rely on the algorithm implemented in the init_genrand().

Please, let me know, if you have more questions or need further help.

Thanks,

Andrey

 

0 Kudos
Ying_H_Intel
Employee
832 Views

Hi Zuch, 

Any news regarding the problem? 

The key is the algorithm of the initial array. 

MKL use vector array, when scalar seed, it take it as 1 size array.  How about you try  the vector function: RandStream.create('',',1);

same result?

Best Regards,
Ying

Example 1

Create a single stream and designate it as the current global stream:

s = RandStream('mt19937ar','Seed',1);
RandStream.setGlobalStream(s);

Example 2

Create three independent streams:

[s1,s2,s3] = RandStream.create('mrg32k3a','NumStreams',3);
r1 = rand(s1,100000,1);

 

0 Kudos
JohnNichols
Valued Contributor III
829 Views

 

I grabbed this code to generate a Weiner Process. The bin file is not in the least bit random, it repeats after about 6 steps, how do I generate a more random seed that actually changes each time. 

 

0 Kudos
Andrey_N_Intel
Employee
832 Views

Hi John,

if you need truly random seed you may want to use non-deterministic generator from Intel MKL to produce it. The generator provides interface to RDRAND instruction available on the latest CPUs. However, if you generate a new seed each time before simulation of the Wiener paths, you will run into non-reproducible results that might not be good for goals of your application. What are the reasons you need a random seed that changes each time for simulation of the Wiener process paths? Why other initialization schemes (like choosing the seed by means of, say,RDRAND and fixing it) that allow getting reproducible results do not work for you?

Also, just in case, this link https://software.intel.com/en-us/articles/initializing-Intel-MKL-BRNGs provides additional details on the initialization of Intel MKL generators.

Thanks,

Andrey 

0 Kudos
JohnNichols
Valued Contributor III
829 Views

Thanks -- the problem is I am trying to get huge data sets that are truly random as this is modelling the flickering thermal input to a timber beam.  The interesting element is the FFT of the actual data shows for constant frequency points a Generalized Poissonian distribution in the results.  I prefer to know the points are not reproduced from run to run as I want to see the underlying stuff well and see the effect of the randomness of the process.

John

0 Kudos
Reply