I have found the answer by reading the source code. This is because the solar positions are updated every 20 days if default settings are used. The solar direction is an averaged direction over the 20 days for the same hour. I set the update frequency to one day, and then the output for Surface Outside Face Beam Solar Incident Angle Cosine Value matched with the calculated ones.

Below is the source code snippet (Solar Shading.cc)

void PerformSolarCalculations(EnergyPlusData &state){

...

    PerDayOfYear = state.dataEnvrn->DayOfYear;
    SumDec = 0.0;
    SumET = 0.0;
    for (Count = 1; Count <= state.dataSolarShading->ShadowingDaysLeft; ++Count) {
        SUN3(PerDayOfYear, SinDec, EqTime);
        SumDec += SinDec;
        SumET += EqTime;
        ++PerDayOfYear;
    }

    //  Compute Period Values
    AvgSinSolarDeclin = SumDec / double(state.dataSolarShading->ShadowingDaysLeft);
    AvgCosSolarDeclin = std::sqrt(1.0 - pow_2(AvgSinSolarDeclin));
    AvgEqOfTime = SumET / double(state.dataSolarShading->ShadowingDaysLeft);
    CalcPerSolarBeam(state, AvgEqOfTime, AvgSinSolarDeclin, AvgCosSolarDeclin);
...

}