hdu-3571:N-dimensional Sphere+高斯消元

#include <iostream>
#include <cstdio>
#include <cstring>
#include <vector>
#include <algorithm>
using namespace std;
typedef long long LL;
const LL MOD = (LL)200000000 * 1000000000 + 3;/* 一个大于2*10^17的质数 */
const LL CORRECT = (LL)100000000 * 1000000000;/* 10^17 */

#define MAXN 50

LL broaden[MAXN][MAXN + 1]; /* 增广矩阵 */
LL ans[MAXN]; /* 用于保存答案 */
LL dataX[MAXN + 1][MAXN]; /* 原数据点 */

/* 计算a的模MOD值 */
LL getMod(LL a)
{
if (a < 0){
a += MOD;
}
else if (a >= MOD){
a -= MOD;
}
return a;
}

/* 二分方法计算a*b的值 */
LL multiMod(LL a, LL b)
{
LL ans = 0, base = a%MOD;
while (b){
if (b & 1){
ans = getMod(base+ans);
}
base = getMod(base << 1);
b >>= 1;
}
return ans;
}

/* 计算第i行的增广矩阵，即第i个式子 Σ[bn(i+1,j)^2-bn(i,j)^2]=2Σ[bn(i+1,j)-bn(i,j)]x(j) */
void calLineBroaden(int i, int n)
{
LL nb = 0;
for (int j = 0; j < n; ++j){
broaden[i][j] = getMod((dataX[i + 1][j] - dataX[i][j]) << 1);
nb += multiMod(dataX[i + 1][j], dataX[i + 1][j]) - multiMod(dataX[i][j], dataX[i][j]);
nb = getMod(nb);
}
broaden[i]
= nb;
}

/* 计算增广矩阵 */
void calBroaden(int n)
{
for (int i = 0; i < n; ++i){
calLineBroaden(i, n);
}
}

/* 扩展欧几里德a*x+b*y=gcd(a,b) */
LL exGCD(LL a, LL b, LL&x, LL&y)
{
LL tmp, d;
if (b == 0){
x = 1, y = 0;
return a;
}
d = exGCD(b, a%b, x, y);
tmp = x, x = y, y = tmp - a / b*y;
return d;
}

/* 高斯消元，求解答案 */
void Gauss(int n)
{
int i, j;
for (i = 0; i < n; ++i){
for (j = 0; j < n; ++j){
if (broaden[i][j])break;
}
if (i < j){
for (int k = 0; k <= n; ++k){
swap(broaden[i][k], broaden[j][k]);
}
}
}
for (i = 0; i < n - 1; ++i){
for (j = i + 1; j < n; ++j){
for (int k = i + 1; k <= n; ++k){
broaden[j][k] = getMod(multiMod(broaden[j][k], broaden[i][i]) - multiMod(broaden[i][k], broaden[j][i]));
}
}
}
LL x, y, g;
for (i = n - 1; i >= 0; --i){
g = exGCD(broaden[i][i], MOD, x, y);
ans[i] = multiMod(x, broaden[i]
);
for (j = 0; j < i; ++j){
broaden[j]
= getMod(broaden[j]
- multiMod(broaden[j][i], ans[i]));
}
}
}

/* 打印结果 */
void printAns(int n)
{
for (int i = 0; i < n; ++i){
printf(i == 0 ? "%lld" : " %lld", (ans[i]-CORRECT)%CORRECT);
}
printf("\n");
}

void gen()
{
int n;
scanf("%d", &n);
for (int i = 0; i <= n; ++i){
for (int j = 0; j < n; ++j){
scanf("%lld", &dataX[i][j]);
dataX[i][j] += CORRECT;
}
}
calBroaden(n);
Gauss(n);
printAns(n);
}

int main()
{
int t;
scanf("%d", &t);
for (int iCase = 1; iCase <= t; ++iCase){
printf("Case %d:\n", iCase);
gen();
}
return 0;
}