lib/initial/channelRand.cc
/**
 ********************************************************************************************************************************************
 * \brief   Function to impose random initial condition for channel flow on the given input velocity field.
 *
 *          The function generates random values for the entire domain (with independent seeds for individual ranks).
 *          The random values are scaled by multiplying with an inverse parabola and then added to the mean velocity field.
 *          The mean field is assumed to have value 1.0
 *
 ********************************************************************************************************************************************
 */
\brief   Function to impose random initial condition for channel flow on the given input velocity field.
The function generates random values for the entire domain (with independent seeds for individual ranks).
The random values are scaled by multiplying with an inverse parabola and then added to the mean velocity field.
The mean field is assumed to have value 1.0
\ brief Function to impose random initial condition for channel flow on the given input velocity field . The function generates random values for the entire domain ( with independent seeds for individual ranks ) . The random values are scaled by multiplying with an inverse parabola and then added to the mean velocity field . The mean field is assumed to have value 1 . 0
void channelRand::initializeField(vfield &uField) {
    if (mesh.rankData.rank == 0) std::cout << "Imposing random initial condition for channel flow" << std::endl << std::endl;
    int randSeed = std::time(0) + mesh.rankData.rank;
    std::srand(randSeed);
    real vScale = 1.0;
#ifdef PLANAR
    for (int i=uField.Vx.F.lbound(0); i <= uField.Vx.F.ubound(0); i++) {
        for (int k=uField.Vx.F.lbound(2); k <= uField.Vx.F.ubound(2); k++) {
            real uP, ndZ;
            real randNum = real(std::rand())/RAND_MAX;
            ndZ = mesh.zStaggr(k)/mesh.zLen;
            uP = randNum*vScale*(4.0*ndZ*(ndZ - 1.0) + 1);
            uField.Vx.F(i, 0, k) = 1.0 + uP;
        }
    }
    for (int i=uField.Vz.F.lbound(0); i <= uField.Vz.F.ubound(0); i++) {
        for (int k=uField.Vz.F.lbound(2); k <= uField.Vz.F.ubound(2); k++) {
            real uP, ndZ;
            real randNum = real(std::rand())/RAND_MAX;
            ndZ = mesh.zStaggr(k)/mesh.zLen;
            uP = randNum*vScale*(4.0*ndZ*(ndZ - 1.0) + 1);
            uField.Vz.F(i, 0, k) = uP;
        }
    }
#else
    for (int i=uField.Vx.F.lbound(0); i <= uField.Vx.F.ubound(0); i++) {
        for (int j=uField.Vx.F.lbound(1); j <= uField.Vx.F.ubound(1); j++) {
            for (int k=uField.Vx.F.lbound(2); k <= uField.Vx.F.ubound(2); k++) {
                real uP, ndZ;
                real randNum = real(std::rand())/RAND_MAX;
                ndZ = mesh.zStaggr(k)/mesh.zLen;
                uP = randNum*vScale*(4.0*ndZ*(ndZ - 1.0) + 1);
                uField.Vx.F(i, j, k) = 1.0 + uP;
            }
        }
    }
    for (int i=uField.Vy.F.lbound(0); i <= uField.Vy.F.ubound(0); i++) {
        for (int j=uField.Vy.F.lbound(1); j <= uField.Vy.F.ubound(1); j++) {
            for (int k=uField.Vy.F.lbound(2); k <= uField.Vy.F.ubound(2); k++) {
                real uP, ndZ;
                real randNum = real(std::rand())/RAND_MAX;
                ndZ = mesh.zStaggr(k)/mesh.zLen;
                uP = randNum*vScale*(4.0*ndZ*(ndZ - 1.0) + 1);
                uField.Vy.F(i, j, k) = uP;
            }
        }
    }
    for (int i=uField.Vz.F.lbound(0); i <= uField.Vz.F.ubound(0); i++) {
        for (int j=uField.Vz.F.lbound(1); j <= uField.Vz.F.ubound(1); j++) {
            for (int k=uField.Vz.F.lbound(2); k <= uField.Vz.F.ubound(2); k++) {
                real uP, ndZ;
                real randNum = real(std::rand())/RAND_MAX;
                ndZ = mesh.zStaggr(k)/mesh.zLen;
                uP = randNum*vScale*(4.0*ndZ*(ndZ - 1.0) + 1);
                uField.Vz.F(i, j, k) = uP;
            }
        }
    }
#endif
}
void channelRand :: initializeField ( vfield & uField ) { if ( mesh . rankData . rank == 0 ) std :: cout << " " << std :: endl << std :: endl ; int randSeed = std :: time ( 0 ) + mesh . rankData . rank ; std :: srand ( randSeed ) ; real vScale = 1.0 ; #ifdef PLANAR for ( int i = uField . Vx . F . lbound ( 0 ) ; i <= uField . Vx . F . ubound ( 0 ) ; i ++ ) { for ( int k = uField . Vx . F . lbound ( 2 ) ; k <= uField . Vx . F . ubound ( 2 ) ; k ++ ) { real uP , ndZ ; real randNum = real ( std :: rand ( ) ) / RAND_MAX ; ndZ = mesh . zStaggr ( k ) / mesh . zLen ; uP = randNum * vScale * ( 4.0 * ndZ * ( ndZ - 1.0 ) + 1 ) ; uField . Vx . F ( i , 0 , k ) = 1.0 + uP ; } } for ( int i = uField . Vz . F . lbound ( 0 ) ; i <= uField . Vz . F . ubound ( 0 ) ; i ++ ) { for ( int k = uField . Vz . F . lbound ( 2 ) ; k <= uField . Vz . F . ubound ( 2 ) ; k ++ ) { real uP , ndZ ; real randNum = real ( std :: rand ( ) ) / RAND_MAX ; ndZ = mesh . zStaggr ( k ) / mesh . zLen ; uP = randNum * vScale * ( 4.0 * ndZ * ( ndZ - 1.0 ) + 1 ) ; uField . Vz . F ( i , 0 , k ) = uP ; } } #else for ( int i = uField . Vx . F . lbound ( 0 ) ; i <= uField . Vx . F . ubound ( 0 ) ; i ++ ) { for ( int j = uField . Vx . F . lbound ( 1 ) ; j <= uField . Vx . F . ubound ( 1 ) ; j ++ ) { for ( int k = uField . Vx . F . lbound ( 2 ) ; k <= uField . Vx . F . ubound ( 2 ) ; k ++ ) { real uP , ndZ ; real randNum = real ( std :: rand ( ) ) / RAND_MAX ; ndZ = mesh . zStaggr ( k ) / mesh . zLen ; uP = randNum * vScale * ( 4.0 * ndZ * ( ndZ - 1.0 ) + 1 ) ; uField . Vx . F ( i , j , k ) = 1.0 + uP ; } } } for ( int i = uField . Vy . F . lbound ( 0 ) ; i <= uField . Vy . F . ubound ( 0 ) ; i ++ ) { for ( int j = uField . Vy . F . lbound ( 1 ) ; j <= uField . Vy . F . ubound ( 1 ) ; j ++ ) { for ( int k = uField . Vy . F . lbound ( 2 ) ; k <= uField . Vy . F . ubound ( 2 ) ; k ++ ) { real uP , ndZ ; real randNum = real ( std :: rand ( ) ) / RAND_MAX ; ndZ = mesh . zStaggr ( k ) / mesh . zLen ; uP = randNum * vScale * ( 4.0 * ndZ * ( ndZ - 1.0 ) + 1 ) ; uField . Vy . F ( i , j , k ) = uP ; } } } for ( int i = uField . Vz . F . lbound ( 0 ) ; i <= uField . Vz . F . ubound ( 0 ) ; i ++ ) { for ( int j = uField . Vz . F . lbound ( 1 ) ; j <= uField . Vz . F . ubound ( 1 ) ; j ++ ) { for ( int k = uField . Vz . F . lbound ( 2 ) ; k <= uField . Vz . F . ubound ( 2 ) ; k ++ ) { real uP , ndZ ; real randNum = real ( std :: rand ( ) ) / RAND_MAX ; ndZ = mesh . zStaggr ( k ) / mesh . zLen ; uP = randNum * vScale * ( 4.0 * ndZ * ( ndZ - 1.0 ) + 1 ) ; uField . Vz . F ( i , j , k ) = uP ; } } } #endif }
\brief   Function to impose random initial condition for channel flow on the given input velocity field.
\ brief Function to impose random initial condition for channel flow on the given input velocity field .
// Seed the random number generator with both time and rank to get different random numbers in different MPI sub-domains // The below factor was recommended to be set to 26.0, but that number appears when V is scaled with friction velocity, U_tau // VELOCITY PERTURBATION FOR PERIODIC CHANNEL FLOW // X-VELOCITY // Z-VELOCITY // VELOCITY PERTURBATION FOR PERIODIC CHANNEL FLOW // X-VELOCITY // ALONG X, THERE IS A MEAN FLOW OF 1.0 UNIT VELOCITY // Y-VELOCITY // Z-VELOCITY