Mostrando entradas con la etiqueta metodos. Mostrar todas las entradas
Mostrando entradas con la etiqueta metodos. Mostrar todas las entradas

sábado, 21 de mayo de 2016

Métodos de generación de números aleatorios en lenguaje C

Aquí muestro la implementación, en lenguaje C, de varios tipos de generadores aleatorios mostrados en el siguiente paper

// Method-1: Sum of Uniform Random Variables
double metodo1_Sum_Of_Uniform_Random_Variables() {
    double X[12];
    double z= 0.0;

    /* Generar 12 numeros con distribucion U(0,1)  */
    for (int i = 0; i < 12; i++ ){
        X[i] = ( (double)rand() / (double)RAND_MAX );
        z+= X[i];
    }

    return z-6;
}

// Box-Muller Method
double metodo2_Box_Muller()
{
    //Generar numero aleatorio uniforme U(0,1)
    double u1 = (double)rand() / (double)RAND_MAX;
    double u2 = (double)rand() / (double)RAND_MAX;

    double z1 = sqrt(-2.0*log(u1))*cos(2*PI*u2); // usa la box-muller transformation
    double z2 = sqrt(-2.0*log(u1))*sin(2*PI*u2); // usa la box-muller transformation

    // Aqui promediamos el valor de las 2 Z para evitar tener que reiterar
    // Este paso es opcional si se retorna el valor de z1 o z2
    double zt= (z1+z2)/2;
    return zt;
}

// Polar Box Muller
double metodo3_Polar_Box_Muller() {

    double S, Z1, Z2, U1, U2, V1,V2;

    do
    {
        // #1 Genera dos números entre -1 y 1
        U1 = (double)rand() / RAND_MAX;
        U2 = (double)rand() / RAND_MAX;

        // #2 V1(-1,1)
        V1 = 2. * U1 - 1.; // V1(-1,1)
        V2 = 2. * U2 - 1.; // V2(-1,1)
        S = V1 * V1 + V2 * V2;

        // #3
    } while(S > 1);

    // #4
    Z1 = sqrt (-2. * log(S) / S)*V1;
    Z2 = sqrt (-2. * log(S) / S)*V2;

    return (Z1+Z2)/2;
}

// Inversion Method
double metodo4_Inversion_Method()
{
    // #1
    double U=(double)rand() / RAND_MAX;

    if(U == 0.5)
        return 0.0;

    // #2 3
    double w=sqrt(-2*log(U < 0.5 ? U : 1-U));

    double a[]={ 2.515517, 0.802853, 0.010328 };
    double b[]={ 1.000000, 1.432788, 0.189269, 0.001308 };

    double Z = -w +      ( (a[0] + a[1]*w + a[2] * w*w) /
            (b[0] + b[1]*w + b[2] * w*w + b[3]*w*w*w));

    // #3
    return U > 0.5 ? -Z : Z;
}

// Acceptance Rejection Method
double metodo5_Acceptance_Rejection_Method()
{
    double Y1=0.0;
    double Y2=0.0;

    // #1 Genera dos Ui ~ U(0,1).
    do{
        double U1=(double)rand()/RAND_MAX;
        double U2=(double)rand()/RAND_MAX;

        Y1= -log(U1);
        Y2= -log(U2);

        // #2
    }while(Y2 > ((Y1-1)*(Y1-1))/2.0 );

    double Z=Y1;

    // #3
    double U3=(double)rand()/RAND_MAX;

    return U3<=0.5 ? fabs(Z) : -fabs(Z);
}

// Using Generalized Exponential Distribution
double metodo6_Using_Generalized_Exponential_Distribution() {

    //Generar numero aleatorio uniforme U(0,1)
    double u = (double)rand() / (double)RAND_MAX;
    double x = -log(1-pow(u,0.0775));

    double z=(log(x)-1.0821)/0.3807;

    //retornar valor de z
    return z;
}

// Bol’shev Formula
double metodo7_Bolshev_Formula() {

    // Array de 12 numeros entre 0 y 1
    double U[12];

    // acumulador
    double X=0.0;
    /* Generar 5 numeros con distribucion U(0,1)  */
    for (int j = 0; j < 12; j++ ){
        U[j] = ( (double)rand() / (double)RAND_MAX );

        X+=sqrt(3)*( 2 * U[j] -1 );
    }
    X*=(1/sqrt(5));
    double z= X-0.01*(3*X-pow(X,3));
    return z;
}

// Inversion Method
double metodo8_Inversion_Method() {
    //Generar numero aleatorio uniforme U(0,1)
    double u = (double)rand() / (double)RAND_MAX;

    //retornar valor de z
    return (-log(1)+ log(u)+ 1);
}

// Proposed Method
double metodo9_Proposed_Method(){
    double U  = (double)rand() / (double)RAND_MAX;
    double X1 = tanh(-31.35694+28.77154*U);
    double X2 = tanh(-2.57136-31.16364*U);
    double X3 = tanh(3.94963-1.66888*U);
    double X4 = tanh(2.31229+1.84289*U);
    double Z=0.46615+90.72192*X1-89.36967*X2-96.55499*X3+97.36346*X4;
    return Z;
}

double calcValorReal(double M, double a, double Z)
{
    //// X=Z*a+M
    return Z*a+M;
}