tozangezan's diary

勝手にソースコードをコピペして利用しないでください。

AOJ 2309: Vector Compression

最小有向全域木。
勉強がてら書いてみました。変な幾何パートもあるし、初書きはこういう問題ではやるべきでないと確信。
あと、こういう問題で英字配列の練習をするべきでないというのも確信。

#include<stdio.h>
#include<algorithm>
#include<vector>
#include<math.h>
using namespace std;
double c[110][110];
double g[210][210];
int rev[210];
int v[210];
int w[210];
int u[210];
double eps=1e-9;
int main(){
	int a,b;
	scanf("%d%d",&a,&b);
	for(int i=0;i<b;i++){
		for(int j=0;j<a;j++){
			scanf("%lf",&c[i][j]);
		}
	}
	for(int i=0;i<210;i++)
		for(int j=0;j<210;j++)
			g[i][j]=999999999;
	int s=b;
	for(int i=0;i<b;i++){
		double dist=0;
		for(int j=0;j<a;j++)dist+=c[i][j]*c[i][j];
		g[s][i]=dist;
		for(int j=0;j<b;j++){
			if(i==j)continue;
			double L=-99999999;
			double R=99999999;
			for(int k=0;k<100;k++){
				double m1=(L*2+R)/3;
				double m2=(L+R*2)/3;
				double d1=0;
				double d2=0;
				for(int l=0;l<a;l++){
					d1+=(c[i][l]*m1-c[j][l])*(c[i][l]*m1-c[j][l]);
					d2+=(c[i][l]*m2-c[j][l])*(c[i][l]*m2-c[j][l]);
				}
				g[i][j]=min(g[i][j],min(d1,d2));
				if(d1<d2)R=m2;
				else L=m1;
			}
		}
	}
	double ret=0;
	int sz=b+1;
	int rt=b;
	while(1){
		for(int i=0;i<sz;i++){
			if(v[i])continue;
			if(rt==i)continue;
			rev[i]=-1;
			for(int j=0;j<sz;j++){
				if(v[j])continue;
				if(i==j)continue;
				if(!~rev[i]||g[rev[i]][i]>g[j][i])rev[i]=j;
			}
			
			if(!u[i])ret+=g[rev[i]][i];
			u[i]=1;
		}
		for(int i=0;i<210;i++)w[i]=0;
		int cnt=0;
		int SZ=sz;
		for(int i=0;i<SZ;i++){
			if(v[i])continue;
			if(w[i])continue;
			if(i==rt)continue;
			int at=i;
			while(1){
				if(w[at]==2){
					cnt++;
					vector<int>vs;
					int a2=at;
					if(w[rt]==2)rt=sz;
					while(1){
						if(w[a2]==1)break;
						vs.push_back(a2);
						w[a2]=1;v[a2]=1;
						a2=rev[a2];
					}
					for(int j=0;j<vs.size();j++){
						for(int k=0;k<sz;k++){
							if(g[vs[j]][k]<999999)g[sz][k]=min(g[sz][k],g[vs[j]][k]);
							if(g[k][vs[j]]<999999)g[k][sz]=min(g[k][sz],g[k][vs[j]]-g[rev[vs[j]]][vs[j]]);
						}
					}
					sz++;
					break;
				}
				if(w[at]==1)break;
				if(at==rt)break;
				w[at]=2;
				at=rev[at];
			}
			at=i;
			while(1){
				if(w[at]==1)break;
				if(at==rt)break;
				w[at]=1;
				at=rev[at];
			}
		}
		if(!cnt)break;
	}
	printf("%f\n",ret);
}