That's correct, the earth is an oblate spheroid (sort of), so any calculation based on a sphere (as mine is) has some amount of error, that varies from place to place on the planet.

In my case, because I have used the equatorial radius (the maximum) the error varies from about 0% at the equator to about +2% at the poles (where we rarely drive or fly). If you bump the radius down to the mean (3959) the error shifts to about +/-1%. And yes, there is a more complicated formula that corrects for this.

However, even that formula incorporates some (smaller) error, because the earth is also not strictly an oblate spheriod. Rather it is slightly asymmetric because the northern hemisphere and the southern hemisphere are not the same in their "oblateness" (one is wider, on the average than the other, but I cannot remember which way it goes). Thus the earth is actually slightly "pear-shaped", and yes, there is a still more complicated formula to correct for this.

But it doesn't stop there, there are significant regional variations in the earth too. Then there is the fact that the Earth's deformations actually change over time... And yes, there are functions/data tables that can correct for these as well, but at some point along the way here, we've stopped being an easy function for everyday SQL Server use and started to be NASA.

The point here is that there is no absolute best formula or approach to use, because you have to trade off your time and performance against what accuracy you actually need. I believe that my function is appropiate for most general uses.