Why doesn't the Cardano Formula yield the result expected by the inverse equation?

Question

My application is about indoor bicycles. Sometimes I need to estimate the virtual speed for certain amount of Watts on the pedal or the other way around (how much watts would a certain speed require on the pedal). I modeled my equations based on Steve Gribble's page: https://www.gribble.org/cycling/power_v_speed.html

The main equation is this:

In my situation the virtual track the cyclist is going is completely flat so the portions of the equation which are driven by the ascent/descent grade can be eliminated.

I'm working with the following constants: Loss_dt is 2%, temperature is 15 Celsius so Rho = 1.226. W = 83 + 9 = 92 kg. C_rr = 0.005, C_d = 0.63, A = 0.509 m^2. When I calculate P_legs for V = 33.21 km/h then the end result is roughly 200 Watts. G = 0 so the sinus portion is eliminated and the cosine comes out to 1.

However I also need the reverse computation. If you look at the source code of the page it performs a successive approximation (with max iterations cutoff at 100) which hones into the desired speed. During the approximation it repeatedly calls the reverse function so the end result will be expected, and it converges.

However this is a cubic polynomial equation where several components are eliminated and I think Cardano's Formula can be used. I was looking at https://proofwiki.org/wiki/Cardano%27s_Formula and https://brilliant.org/wiki/cardano-method/ for a series of implementable steps.

In our case:

• loss L = 1/(1 - 2/100) = 1/0.98 = 1.020408163
• a = L * 0.5 * C_d * A * Rho = 1.020408163 * 0.5 * 0.63 * 0.509 * 1.226 = 0.200582357
• b = 0
• c = L * g * W * C_rr = 1.020408163 * 9.8067 * 92 * 0.005 = 4.603144897
• d = -P_legs = -200
• Q = (3ac - b^2) / 9a^3 = 3ac/9a^3 = c/3a^2 = 4.603144897 / 3*0.200582357^2 = 38.137123256
• R = (9abc - 27a^2d - 2b^3) / 54a^3 = -27a^2d / 54a^3 = -d / 2a = 0.5 * -P / a = 0.5 * 200 / 0.200582357 = 498.548334438
• D = Q^3 + R^2 > 0 => one root is real and the other two roots are complex
• intermediate variable for S and T: e = sqrt(Q^3 + R^2) = sqrt(55468.164041048 + 248550.441770912) = sqrt(304018.60581196) = 551.378822419
• S = cubert(R + sqrt(Q^3 + R^2)) = cubert(R + e) = cubert(498.548334438 + 551.378822419) = 1049.927156857^(1/3) = 10.163728523
• T = cubert(R - sqrt(Q^3 + R^2)) = cubert(R - e) = cubert(498.548334438 - 551.378822419) = −52.830487981^(1/3) = −3.752276851
• V = x1 = S + T - b/3a = S + T = 10.163728523 −3.752276851 = 6.411451672 m/s
• Multiply by 3.6 to get 23.081226019 km/h which is far from the expected 33.21 km/h

I don't know if I'm misinterpreting some computation during the Cardano Formula, if I am making some mistakes when translating the formula into pseudo code, if I am making some calculation errors with the numbers, or the Cardano is not even applicable? I think it'd be much better to have a closed formula than an iterative estimation algorithm.

---

Dart source code I wrote based on the above:

  double velocityForPowerCardano(int power) {
    final driveTrainFraction = 1 / (1.0 - (driveTrainLoss / 100.0));
    final a = driveTrainFraction * 0.5 * dragCoefficient * frontalArea * airDensity;
    final w = athleteWeight + bikeWeight;
    final c = driveTrainFraction * gConst * w * rollingResistanceCoefficient;
    final q = c / (3 * a * a);
    final r = 0.5 * power / a;
    final e = sqrt(q * q * q + r * r);
    const third = 1 / 3;
    final s = pow(r + e, third);
    final t = pow(r - e, third);
    return (s + t).toDouble();
  }

Answer 1

I contacted Steve Gribble, he updated the page and I compared my code semantically to his JavaScript. My mistake was q = c / (3 * a * a); instead of q = c / (3 * a);. Claude Leibovici properly pointed out that it's not necessarily all about Cardano, but I can use hyperbolic functions too.

As for Dart goes it doesn't like negative numbers to be raised to the power. So I need to work that around and change the signs, currently I'm only doing it for t:

    final driveTrainFraction = 1.0 - (driveTrainLoss / 100.0);
    final a = 0.5 * dragCoefficient * frontalArea * airDensity;
    final c = gConst * (athleteWeight + bikeWeight) * rollingResistanceCoefficient;
    final dNeg = driveTrainFraction * power;
    final q = c / (3 * a);
    final r = dNeg / (2 * a);
    final e = sqrt(q * q * q + r * r);
    const third = 1 / 3;
    final s = pow(r + e, third);
    final rSubE = r - e;
    final negate = rSubE < 0;
    final t = pow(rSubE.abs(), third) * (negate ? -1 : 1);
    return (s + t).toDouble();

So I'm still contemplating about the hyperbolic solution, since my application is on mobile, I'll look if trigonometric functions are speedier or square root / cube root.