Slow projection maths in C++ code I'm optimising our codebase and currently, the slowest part of the code is a method called calculateResourceZ which is used in projections. I'm looking for suggestions on areas of improvement. It needs the accuracy so it has to use doubles.

double calculateResourceZ(const double resource0, const double fwAtWicket, const double decimalOversRemaining, const double lambdaLocal, const double power) const {
    if (fwAtWicket == 0) {
        return 0;
    }

    const double resourceZPart1 = resource0 * fwAtWicket * pow(lambdaLocal, (power + 1.0));
    const double fg = calculateFunctionG(lambdaLocal, decimalOversRemaining);

    const double resourceZPart2 = (1.0 - exp(-((decimalOversRemaining * expDecay * fg * pow(lambdaLocal, -power)) / fwAtWicket)));
    return resourceZPart1 * resourceZPart2;
}

double calculateFunctionG(const double lambda, const double oversRemaining) const {
    if (version &lt; Version::OCT2014) {
        return 1;
    }

    if (oversRemaining == 0 || lambda == 0) {
        return 0;
    }

    const double alphaLambda = calculateAlphaLambda(lambda);
    const double betaLambda = calculateBetaLambda(lambda);

    return pow((oversRemaining / 50.0), -(1.0 + alphaLambda + betaLambda * oversRemaining));
}

double calculateAlphaLambda(const double lambda) const {
    return -1.0 / (1.0 + c1 * (lambda - 1) * exp(-c2 * (lambda - 1.0)));
}

double calculateBetaLambda(const double lambda) const {
    return -c3 * (lambda - 1.0) * exp(-c4 * (lambda - 1.0));
}
This doesn't do much, so if it's the bottleneck it must be because you're calling it a lot. The first thing to do is to look at the calling pattern. If there are lots of calls with the same value of lambda then you may be able to calculate alpha and beta once rather than each time; you may even be able to cache fg and share it between lots of calls with the same lambda and oversRemaining.

On the level of microoptimisation, pow(lambdaLocal, power + 1.0) is mathematically equal to lambdaLocal * pow(lambdaLocal, power) and * pow(lambdaLocal, -power) is mathematically equal to / pow(lambdaLocal, power); I don't see a numerical reason that pulling that pow out would hurt, and it might shave a couple of cycles off. Then (... / lambdaPower) / fwAtWicket is (...) / (fwAtWicket * lambdaPower) which is another common subexpression.

There's potentially another algebraic tweak: decimalOversRemaining * fg could be refactored by including that extra product in the method:

double calculateFunctionGPrime(const double lambda, const double oversRemaining) const {
    if (version &lt; Version::OCT2014) {
        return oversRemaining;
    }

    if (oversRemaining == 0 || lambda == 0) {
        return 0;
    }

    const double alphaLambda = calculateAlphaLambda(lambda);
    const double betaLambda = calculateBetaLambda(lambda);

    return 50 * pow((oversRemaining / 50.0), -(alphaLambda + betaLambda * oversRemaining));
}


which saves a multiplication in the exceptional cases and an addition in the primary case. Then since the calculations of both alpha and beta have a leading negation those can be removed along with the negation in that last pow:

double calculateFunctionGPrime(const double lambda, const double oversRemaining) const {
    if (version &lt; Version::OCT2014) {
        return oversRemaining;
    }

    if (oversRemaining == 0 || lambda == 0) {
        return 0;
    }

    const double alphaLambda = calculateAlphaLambda(lambda);
    const double betaLambda = calculateBetaLambda(lambda);

    return 50 * pow((oversRemaining / 50.0), alphaLambda + betaLambda * oversRemaining);
}

double calculateAlphaLambda(const double lambda) const {
    return 1.0 / (1.0 + c1 * (lambda - 1) * exp(-c2 * (lambda - 1.0)));
}

double calculateBetaLambda(const double lambda) const {
    return c3 * (lambda - 1.0) * exp(-c4 * (lambda - 1.0));
}