Kickstart Round B 2017——Problem B. Center（及一點延伸）
阿新 • • 發佈：2019-01-06
題目介紹
相關知識
先介紹一些概念性的東西（機器學習的東西不是很懂，如有錯誤，歡迎指出，也請見諒）：
- 機器學習中常用這些距離公式估算不同樣本之間的相似性，可以依此找到聚類質心點，進行聚類。之前在師兄的論文中也看到，使用閔氏距離估算某一個特徵在正負樣本上的區分度，距離值越大，表明區分度越好，特徵越重要。
- 但是在實際問題中，很多樣本不是單純地距離關係，也就不能只依據距離大小進行聚類，比如本題中每個點還有權重的影響。於是，才理解本題設計的初衷：要求使用切比雪夫距離公式寫一個聚類演算法，找到聚類的質心
切比雪夫距離具有這樣一種性質：
- 對於平面中原座標系中兩點間的 Chebyshev 距離，是將座標軸順(逆)時針旋轉45度並將所有點的座標值放大sqrt(2)倍所得到的新座標系中的Manhattan距離的二分之一。通過畫圖可以發現，在切比雪夫座標系下的點A（x,y），假設x>y，在曼哈頓座標系下對應座標為（(x+y)/sqrt(2)， (x-y)/sqrt(2)）,則在原座標系下點O(0,0)到A點的切比雪夫距離等於x，等於旋轉後的座標擴大sqrt(2)倍後的曼哈頓距離的二分之一，即座標為（(x+y)/2， (x-y)/2）到O點的曼哈頓距離。（剛手畫了一個圖，太醜了，就不放上了，稍後補上）
曼哈頓距離和切比雪夫距離的兩個基礎題目
- 相關題目：HDU 4311為典型的曼哈頓距離題目，HDU 4312則求切比雪夫距離，建議可以先拿這兩個比較直觀的題目練練手。此類題目有兩種解法。因為兩個題目的解題方法大同小異，只是4312需要將座標轉化一下，然後便可以完全按照4311的做法進行求解，這裡以4311為例，題意為在n個點鐘中選擇某一個點（X,Y），使得到其他所有點的曼哈頓距離和最小，即sigma(|X-Xi|+|Y-Yi|)最小：
http://acm.hdu.edu.cn/showproblem.php?pid=4311
http://acm.hdu.edu.cn/showproblem.php?pid=4312 - 解法一：我們可以發現在曼哈頓距離中，X軸方向的距離和Y軸方向的距離無關，且相互沒有影響，故可以單獨計算。即sigma(|X-Xi|+|Y-Yi|) = sigma(|X-Xi|) + sigma(|Y-Yi|)，據此我們可以先按照x排序，計算以每個x值為中心時X軸方向上的距離和，然後按照y排序，計算以每個y值為中心時Y軸方向上的距離和。
– 複雜度分析：排序的時間複雜度為O(nlogn)。列舉n個點，暴力計算以每個點為中心時的距離和為O(n^2)，在這一步我們可以使用字首和的方式將時間優化到O(n)，故總時間複雜度為O(nlogn)。
– 字首和公式：稍微推算一下就會發現，X軸方向上排序後以第j個點作為中心點的距離和sumX[j] = sumX[j-1] + (2*j - n) * (X[j] - X[j-1])，Y軸方向上同理。
– HDU 4311程式碼：
- 相關題目：HDU 4311為典型的曼哈頓距離題目，HDU 4312則求切比雪夫距離，建議可以先拿這兩個比較直觀的題目練練手。此類題目有兩種解法。因為兩個題目的解題方法大同小異，只是4312需要將座標轉化一下，然後便可以完全按照4311的做法進行求解，這裡以4311為例，題意為在n個點鐘中選擇某一個點（X,Y），使得到其他所有點的曼哈頓距離和最小，即sigma(|X-Xi|+|Y-Yi|)最小：
#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const int N = 100010;
const ll INF = LLONG_MAX;
struct Node {
int x;
int y;
int id;
}node[N];
bool CmpX(Node a, Node b) {
return a.x < b.x;
}
bool CmpY(Node a, Node b) {
return a.y < b.y;
}
ll sumX[N];
ll sumY[N];
int main() {
int t;
cin >> t;
while(t--) {
int n;
cin >> n;
for(int i = 0; i< n; i++) {
cin>> node[i].x >> node[i].y;
node[i].id = i;
}
sort(node, node + n, CmpX);
memset(sumX, 0, sizeof(sumX));
for(int i = 1; i< n; i++) {
sumX[node[0].id] += (node[i].x-node[0].x);
}
for(int i = 1; i < n; i++) {
sumX[node[i].id] = sumX[node[i-1].id] + (ll)(2*i-n) * (node[i].x - node[i-1].x);
}
sort(node, node + n, CmpY);
memset(sumY, 0, sizeof(sumY));
for(int i = 1; i< n; i++) {
sumY[node[0].id] += (node[i].y-node[0].y);
}
for(int i = 1; i < n; i++) {
sumY[node[i].id] = sumY[node[i-1].id] + (ll)(2*i-n) * (node[i].y - node[i-1].y);
}
ll ans = INF;
for(int i = 0; i< n; ++i) {
// cout<< sumX[i] << " "<<sumY[i]<<endl;
ans = min(ans, sumX[i]+sumY[i]);
}
cout<<ans<<endl;
}
return 0;
}– HDU 4312程式碼：
#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const int N = 100010;
const ll INF = LLONG_MAX;
struct Node {
ll x;
ll y;
ll id;
Node(int a=0, int b=0, int id = 0):x(a), y(b), id(id){}
}node[N];
bool CmpX(Node a, Node b) {
return a.x < b.x;
}
bool CmpY(Node a, Node b) {
return a.y < b.y;
}
ll sumX[N];
ll sumY[N];
int main() {
int t;
cin >> t;
while(t--) {
int n;
cin >> n;
int x, y;
for(int i = 0; i< n; i++) {
//cin>> node[i].x >> node[i].y;
//node[i].id = i;
cin>>x>>y;
node[i] = Node(x-y, x+y, i);
}
sort(node, node + n, CmpX);
memset(sumX, 0, sizeof(sumX));
for(int i = 1; i< n; i++) {
sumX[node[0].id] += (node[i].x-node[0].x);
}
for(int i = 1; i < n; i++) {
sumX[node[i].id] = sumX[node[i-1].id] + (ll)(2*i-n) * (node[i].x - node[i-1].x);
}
sort(node, node + n, CmpY);
memset(sumY, 0, sizeof(sumY));
for(int i = 1; i< n; i++) {
sumY[node[0].id] += (node[i].y-node[0].y);
}
for(int i = 1; i < n; i++) {
sumY[node[i].id] = sumY[node[i-1].id] + (ll)(2*i-n) * (node[i].y - node[i-1].y);
}
ll ans = INF;
for(int i = 0; i< n; ++i) {
// cout<< sumX[i] << " "<<sumY[i]<<endl;
ans = min(ans, sumX[i]+sumY[i]);
}
cout<<ans/2<<endl;
}
return 0;
}–
- 解法二：二分法逐漸逼近最優解（未完待續……）
本題解題思路
本題與HDU 4312的區別在於：（1）選中的中心點可以是平面中的任一點，意味著可能不是N個點中的某一個點；（2）每個點有一個權重Wi。
**本題兩種解法：
（1）解法一： 字首和，做法同HDU 4311，但是有點小區別：我們選擇最小的sumX[i]值，然後選擇最小的sumY[j]值，最後選擇（min（sumX[i]）+min（sumY[j]））的最小值即為最終結果。而HDU 4311中的（min（sumX））最終的中心點為（X,Y）
（2）解法二：三分搜尋，逼近最優解- 回顧一下本題的題意，max（|X-Xi|， |Y-Yi|）*Wi
- 二分適用於單調函式，或者導數為單調函式的函式，求導後，使用二分。三分適用於單峰函式，可以求得最值。
http://blog.csdn.net/fjsd155/article/details/6918873 二分法作為分治中最常見的方法，適用於單調函式，逼近求解某點的值。但當函式是凸性函式時，二分法就無法適用，這時三分法就可以“大顯身手”～～ 在求解函式f的極大極小問題時，通常使用三分。但也可以轉化為求解函式f的導數g = 0的問題，此時如果g在解所在的區間內是單調的，則可以使用二分求解。
在Kickstart Round A 2017的Problem C Jane’s Flower Shop，使用二分法求解。我沒有證明出來是單調的，哪位大神如果知道，還望不吝賜教。師兄提供了另一種解法，牛頓法。P.S. 想說谷歌對數學的要求還真是挺高的，很多人即使題做出來了，也是知其然不知其所以然。當然，我連知其然還沒做到，加油咯。
程式碼
解法一：字首和。
#define _CRT_SECURE_NO_WARNINGS
//#include <bits/stdc++.h>
#include<iostream>
#include<algorithm>
#include<cstring>
#include<ctime>
using namespace std;
const int N = 50010;
double sumx[N], sumy[N];
struct Node{
double x;
double y;
double w;
int id;
Node(double x = 0.0, double y = 0.0, double w = 0.0, int id = 0) :x(x), y(y), w(w), id(id){
}
}node[N];
bool Cmpx(Node a, Node b) {
return a.x < b.x;
}
bool Cmpy(Node a, Node b) {
return a.y < b.y;
}
int main() {
freopen("B-small-practice (2).in", "r", stdin);
freopen("B-small-practice-my.out", "w", stdout);
int t;
cin >> t;
for (int k = 1; k <= t; k++) {
int n;
cin >> n;
double x, y, w;
double rightw = 0.0;
for (int j = 0; j<n; j++) {
cin >> x >> y >> w;
node[j] = Node((x - y)*0.5, (x + y) *0.5, w, j);
//cout<<node[j].x <<" " <<node[j].y<<endl;
rightw += w;
}
double righty = rightw;
sort(node, node + n, Cmpx);
for (int i = 0; i < n; i++) {
sumx[i] = 0;
sumy[i] = 0;
}
for (int i = 1; i< n; i++) {
sumx[node[0].id] = sumx[node[0].id]+ (node[i].x - node[0].x) * node[i].w;
// cout<<sumx[node[0].id]<<endl;
}
//cout<<sumx[node[0].id]<<endl;
double leftw = 0.0;
double temp = DBL_MAX;
for (int i = 1; i < n; i++) {
leftw += node[i - 1].w;
rightw -= node[i - 1].w;
sumx[node[i].id] = sumx[node[i - 1].id] + leftw * (node[i].x - node[i - 1].x) - rightw * (node[i].x - node[i - 1].x);
//cout<<sumx[node[i].id]<<" ";
temp = min(temp, sumx[node[i].id]);
}
//cout<<endl;
sort(node, node + n, Cmpy);
//memset(sumy, 0.0000, sizeof(sumy));
for (int i = 1; i< n; i++) {
sumy[node[0].id] = sumy[node[0].id]+(node[i].y - node[0].y) * node[i].w;
}
double lefty = 0.0;
double tempy = DBL_MAX;
for (int i = 1; i < n; i++) {
lefty += node[i - 1].w;
righty -= node[i - 1].w;
sumy[node[i].id] = sumy[node[i - 1].id] + lefty * (node[i].y - node[i - 1].y) - righty * (node[i].y - node[i - 1].y);
tempy = min(tempy, sumy[node[i].id]);
}
/*cout << temp << " && " << tempy << endl;
double ans = sumx[0]+sumy[0];
for (int i = 0; i < n; i++) {
ans = min(ans, sumx[i] + sumy[i]);
}*/
double ans = temp + tempy;
printf("Case #%d: %.7lf\n", k, ans);
//cout << "Case #" << k << ": " << ans << endl;
}
//system("pause");
return 0;
}解法二：將切比雪夫距離轉化成曼哈頓距離，然後分別對X，Y進行三分，X方向距離和的最小值加上Y方向距離和的最小值，即為最終解。
#define _CRT_SECURE_NO_WARNINGS
#include <stdio.h>
#include <stdlib.h>
#include <algorithm>
#include <string.h>
#include <iostream>
#include <map>
#include <vector>
#include <queue>
#include <set>
#include <string>
#include <math.h>
using namespace std;
const int N = 10010;
const int LG = 100;
int n;
struct Node {
double x;
double y;
double w;
Node(double a = 0.0, double b = 0.0, double c = 0.0) :x(a), y(b), w(c){
}
}node[N];
double Calc(double p, int flag) {
double ans = 0.0;
if (flag == 1){
for (int i = 0; i< n; i++) {
ans += abs((node[i].x - p))*node[i].w;
}
}
else {
for (int i = 0; i< n; i++) {
ans += abs((node[i].y - p))*node[i].w;
}
}
return ans;
}
double Search(double l, double r, int flag) {
double ans = min(Calc(l, flag), Calc(r, flag));
for (int i = 0; i< LG; i++) {
double mid1 = (l + l + r) / 3;
double mid2 = (l + r + r) / 3;
double res1 = Calc(mid1, flag);
double res2 = Calc(mid2, flag);
if (res1 > res2) {
l = mid1;
}
else{
r = mid2;
}
ans = min(ans, min(res1, res2));
}
return ans;
}
int main() {
freopen("B-large-practice.in", "r", stdin);
freopen("B-large-practice-my.out", "w", stdout);
int t;
cin >> t;
for (int i = 1; i <= t; i++) {
cin >> n;
double x, y, w;
double minx = DBL_MAX;
double maxx = 0.0;
double miny = DBL_MAX;
double maxy = 0.0;
for (int j = 0; j < n; j++) {
cin >> x >> y >> w;
node[j] = Node((x + y) / 2, (x - y) / 2, w);
minx = min(minx, node[j].x);
miny = min(miny, node[j].y);
maxx = max(maxx, node[j].x);
maxy = max(maxy, node[j].y);
}
double result = Search(minx, maxx, 1) + Search(miny, maxy, 0);
//cout<<"Case #"<<i<<": " <<result<<endl;
printf("Case #%d: %.6lf\n", i, result);
}
return 0;
}