programming - (C++) Monte Carlo pricer for SABR model to test Hagan / Paulot formulas

I'm trying to test the so-called Hagan formula (p.6 of this paper) and the Paulot formula, order 1 only (eq. (43) p.19 of this paper. For this, i'm trying to use both Euler and Milstein scheme described here (p.9, eq. (3.1) and (3.2)) to price a call option, but the results seem not very consistent, so i'm asking myself if my code's right...

That is my C++ function:

double MC_SABR_price(const int& num_sims, const int& num_intervals, const double& F_0, const double& K, const double& alpha, const double& beta, const double& rho, const double& nu, const double& r, const double& T)
{
    double dt = T / num_intervals;
    double F[num_intervals];
    double V[num_intervals];
    F[0] = F_0;
    V[0] = alpha;

    double payoff_sum = 0.0;
    for (int i=1; i
    {
        for (int j=1; j        {
            double Z1 = NormalSimulation();
            double Z2 = NormalSimulation();
            V[j] = V[j-1] * exp((nu * sqrt(dt) * Z1) - (0.5 * nu * nu * dt));
            F[j] = F[j-1] + (V[j-1] * pow(F[j-1], beta) * sqrt(dt) * ((rho * Z1) + (sqrt(1 - (rho * rho)) * Z2)));
            F[j] = max(F[j], 0.0);
        }
        payoff_sum += max(F[num_intervals-1] - K, 0.0);

    }
    return (payoff_sum / num_sims) * exp(-r*T);
}

With those parameters:

double num_sims = 100000;   // Number of simulated asset paths
double num_intervals = 1000;  // Number of intervals for the asset path to be sampled

double F_0 = 5.0;       // Initial forward price
vector K(10);

for (int i=0; i
double r = 0.0;         // Risk-free rate
double T = 2.5;         // One year until expiry

double alpha = 0.3;     // Initial volatility
double beta = 0.7;      // Elasticity
double rho = -0.5;      // Correlation of asset and volatility
double nu = 0.4;        // "Vol of vol"

These are the results I get:

As you can see nothing coincide...