Go: // Great-circular distance
// Reference: https://ja.wikipedia.org/wiki/大円距離

package main

import (
	"bufio"
	"os"
	"strings"
	"strconv"
	"math"
	"fmt"
)

var scanner = bufio.NewScanner(os.Stdin)
const (
	R = 6378.1
	M = math.Pi/180
)

func main() {
	for scanner.Scan() {
		buf := strings.Split(scanner.Text(), " ")
		input := make([]float64, len(buf))
		for i := range input {
			input[i], _ = strconv.ParseFloat(buf[i], 64)
		}
		if input[0] == -1 && input[1]== -1&&
			input[2] == -1 && input[3] == -1 {
			break
		}
		y1, x1, y2, x2 := input[0], input[1], input[2], input[3]
		a := math.Cos(y1*M)*math.Cos(y2*M)*math.Cos((x1-x2)*M) +
			math.Sin(y1*M)*math.Sin(y2*M)
		fmt.Println(int(R*math.Acos(a)+0.5))
	}
}

C++: #include <cstdio>
#include <cstdlib>
#include <cmath>
#include <climits>
#include <cfloat>
#include <map>
#include <utility>
#include <set>
#include <iostream>
#include <memory>
#include <string>
#include <vector>
#include <algorithm>
#include <functional>
#include <sstream>
#include <complex>
#include <stack>
#include <queue>
#include <cstring>
#include <sstream>
#include <cassert>
#include <list>
using namespace std;
static const double EPS = 1e-10;
typedef long long ll;
typedef pair<int,int> PI;
typedef vector<int> vi;
#ifndef M_PI
#define M_PI 3.14159265358979323846
#endif
#define rep(i,n) for(int i=0;i<(int)n;++i)
#define FOR(i,c) for(__typeof((c).begin())i=(c).begin();i!=(c).end();++i)
#define ALL(c) (c).begin(), (c).end()
#define mp(a,b) make_pair(a,b)
#define pb(a) push_back(a)
#define F first
#define S second

int dx[]={0,1,0,-1},dy[]={1,0,-1,0};

main(){
  long double a,b,c,d;
  while(cin>>a>>b>>c>>d,b>0){
    a*=M_PI/180;
    b*=M_PI/180;
    c*=M_PI/180;
    d*=M_PI/180;
    long double the=acos(cos(M_PI/2-a)*cos(M_PI/2-c)+sin(M_PI/2-a)*sin(M_PI/2-c)*cos(d-b+2*M_PI));
    printf("%.0Lf\n",the*6378.1);
  }
}

C++: #include <iostream>
#include <vector>
#include <algorithm>
#include <cmath>

using namespace std;

int main(){

	const double PI = 2*acos(0.0);
	double a,b,c,d;
	while(cin>>a>>b>>c>>d&&!(a==-1&&b==-1&&c==-1&&d==-1)){
		pair<double,double> angle[2];
		vector<pair<double,double> > v;
		v.push_back(make_pair(a,b));
		v.push_back(make_pair(c,d));
		for(int i = 0; i < 2; i++){
			angle[i].first=(v[i].first/180)*PI;
			angle[i].second=(v[i].second/180)*PI;
		}
		const double r=6378.1;
		vector<double> points[2];
		for(int i = 0; i < 2; i++){
			points[i].push_back(r*cos(angle[i].first)*cos(angle[i].second));
			points[i].push_back(r*cos(angle[i].first)*sin(angle[i].second));
			points[i].push_back(r*sin(angle[i].first));
		}
		double dist = sqrt((points[0][0]-points[1][0])*(points[0][0]-points[1][0])
			+(points[0][1]-points[1][1])*(points[0][1]-points[1][1])
			+(points[0][2]-points[1][2])*(points[0][2]-points[1][2]));
		double diffAngle=(r*r+r*r-dist*dist)/(2*r*r);


		double anp = acos(diffAngle);
		double comp=PI*r*(anp/PI);
		cout<<(int)(comp+0.5)<<endl;
	}

	return 0;
}

C++: #include<iostream>
#include<cmath>
#define LEN 6378.1
#define PI 3.141592653589793238
using namespace std;
double a, b, c, d, BC, CA, C;
int main() {
    while (true) {
        cin >> a >> b >> c >> d;
        if (a == -1) { break; }
        BC = PI / 180 * (90.0 - a);
        CA = PI / 180 * (90.0 - c);
        C = PI*(d - b) / 180;
        cout << (int)(0.5 + LEN*acos(cos(BC)*cos(CA) + sin(BC)*sin(CA)*cos(C)))<< endl;
    }
}

C++: #include <stdio.h>
#include <string.h>
#include <math.h>
#include <bits/stdc++.h>
using namespace std;

#define EPS 1e-10
#define EQ(a, b) (fabs(a - b) < EPS)
 
int main(void)
{
    const long double radius = 6378.1;
    const long double pi = acos(-1);
    long double x[2], y[2];
    long double dx;
    long double d;
     
    while (1){
        int i;
        
		bool f = true;
		
        for (i = 0; i < 2; i++){
            cin >> y[i] >> x[i];
            if (y[i] != -1 || x[i] != -1) f = false;
			 
            x[i] = x[i] * pi / 180;
            y[i] = y[i] * pi / 180;
        }
		
		if (f) return 0;
         
        dx = x[1] - x[0];
         
        d = radius * acos( sin(y[0]) * sin(y[1]) + cos(y[0]) * cos(y[1]) * cos(dx) );
         
        cout << round(d) << endl;
    }
     
    return 0;
}

C++: #include <iostream>
#include <vector>
#include <string>
#include <sstream>
#include <algorithm>
#include <map>
#include <set>
#include <cstdio>
#include <cmath>
#define rep(i,l,n) for(lint i=l;i<n;i++)
#define rer(i,l,n) for(lint i=l;i<=n;i++)
#define all(a) a.begin(),a.end()
#define o(a) cout<<a<<endl
#define pb(a) push_back(a)
#define mk(a,b) make_pair(a,b)
#define fi first
#define se second
using namespace std;
typedef long long lint;
typedef vector<int> vi;
typedef vector<lint> vli;
typedef vector<vi> vvi;
typedef pair<int,int> pii;

double M_PI=acos(-1);

int main(){
	double a[2],b[2],r=6378.1;
	while(1){
		rep(i,0,2) cin>>a[i]>>b[i];
		if(a[0]==-1 && b[0]==-1 && a[1]==-1 && b[1]==-1) break;
		double x[2],y[2],z[2],d2,l,c=M_PI/180.0;
		/*rep(i,0,2) x[i]=r*cos(a[i]*c)*cos(b[i]*c);
		rep(i,0,2) y[i]=r*sin(a[i]*c)*cos(b[i]*c);
		rep(i,0,2) z[i]=r*sin(b[i]*c);
		d2=(x[1]-x[0])*(x[1]-x[0])+(y[1]-y[0])*(y[1]-y[0])+(z[1]-z[0])*(z[1]-z[0]);
		o(" "<<d2);
		double t=acos(1-d2/(2.0*r*r));
		o(t);
		l=r*t;
		l=int(l+0.5);
		o(l);*/
		double t=acos(cos(a[0]*c)*cos(a[1]*c)*cos(b[0]*c-b[1]*c)+sin(a[0]*c)*sin(a[1]*c));
		double ans=r*t;
		ans=(int)(ans+0.5);
		o(ans);
	}
}

C++: #include <bits/stdc++.h>
#define FOR(i,a,b) for(int i=(a);i<(b);i++)
#define RFOR(i,a,b) for(int i=(b) - 1;i>=(a);i--)
#define REP(i,n) for(int i=0;i<(n);i++)
#define RREP(i,n) for(int i=n-1;i>=0;i--)

#define PB push_back
#define INF INT_MAX/3
#define ALL(a) (a).begin(),(a).end()
#define CLR(a) memset(a,0,sizeof(a))

typedef long long int ll;

using namespace std;

int main(){
    while(true){
	double a,b,c,d;
	const double r = 6378.1;
	const double p = 3.14159265;
	cin >> a >> b >> c >> d;
	if(a==-1&&b==-1&&c==-1&&d==-1)break;
	a *= p / 180.0;
	b *= p / 180.0;
	c *= p / 180.0;
	d *= p / 180.0;
	/*
	if(b > 180.0) b -= 360.0;
	if(d > 180.0) d -= 360.0;
	if(a < c) swap(a,c);
	if(b < d) swap(b,d);
	*/
	double ans;
	double s = acos( cos(a) * cos(c) * cos(b - d) + sin(a)*sin(c) );
	ans = r * s;
	printf("%.f\n",ans);

    }
    return 0;
}

C++: #include <cstdio>
#include <iostream>
#include <sstream>
#include <iomanip>
#include <algorithm>
#include <cmath>
#include <string>
#include <vector>
#include <list>
#include <queue>
#include <stack>
#include <set>
#include <map>
#include <bitset>
#include <numeric>
#include <climits>
#include <cfloat>
using namespace std;

const double L = 6378.1;
const double PI = acos(-1.0);

void solve(double n0, double e0, double& x, double& y, double& z)
{
    double n = n0 / 180 * PI;
    double e = e0 / 180 * PI;

    x = cos(n) * cos(e);
    y = cos(n) * sin(e);
    z = sin(n);
}

int main()
{
    for(;;){
        double a, b, c, d;
        cin >> a >> b >> c >> d;
        if(b < -0.5)
            return 0;

        double x1, y1, z1, x2, y2, z2;
        solve(a, b, x1, y1, z1);
        solve(c, d, x2, y2, z2);

        double dist = sqrt(pow(x1-x2, 2) + pow(y1-y2, 2) + pow(z1-z2, 2));
        double theta = asin(dist / 2) * 2;

        int ret = (int)(L * theta + 0.5);
        cout << ret << endl;
    }
}

C++: #include <cstdio>
#include <cmath>
 
#define r 6378.1
 
int main()
{
    double a, b, c, d;
 
    while (true) {
	scanf("%lf %lf %lf %lf" ,&a ,&b ,&c ,&d);
	if (a == -1 && b == -1 && c == -1 && d == -1) {
	    break;
	}
 
	a = a * M_PI / 180;
	b = b * M_PI / 180;
	c = c * M_PI / 180;
	d = d * M_PI / 180;
 
	double e = cos(a)*cos(b);
	double f = sin(b)*cos(a);
	double g = sin(a);
	double h = cos(c)*cos(d);
	double i = sin(d)*cos(c);
	double j = sin(c);
 
	double dist = acos(e*h + f*i + g*j);  
	printf("%d\n", (int)(r*dist+0.5));
    }
    return 0;
}