JSOI球形空间产生器 (高斯消元)
按照朴素的列方程,可以列出n+1个n元2次方程。
将相邻的两个方程相减就可以得到n个n元1次方程,进行高斯消元就可以了。
1 var a,b:array[0..100,0..100] of extended; 2 temp,ans:array[0..100] of extended; 3 i,j,k,n:longint; 4 cnt:extended; 5 begin 6 readln(n); 7 for i:=1 to n+1 do 8 for j:=1 to n do 9 read(b[i,j]); 10 for i:=1 to n do 11 for j:=1 to n do 12 begin 13 a[i,j]:=2*(b[i+1,j]-b[i,j]); 14 a[i,n+1]:=a[i,n+1]+b[i+1,j]*b[i+1,j]-b[i,j]*b[i,j]; 15 end; 16 //--------------------------------rpCardinal Orz-------------------------------------- 17 for i:=1 to n-1 do 18 begin 19 for j:=i+1 to n do 20 if a[j,i]>a[i,i] then 21 begin 22 temp:=a[i]; 23 a[i]:=a[j]; 24 a[j]:=temp; 25 end; 26 for j:=i+1 to n do 27 begin 28 cnt:=a[j,i]/a[i,i]; a[j,i]:=0; 29 for k:=i+1 to n+1 do 30 a[j,k]:=a[i,k]*cnt-a[j,k]; 31 end; 32 end; 33 ans[n]:=a[n,n+1]/a[n,n]; 34 for i:=n-1 downto 1 do 35 begin 36 for j:=i+1 to n do a[i,n+1]:=a[i,n+1]-ans[j]*a[i,j]; 37 ans[i]:=a[i,n+1]/a[i,i]; 38 end; 39 for i:=1 to n do write(ans[i]:0:3,' '); 40 end.