4.5
Your rating: None Average: 4.5 (2 votes)

CodeChef submission 20886 (C++ 4.0.0-8)

CodeChef submission 20886 (C++ 4.0.0-8) plaintext list. Status: AC, problem BX, contest APRIL09. By anshuman_singh (Anshuman Singh), 2009-04-15 14:44:31.
  1. #include<iostream>
  2. #include<cstdio>
  3. #include<cstring>
  4. #include<cstdlib>
  5. #include<algorithm>
  6. #include<iterator>
  7. #include<map>
  8. #include<vector>
  9. #include<list>
  10. #include<set>
  11. #include<queue>
  12. #include<cassert>
  13. #include<deque>
  14. #include<stack>
  15. #include<numeric>
  16. #include<sstream>
  17. #include<string>
  18. #include<cmath>
  19. using namespace std;
  20.  
  21. #define LET(x,a) 	__typeof(a) x(a)
  22. #define IFOR(i,a,b) 	for(LET(i,a);i!=(b);++i)
  23. #define EACH(it,v)  	IFOR(it,v.begin(),v.end())
  24. #define FOR(i,a,b)  	for(int i=(int)(a) ; i < (int)(b);++i)
  25. #define	REP(i,n)    	FOR(i,0,n)
  26. #define SZ		size()
  27. #define PB		push_back
  28. #define EPS		1e-9
  29. #define D_EQ(a,b) 	(a>b-EPS && a<b+EPS)
  30. #define D_LT(a,b) 	((a)<((b)-EPS))
  31. #define D_LTEQ(a,b) 	((a)<((b)+EPS))
  32. #define D_GT(a,b) 	((a)>((b)+EPS))
  33. #define D_GTEQ(a,b) 	((a)>((b)-EPS))
  34. #define ABS(a) 		((a)>=0?(a):(-1.0*(a)))
  35. #define V(x) vector< x >
  36. #define INF 100000000
  37. #define SQ(a) ((a)*(a))
  38.  
  39. typedef vector<int>	VI;
  40. typedef long long	LL;
  41. typedef pair<int,int>	PI;
  42. double x[2000],y[2000],dist[2000];
  43. int n,connectedTo[2000],edgeCnt;
  44. bool done[2000];
  45. PI TreeEdges[10000],tempEdges[10000];
  46.  
  47. /* Basic definitions of point class and line class */
  48. struct Vector {
  49. 	double x, y;
  50. 	inline Vector operator-(const Vector &v) const{
  51. 		return Vector(x-v.x, y-v.y);
  52. 	}
  53. 	inline Vector operator+(const Vector &v) const{
  54. 		return Vector(x+v.x, y+v.y);
  55. 	}
  56. 	inline Vector operator/(double fact) const{
  57. 		return Vector(x/fact, y/fact);
  58. 	}
  59. 	inline Vector operator*(double fact) const{
  60. 		return Vector(x*fact, y*fact);
  61. 	}
  62. 	inline bool operator==(const Vector &v) const{
  63. 		return ABS(x-v.x) < EPS && ABS(y-v.y) < EPS;
  64. 	}
  65. 	inline bool operator<(const Vector &v)const{
  66. 		if(D_EQ(x,v.x)){
  67. 			return (D_LT(y,v.y));
  68. 		}
  69. 		return D_LT(x,v.x);
  70. 	}
  71. 	inline Vector(){};
  72. 	inline Vector(double x, double y):x(x),y(y){};
  73. };
  74. typedef Vector Point; 
  75.  
  76. struct Segment{
  77. 	Point c1,c2;
  78. 	inline bool operator==(const Segment &s) const {return c1==s.c1 && c2==s.c2;}
  79. 	inline Segment(){};
  80. 	inline Segment(const Point &c1, const Point &c2) : c1(c1), c2(c2) {};
  81. }; 
  82.  
  83. struct Line{
  84. 	Point c1, c2;
  85. 	inline bool operator==(const Line &l) const {return c1==l.c1 && c2==l.c2;}
  86. 	inline Line(){};
  87. 	inline Line(const Point &c1, const Point &c2) : c1(c1), c2(c2) {};
  88. 	inline Line(const Segment &s):c1(s.c1),c2(s.c2){};
  89. };
  90.  
  91.  
  92. inline double lensqr(const Vector &v){
  93. 	return v.x*v.x + v.y*v.y;
  94. }
  95.  
  96. inline double len(const Vector &v){
  97. 	return sqrt(lensqr(v));
  98. }
  99.  
  100. inline double dotp(const Vector &v1, const Vector &v2){
  101. 	return v1.x*v2.x + v1.y*v2.y;
  102. }
  103.  
  104. inline double crossp(const Vector &v1, const Vector &v2){
  105. 	return v1.x*v2.y - v1.y*v2.x;
  106. }
  107.  
  108. inline double dotp(const Line &l, const Point &p){
  109. 	return dotp(l.c2-l.c1, p-l.c1);
  110. }
  111.  
  112. inline double crossp(const Line &l, const Point &p){
  113. 	return crossp(l.c2-l.c1, p-l.c1);
  114. } 
  115.  
  116. inline double dotp(const Segment &l, const Point &p){
  117. 	return dotp(l.c2-l.c1, p-l.c1);
  118. }
  119.  
  120. inline double crossp(const Segment &l, const Point &p){
  121. 	return crossp(l.c2-l.c1, p-l.c1);
  122. } 
  123.  
  124. inline double operator*(const Vector &v1,const Vector &v2){
  125. 	return dotp(v1,v2);
  126. }
  127.  
  128. inline double operator^(const Vector &v1,const Vector &v2){
  129. 	return crossp(v1,v2);
  130. }
  131.  
  132. inline double pointLineDistance(const Line &l,const Point &p){
  133. 	double dist = ((l.c2-l.c1)^(p-l.c1))/len(l.c2-l.c1);
  134. 	return ABS(dist);
  135. }
  136.  
  137. inline double pointSegmentDistance(const Segment &s,const Point &p){
  138. 	double dist = ((s.c2-s.c1)^(p-s.c1))/len(s.c2-s.c1),dot1,dot2;
  139. 	dot1 = (p-s.c2)*(s.c2-s.c1);if(dot1 > EPS) return dot1;
  140. 	dot2 = (p-s.c1)*(s.c1-s.c2);if(dot2 > EPS) return dot2;
  141. 	return ABS(dist);
  142. }
  143.  
  144. inline double area(vector<Point> &p){
  145. 	double res=0.0;
  146. 	REP(i,p.SZ) res += p[i].x*p[(i+1)%p.SZ].y - p[(i+1)%p.SZ].x*p[i].y;
  147. 	return ABS(res)/2.0;
  148. }
  149.  
  150. inline Point rotate(Point p,double theta){
  151. 	theta = theta*acos(0)/90;
  152. 	return Point(p.x*cos(theta) - p.y*sin(theta),p.x*sin(theta)+p.y*cos(theta));
  153. }
  154.  
  155.  
  156. Point intersectionPoint( const Line &l1, const Line &l2 ){
  157. 	Point ret(INF, INF);
  158. 	double d1 = crossp(l1, l2.c1);
  159. 	double d2 = crossp(l1, l2.c2);
  160. 	if(ABS(d1-d2) < EPS ) return ret;
  161. 	ret.x = l2.c1.x + (l2.c2.x-l2.c1.x)*d1/(d1-d2);
  162. 	ret.y = l2.c1.y + (l2.c2.y-l2.c1.y)*d1/(d1-d2);
  163. 	return ret;
  164. } 
  165.  
  166. Point intersectionPoint(const Segment &l1, const Segment &l2){
  167. 	Point ret(INF, INF);
  168. 	double d1 = crossp(l1, l2.c1);
  169. 	double d2 = crossp(l1, l2.c2);
  170. 	if(ABS(d1-d2) < EPS ) return ret;
  171. 	ret.x = l2.c1.x + (l2.c2.x-l2.c1.x)*d1/(d1-d2);
  172. 	ret.y = l2.c1.y + (l2.c2.y-l2.c1.y)*d1/(d1-d2);
  173. 	if(ret.x >= (l1.c1.x<?l1.c2.x)-EPS && ret.x <= (l1.c1.x>?l1.c2.x)+EPS){
  174. 		if(ret.y >= (l1.c1.y<?l1.c2.y)-EPS && ret.y <= (l1.c1.y>?l1.c2.y)+EPS){
  175. 			if(ret.x >= (l2.c1.x<?l2.c2.x)-EPS && ret.x <= (l2.c1.x>?l2.c2.x)+EPS){
  176. 				if(ret.y >= (l2.c1.y<?l2.c2.y)-EPS && ret.y <= (l2.c1.y>?l2.c2.y)+EPS){
  177. 					return ret;
  178. 				}
  179. 			}
  180. 		}
  181. 	}
  182. 	return Point(INF,INF);
  183. }
  184.  
  185. double highestAngle(Point p1,Point p2,Point p3){
  186. 	double s[3];
  187. 	s[0]=len(p1-p2),s[1]=len(p2-p3),s[2]=len(p3-p1);
  188. 	sort(s,s+3);
  189. 	REP(i,3)if(D_EQ(s[i],0))return 180;
  190. 	double ang = (double)acosl((s[0]*s[0]+s[1]*s[1]-s[2]*s[2])/(2.0*s[0]*s[1]));
  191. 	return (ang*90.0)/acos(0);
  192. }
  193.  
  194. Line GetLine(Point p1,Point p2,Point p3){ //getting the point on the equilateral triangle of p1,p2,m and forming line p3,m
  195. 	Line L(p1,p2);
  196. 	Point mid = (p1 + p2);
  197. 	mid.x/=2.0;mid.y/=2.0;
  198. 	Vector tmp = mid - p1; 
  199. 	tmp=tmp/len(tmp);
  200. 	tmp = rotate(tmp,90);
  201. 	double SideLen = len(p1-p2);
  202. 	SideLen*=(sqrtl(3.0)/2.0);
  203. 	Point m=tmp*SideLen;
  204. 	m=m+mid;
  205. 	if(crossp(L,m)*crossp(L,p3)<0){
  206. 		return Line(m,p3);
  207. 	}
  208. 	else {
  209. 		m=m-mid;
  210. 		m=mid-m;
  211. 		return Line(m,p3);
  212. 	}
  213. }
  214. Point ToricelliPoint(Point p1,Point p2,Point p3){
  215. 	Line l1=GetLine(p1,p2,p3),l2=GetLine(p1,p3,p2);
  216. 	return intersectionPoint(l1,l2);
  217. }
  218.  
  219. Point locallyOptimize(Point t,Point p1,Point p2,Point p3){
  220. 	bool bettered=1;
  221. 	double d1=len(Point(t-p1));d1/=max(1.0,log(d1));
  222. 	double d2=len(Point(t-p2));d2/=max(1.0,log(d2));
  223. 	double d3=len(Point(t-p3));d3/=max(1.0,log(d3));
  224. 	double curDist = d1+d2+d3;
  225. 	vector<pair<int,double> > v;
  226. 	REP(i,edgeCnt){
  227. 		int k1=TreeEdges[i].first,k2=TreeEdges[i].second;
  228. 		bool poss1=0,poss2=0;
  229. 		Point n1(x[k1],y[k1]),n2(x[k2],y[k2]);
  230. 		if(n1==p1||n1==p2||n1==p3)poss1=1;
  231. 		if(n2==p1||n2==p2||n2==p3)poss2=1;
  232. 		double l=len(n1-n2);l/=max(1.0,log(l));
  233. 		if(!poss1 && poss2){
  234. 			v.PB(make_pair(k1,l));
  235. 		}
  236. 		else if(poss1 && !poss2)v.PB(make_pair(k2,l));
  237. 	}
  238.  
  239. 	while(bettered){
  240. 		bettered=0;
  241. 		for(double x1=t.x-5;x1<=t.x+5;x1+=0.05){
  242. 			for(double y1=t.y-5;y1<=t.y+5;y1+=0.05){
  243. 				Point t1(x1,y1);
  244. 				double d11=len(Point(t1-p1));d11/=max(1.0,log(d11));
  245. 				double d12=len(Point(t1-p2));d12/=max(1.0,log(d12));
  246. 				double d13=len(Point(t1-p3));d13/=max(1.0,log(d13));
  247. 				double sub=0;
  248. 				REP(i,(int)v.size()){
  249. 					Point n1(x[v[i].first],y[v[i].first]);
  250. 					double dist=len(n1-t1);dist/=max(1.0,log(dist));
  251. 					sub+=max(0.0,v[i].second-dist);
  252. 				}
  253. 				if(d11+d12+d13-sub<curDist){
  254. 					curDist=d11+d12+d13-sub;
  255. 					t=t1;
  256. 					bettered=1;
  257. 				}
  258. 			}
  259. 		}
  260. 	}
  261. 	return t;
  262. }
  263. int graph[1501][1501];
  264. int main(){
  265. 	scanf("%d",&n);
  266. 	REP(i,n)scanf("%lf%lf",&x[i],&y[i]);
  267. 	//prim MST
  268. 	memset(done,0,sizeof(done));
  269. 	done[0]=1;
  270. 	REP(i,2000){
  271. 		dist[i]=1e20;
  272. 	}
  273. 	REP(i,n)if(i){
  274. 		dist[i]=SQ(x[i]-x[0])+SQ(y[i]-y[0]);
  275. 		connectedTo[i]=0;
  276. 	}
  277. 	edgeCnt=0;
  278. 	REP(i,n-1){
  279. 		double min1=1e20;
  280. 		int pt=-1;
  281. 		REP(i,n)if(!done[i] && dist[i]<min1){
  282. 			min1=dist[i];
  283. 			pt=i;
  284. 		}
  285. 		assert(pt!=-1);
  286. 		done[pt]=1;
  287. 		TreeEdges[edgeCnt++]=make_pair(pt,connectedTo[pt]);
  288. 		REP(i,n)if(!done[i]){
  289. 			double newDist=(SQ(x[i]-x[pt])+SQ(y[i]-y[pt]));
  290. 			if(newDist<dist[i]){
  291. 				dist[i]=newDist;
  292. 				connectedTo[i]=pt;
  293. 			}
  294. 		}
  295. 	}
  296. 	int added=0,cur=n,oldEdgeCnt=edgeCnt;
  297. 	long double profit=0;
  298. 	PI TempTreeEdge[10000];
  299. 	REP(i,edgeCnt)TempTreeEdge[i]=TreeEdges[i];
  300. 	REP(i,edgeCnt){
  301. 		graph[TreeEdges[i].first][TreeEdges[i].second]=1;
  302. 		graph[TreeEdges[i].second][TreeEdges[i].first]=1;
  303. 	}
  304. 	double best = 0;
  305. 	pair<pair<int,int>,int> bestOption;
  306. 	REP(i,edgeCnt+1){
  307. 		if(i==edgeCnt && best>1e-7){
  308. 			int index1=bestOption.first.first,index2=bestOption.first.second,j=bestOption.second;
  309. 			Point p3=Point(x[j],y[j]),p1=Point(x[index1],y[index1]),p2=Point(x[index2],y[index2]);
  310. 			if(highestAngle(p1,p2,p3)+EPS<120.0){
  311. 				Point t1 = ToricelliPoint(p1,p2,p3);
  312. 				Point t=locallyOptimize(t1,p1,p2,p3); //try including other n-1 points in locally optimize too
  313. 				double d1=len(p1-p2),d2=len(p1-p3);
  314. 				d1/=max(1.0,log(d1));d2/=max(1.0,log(d2));
  315. 				double d3=len(p1-t),d4=len(p2-t),d5=len(p3-t);
  316. 				d3/=max(1.0,log(d3));d4/=max(1.0,log(d4));d5/=max(1.0,log(d5));
  317. 				if(d3+d4+d5+log(added+2)-log(added+1)+EPS<d1+d2){
  318. 					int ind = -1,ind2=-1;
  319. 					REP(tmp,edgeCnt){
  320. 						if((j==TreeEdges[tmp].first && index1==TreeEdges[tmp].second)||(j==TreeEdges[tmp].second && index1==TreeEdges[tmp].first)){
  321. 							ind=tmp;
  322. 							break;
  323. 						}
  324. 					}
  325. 					REP(tmp,edgeCnt){
  326. 						if((index2==TreeEdges[tmp].first && index1==TreeEdges[tmp].second)||(index2==TreeEdges[tmp].second && index1==TreeEdges[tmp].first)){
  327. 							ind2=tmp;
  328. 							break;
  329. 						}
  330. 					}
  331. 					graph[index1][index2]=0;
  332. 					graph[index2][index1]=0;
  333. 					graph[index1][j]=0;
  334. 					graph[j][index1]=0;
  335. 					graph[index1][cur]=1;
  336. 					graph[index2][cur]=1;
  337. 					graph[j][cur]=1;
  338. 					graph[cur][index1]=1;
  339. 					graph[cur][index2]=1;
  340. 					graph[cur][j]=1;
  341.  
  342. 					REP(i1,edgeCnt){
  343. 						int k1=TreeEdges[i1].first,k2=TreeEdges[i1].second;
  344. 						bool poss1=0,poss2=0;
  345. 						Point n1(x[k1],y[k1]),n2(x[k2],y[k2]);
  346. 						if(k1==index1||k1==j||k1==index2)poss1=1;
  347. 						if(k2==index1||k2==j||k2==index2)poss2=1;
  348. 						double l=len(n1-n2);l/=max(1.0,log(l));
  349. 						if(!poss1 && poss2){
  350. 							double sub = len(t-n1);sub/=max(1.0,log(sub));
  351. 							sub=l-sub;
  352. 							if(sub>EPS){
  353. 								graph[k1][k2]=graph[k2][k1]=0;
  354. 								graph[k1][cur]=graph[cur][k1]=1;
  355. 								TreeEdges[i1]=make_pair(k1,cur);
  356. 							}
  357. 						}
  358. 						else if(poss1 && !poss2){
  359. 							double sub = len(t-n2);sub/=max(1.0,log(sub));
  360. 							sub=l-sub;
  361. 							if(sub>1e-5){
  362. 								graph[k1][k2]=graph[k2][k1]=0;
  363. 								graph[k2][cur]=graph[cur][k2]=1;
  364. 								TreeEdges[i1]=make_pair(k2,cur);
  365. 							}
  366. 						}
  367. 					}
  368. 					
  369.  
  370. 					TreeEdges[ind2]=make_pair(index1,cur);
  371. 					TreeEdges[ind]=make_pair(j,cur);
  372. 					TreeEdges[edgeCnt++]=make_pair(index2,cur);
  373. 					x[cur]=t.x;y[cur++]=t.y;
  374. 					added++;
  375. 					i=-1;
  376. 					best=0;
  377. 					continue;
  378. 				}
  379. 			}
  380. 		}
  381. 		int index1=TreeEdges[i].first,index2=TreeEdges[i].second;
  382. 		Point p1=Point(x[index1],y[index1]),p2=Point(x[index2],y[index2]);
  383. 		REP(j,cur)if(graph[index1][j] && j!=index1 && j!=index2){
  384. 			Point p3=Point(x[j],y[j]);
  385. 			if(highestAngle(p1,p2,p3)+EPS<120.0){
  386. 				Point t = ToricelliPoint(p1,p2,p3);
  387. 				double d1=len(p1-p2),d2=len(p1-p3);
  388. 				d1/=max(1.0,log(d1));d2/=max(1.0,log(d2));
  389. 				double d3=len(p1-t),d4=len(p2-t),d5=len(p3-t);
  390. 				d3/=max(1.0,log(d3));d4/=max(1.0,log(d4));d5/=max(1.0,log(d5));
  391. 				if(d3+d4+d5+EPS+log(added+2)-log(added+1)<d1+d2 && best<d1+d2-d3-d4-d5-log(added+2)+log(added+1)){
  392. 					best = d1+d2-d3-d4-d5-log(added+2)+log(added+1);
  393. 					bestOption = make_pair(make_pair(index1,index2),j);
  394. 				}
  395. 			}
  396. 		}
  397. 		REP(j,cur)if(graph[index2][j] && j!=index2 && j!=index1){
  398. 			Point p3=Point(x[j],y[j]);
  399. 			if(highestAngle(p1,p2,p3)+EPS<120.0){
  400. 				Point t = ToricelliPoint(p1,p2,p3);
  401. 				double d1=len(p1-p2),d2=len(p2-p3);
  402. 				d1/=max(1.0,log(d1));d2/=max(1.0,log(d2));
  403. 				double d3=len(p1-t),d4=len(p2-t),d5=len(p3-t);
  404. 				d3/=max(1.0,log(d3));d4/=max(1.0,log(d4));d5/=max(1.0,log(d5));
  405. 				if(d3+d4+d5+log(added+2)-log(added+1)+EPS<d1+d2 && best<d1+d2-d3-d4-d5-log(added+2)+log(added+1)){
  406. 					best = d1+d2-d3-d4-d5-log(added+2)+log(added+1);
  407. 					bestOption = make_pair(make_pair(index2,index1),j);
  408. 				}
  409. 			}
  410. 		}
  411. 	}
  412. 	long double cost1=0,cost2=0;
  413. 	cost1 = log(added+1);
  414. 	REP(i,edgeCnt){
  415. 		Point p1=Point(x[TreeEdges[i].first],y[TreeEdges[i].first]),p2=Point(x[TreeEdges[i].second],y[TreeEdges[i].second]);
  416. 		double d=len((p1-p2));
  417. 		d/=max(1.0,log(d));
  418. 		cost1+=d;
  419. 	}
  420. 	REP(i,oldEdgeCnt){
  421. 		Point p1=Point(x[TempTreeEdge[i].first],y[TempTreeEdge[i].first]),p2=Point(x[TempTreeEdge[i].second],y[TempTreeEdge[i].second]);
  422. 		double d=len((p1-p2));
  423. 		d/=max(1.0,log(d));
  424. 		cost2+=d;
  425. 	}
  426. 	if(cost1<cost2){
  427. 		//output the answer
  428. 		printf("%d\n",added);
  429. 		REP(i,added)printf("%lf %lf\n",x[i+n],y[i+n]);
  430. 		printf("%d\n",edgeCnt);
  431. 		REP(i,edgeCnt)printf("%d %d\n",TreeEdges[i].first,TreeEdges[i].second);
  432. 	}
  433. 	else {
  434. 		puts("0");
  435. 		printf("%d\n",oldEdgeCnt);
  436. 		REP(i,oldEdgeCnt)printf("%d %d\n",TempTreeEdge[i].first,TempTreeEdge[i].second);
  437. 	}
  438. 	return 0;	
  439. }
  440.

Comments