Notes for Intel® oneAPI Math Kernel Library Vector Statistics

ID 772987
Date 12/04/2020
Public
Document Table of Contents

Example of VS Use

A typical algorithm for VS generators is as follows:

  1. Create and initialize a stream/streams. You can use functions vslNewStream, vslNewStreamEx, vslCopyStream, vslCopyStreamState, vslLeapfrogStream, vslSkipAheadStream, vslSkipAheadStreamEx.

  2. Call one or more RNGs.

  3. Process the output.

  4. Delete the stream/streams. Use function vslDeleteStream.

You may reiterate steps 2-3. Random number streams may be generated for different threads.

The following example demonstrates generation of two random streams:

  1. The first stream is the output of the basic generator MCG31m1. This stream is used to generate 1,000 normally distributed random numbers in blocks of 100 random numbers with parameters a = 5 and sigma = 2.

  2. The second stream is the output of the basic generator R250. This stream is used to produce 1,000 exponentially distributed random numbers in blocks of 100 random numbers with parameters a = -3 and beta = 2.

For each of the streams, the seeds are equal to 1. Delete the streams after completing the generation. The purpose is to calculate the sample mean for normal and exponential distributions with the given parameters.

#include <stdio.h>
#include "mkl.h"
float rn[100], re[100];        /* buffers 
 for random numbers */
float sn, se;                  /* 
 averages */
VSLStreamStatePtr streamn, streame;
int i, j;
/* Initializing the streams*/
sn = 0.0f;
se = 0.0f;
vslNewStream( &streamn, VSL_BRNG_MCG31, 1 );
vslNewStream( &streame, VSL_BRNG_R250, 1 );
/* Generating */
for ( i=0; i<10; i++ )
{
  vsRngGaussian( VSL_METHOD_SGAUSSIAN_BOXMULLER2,
                 streamn, 
 100, rn, 5.0f, 2.0f );
  vsRngExponential(VSL_RNG_METHOD_EXPONENTIAL_ICDF,
                        streame, 
 100, re, -3.0f, 4.0f );
  for ( j=0; j<100; j++ )
  {
    sn += rn[j];
    se += re[j];
  }
}
sn /= 1000.0f;
se /= 1000.0f;
/* Deleting the streams */
vslDeleteStream( &streamn );
vslDeleteStream( &streame );
/* Printing the results */
 
printf( "Sample mean of normal distribution = %f\n", sn );
printf( "Sample mean of exponential distribution = %f\n", se 
 );

When you call a generator of random numbers of normal (Gaussian) distribution, use the named constant VSL_RNG_METHOD_GAUSSIAN_BOXMULLER2 to invoke the Box-Muller2 generation method. In the case of a generator of exponential distribution, assign the method by the named constant VSL_RNG_METHOD_EXPONENTIAL_ICDF.

The following example generates 100 three-dimensional quasi-random vectors in the (2,3)3 hypercube using SOBOL BRNG.

#include <stdio.h>
#include "mkl.h"
float r[100][3]; /* buffer for quasi-random numbers */
VSLStreamStatePtr stream;
 
/* Initializing the streams*/
vslNewStream( &stream, VSL_BRNG_SOBOL, 3 );
/* Generating */
vsRngUniform( VSL_RNG_METHOD_UNIFORM_STD,
              stream, 
 100*3, (float*)r, 2.0f, 3.0f );
/* Deleting the streams */
vslDeleteStream( &stream );

The example below demonstrates the use of a non-deterministic random number generator for initialization of Mersenne Twister generator in parallel computations, so that each run of the application produces different random number sequences:

#define N 6024
#define METHOD1 VSL_RNG_METHOD_UNIFORMBITS_STD
#define METHOD2 VSL_RNG_METHOD_GAUSSIAN_ICDF
#define M 1024
int main()
{
    int i;
    VSLStreamStatePtr istream;
    VSLStreamStatePtr stream[N];
    unsigned seed[N];
    double buf[M];
    /* Initialize non-deterministic random number generator 
 */
    errcode = vslNewStream( &stream, VSL_BRNG_NONDETERM, 
 VSL_BRNG_RDRAND );
    if (errcode == VSL_RNG_ERROR_NONDETERM_NOT_SUPPORTED 
 )
    {
         printf(“CPU does 
 not support non-deterministic generator”);
         exit( 1 );
    }
 
    /* Generate seeds for MT2203 generator */
    errcode = viRngUniformBits( METHOD1, istream, N, 
 seed );
    /* De-initialize */
    errcode = vslDeleteStream( &istream );
    for ( i = 0; i < N; i++ )
    {
        errcode = vslNewStream( 
 &stream[i], VSL_BRNG_MT2203, seed[i] );
        errcode = vdRngGaussian( 
 METHOD2, stream[i], M, buf, 0.0, 1.0 );
        /* Process Gaussian numbers 
 in the buffer buf */    
        ...
         errcode = vslDeleteStream( 
 &stream[i] );
    }
    return 0;
}

The example below demonstrates the use of a non-deterministic random number generator in parallel Monte Carlo simulations:

#define NCORE 4
#define METHOD VSL_RNG_METHOD_GAUSSIAN_ICDF
#define M 1024
int main()
{
    int i;
    VSLStreamStatePtr stream[NCORE];
    for ( i = 0; i < NCORE; i++ )
    {
        /* Initialize non-deterministic 
 random number generator */
        errcode = vslNewStream( 
 &stream[i], VSL_BRNG_NONDETERM, VSL_BRNG_RDRAND );
        if (errcode == VSL_RNG_ERROR_NONDETERM_NOT_SUPPORTED 
 )
        {
            printf(“CPU 
 does not support non-deterministic generator”);
            exit( 
 1 );
        }
    }
 
    #pragma omp parallel for
    for ( i = 0; i < NCORE; i++ )
    {
        double buf[M];
        errcode = vdRngGaussian( 
 METHOD, stream[i], M, buf, 0.0, 1.0 );
        /* Process Gaussian numbers 
 in the buffer buf */    
        ...
    }
    /* De-initialize */
    for ( i = 0; i < NCORE; i++ )
    {
        errcode = vslDeleteStream( 
 &stream[i] );
    }
    return 0;
}