http://acm.hdu.edu.cn/showproblem.php?pid=3571

比赛前一天晚上刚跟haozi研究了最小包围球,其中要根据四个点求出球心坐标,然后我提出了这么拓展到N维的情况。今天比赛的时候居然就出了这么一题,虽然知道可以用高斯消元,但是这题的数据好大,有1017那么大,而且要是整数解,用double精度根本不够。我想利用大整数,要队友写了一个java的,每次消元的时候取最小公倍数,结果TLE了。。。比赛的时候只有haozi他们队和电子科大的一个队伍过了,Orz!!!

赛后,问了一下haozi,因为最后的答案在[ -1017  , 1017] 的范围内,可以取一个大于2×1017的素数p,计算的过程中都 mod p,由于坐标有正有负,可以先将球心坐标都加上1017  这么一个偏移量,最后求出答案之后在剪掉。
PS: 高斯消元改自浙大模板

  1 #include <cstdio>
  2 #include <iostream>
  3 #include <map>
  4 #include <queue>
  5 #include <complex>
  6 #include <algorithm>
  7 #include <cmath>
  8 #include <cstring>
  9 using namespace std;
 10 typedef __int64 ll;
 11 
 12 const ll p =   1000000000000000003LL;
 13 const int maxn = 60;
 14 const ll inf = 100000000000000000LL;
 15 
 16 ll mod( ll a, ll n ) { return ( a % n + n ) % n; }
 17 
 18 ll mul_mod ( ll a, ll b )
 19 
 20     ll ret = 0, tmp = a % p;
 21     while( b )
 22     {     
 23         if( b & 1 )
 24             if( ( ret += tmp ) >= p ) 
 25                 ret -= p;     
 26         if( ( tmp <<= 1 ) >= p ) tmp -= p;  
 27         b >>= 1;   
 28     }  
 29     return ret;
 30 }
 31 
 32 void gcd( ll a, ll b, ll & d, ll & x, ll & y )
 33 {
 34     if! b ) { d = a; x = 1; y = 0; }
 35     else
 36     {
 37         gcd( b, a % b, d, y, x );
 38         y -= x * ( a / b );
 39     }
 40 }
 41 
 42 ll inv( ll a, ll n ) 
 43 {
 44     ll d, x, y;
 45     gcd( a, p, d, x, y );
 46     return d == 1 ? mod( x, p ) : -1;
 47 }
 48 
 49 ll ABS( ll a ) { return ( a < 0 ? -a : a ); }
 50 
 51 int Gauss( int n, ll a[][maxn], ll b[] )
 52 {
 53     int i, j, k, row;
 54     ll maxp, t;
 55     for( k = 0; k < n; k++ )
 56     {
 57         for( maxp = 0, i = k; i < n; i++ )
 58             if( ABS( a[i][k] ) > ABS( maxp ) )
 59                 maxp = a[row=i][k];
 60         if( maxp == 0 ) return 0;
 61         if( row != k )
 62         {
 63             for( j = k; j < n; j++ )
 64                 t = a[k][j], a[k][j] = a[row][j], a[row][j] = t;
 65             t = b[k], b[k] = b[row], b[row] = t;
 66         }
 67         ll INV = inv( maxp, p );
 68         for( j = k + 1; j < n; j++ )
 69         {
 70             a[k][j] = mul_mod( a[k][j], INV );
 71             for( i = k + 1; i < n; i++ )
 72                 a[i][j] = mod( a[i][j] - mul_mod( a[i][k], a[k][j] ), p );
 73         }
 74         b[k] = mul_mod( b[k], INV );
 75         for( i = k + 1; i < n; i++ )
 76             b[i] = mod( b[i] - mul_mod( b[k], a[i][k] ), p );
 77     }
 78     for( i = n - 1; i >= 0; i-- )
 79         for( j = i + 1; j < n; j++ )
 80             b[i] = mod( b[i] - mul_mod( a[i][j], b[j] ), p );
 81     return 1;
 82 }
 83 
 84 int main(int argc, char *argv[])
 85 {
 86     int cas, n;
 87     ll a[maxn][maxn], b[maxn], c[maxn][maxn], d[maxn];
 88     scanf("%d",&cas);
 89     forint t = 1; t <= cas; t++ )
 90     {
 91         scanf("%d",&n);
 92         forint i = 0; i <= n; i++ )
 93         {
 94             b[i] = 0;
 95             forint j = 0; j < n ; j++ )
 96             {
 97                 scanf("%I64d",&a[i][j]);
 98                 a[i][j] += inf;
 99                 b[i] = ( b[i] + mul_mod( a[i][j], a[i][j] ) ) % p;
100                 a[i][j] = ( a[i][j] + a[i][j] ) % p;
101             }
102         }
103         forint i = 0; i < n; i++ )
104         {
105             forint j = 0; j < n; j++ )
106                 c[i][j] = mod( a[i][j] - a[n][j], p );
107             d[i] = mod( b[i] - b[n], p );
108         }
109         Gauss( n, c, d );
110         //gauss_cpivot( n, c, d );
111         printf("Case %d:\n",t);
112         printf("%I64d",d[0]-inf);
113         forint i = 1; i < n; i++ )
114             printf(" %I64d",d[i]-inf);
115         printf("\n");
116     }
117     return 0;
118 }