根据地球上两个地点的经度和纬度，如何获得这两点的球面距离或直线距离？有无定理公式？

Question

根据地球上两个地点的经度和纬度，如何获得这两点的球面距离或直线距离？
有无定理公式？
假设前提是两个地点均在地表面的零海拔，且地球为理想球体。
假设A点的经度、纬度分别为λA和ΦA， B点的经度、纬度分别为λB和ΦB，d为距离。
D = arc cos((sin北纬A×sin北纬B)＋(cos北纬A×cos北纬B×cosAB两地经度差绝对值) )×地球平均半径
= 6371.004×cos-1[sinΦAsinΦB十cosΦAcosΦBcos(λB—λA)]
其中地球平均半径为6371.004 km，D的单位为km
-------------------------------------------------------

至于网上流传的以下公式，经推导验证都是错误的。
D＝111.12×cos{1/[sinΦAsinΦB十cosΦAcosΦBcos(λB—λA)]}
D＝111.12×cos-1[sinΦAsinΦB十cosΦAcosΦBcos(λB—λA)]
D = arc cos(sin北纬A×sin北纬B＋cos北纬A×cos北纬B×cosAB两地经度差绝对值)÷360×2PI×6371
验证条件：
纬度只差1度时的距离为
D1 = 地球经线或赤道周长÷360
= 6371.004×2×3.1415926536÷360
= 111.19499645809008 km
约111.2千米。

Answer 1

假设前提是两个地点均在地表面的零海拔，且地球为理想球体。
假设A点的经度、纬度分别为λA和ΦA， B点的经度、纬度分别为λB和ΦB，d为距离。
D = arc cos((sin北纬A×sin北纬B)＋(cos北纬A×cos北纬B×cosAB两地经度差绝对值) )×地球平均半径
= 6371.004×cos-1[sinΦAsinΦB十cosΦAcosΦBcos(λB—λA)]
其中地球平均半径为6371.004 km，D的单位为km
-------------------------------------------------------

至于网上流传的以下公式，经推导验证都是错误的。
D＝111.12×cos{1/[sinΦAsinΦB十cosΦAcosΦBcos(λB—λA)]}
D＝111.12×cos-1[sinΦAsinΦB十cosΦAcosΦBcos(λB—λA)]
D = arc cos(sin北纬A×sin北纬B＋cos北纬A×cos北纬B×cosAB两地经度差绝对值)÷360×2PI×6371
验证条件：
纬度只差1度时的距离为
D1 = 地球经线或赤道周长÷360
= 6371.004×2×3.1415926536÷360
= 111.19499645809008 km
约111.2千米。

Answer 2

因为地球是椭球体，所以计算出准确的距离需要更复杂的公式。如果假设地球是球体就容易的多了，
设地球上某点的经度为lambda,纬度为phi，
则这点的空间坐标是
x=cos(phi)*cos(lambda)
y=cos(phi)*sin(lambda)
z=sin(phi)
设地球上两点的空间坐标分别为(x1,y1,z1),(x2,y2,z2)
直线距离即为R*sqrt((x2-x1)*(x2-x1)+(y2-y1)*(y2-y1)+(z2-z1)*(z2-z1)),不过没什么意义
则它们的夹角为
A=acos(x1*x2+y1*y2+z1*z2)，A是角度
则两地距离为
A/180*pi*R，其中R为地球平均半径6371
误差不超过1%

Answer 3

设地球上某点的经度为lambda,纬度为phi，
则这点的空间坐标是
x=cos(phi)*cos(lambda)
y=cos(phi)*sin(lambda)
z=sin(phi)
设地球上两点的空间坐标分别为(x1,y1,z1),(x2,y2,z2)
直线距离即为R*sqrt((x2-x1)*(x2-x1)+(y2-y1)*(y2-y1)+(z2-z1)*(z2-z1)),不过没什么意义
则它们的夹角为
A=acos(x1*x2+y1*y2+z1*z2)，A是角度
则两地距离为
A/180*pi*R，其中R为地球平均半径6371
误差不超过1%

Answer 4

同纬度：#include <stdio.h>/*其中A代表纬度，B和C代表经度，S为所求球面距离*/
#include <math.h>
#define PI 3.1415926
main()
{
float a,b,c,s,m;
scanf("%f%f%f",&a,&b,&c);
m=fabs(b-c);
if(m<180)
s=2*6371*asin(cos(PI/180*a)*sin(PI/360*m));
else
s=2*6371*asin(cos(PI/180*a)*sin(PI/360*(360-m)));
printf("距离为%f\n",s);
}
不同纬度不同经度：/*A(a1,a2)B(b1,b2)代表两点的经纬度*/
#include <stdio.h>
#include <math.h>
#define PI 3.1415926
main()
{
float a1,a2,b1,b2,m,s;
scanf("%f%f%f%f",&a1,&a2,&b1,&b2);
m=fabs(a1-b1);
if(m<180)
s=6371*2*asin(sqrt(0.5*(1-sin(PI/180*a2)*sin(PI/180*b2)-cos(PI/180*a2)*cos(PI/180*b2)*
cos(PI/180*m))));
else
s=6371*2*asin(sqrt(0.5*(1-sin(PI/180*a2)*sin(PI/180*b2)-cos(PI/180*a2)*cos(PI/180*b2)*
cos(PI/180*(360-m)))));
printf("%f\n",s);
}

Answer 5

下面是所写程序，A（A1，A2），B（B1，B2）代表经纬度。
#include <stdio.h>
#include <math.h>
#define PI 3.1415926
main()
{
float a1,a2,b1,b2,m,s;
scanf("%f%f%f%f",&a1,&a2,&b1,&b2);
m=fabs(a1-b1);
if(m<180)
s=6371*2*asin(sqrt(0.5*(1-sin(PI/180*a2)*sin(PI/180*b2)-cos(PI/180*a2)*cos(PI/180*b2)*
cos(PI/180*m))));
else
s=6371*2*asin(sqrt(0.5*(1-sin(PI/180*a2)*sin(PI/180*b2)-cos(PI/180*a2)*cos(PI/180*b2)*
cos(PI/180*(360-m)))));
printf("%f\n",s);