[动态规划] pku1185 状态压缩
发表评论45次阅读2007.09.30 22:09 作者：Felicia 编辑
[阅读更多]
经典的状态压缩DP，状态是f[i][j]，表示第i行，以3进制j为状态。j的位代表一个格子，只能是：0表示第i行和第i – 1行都没有炮兵，1表示第i行没有炮兵而第i-1行有炮兵，2表示第i行有炮兵。然后用DFS进行状态转移。一开始我做了超时，后来预处理了一下合法状态，快了不少，才AC。
下面是我的代码
下载: pku1185.cpp
/**********************************************************************
Author: WHU_GCC
Created Time: 2007-9-29 20:53:51
File Name: pku1185.cpp
Description:
**********************************************************************/
#include <iostream>
using namespace std;
#define out(x) (cout << #x << ": " << x << endl)
typedef long long int64;
const int maxint = 0x7FFFFFFF;
const int64 maxint64 = 0x7FFFFFFFFFFFFFFFLL;
template <class T> void show(T a, int n) { for (int i = 0; i < n; ++i) cout << a[i] << ' '; cout << endl; }
template <class T> void show(T a, int r, int l) { for (int i = 0; i < r; ++i) show(a[i], l); cout << endl; }
const int maxn = 110;
const int maxm = 11;
const int mask[] = {1, 3, 9, 27, 81, 243, 729, 2187, 6561, 19683, 59049, 177147, 531441, 1594323};
int n, m;
int g[maxn][maxm];
int f[maxn][59049];
int ok_state[10000];
int len;
int bits[14];
inline int get_bit(int x, int k)
{
return x % mask[k + 1] / mask[k];
}
void dfs(int origin_state, int origin_line, int now_col, int delta)
{
if (now_col >= m)
{
int new_state = 0;
for (int i = 0; i < m; i++)
if (bits[i] > 0)
new_state += bits[i] * mask[i];
f[origin_line + 1][new_state] >?= f[origin_line][origin_state] + delta;
return;
}
if (bits[now_col] == -1 && g[origin_line + 1][now_col] == 0)
{
bits[now_col] = 2;
dfs(origin_state, origin_line, now_col + 3, delta + 1);
bits[now_col] = -1;
}
dfs(origin_state, origin_line, now_col + 1, delta);
}
void dfs1(int now_col, int new_state, int last1, int last2)
{
if (now_col >= m)
{
ok_state[len++] = new_state;
return;
}
if (last1 != 2 && last2 != 2)
dfs1(now_col + 1, new_state + 2 * mask[now_col], 2, last1);
if (last1 != 1 && last2 != 1)
dfs1(now_col + 1, new_state + mask[now_col], 1, last1);
dfs1(now_col + 1, new_state, 0, last1);
}
void pre_process()
{
len = 0;
dfs1(0, 0, -1, -1);
}
int dp()
{
memset(f, -1, sizeof(f));
f[0][0] = 0;
for (int i = 0; i < n; i++)
{
for (int j = 0; j < len; j++)
if (f[i][ok_state[j]] != -1)
{
for (int k = 0; k < m; k++)
bits[k] = get_bit(ok_state[j], k) - 1;
dfs(ok_state[j], i, 0, 0);
}
}
int ret = 0;
for (int i = 0; i < mask[m]; i++)
ret >?= f[n][i];
return ret;
}
int main()
{
char s[22];
scanf("%d%d", &n, &m);
for (int i = 1; i <= n; i++)
{
scanf("%s", s);
for (int j = 0; j < m; j++)
if (s[j] == 'P')
g[i][j] = 0;
else
g[i][j] = 1;
}
pre_process();
printf("%d\n", dp());
return 0;
}
Author: WHU_GCC
Created Time: 2007-9-29 20:53:51
File Name: pku1185.cpp
Description:
**********************************************************************/
#include <iostream>
using namespace std;
#define out(x) (cout << #x << ": " << x << endl)
typedef long long int64;
const int maxint = 0x7FFFFFFF;
const int64 maxint64 = 0x7FFFFFFFFFFFFFFFLL;
template <class T> void show(T a, int n) { for (int i = 0; i < n; ++i) cout << a[i] << ' '; cout << endl; }
template <class T> void show(T a, int r, int l) { for (int i = 0; i < r; ++i) show(a[i], l); cout << endl; }
const int maxn = 110;
const int maxm = 11;
const int mask[] = {1, 3, 9, 27, 81, 243, 729, 2187, 6561, 19683, 59049, 177147, 531441, 1594323};
int n, m;
int g[maxn][maxm];
int f[maxn][59049];
int ok_state[10000];
int len;
int bits[14];
inline int get_bit(int x, int k)
{
return x % mask[k + 1] / mask[k];
}
void dfs(int origin_state, int origin_line, int now_col, int delta)
{
if (now_col >= m)
{
int new_state = 0;
for (int i = 0; i < m; i++)
if (bits[i] > 0)
new_state += bits[i] * mask[i];
f[origin_line + 1][new_state] >?= f[origin_line][origin_state] + delta;
return;
}
if (bits[now_col] == -1 && g[origin_line + 1][now_col] == 0)
{
bits[now_col] = 2;
dfs(origin_state, origin_line, now_col + 3, delta + 1);
bits[now_col] = -1;
}
dfs(origin_state, origin_line, now_col + 1, delta);
}
void dfs1(int now_col, int new_state, int last1, int last2)
{
if (now_col >= m)
{
ok_state[len++] = new_state;
return;
}
if (last1 != 2 && last2 != 2)
dfs1(now_col + 1, new_state + 2 * mask[now_col], 2, last1);
if (last1 != 1 && last2 != 1)
dfs1(now_col + 1, new_state + mask[now_col], 1, last1);
dfs1(now_col + 1, new_state, 0, last1);
}
void pre_process()
{
len = 0;
dfs1(0, 0, -1, -1);
}
int dp()
{
memset(f, -1, sizeof(f));
f[0][0] = 0;
for (int i = 0; i < n; i++)
{
for (int j = 0; j < len; j++)
if (f[i][ok_state[j]] != -1)
{
for (int k = 0; k < m; k++)
bits[k] = get_bit(ok_state[j], k) - 1;
dfs(ok_state[j], i, 0, 0);
}
}
int ret = 0;
for (int i = 0; i < mask[m]; i++)
ret >?= f[n][i];
return ret;
}
int main()
{
char s[22];
scanf("%d%d", &n, &m);
for (int i = 1; i <= n; i++)
{
scanf("%s", s);
for (int j = 0; j < m; j++)
if (s[j] == 'P')
g[i][j] = 0;
else
g[i][j] = 1;
}
pre_process();
printf("%d\n", dp());
return 0;
}
[动态规划] pku1163 二维DP
发表评论12次阅读2007.09.29 22:43 作者：Felicia 编辑
[阅读更多]
今天郁闷了，贴个小代码。经典的DP入门题。
下面是我的代码
下载: pku1163.cpp
#include <stdio.h>
#include <string.h>
int max,n,i,j,a[110][110],f[110][110];
int main() {
memset(f,0,sizeof(f));
scanf("%d",&n);
for (i=1;i<=n;i++) for (j=1;j<=i;j++) scanf("%d",&a[i][j]);
f[1][1]=a[1][1];
for (i=1;i<n;i++) for (j=1;j<=i;j++) {
f[i+1][j]>?=f[i][j]+a[i+1][j];
f[i+1][j+1]>?=f[i][j]+a[i+1][j+1];
}
for (max=0,i=1;i<=n;i++) max>?=f[n][i];
printf("%d\n",max);
return 0;
}
#include <string.h>
int max,n,i,j,a[110][110],f[110][110];
int main() {
memset(f,0,sizeof(f));
scanf("%d",&n);
for (i=1;i<=n;i++) for (j=1;j<=i;j++) scanf("%d",&a[i][j]);
f[1][1]=a[1][1];
for (i=1;i<n;i++) for (j=1;j<=i;j++) {
f[i+1][j]>?=f[i][j]+a[i+1][j];
f[i+1][j+1]>?=f[i][j]+a[i+1][j+1];
}
for (max=0,i=1;i<=n;i++) max>?=f[n][i];
printf("%d\n",max);
return 0;
}
[计算几何] pku3391 平面欧几里得最小生成树
发表评论116次阅读2007.09.28 20:17 作者：Felicia 编辑
[阅读更多]
问题是求平面欧几里得最小生成树的第n – k小边。
平面欧几里德最小生成树是经典问题，可以做到O(nlogn)。具体做法是先对平面点进行三角剖分，时间复杂度是O(nlogn)，三角剖分的边就是可能的在最小生成树的边。因为是平面图，所以有O(n)条边，在其上应用 Kruscal 算法即可。
下面是我的代码
下载: pku3391.cpp
#include <cstdio>
#include <cmath>
#include <cstring>
#include <cstdlib>
#include <algorithm>
using namespace std;
#define TRUE 1
#define FALSE 0
#define Vector(p1, p2, u, v) (u = p2->x - p1->x, v = p2->y - p1->y)
#define Cross_product_2v(u1, v1, u2, v2) (u1 * v2 - v1 * u2)
#define Cross_product_3p(p1, p2, p3) ((p2->x - p1->x) * (p3->y - p1->y) - (p2->y - p1->y) * (p3->x - p1->x))
#define Dot_product_2v(u1, v1, u2, v2) (u1 * u2 + v1 * v2)
#define Org(e) ((e)->org)
#define Dest(e) ((e)->dest)
#define Onext(e) ((e)->onext)
#define Oprev(e) ((e)->oprev)
#define Dnext(e) ((e)->dnext)
#define Dprev(e) ((e)->dprev)
#define Other_point(e, p) ((e)->org == p ? (e)->dest : (e)->org)
#define Next(e, p) ((e)->org == p ? (e)->onext : (e)->dnext)
#define Prev(e, p) ((e)->org == p ? (e)->oprev : (e)->dprev)
#define Visited(p) ((p)->f)
#define Identical_refs(e1, e2) (e1 == e2)
typedef enum {right, left} side;
typedef double Real;
typedef int boolean;
typedef struct point point;
typedef struct edge edge;
struct point
{
Real x, y;
edge *entry_pt;
};
struct edge
{
point *org, *dest;
edge *onext, *oprev, *dnext, *dprev;
};
#define MAXV 1000001
#define MAXE 3000001
struct EDGE
{
int u, v;
double w;
};
EDGE E[MAXE], mst[MAXE];
int n, en, k;
double dist[MAXV];
int p[MAXV], d[MAXV];
void UFinit(int m)
{
int i;
for(i = 0; i < m; i++)
{
p[i] = i, d[i] = 0;
}
}
int UFfind(int x)
{
int i = x, j, k;
while(p[i] != i) i = p[i];
j = x;
while(p[j] != i)
{
k = p[j];
p[j] = i;
j = k;
}
return i;
}
void UFunion(int x, int y)
{
int i = UFfind(x), j = UFfind(y);
if(d[i] > d[j]) p[j] = i;
else
{
p[i] = j;
if(d[i] == d[j]) d[j]++;
}
}
bool cmp(const EDGE &e1, const EDGE &e2)
{
return e1.w < e2.w;
}
void kruskal(void)
{
if(k > n) k = n;
int i, kk;
sort(E, E + en, cmp);
UFinit(n);
for(i = 0, kk = 0; i < en && kk < n - 1; i++)
if(UFfind(E[i].u) != UFfind(E[i].v))
{
UFunion(E[i].u, E[i].v);
mst[kk++] = E[i];
}
printf("%d\n", int(ceil(mst[n - k - 1].w) + 0.00000001));
}
point *p_array;
static edge *e_array;
static edge **free_list_e;
static int n_free_e;
void alloc_memory(int n)
{
edge *e;
int i;
p_array = (point *)calloc(n, sizeof(point));
n_free_e = 3 * n;
e_array = e = (edge *)calloc(n_free_e, sizeof(edge));
free_list_e = (edge **)calloc(n_free_e, sizeof(edge *));
for(i = 0; i < n_free_e; i++, e++)
free_list_e[i] = e;
}
void free_memory(void)
{
free(p_array);
free(e_array);
free(free_list_e);
}
edge *get_edge(void)
{
if(n_free_e == 0)
printf("Out of memory for edges\n");
return (free_list_e[--n_free_e]);
}
void free_edge(edge *e)
{
free_list_e[n_free_e++] = e;
}
void splice(edge *a, edge *b, point *v)
{
edge *next;
if(Org(a) == v)
{
next = Onext(a);
Onext(a) = b;
}
else
{
next = Dnext(a);
Dnext(a) = b;
}
if(Org(next) == v)
Oprev(next) = b;
else
Dprev(next) = b;
if(Org(b) == v)
{
Onext(b) = next;
Oprev(b) = a;
}
else
{
Dnext(b) = next;
Dprev(b) = a;
}
}
edge *make_edge(point *u, point *v)
{
edge *e;
e = get_edge();
e->onext = e->oprev = e->dnext = e->dprev = e;
e->org = u;
e->dest = v;
if(u->entry_pt == NULL)
u->entry_pt = e;
if(v->entry_pt == NULL)
v->entry_pt = e;
return e;
}
edge *join(edge *a, point *u, edge *b, point *v, side s)
{
edge *e;
e = make_edge(u, v);
if(s == left)
{
if (Org(a) == u)
splice(Oprev(a), e, u);
else
splice(Dprev(a), e, u);
splice(b, e, v);
}
else
{
splice(a, e, u);
if(Org(b) == v)
splice(Oprev(b), e, v);
else
splice(Dprev(b), e, v);
}
return e;
}
void delete_edge(edge *e)
{
point *u, *v;
u = Org(e);
v = Dest(e);
if(u->entry_pt == e)
u->entry_pt = e->onext;
if(v->entry_pt == e)
v->entry_pt = e->dnext;
if(Org(e->onext) == u)
e->onext->oprev = e->oprev;
else
e->onext->dprev = e->oprev;
if(Org(e->oprev) == u)
e->oprev->onext = e->onext;
else
e->oprev->dnext = e->onext;
if(Org(e->dnext) == v)
e->dnext->oprev = e->dprev;
else
e->dnext->dprev = e->dprev;
if(Org(e->dprev) == v)
e->dprev->onext = e->dnext;
else
e->dprev->dnext = e->dnext;
free_edge(e);
}
double px[MAXV], py[MAXV];
void read_points()
{
int i = 0;
while (1)
{
int t1, t2;
scanf("%d", &t1);
px[i] = t1;
if (t1 == -1)
break;
scanf("%d", &t2);
py[i] = t2;
i++;
}
n = i;
}
static void print_edges(int n)
{
edge *e_start, *e;
point *u, *v;
int i;
for(i = 0; i < n; i++)
{
u = &p_array[i];
e_start = e = u->entry_pt;
do
{
v = Other_point(e, u);
if(u < v)
{
E[en].u = u - p_array, E[en].v = v - p_array;
E[en++].w = hypot(u->x - v->x, u->y - v->y);
}
e = Next(e, u);
}while(!Identical_refs(e, e_start));
}
}
void merge_sort(point *p[], point *p_temp[], int l, int r)
{
int i, j, k, m;
if(r - l > 0)
{
m = (r + l) / 2;
merge_sort(p, p_temp, l, m);
merge_sort(p, p_temp, m + 1, r);
for(i = m + 1; i > l; i--)
p_temp[i - 1] = p[i - 1];
for(j = m; j < r; j++)
p_temp[r + m - j] = p[j + 1];
for(k = l; k <= r; k++)
if(p_temp[i]->x < p_temp[j]->x)
{
p[k] = p_temp[i];
i = i + 1;
}
else if(p_temp[i]->x == p_temp[j]->x && p_temp[i]->y < p_temp[j]->y)
{
p[k] = p_temp[i];
i = i + 1;
}
else
{
p[k] = p_temp[j];
j = j - 1;
}
}
}
static void lower_tangent(edge *r_cw_l, point *s, edge *l_ccw_r, point *u, edge **l_lower, point **org_l_lower, edge **r_lower, point **org_r_lower)
{
edge *l, *r;
point *o_l, *o_r, *d_l, *d_r;
boolean finished;
l = r_cw_l;
r = l_ccw_r;
o_l = s;
d_l = Other_point(l, s);
o_r = u;
d_r = Other_point(r, u);
finished = FALSE;
while(!finished)
if(Cross_product_3p(o_l, d_l, o_r) > 0.0)
{
l = Prev(l, d_l);
o_l = d_l;
d_l = Other_point(l, o_l);
}
else if(Cross_product_3p(o_r, d_r, o_l) < 0.0)
{
r = Next(r, d_r);
o_r = d_r;
d_r = Other_point(r, o_r);
}
else
finished = TRUE;
*l_lower = l;
*r_lower = r;
*org_l_lower = o_l;
*org_r_lower = o_r;
}
static void merge(edge *r_cw_l, point *s, edge *l_ccw_r, point *u, edge **l_tangent)
{
edge *base, *l_cand, *r_cand;
point *org_base, *dest_base;
Real u_l_c_o_b, v_l_c_o_b, u_l_c_d_b, v_l_c_d_b;
Real u_r_c_o_b, v_r_c_o_b, u_r_c_d_b, v_r_c_d_b;
Real c_p_l_cand, c_p_r_cand;
Real d_p_l_cand, d_p_r_cand;
boolean above_l_cand, above_r_cand, above_next, above_prev;
point *dest_l_cand, *dest_r_cand;
Real cot_l_cand, cot_r_cand;
edge *l_lower, *r_lower;
point *org_r_lower, *org_l_lower;
lower_tangent(r_cw_l, s, l_ccw_r, u, &l_lower, &org_l_lower, &r_lower, &org_r_lower);
base = join(l_lower, org_l_lower, r_lower, org_r_lower, right);
org_base = org_l_lower;
dest_base = org_r_lower;
*l_tangent = base;
do
{
l_cand = Next(base, org_base);
r_cand = Prev(base, dest_base);
dest_l_cand = Other_point(l_cand, org_base);
dest_r_cand = Other_point(r_cand, dest_base);
Vector(dest_l_cand, org_base, u_l_c_o_b, v_l_c_o_b);
Vector(dest_l_cand, dest_base, u_l_c_d_b, v_l_c_d_b);
Vector(dest_r_cand, org_base, u_r_c_o_b, v_r_c_o_b);
Vector(dest_r_cand, dest_base, u_r_c_d_b, v_r_c_d_b);
c_p_l_cand = Cross_product_2v(u_l_c_o_b, v_l_c_o_b, u_l_c_d_b, v_l_c_d_b);
c_p_r_cand = Cross_product_2v(u_r_c_o_b, v_r_c_o_b, u_r_c_d_b, v_r_c_d_b);
above_l_cand = c_p_l_cand > 0.0;
above_r_cand = c_p_r_cand > 0.0;
if(!above_l_cand && !above_r_cand)
break;
if(above_l_cand)
{
Real u_n_o_b, v_n_o_b, u_n_d_b, v_n_d_b;
Real c_p_next, d_p_next, cot_next;
edge *next;
point *dest_next;
d_p_l_cand = Dot_product_2v(u_l_c_o_b, v_l_c_o_b, u_l_c_d_b, v_l_c_d_b);
cot_l_cand = d_p_l_cand / c_p_l_cand;
do
{
next = Next(l_cand, org_base);
dest_next = Other_point(next, org_base);
Vector(dest_next, org_base, u_n_o_b, v_n_o_b);
Vector(dest_next, dest_base, u_n_d_b, v_n_d_b);
c_p_next = Cross_product_2v(u_n_o_b, v_n_o_b, u_n_d_b, v_n_d_b);
above_next = c_p_next > 0.0;
if(!above_next)
break;
d_p_next = Dot_product_2v(u_n_o_b, v_n_o_b, u_n_d_b, v_n_d_b);
cot_next = d_p_next / c_p_next;
if(cot_next > cot_l_cand)
break;
delete_edge(l_cand);
l_cand = next;
cot_l_cand = cot_next;
}while(TRUE);
}
if(above_r_cand)
{
Real u_p_o_b, v_p_o_b, u_p_d_b, v_p_d_b;
Real c_p_prev, d_p_prev, cot_prev;
edge *prev;
point *dest_prev;
d_p_r_cand = Dot_product_2v(u_r_c_o_b, v_r_c_o_b, u_r_c_d_b, v_r_c_d_b);
cot_r_cand = d_p_r_cand / c_p_r_cand;
do
{
prev = Prev(r_cand, dest_base);
dest_prev = Other_point(prev, dest_base);
Vector(dest_prev, org_base, u_p_o_b, v_p_o_b);
Vector(dest_prev, dest_base, u_p_d_b, v_p_d_b);
c_p_prev = Cross_product_2v(u_p_o_b, v_p_o_b, u_p_d_b, v_p_d_b);
above_prev = c_p_prev > 0.0;
if(!above_prev)
break;
d_p_prev = Dot_product_2v(u_p_o_b, v_p_o_b, u_p_d_b, v_p_d_b);
cot_prev = d_p_prev / c_p_prev;
if(cot_prev > cot_r_cand)
break;
delete_edge(r_cand);
r_cand = prev;
cot_r_cand = cot_prev;
}while(TRUE);
}
dest_l_cand = Other_point(l_cand, org_base);
dest_r_cand = Other_point(r_cand, dest_base);
if(!above_l_cand || (above_l_cand && above_r_cand && cot_r_cand < cot_l_cand))
{
base = join(base, org_base, r_cand, dest_r_cand, right);
dest_base = dest_r_cand;
}
else
{
base = join(l_cand, dest_l_cand, base, dest_base, right);
org_base = dest_l_cand;
}
}while(TRUE);
}
void divide(point *p_sorted[], int l, int r, edge **l_ccw, edge **r_cw)
{
int n;
int split;
edge *l_ccw_l, *r_cw_l, *l_ccw_r, *r_cw_r, *l_tangent;
edge *a, *b, *c;
Real c_p;
n = r - l + 1;
if(n == 2)
{
*l_ccw = *r_cw = make_edge(p_sorted[l], p_sorted[r]);
}
else if(n == 3)
{
a = make_edge(p_sorted[l], p_sorted[l + 1]);
b = make_edge(p_sorted[l + 1], p_sorted[r]);
splice(a, b, p_sorted[l + 1]);
c_p = Cross_product_3p(p_sorted[l], p_sorted[l + 1], p_sorted[r]);
if(c_p > 0.0)
{
c = join(a, p_sorted[l], b, p_sorted[r], right);
*l_ccw = a;
*r_cw = b;
}
else if(c_p < 0.0)
{
c = join(a, p_sorted[l], b, p_sorted[r], left);
*l_ccw = c;
*r_cw = c;
}
else
{
*l_ccw = a;
*r_cw = b;
}
}
else if(n > 3)
{
split = (l + r) / 2;
divide(p_sorted, l, split, &l_ccw_l, &r_cw_l);
divide(p_sorted, split+1, r, &l_ccw_r, &r_cw_r);
merge(r_cw_l, p_sorted[split], l_ccw_r, p_sorted[split + 1], &l_tangent);
if(Org(l_tangent) == p_sorted[l])
l_ccw_l = l_tangent;
if(Dest(l_tangent) == p_sorted[r])
r_cw_r = l_tangent;
*l_ccw = l_ccw_l;
*r_cw = r_cw_r;
}
}
int main(void)
{
int ca;
for (scanf("%d", &ca); ca--;)
{
scanf("%d", &k);
edge *l_cw, *r_ccw;
int i;
point **p_sorted, **p_temp;
read_points();
if (k > n) k = n;
if (n == 1)
{
printf("0\n");
continue;
}
alloc_memory(n);
for (i = 0; i < n; i++)
{
p_array[i].x = px[i];
p_array[i].y = py[i];
}
for(i = 0; i < n; i++)
p_array[i].entry_pt = NULL;
p_sorted = (point **)malloc((unsigned)n * sizeof(point *));
p_temp = (point **)malloc((unsigned)n * sizeof(point *));
for(i = 0; i < n; i++)
p_sorted[i] = p_array + i;
merge_sort(p_sorted, p_temp, 0, n - 1);
free((char *)p_temp);
divide(p_sorted, 0, n - 1, &l_cw, &r_ccw);
free((char *)p_sorted);
en = 0;
print_edges(n);
kruskal();
free_memory();
}
return 0;
}
#include <cmath>
#include <cstring>
#include <cstdlib>
#include <algorithm>
using namespace std;
#define TRUE 1
#define FALSE 0
#define Vector(p1, p2, u, v) (u = p2->x - p1->x, v = p2->y - p1->y)
#define Cross_product_2v(u1, v1, u2, v2) (u1 * v2 - v1 * u2)
#define Cross_product_3p(p1, p2, p3) ((p2->x - p1->x) * (p3->y - p1->y) - (p2->y - p1->y) * (p3->x - p1->x))
#define Dot_product_2v(u1, v1, u2, v2) (u1 * u2 + v1 * v2)
#define Org(e) ((e)->org)
#define Dest(e) ((e)->dest)
#define Onext(e) ((e)->onext)
#define Oprev(e) ((e)->oprev)
#define Dnext(e) ((e)->dnext)
#define Dprev(e) ((e)->dprev)
#define Other_point(e, p) ((e)->org == p ? (e)->dest : (e)->org)
#define Next(e, p) ((e)->org == p ? (e)->onext : (e)->dnext)
#define Prev(e, p) ((e)->org == p ? (e)->oprev : (e)->dprev)
#define Visited(p) ((p)->f)
#define Identical_refs(e1, e2) (e1 == e2)
typedef enum {right, left} side;
typedef double Real;
typedef int boolean;
typedef struct point point;
typedef struct edge edge;
struct point
{
Real x, y;
edge *entry_pt;
};
struct edge
{
point *org, *dest;
edge *onext, *oprev, *dnext, *dprev;
};
#define MAXV 1000001
#define MAXE 3000001
struct EDGE
{
int u, v;
double w;
};
EDGE E[MAXE], mst[MAXE];
int n, en, k;
double dist[MAXV];
int p[MAXV], d[MAXV];
void UFinit(int m)
{
int i;
for(i = 0; i < m; i++)
{
p[i] = i, d[i] = 0;
}
}
int UFfind(int x)
{
int i = x, j, k;
while(p[i] != i) i = p[i];
j = x;
while(p[j] != i)
{
k = p[j];
p[j] = i;
j = k;
}
return i;
}
void UFunion(int x, int y)
{
int i = UFfind(x), j = UFfind(y);
if(d[i] > d[j]) p[j] = i;
else
{
p[i] = j;
if(d[i] == d[j]) d[j]++;
}
}
bool cmp(const EDGE &e1, const EDGE &e2)
{
return e1.w < e2.w;
}
void kruskal(void)
{
if(k > n) k = n;
int i, kk;
sort(E, E + en, cmp);
UFinit(n);
for(i = 0, kk = 0; i < en && kk < n - 1; i++)
if(UFfind(E[i].u) != UFfind(E[i].v))
{
UFunion(E[i].u, E[i].v);
mst[kk++] = E[i];
}
printf("%d\n", int(ceil(mst[n - k - 1].w) + 0.00000001));
}
point *p_array;
static edge *e_array;
static edge **free_list_e;
static int n_free_e;
void alloc_memory(int n)
{
edge *e;
int i;
p_array = (point *)calloc(n, sizeof(point));
n_free_e = 3 * n;
e_array = e = (edge *)calloc(n_free_e, sizeof(edge));
free_list_e = (edge **)calloc(n_free_e, sizeof(edge *));
for(i = 0; i < n_free_e; i++, e++)
free_list_e[i] = e;
}
void free_memory(void)
{
free(p_array);
free(e_array);
free(free_list_e);
}
edge *get_edge(void)
{
if(n_free_e == 0)
printf("Out of memory for edges\n");
return (free_list_e[--n_free_e]);
}
void free_edge(edge *e)
{
free_list_e[n_free_e++] = e;
}
void splice(edge *a, edge *b, point *v)
{
edge *next;
if(Org(a) == v)
{
next = Onext(a);
Onext(a) = b;
}
else
{
next = Dnext(a);
Dnext(a) = b;
}
if(Org(next) == v)
Oprev(next) = b;
else
Dprev(next) = b;
if(Org(b) == v)
{
Onext(b) = next;
Oprev(b) = a;
}
else
{
Dnext(b) = next;
Dprev(b) = a;
}
}
edge *make_edge(point *u, point *v)
{
edge *e;
e = get_edge();
e->onext = e->oprev = e->dnext = e->dprev = e;
e->org = u;
e->dest = v;
if(u->entry_pt == NULL)
u->entry_pt = e;
if(v->entry_pt == NULL)
v->entry_pt = e;
return e;
}
edge *join(edge *a, point *u, edge *b, point *v, side s)
{
edge *e;
e = make_edge(u, v);
if(s == left)
{
if (Org(a) == u)
splice(Oprev(a), e, u);
else
splice(Dprev(a), e, u);
splice(b, e, v);
}
else
{
splice(a, e, u);
if(Org(b) == v)
splice(Oprev(b), e, v);
else
splice(Dprev(b), e, v);
}
return e;
}
void delete_edge(edge *e)
{
point *u, *v;
u = Org(e);
v = Dest(e);
if(u->entry_pt == e)
u->entry_pt = e->onext;
if(v->entry_pt == e)
v->entry_pt = e->dnext;
if(Org(e->onext) == u)
e->onext->oprev = e->oprev;
else
e->onext->dprev = e->oprev;
if(Org(e->oprev) == u)
e->oprev->onext = e->onext;
else
e->oprev->dnext = e->onext;
if(Org(e->dnext) == v)
e->dnext->oprev = e->dprev;
else
e->dnext->dprev = e->dprev;
if(Org(e->dprev) == v)
e->dprev->onext = e->dnext;
else
e->dprev->dnext = e->dnext;
free_edge(e);
}
double px[MAXV], py[MAXV];
void read_points()
{
int i = 0;
while (1)
{
int t1, t2;
scanf("%d", &t1);
px[i] = t1;
if (t1 == -1)
break;
scanf("%d", &t2);
py[i] = t2;
i++;
}
n = i;
}
static void print_edges(int n)
{
edge *e_start, *e;
point *u, *v;
int i;
for(i = 0; i < n; i++)
{
u = &p_array[i];
e_start = e = u->entry_pt;
do
{
v = Other_point(e, u);
if(u < v)
{
E[en].u = u - p_array, E[en].v = v - p_array;
E[en++].w = hypot(u->x - v->x, u->y - v->y);
}
e = Next(e, u);
}while(!Identical_refs(e, e_start));
}
}
void merge_sort(point *p[], point *p_temp[], int l, int r)
{
int i, j, k, m;
if(r - l > 0)
{
m = (r + l) / 2;
merge_sort(p, p_temp, l, m);
merge_sort(p, p_temp, m + 1, r);
for(i = m + 1; i > l; i--)
p_temp[i - 1] = p[i - 1];
for(j = m; j < r; j++)
p_temp[r + m - j] = p[j + 1];
for(k = l; k <= r; k++)
if(p_temp[i]->x < p_temp[j]->x)
{
p[k] = p_temp[i];
i = i + 1;
}
else if(p_temp[i]->x == p_temp[j]->x && p_temp[i]->y < p_temp[j]->y)
{
p[k] = p_temp[i];
i = i + 1;
}
else
{
p[k] = p_temp[j];
j = j - 1;
}
}
}
static void lower_tangent(edge *r_cw_l, point *s, edge *l_ccw_r, point *u, edge **l_lower, point **org_l_lower, edge **r_lower, point **org_r_lower)
{
edge *l, *r;
point *o_l, *o_r, *d_l, *d_r;
boolean finished;
l = r_cw_l;
r = l_ccw_r;
o_l = s;
d_l = Other_point(l, s);
o_r = u;
d_r = Other_point(r, u);
finished = FALSE;
while(!finished)
if(Cross_product_3p(o_l, d_l, o_r) > 0.0)
{
l = Prev(l, d_l);
o_l = d_l;
d_l = Other_point(l, o_l);
}
else if(Cross_product_3p(o_r, d_r, o_l) < 0.0)
{
r = Next(r, d_r);
o_r = d_r;
d_r = Other_point(r, o_r);
}
else
finished = TRUE;
*l_lower = l;
*r_lower = r;
*org_l_lower = o_l;
*org_r_lower = o_r;
}
static void merge(edge *r_cw_l, point *s, edge *l_ccw_r, point *u, edge **l_tangent)
{
edge *base, *l_cand, *r_cand;
point *org_base, *dest_base;
Real u_l_c_o_b, v_l_c_o_b, u_l_c_d_b, v_l_c_d_b;
Real u_r_c_o_b, v_r_c_o_b, u_r_c_d_b, v_r_c_d_b;
Real c_p_l_cand, c_p_r_cand;
Real d_p_l_cand, d_p_r_cand;
boolean above_l_cand, above_r_cand, above_next, above_prev;
point *dest_l_cand, *dest_r_cand;
Real cot_l_cand, cot_r_cand;
edge *l_lower, *r_lower;
point *org_r_lower, *org_l_lower;
lower_tangent(r_cw_l, s, l_ccw_r, u, &l_lower, &org_l_lower, &r_lower, &org_r_lower);
base = join(l_lower, org_l_lower, r_lower, org_r_lower, right);
org_base = org_l_lower;
dest_base = org_r_lower;
*l_tangent = base;
do
{
l_cand = Next(base, org_base);
r_cand = Prev(base, dest_base);
dest_l_cand = Other_point(l_cand, org_base);
dest_r_cand = Other_point(r_cand, dest_base);
Vector(dest_l_cand, org_base, u_l_c_o_b, v_l_c_o_b);
Vector(dest_l_cand, dest_base, u_l_c_d_b, v_l_c_d_b);
Vector(dest_r_cand, org_base, u_r_c_o_b, v_r_c_o_b);
Vector(dest_r_cand, dest_base, u_r_c_d_b, v_r_c_d_b);
c_p_l_cand = Cross_product_2v(u_l_c_o_b, v_l_c_o_b, u_l_c_d_b, v_l_c_d_b);
c_p_r_cand = Cross_product_2v(u_r_c_o_b, v_r_c_o_b, u_r_c_d_b, v_r_c_d_b);
above_l_cand = c_p_l_cand > 0.0;
above_r_cand = c_p_r_cand > 0.0;
if(!above_l_cand && !above_r_cand)
break;
if(above_l_cand)
{
Real u_n_o_b, v_n_o_b, u_n_d_b, v_n_d_b;
Real c_p_next, d_p_next, cot_next;
edge *next;
point *dest_next;
d_p_l_cand = Dot_product_2v(u_l_c_o_b, v_l_c_o_b, u_l_c_d_b, v_l_c_d_b);
cot_l_cand = d_p_l_cand / c_p_l_cand;
do
{
next = Next(l_cand, org_base);
dest_next = Other_point(next, org_base);
Vector(dest_next, org_base, u_n_o_b, v_n_o_b);
Vector(dest_next, dest_base, u_n_d_b, v_n_d_b);
c_p_next = Cross_product_2v(u_n_o_b, v_n_o_b, u_n_d_b, v_n_d_b);
above_next = c_p_next > 0.0;
if(!above_next)
break;
d_p_next = Dot_product_2v(u_n_o_b, v_n_o_b, u_n_d_b, v_n_d_b);
cot_next = d_p_next / c_p_next;
if(cot_next > cot_l_cand)
break;
delete_edge(l_cand);
l_cand = next;
cot_l_cand = cot_next;
}while(TRUE);
}
if(above_r_cand)
{
Real u_p_o_b, v_p_o_b, u_p_d_b, v_p_d_b;
Real c_p_prev, d_p_prev, cot_prev;
edge *prev;
point *dest_prev;
d_p_r_cand = Dot_product_2v(u_r_c_o_b, v_r_c_o_b, u_r_c_d_b, v_r_c_d_b);
cot_r_cand = d_p_r_cand / c_p_r_cand;
do
{
prev = Prev(r_cand, dest_base);
dest_prev = Other_point(prev, dest_base);
Vector(dest_prev, org_base, u_p_o_b, v_p_o_b);
Vector(dest_prev, dest_base, u_p_d_b, v_p_d_b);
c_p_prev = Cross_product_2v(u_p_o_b, v_p_o_b, u_p_d_b, v_p_d_b);
above_prev = c_p_prev > 0.0;
if(!above_prev)
break;
d_p_prev = Dot_product_2v(u_p_o_b, v_p_o_b, u_p_d_b, v_p_d_b);
cot_prev = d_p_prev / c_p_prev;
if(cot_prev > cot_r_cand)
break;
delete_edge(r_cand);
r_cand = prev;
cot_r_cand = cot_prev;
}while(TRUE);
}
dest_l_cand = Other_point(l_cand, org_base);
dest_r_cand = Other_point(r_cand, dest_base);
if(!above_l_cand || (above_l_cand && above_r_cand && cot_r_cand < cot_l_cand))
{
base = join(base, org_base, r_cand, dest_r_cand, right);
dest_base = dest_r_cand;
}
else
{
base = join(l_cand, dest_l_cand, base, dest_base, right);
org_base = dest_l_cand;
}
}while(TRUE);
}
void divide(point *p_sorted[], int l, int r, edge **l_ccw, edge **r_cw)
{
int n;
int split;
edge *l_ccw_l, *r_cw_l, *l_ccw_r, *r_cw_r, *l_tangent;
edge *a, *b, *c;
Real c_p;
n = r - l + 1;
if(n == 2)
{
*l_ccw = *r_cw = make_edge(p_sorted[l], p_sorted[r]);
}
else if(n == 3)
{
a = make_edge(p_sorted[l], p_sorted[l + 1]);
b = make_edge(p_sorted[l + 1], p_sorted[r]);
splice(a, b, p_sorted[l + 1]);
c_p = Cross_product_3p(p_sorted[l], p_sorted[l + 1], p_sorted[r]);
if(c_p > 0.0)
{
c = join(a, p_sorted[l], b, p_sorted[r], right);
*l_ccw = a;
*r_cw = b;
}
else if(c_p < 0.0)
{
c = join(a, p_sorted[l], b, p_sorted[r], left);
*l_ccw = c;
*r_cw = c;
}
else
{
*l_ccw = a;
*r_cw = b;
}
}
else if(n > 3)
{
split = (l + r) / 2;
divide(p_sorted, l, split, &l_ccw_l, &r_cw_l);
divide(p_sorted, split+1, r, &l_ccw_r, &r_cw_r);
merge(r_cw_l, p_sorted[split], l_ccw_r, p_sorted[split + 1], &l_tangent);
if(Org(l_tangent) == p_sorted[l])
l_ccw_l = l_tangent;
if(Dest(l_tangent) == p_sorted[r])
r_cw_r = l_tangent;
*l_ccw = l_ccw_l;
*r_cw = r_cw_r;
}
}
int main(void)
{
int ca;
for (scanf("%d", &ca); ca--;)
{
scanf("%d", &k);
edge *l_cw, *r_ccw;
int i;
point **p_sorted, **p_temp;
read_points();
if (k > n) k = n;
if (n == 1)
{
printf("0\n");
continue;
}
alloc_memory(n);
for (i = 0; i < n; i++)
{
p_array[i].x = px[i];
p_array[i].y = py[i];
}
for(i = 0; i < n; i++)
p_array[i].entry_pt = NULL;
p_sorted = (point **)malloc((unsigned)n * sizeof(point *));
p_temp = (point **)malloc((unsigned)n * sizeof(point *));
for(i = 0; i < n; i++)
p_sorted[i] = p_array + i;
merge_sort(p_sorted, p_temp, 0, n - 1);
free((char *)p_temp);
divide(p_sorted, 0, n - 1, &l_cw, &r_ccw);
free((char *)p_sorted);
en = 0;
print_edges(n);
kruskal();
free_memory();
}
return 0;
}
[计算几何] pku1673 求三角形垂心
发表评论40次阅读2007.09.27 17:18 作者：Felicia 编辑
[阅读更多]
此问题可转化为求三角形垂心。我的做法是设垂心坐标为(x, y)，然后利用垂直关系解方程。
下面是我的代码
下载: pku1673.cpp
/**********************************************************************
Author: WHU_GCC
Created Time: 2007-9-27 16:49:56
File Name: pku1673.cpp
Description:
**********************************************************************/
#include <iostream>
using namespace std;
#define out(x) (cout << #x << ": " << x << endl)
typedef long long int64;
const int maxint = 0x7FFFFFFF;
const int64 maxint64 = 0x7FFFFFFFFFFFFFFFLL;
template <class T> void show(T a, int n) { for (int i = 0; i < n; ++i) cout << a[i] << ' '; cout << endl; }
template <class T> void show(T a, int r, int l) { for (int i = 0; i < r; ++i) show(a[i], l); cout << endl; }
double cross(double a, double b, double c, double d)
{
return a * d - b * c;
}
int main()
{
double x1, x2, x3, y1, y2, y3;
double A1, B1, C1, A2, B2, C2;
int ca;
for (scanf("%d", &ca); ca--;)
{
scanf("%lf%lf%lf%lf%lf%lf", &x1, &y1, &x2, &y2, &x3, &y3);
A1 = x2 - x3;
B1 = y2 - y3;
C1 = x1 * x2 - x1 * x3 + y1 * y2 - y1 * y3;
A2 = x1 - x3;
B2 = y1 - y3;
C2 = x1 * x2 - x2 * x3 + y1 * y2 - y2 * y3;
printf("%.4lf %.4lf\n", cross(C1, B1, C2, B2) / cross(A1, B1, A2, B2), cross(A1, C1, A2, C2) / cross(A1, B1, A2, B2));
}
return 0;
}
Author: WHU_GCC
Created Time: 2007-9-27 16:49:56
File Name: pku1673.cpp
Description:
**********************************************************************/
#include <iostream>
using namespace std;
#define out(x) (cout << #x << ": " << x << endl)
typedef long long int64;
const int maxint = 0x7FFFFFFF;
const int64 maxint64 = 0x7FFFFFFFFFFFFFFFLL;
template <class T> void show(T a, int n) { for (int i = 0; i < n; ++i) cout << a[i] << ' '; cout << endl; }
template <class T> void show(T a, int r, int l) { for (int i = 0; i < r; ++i) show(a[i], l); cout << endl; }
double cross(double a, double b, double c, double d)
{
return a * d - b * c;
}
int main()
{
double x1, x2, x3, y1, y2, y3;
double A1, B1, C1, A2, B2, C2;
int ca;
for (scanf("%d", &ca); ca--;)
{
scanf("%lf%lf%lf%lf%lf%lf", &x1, &y1, &x2, &y2, &x3, &y3);
A1 = x2 - x3;
B1 = y2 - y3;
C1 = x1 * x2 - x1 * x3 + y1 * y2 - y1 * y3;
A2 = x1 - x3;
B2 = y1 - y3;
C2 = x1 * x2 - x2 * x3 + y1 * y2 - y2 * y3;
printf("%.4lf %.4lf\n", cross(C1, B1, C2, B2) / cross(A1, B1, A2, B2), cross(A1, C1, A2, C2) / cross(A1, B1, A2, B2));
}
return 0;
}
[计算几何] pku1624 枚举
发表评论5次阅读2007.09.26 20:46 作者：Felicia 编辑
[阅读更多]
简单计算几何，我的做法是列出所有可能的切法（一共18种），求最优值。
下面是我的代码
下载: pku1624.cpp
/**********************************************************************
Author: WHU_GCC
Created Time: 2007-9-26 17:48:45
File Name: pku1624.cpp
Description:
**********************************************************************/
// 平面计算几何
// Author: Felicia @ GCC
// Last Modified: 2007.09.17
#include <iostream>
#include <cmath>
using namespace std;
typedef long long int64;
const int maxint = 0x7FFFFFFF;
const int64 maxint64 = 0x7FFFFFFFFFFFFFFFLL;
const double eps = 1E-9;
const double pi = acos(-1.0);
const double inf = 1E200;
#define out(x) (cout<<#x<<": "<<x<<endl)
template <class T> void show(T a, int n) {for (int i = 0; i < n; ++i) cout << a[i] << ' '; cout << endl;}
template <class T> void show(T a, int r, int l) {for (int i = 0; i < r; ++i) show(a[i],l); cout << endl;}
/********************
基本几何结构
点
线段
直线
多边形
圆
********************/
// 点, 同时也可以看成向量
struct point_t
{
double x, y;
point_t(double a = 0, double b = 0)
{
x = a;
y = b;
}
};
// 线段
struct lineseg_t
{
point_t s, e;
lineseg_t()
{
}
lineseg_t(point_t a, point_t b)
{
s = a;
e = b;
}
};
// 直线
// 解析方程 ax + by + c = 0 为统一表示，约定 a >= 0
struct line_t
{
double a, b, c;
line_t(double d1 = 1, double d2 = -1, double d3 = 0)
{
a = d1;
b = d2;
c = d3;
}
};
// 这里定义多边形的最大点数
const int max_polygon_size = 1000;
// 多边形, 规定逆时针为正方向
struct polygon_t
{
int n;
point_t p[max_polygon_size];
};
// 圆
struct circle_t
{
point_t center;
double r;
};
/********************
常用小函数与算符
浮点数比较
平方
点到原点距离
两点距离
两点重合
********************/
// 浮点数与0比较. x == 0 返回 0, x > 0 返回 1, x < 0 返回 -1
int dcmp(double x)
{
if (-eps < x && x < eps)
return 0;
else if (x > 0)
return 1;
else
return -1;
}
// 判断两个点是否重合
bool operator ==(const point_t &a, const point_t &b)
{
return dcmp(a.x - b.x) == 0 && dcmp(a.y - b.y) == 0;
}
bool operator !=(const point_t &a, const point_t &b)
{
return !(a == b);
}
// 向量加法
point_t operator +(const point_t &a, const point_t &b)
{
point_t ret(a.x + b.x, a.y + b.y);
return ret;
}
// 向量减法
point_t operator -(const point_t &a, const point_t &b)
{
point_t ret(a.x - b.x, a.y - b.y);
return ret;
}
// 平方
inline double sqr(double x)
{
return x * x;
}
// 点到原点距离
double dist(const point_t &p)
{
return sqrt(sqr(p.x) + sqr(p.y));
}
// 两点距离
double dist(const point_t &a, const point_t &b)
{
return sqrt(dist(a - b));
}
/********************\
* *
* 点的基本运算 *
* *
\********************/
/******************************************************************************
r=multiply(sp,ep,op),得到(sp-op)*(ep-op)的叉积
r>0:ep在矢量opsp的逆时针方向；
r=0：opspep三点共线；
r<0:ep在矢量opsp的顺时针方向
*******************************************************************************/
double cross_mul(const point_t &a, const point_t &b)
{
return a.x * b.y - a.y * b.x;
}
double dot_mul(const point_t &a, const point_t &b)
{
return a.x * b.x + a.y * b.y;
}
double cross_mul(const point_t &a, const point_t &b, const point_t &c)
{
return cross_mul(a - c, b - c);
}
/*******************************************************************************
r=dotmultiply(p1,p2,op),得到矢量(p1-op)和(p2-op)的点积，如果两个矢量都非零矢量
r<0:两矢量夹角为锐角；r=0：两矢量夹角为直角；r>0:两矢量夹角为钝角
*******************************************************************************/
double dot_mul(const point_t &a, const point_t &b, const point_t &c)
{
return dot_mul(a - c, b - c);
}
// 判断点p是否在线段l上
// 条件：p在线段l所在的直线上 && 点p在以线段l为对角线的矩形内
bool online(const lineseg_t l, const point_t p)
{
return dcmp(cross_mul(l.e, p, l.s)) == 0 && (p.x - l.s.x) * (p.x - l.e.x) <= 0 && (p.y - l.s.y) * (p.y - l.e.y) <= 0;
}
// 返回点p以点o为圆心逆时针旋转alpha(单位：弧度)后所在的位置
point_t rotate(const point_t &o, double alpha, point_t p)
{
point_t ret;
p.x -= o.x;
p.y -= o.y;
ret.x = p.x * cos(alpha) - p.y * sin(alpha) + o.x;
ret.y = p.y * cos(alpha) + p.x * sin(alpha) + o.y;
return ret;
}
// 返回向量a按逆时针方向,旋转到向量b的角度
// 角度小于pi，返回正值
// 角度大于pi，返回负值
double angle(const point_t &a, const point_t &b)
{
double ret = acos(dot_mul(a, b) / (dist(a) * dist(b)));
if (cross_mul(a, b) < 0)
return ret;
else
return -ret;
}
// 返回顶角在o点,起始边为os,终止边为oe的夹角(单位:弧度),规定逆时针为正方向
// 角度小于pi，返回正值
// 角度大于pi，返回负值
// 可以用于求线段之间的夹角
double angle(const point_t &o, const point_t &s, const point_t &e)
{
return angle(s - o, e - o);
}
// 线段及直线的基本运算
/* 判断点与线段的关系,用途很广泛
本函数是根据下面的公式写的,P是点C到线段AB所在直线的垂足
AC dot AB
r = ---------
||AB||^2
(Cx-Ax)(Bx-Ax) + (Cy-Ay)(By-Ay)
= -------------------------------
L^2
r has the following meaning:
r=0 P = A
r=1 P = B
r<0 P is on the backward extension of AB
r>1 P is on the forward extension of AB
0<r<1 P is interior to AB
*/
double relation(const point_t &p, const lineseg_t &l)
{
return dot_mul(p, l.e, l.s) / (dist(l.s, l.e) * dist(l.s, l.e));
}
// 求点p到线段l所在直线的垂足
point_t foot(const point_t &p, const lineseg_t &l)
{
double r = relation(p, l);
point_t ret;
ret.x = l.s.x + r * (l.e.x - l.s.x);
ret.y = l.s.y + r * (l.e.y - l.s.y);
return ret;
}
// 求点p到线段l的最短距离,并返回线段上距该点最近的点ret
// 注意: ret是线段l上到点p最近的点,不一定是垂足
double dist(const point_t p, const lineseg_t &l, point_t &ret)
{
double r = relation(p,l);
if (r < 0)
ret = l.s;
else if (r > 1)
ret = l.e;
else
ret = foot(p,l);
return dist(p, ret);
}
// 求点p到线段l所在直线的距离,请注意本函数与上个函数的区别
double dist(const point_t p, const lineseg_t l)
{
return abs(cross_mul(p, l.e, l.s)) / dist(l.s, l.e);
}
// 计算点到折线集的最近距离,并返回最近点
double dist(int cnt_v, point_t point_set[], const point_t &p, point_t &ret_p)
{
double ret = inf;
for (int i = 0; i < cnt_v - 1; i++)
{
lineseg_t l(point_set[i], point_set[i + 1]);
point_t tmp_p;
double tmp = dist(p, l, tmp_p);
if (tmp < ret)
{
ret = tmp;
ret_p = tmp_p;
}
}
return ret;
}
// 判断圆是否在多边形内
bool inside(const circle_t c, int cnt_v, point_t poly[])
{
point_t q;
double d = dist(cnt_v, poly, c.center, q);
return d < c.r || fabs(d - c.r) < eps;
}
// 如果线段a和b相交(包括相交在端点处)时返回true
bool intersect_e(const lineseg_t &a, const lineseg_t &b)
{
return
//排斥实验
max(a.s.x, a.e.x) >= min(b.s.x, b.e.x) &&
max(b.s.x, b.e.x) >= min(a.s.x, a.e.x) &&
max(a.s.y, a.e.y) >= min(b.s.y, b.e.y) &&
max(b.s.y, b.e.y) >= min(a.s.y, a.e.y) &&
//跨立实验
cross_mul(b.s, a.e, a.s) * cross_mul(a.e, b.e, a.s) >= 0 &&
cross_mul(a.s, b.e, b.s) * cross_mul(b.e, a.e, b.s) >= 0;
}
// 线段a和b相交 && 交点不是双方的端点时返回true
bool intersect_ne(const lineseg_t &a, const lineseg_t &b)
{
return
intersect_e(a, b) &&
!online(a, b.s) &&
!online(a, b.e) &&
!online(b, a.e) &&
!online(b, a.s);
}
// 线段l所在直线与线段a相交(包括相交在端点处)时返回true
// 方法:判断线段a是否跨立线段l
bool intersect_l(const lineseg_t &a, const lineseg_t &l)
{
return cross_mul(a.s, l.e, l.s) * cross_mul(l.e, a.e, l.s) >= 0;
}
// 根据已知两点坐标,求过这两点的直线解析方程: ax + by + c = 0 (a >= 0)
// 若两点不重合,返回true,否则返回false
bool make_line(const point_t &a, const point_t &b, line_t &ret)
{
int sign = 1;
ret.a = b.y - a.y;
if (dcmp(ret.a) == 0 && dcmp(b.x - a.x) == 0)
return false;
if (dcmp(ret.a) == 0)
{
ret.a = 0;
ret.b = 1;
ret.c = -a.y;
return true;
}
if (ret.a < 0)
{
sign = -1;
ret.a = -ret.a;
}
ret.b = sign * (a.x - b.x);
ret.c = sign * (a.y * b.x - a.x * b.y);
return true;
}
// 根据直线解析方程返回直线的斜率k,水平线返回0,竖直线返回inf
double slope(const line_t &l)
{
if (dcmp(l.a) == 0)
return 0;
if (dcmp(l.b) == 0)
return inf;
return -(l.a / l.b);
}
// 返回直线的倾斜角alpha (0 - pi)
double alpha(const line_t &l)
{
if (dcmp(l.a) == 0)
return 0;
if (dcmp(l.b) == 0)
return pi / 2;
double k = slope(l);
return k > 0 ? atan(k) : pi + atan(k);
}
// 求点p关于直线l的对称点
point_t symmetry(const line_t &l, const point_t &p)
{
point_t ret;
double sla = sqr(l.a), slb = sqr(l.b);
ret.x = ((slb - sla) * p.x - 2 * l.a * l.b * p.y - 2 * l.a * l.c) / (sla + slb);
ret.y = ((sla - slb) * p.y - 2 * l.a * l.b * p.x - 2 * l.b * l.c) / (sla + slb);
return ret;
}
// 如果两条直线 l1(a1x + b1y + c1 = 0), l2(a2x + b2y + c2 = 0)相交, 返回true, 且返回交点p
bool intersect(const line_t &l1, const line_t &l2, point_t &p)
{
double d = l1.a * l2.b - l2.a * l1.b;
if (dcmp(d) == 0)
return false;
p.x = (l2.c * l1.b - l1.c * l2.b) / d;
p.y = (l2.a * l1.c - l1.a * l2.c) / d;
return true;
}
// 如果线段l1和l2相交，返回true且交点由(inter)返回，否则返回false
bool intersect(const lineseg_t &l1, const lineseg_t &l2, point_t &ret)
{
line_t t1, t2;
if (!make_line(l1.s, l1.e, t1))
return false;
if (!make_line(l2.s, l2.e, t2))
return false;
if (intersect(t1, t2, ret))
return online(l1, ret) && online(l2, ret);
else
return false;
}
//点到直线距离
double dist(const point_t &a, const line_t &b)
{
return abs(b.a * a.x + b.b * a.y + b.c) / sqrt(sqr(b.a) + sqr(b.b));
}
const int size[] =
{
3, 3, 3, 3, 4, 4, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3
};
const int poly[18][4] =
{
{0, 1, 7},
{1, 2, 3},
{3, 4, 5},
{5, 6, 7},
{1, 2, 4, 5},
{0, 2, 3, 7},
{0, 2, 3},
{0, 2, 4},
{0, 2, 5},
{2, 4, 5},
{2, 4, 6},
{2, 7, 0},
{4, 6, 7},
{4, 6, 0},
{4, 1, 2},
{6, 0, 1},
{6, 0, 2},
{6, 3, 4},
};
point_t p[8];
int main()
{
int t1, t2, t3, t4, t5, t6, t7, t8;
int ca = 1;
while (scanf("%d%d%d%d%d%d%d%d", &t1, &t2, &t3, &t4, &t5, &t6, &t7, &t8),
t1 != 0 || t2 != 0 || t3 != 0 || t4 != 0 || t5 != 0 || t6 != 0 || t7 != 0 || t8 != 0)
{
p[0].x = t1;
p[0].y = t2;
p[2].x = t3;
p[2].y = t4;
p[4].x = t5;
p[4].y = t6;
p[6].x = t7;
p[6].y = t8;
p[1].x = (p[0].x + p[2].x) / 2;
p[1].y = (p[0].y + p[2].y) / 2;
p[3].x = (p[2].x + p[4].x) / 2;
p[3].y = (p[2].y + p[4].y) / 2;
p[5].x = (p[4].x + p[6].x) / 2;
p[5].y = (p[4].y + p[6].y) / 2;
p[7].x = (p[6].x + p[0].x) / 2;
p[7].y = (p[6].y + p[0].y) / 2;
double area = cross_mul(p[0], p[2]) + cross_mul(p[2], p[4]) + cross_mul(p[4], p[6]) + cross_mul(p[6], p[0]);
double ans = 1e200, ans_a, ans_b;
for (int i = 0; i < 18; i++)
{
double tmp = 0;
for (int j = 0; j < size[i]; j++)
tmp += cross_mul(p[poly[i][j]], p[poly[i][(j + 1) % size[i]]]);
if (fabs(area - tmp - tmp) < ans)
{
ans = fabs(area - tmp - tmp);
ans_a = tmp;
ans_b = area - tmp;
}
}
if (ans_a > ans_b)
swap(ans_a, ans_b);
printf("Cake %d: %.3lf %.3lf\n", ca++, ans_a / 2, ans_b / 2);
}
return 0;
}
Author: WHU_GCC
Created Time: 2007-9-26 17:48:45
File Name: pku1624.cpp
Description:
**********************************************************************/
// 平面计算几何
// Author: Felicia @ GCC
// Last Modified: 2007.09.17
#include <iostream>
#include <cmath>
using namespace std;
typedef long long int64;
const int maxint = 0x7FFFFFFF;
const int64 maxint64 = 0x7FFFFFFFFFFFFFFFLL;
const double eps = 1E-9;
const double pi = acos(-1.0);
const double inf = 1E200;
#define out(x) (cout<<#x<<": "<<x<<endl)
template <class T> void show(T a, int n) {for (int i = 0; i < n; ++i) cout << a[i] << ' '; cout << endl;}
template <class T> void show(T a, int r, int l) {for (int i = 0; i < r; ++i) show(a[i],l); cout << endl;}
/********************
基本几何结构
点
线段
直线
多边形
圆
********************/
// 点, 同时也可以看成向量
struct point_t
{
double x, y;
point_t(double a = 0, double b = 0)
{
x = a;
y = b;
}
};
// 线段
struct lineseg_t
{
point_t s, e;
lineseg_t()
{
}
lineseg_t(point_t a, point_t b)
{
s = a;
e = b;
}
};
// 直线
// 解析方程 ax + by + c = 0 为统一表示，约定 a >= 0
struct line_t
{
double a, b, c;
line_t(double d1 = 1, double d2 = -1, double d3 = 0)
{
a = d1;
b = d2;
c = d3;
}
};
// 这里定义多边形的最大点数
const int max_polygon_size = 1000;
// 多边形, 规定逆时针为正方向
struct polygon_t
{
int n;
point_t p[max_polygon_size];
};
// 圆
struct circle_t
{
point_t center;
double r;
};
/********************
常用小函数与算符
浮点数比较
平方
点到原点距离
两点距离
两点重合
********************/
// 浮点数与0比较. x == 0 返回 0, x > 0 返回 1, x < 0 返回 -1
int dcmp(double x)
{
if (-eps < x && x < eps)
return 0;
else if (x > 0)
return 1;
else
return -1;
}
// 判断两个点是否重合
bool operator ==(const point_t &a, const point_t &b)
{
return dcmp(a.x - b.x) == 0 && dcmp(a.y - b.y) == 0;
}
bool operator !=(const point_t &a, const point_t &b)
{
return !(a == b);
}
// 向量加法
point_t operator +(const point_t &a, const point_t &b)
{
point_t ret(a.x + b.x, a.y + b.y);
return ret;
}
// 向量减法
point_t operator -(const point_t &a, const point_t &b)
{
point_t ret(a.x - b.x, a.y - b.y);
return ret;
}
// 平方
inline double sqr(double x)
{
return x * x;
}
// 点到原点距离
double dist(const point_t &p)
{
return sqrt(sqr(p.x) + sqr(p.y));
}
// 两点距离
double dist(const point_t &a, const point_t &b)
{
return sqrt(dist(a - b));
}
/********************\
* *
* 点的基本运算 *
* *
\********************/
/******************************************************************************
r=multiply(sp,ep,op),得到(sp-op)*(ep-op)的叉积
r>0:ep在矢量opsp的逆时针方向；
r=0：opspep三点共线；
r<0:ep在矢量opsp的顺时针方向
*******************************************************************************/
double cross_mul(const point_t &a, const point_t &b)
{
return a.x * b.y - a.y * b.x;
}
double dot_mul(const point_t &a, const point_t &b)
{
return a.x * b.x + a.y * b.y;
}
double cross_mul(const point_t &a, const point_t &b, const point_t &c)
{
return cross_mul(a - c, b - c);
}
/*******************************************************************************
r=dotmultiply(p1,p2,op),得到矢量(p1-op)和(p2-op)的点积，如果两个矢量都非零矢量
r<0:两矢量夹角为锐角；r=0：两矢量夹角为直角；r>0:两矢量夹角为钝角
*******************************************************************************/
double dot_mul(const point_t &a, const point_t &b, const point_t &c)
{
return dot_mul(a - c, b - c);
}
// 判断点p是否在线段l上
// 条件：p在线段l所在的直线上 && 点p在以线段l为对角线的矩形内
bool online(const lineseg_t l, const point_t p)
{
return dcmp(cross_mul(l.e, p, l.s)) == 0 && (p.x - l.s.x) * (p.x - l.e.x) <= 0 && (p.y - l.s.y) * (p.y - l.e.y) <= 0;
}
// 返回点p以点o为圆心逆时针旋转alpha(单位：弧度)后所在的位置
point_t rotate(const point_t &o, double alpha, point_t p)
{
point_t ret;
p.x -= o.x;
p.y -= o.y;
ret.x = p.x * cos(alpha) - p.y * sin(alpha) + o.x;
ret.y = p.y * cos(alpha) + p.x * sin(alpha) + o.y;
return ret;
}
// 返回向量a按逆时针方向,旋转到向量b的角度
// 角度小于pi，返回正值
// 角度大于pi，返回负值
double angle(const point_t &a, const point_t &b)
{
double ret = acos(dot_mul(a, b) / (dist(a) * dist(b)));
if (cross_mul(a, b) < 0)
return ret;
else
return -ret;
}
// 返回顶角在o点,起始边为os,终止边为oe的夹角(单位:弧度),规定逆时针为正方向
// 角度小于pi，返回正值
// 角度大于pi，返回负值
// 可以用于求线段之间的夹角
double angle(const point_t &o, const point_t &s, const point_t &e)
{
return angle(s - o, e - o);
}
// 线段及直线的基本运算
/* 判断点与线段的关系,用途很广泛
本函数是根据下面的公式写的,P是点C到线段AB所在直线的垂足
AC dot AB
r = ---------
||AB||^2
(Cx-Ax)(Bx-Ax) + (Cy-Ay)(By-Ay)
= -------------------------------
L^2
r has the following meaning:
r=0 P = A
r=1 P = B
r<0 P is on the backward extension of AB
r>1 P is on the forward extension of AB
0<r<1 P is interior to AB
*/
double relation(const point_t &p, const lineseg_t &l)
{
return dot_mul(p, l.e, l.s) / (dist(l.s, l.e) * dist(l.s, l.e));
}
// 求点p到线段l所在直线的垂足
point_t foot(const point_t &p, const lineseg_t &l)
{
double r = relation(p, l);
point_t ret;
ret.x = l.s.x + r * (l.e.x - l.s.x);
ret.y = l.s.y + r * (l.e.y - l.s.y);
return ret;
}
// 求点p到线段l的最短距离,并返回线段上距该点最近的点ret
// 注意: ret是线段l上到点p最近的点,不一定是垂足
double dist(const point_t p, const lineseg_t &l, point_t &ret)
{
double r = relation(p,l);
if (r < 0)
ret = l.s;
else if (r > 1)
ret = l.e;
else
ret = foot(p,l);
return dist(p, ret);
}
// 求点p到线段l所在直线的距离,请注意本函数与上个函数的区别
double dist(const point_t p, const lineseg_t l)
{
return abs(cross_mul(p, l.e, l.s)) / dist(l.s, l.e);
}
// 计算点到折线集的最近距离,并返回最近点
double dist(int cnt_v, point_t point_set[], const point_t &p, point_t &ret_p)
{
double ret = inf;
for (int i = 0; i < cnt_v - 1; i++)
{
lineseg_t l(point_set[i], point_set[i + 1]);
point_t tmp_p;
double tmp = dist(p, l, tmp_p);
if (tmp < ret)
{
ret = tmp;
ret_p = tmp_p;
}
}
return ret;
}
// 判断圆是否在多边形内
bool inside(const circle_t c, int cnt_v, point_t poly[])
{
point_t q;
double d = dist(cnt_v, poly, c.center, q);
return d < c.r || fabs(d - c.r) < eps;
}
// 如果线段a和b相交(包括相交在端点处)时返回true
bool intersect_e(const lineseg_t &a, const lineseg_t &b)
{
return
//排斥实验
max(a.s.x, a.e.x) >= min(b.s.x, b.e.x) &&
max(b.s.x, b.e.x) >= min(a.s.x, a.e.x) &&
max(a.s.y, a.e.y) >= min(b.s.y, b.e.y) &&
max(b.s.y, b.e.y) >= min(a.s.y, a.e.y) &&
//跨立实验
cross_mul(b.s, a.e, a.s) * cross_mul(a.e, b.e, a.s) >= 0 &&
cross_mul(a.s, b.e, b.s) * cross_mul(b.e, a.e, b.s) >= 0;
}
// 线段a和b相交 && 交点不是双方的端点时返回true
bool intersect_ne(const lineseg_t &a, const lineseg_t &b)
{
return
intersect_e(a, b) &&
!online(a, b.s) &&
!online(a, b.e) &&
!online(b, a.e) &&
!online(b, a.s);
}
// 线段l所在直线与线段a相交(包括相交在端点处)时返回true
// 方法:判断线段a是否跨立线段l
bool intersect_l(const lineseg_t &a, const lineseg_t &l)
{
return cross_mul(a.s, l.e, l.s) * cross_mul(l.e, a.e, l.s) >= 0;
}
// 根据已知两点坐标,求过这两点的直线解析方程: ax + by + c = 0 (a >= 0)
// 若两点不重合,返回true,否则返回false
bool make_line(const point_t &a, const point_t &b, line_t &ret)
{
int sign = 1;
ret.a = b.y - a.y;
if (dcmp(ret.a) == 0 && dcmp(b.x - a.x) == 0)
return false;
if (dcmp(ret.a) == 0)
{
ret.a = 0;
ret.b = 1;
ret.c = -a.y;
return true;
}
if (ret.a < 0)
{
sign = -1;
ret.a = -ret.a;
}
ret.b = sign * (a.x - b.x);
ret.c = sign * (a.y * b.x - a.x * b.y);
return true;
}
// 根据直线解析方程返回直线的斜率k,水平线返回0,竖直线返回inf
double slope(const line_t &l)
{
if (dcmp(l.a) == 0)
return 0;
if (dcmp(l.b) == 0)
return inf;
return -(l.a / l.b);
}
// 返回直线的倾斜角alpha (0 - pi)
double alpha(const line_t &l)
{
if (dcmp(l.a) == 0)
return 0;
if (dcmp(l.b) == 0)
return pi / 2;
double k = slope(l);
return k > 0 ? atan(k) : pi + atan(k);
}
// 求点p关于直线l的对称点
point_t symmetry(const line_t &l, const point_t &p)
{
point_t ret;
double sla = sqr(l.a), slb = sqr(l.b);
ret.x = ((slb - sla) * p.x - 2 * l.a * l.b * p.y - 2 * l.a * l.c) / (sla + slb);
ret.y = ((sla - slb) * p.y - 2 * l.a * l.b * p.x - 2 * l.b * l.c) / (sla + slb);
return ret;
}
// 如果两条直线 l1(a1x + b1y + c1 = 0), l2(a2x + b2y + c2 = 0)相交, 返回true, 且返回交点p
bool intersect(const line_t &l1, const line_t &l2, point_t &p)
{
double d = l1.a * l2.b - l2.a * l1.b;
if (dcmp(d) == 0)
return false;
p.x = (l2.c * l1.b - l1.c * l2.b) / d;
p.y = (l2.a * l1.c - l1.a * l2.c) / d;
return true;
}
// 如果线段l1和l2相交，返回true且交点由(inter)返回，否则返回false
bool intersect(const lineseg_t &l1, const lineseg_t &l2, point_t &ret)
{
line_t t1, t2;
if (!make_line(l1.s, l1.e, t1))
return false;
if (!make_line(l2.s, l2.e, t2))
return false;
if (intersect(t1, t2, ret))
return online(l1, ret) && online(l2, ret);
else
return false;
}
//点到直线距离
double dist(const point_t &a, const line_t &b)
{
return abs(b.a * a.x + b.b * a.y + b.c) / sqrt(sqr(b.a) + sqr(b.b));
}
const int size[] =
{
3, 3, 3, 3, 4, 4, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3
};
const int poly[18][4] =
{
{0, 1, 7},
{1, 2, 3},
{3, 4, 5},
{5, 6, 7},
{1, 2, 4, 5},
{0, 2, 3, 7},
{0, 2, 3},
{0, 2, 4},
{0, 2, 5},
{2, 4, 5},
{2, 4, 6},
{2, 7, 0},
{4, 6, 7},
{4, 6, 0},
{4, 1, 2},
{6, 0, 1},
{6, 0, 2},
{6, 3, 4},
};
point_t p[8];
int main()
{
int t1, t2, t3, t4, t5, t6, t7, t8;
int ca = 1;
while (scanf("%d%d%d%d%d%d%d%d", &t1, &t2, &t3, &t4, &t5, &t6, &t7, &t8),
t1 != 0 || t2 != 0 || t3 != 0 || t4 != 0 || t5 != 0 || t6 != 0 || t7 != 0 || t8 != 0)
{
p[0].x = t1;
p[0].y = t2;
p[2].x = t3;
p[2].y = t4;
p[4].x = t5;
p[4].y = t6;
p[6].x = t7;
p[6].y = t8;
p[1].x = (p[0].x + p[2].x) / 2;
p[1].y = (p[0].y + p[2].y) / 2;
p[3].x = (p[2].x + p[4].x) / 2;
p[3].y = (p[2].y + p[4].y) / 2;
p[5].x = (p[4].x + p[6].x) / 2;
p[5].y = (p[4].y + p[6].y) / 2;
p[7].x = (p[6].x + p[0].x) / 2;
p[7].y = (p[6].y + p[0].y) / 2;
double area = cross_mul(p[0], p[2]) + cross_mul(p[2], p[4]) + cross_mul(p[4], p[6]) + cross_mul(p[6], p[0]);
double ans = 1e200, ans_a, ans_b;
for (int i = 0; i < 18; i++)
{
double tmp = 0;
for (int j = 0; j < size[i]; j++)
tmp += cross_mul(p[poly[i][j]], p[poly[i][(j + 1) % size[i]]]);
if (fabs(area - tmp - tmp) < ans)
{
ans = fabs(area - tmp - tmp);
ans_a = tmp;
ans_b = area - tmp;
}
}
if (ans_a > ans_b)
swap(ans_a, ans_b);
printf("Cake %d: %.3lf %.3lf\n", ca++, ans_a / 2, ans_b / 2);
}
return 0;
}
