題目：NxM矩陣中散落K點，求一個點，該點到各個的距離之和最短，兩點之間距離的定義：折線距離（即只能橫著走或豎著走）。
解法：求K個點的重心，從矩陣左上角，先橫后豎遍歷，找到第K/2個點所在的行；先豎后橫找到第K/2個點所在的列；用這個行和這個列定位的點即為所求。
原理：假定P點為所求，距離和為S，P點到任何一點的距離都由水平距離＋垂直距離構成；若水平移動，豎直方向的距離之和保持不變；同理垂直方向。因此，該問題可以歸約為一個一維問題：數軸上散落N個點，求一個點到各個點的距離之和最小，下面進行簡單論證：
a b c d e f g七個點，d點即位所求；如果是f點，那么該問題可以歸約為b c d e f五個點的問題（在比較d和f誰更優時，a和g同時存在或不存在對問題的影響是一致的），此時，顯然f不如d。
推廣：K個點，每個點都有一個系數q，兩點之間距離定義為：折線距離 x 系數k；解法依舊：先計算K個點的總系數Q，遍歷方法依舊，只不過這次是找到一個點，截至這個點累加起來的系數之和剛好達到Q/2。
領悟：實則就是求重心。
代碼：和窮舉法做了對比，結果一致。
#include <iostream>
#include <stdlib.h>
#include <string>
#include <vector>
#include <map>
#include <math.h>
#include <set>
using namespace std;
class MatchPointFinder {
public:
double GetDistanceBetweenTwoPoints(int x1, int y1, int x2, int y2, double alpha) const {
return (abs(x1 - x2) + abs(y1 - y2)) * alpha;
//return (x1 - x2) * (x1 - x2) + (y1 - y2) * (y1 - y2);
//two methods will *not* be equivalent under the following way of calculating distance, for the case in "main"
//return sqrt((x1 - x2) * (x1 - x2) + (y1 - y2) * (y1 - y2));
}
double CalcSumDistance(int x, int y) {
double dist = 0.0f;
for (typeof(m_node_pos.begin()) it = m_node_pos.begin(), it4End = m_node_pos.end(); it4End != it; ++it) {
dist += GetDistanceBetweenTwoPoints(it->first.first, it->first.second, x, y, it->second);
}
return dist;
}
MatchPointFinder(): m_width(0), m_height(0), m_match_point_x(0), m_match_point_y(0), m_total_weight(0), m_nodes(NULL) {
}
void ValidateEquivalenceOfTwoMethods(int trial_num = 1000) {
for (int i = 0; i < trial_num; ++i) {
int width = 5 + rand() % 100;
int height = 5 + rand() % 100;
Reset();
RandomInit(width, height);
if (!CompareTwoMethod()) {
Try();
break;
} else {
printf("two methods are equivalent under %d x %d array scale\n", height, width);
}
}
}
void Try() {
FindMatchPointByEnumerating(m_match_point_x, m_match_point_y);
printf("total weight is %f\n", m_total_weight);
PrintMap(false);
PrintMap();
int old_mp_x = m_match_point_x;
int old_mp_y = m_match_point_y;
FindMatchPointByMedianMethod(m_match_point_x, m_match_point_y);
PrintMap(false);
if (old_mp_x == m_match_point_x && old_mp_y == m_match_point_y) {
printf("two methods are equivalent under %d x %d array scale\n", m_height, m_width);
} else {
printf("two methods are *not* equivalent under %d x %d array scale\n", m_height, m_width);
}
}
template <typename T>
T abs(const T& t) const {
return t > 0 ? t : -1 * t;
}
bool CompareTwoMethod() {
FindMatchPointByEnumerating(m_match_point_x, m_match_point_y);
int old_mp_x = m_match_point_x;
int old_mp_y = m_match_point_y;
FindMatchPointByMedianMethod(m_match_point_x, m_match_point_y);
if (old_mp_x == m_match_point_x && old_mp_y == m_match_point_y) {
return true;
} else if (abs(m_node_dist[std::make_pair(old_mp_x, old_mp_y)] - m_node_dist[std::make_pair(m_match_point_x, m_match_point_y)]) < 0.000001) {
return true;
}
return false;
}
~MatchPointFinder() {
Reset();
}
void RandomInit(int width = -1, int height = -1, bool init_random = true) {
if (-1 != width) {
m_width = width;
}
if (-1 != height) {
m_height = height;
}
srand(time(NULL));
m_nodes = new double*[m_height];
memset(m_nodes, 0, sizeof(double*) * m_height);
for (int i = 0; i < m_height; ++i) {
if (!init_random) {
m_nodes[i] = new double[m_width];
memset(m_nodes[i], 0, sizeof(int) * m_width);
} else {
for (int j = 0; j < m_width; ++j) {
if (0 == j) {
m_nodes[i] = new double[m_width];
memset(m_nodes[i], 0, sizeof(int) * m_width);
}
m_nodes[i][j] = ((rand() % 10) >= 5 ? (double)(1 + rand() % 3) / 5 : 0);
if (m_nodes[i][j]) {
m_node_pos[std::make_pair(i, j)] = m_nodes[i][j];
m_total_weight += m_nodes[i][j];
}
}
}
}
}
void AddNode(int x, int y) {
m_nodes[x][y] = 1;
m_node_pos[std::make_pair(x, y)] = 1;
m_total_weight += 1;
}
void FindMatchPointByMedianMethod(int& mp_x, int& mp_y) {
double collected_weight = 0;
mp_x = -1;
for (int i = 0; i < m_height; ++i) {
for (int j = 0; j < m_width; ++j) {
if (m_nodes[i][j] <= 0) {
continue;
}
collected_weight += m_nodes[i][j];
if (collected_weight >= m_total_weight / 2) {
mp_x = i;
break;
}
}
if (-1 != mp_x) {
break;
}
}
collected_weight = 0;
mp_y = -1;
for (int j = 0; j < m_width; ++j) {
for (int i = 0; i < m_height; ++i) {
if (m_nodes[i][j] <= 0) {
continue;
}
collected_weight += m_nodes[i][j];
if (collected_weight >= m_total_weight / 2) {
mp_y = j;
break;
}
}
if (-1 != mp_y) {
break;
}
}
}
void FindMatchPointByEnumerating(int& mp_x, int& mp_y) {
double shortest_sum_dist = std::numeric_limits<double>::max();
for (int i = 0; i < m_height; ++i) {
for (int j = 0; j < m_width; ++j) {
double dist = CalcSumDistance(i, j);
m_node_dist[std::make_pair(i, j)] = dist;
if (dist < shortest_sum_dist) {
shortest_sum_dist = dist;
mp_x = i;
mp_y = j;
}
}
}
}
void PrintMap(bool print_dist = true) {
printf("%4s ", " ");
for (int i = 0; i < m_width; ++i) {
printf(" (%2d) ", i);
}
printf("\n");
for (int i = 0; i < m_height; ++i) {
printf("(%2d) ", i);
for (int j = 0; j < m_width; ++j) {
if (i == m_match_point_x && j == m_match_point_y) {
if (print_dist) {
printf("*%.1f[%2.1f] ", m_nodes[i][j], m_node_dist[std::make_pair(i, j)]);
} else {
printf(" *%.1f ", m_nodes[i][j]);
}
} else {
if (print_dist) {
printf("%.1f[%2.1f] ", m_nodes[i][j], m_node_dist[std::make_pair(i, j)]);
} else {
printf(" %.1f ", m_nodes[i][j]);
}
}
}
printf("\n");
}
printf("\n");
}
void Reset() {
if (m_nodes) {
for (int i = 0; i < m_height; ++i) {
delete [] m_nodes[i];
}
delete [] m_nodes;
}
m_nodes = NULL;
m_width = 0;
m_height = 0;
m_match_point_x = 0;
m_match_point_y = 0;
m_total_weight = 0;
m_node_pos.clear();
m_node_dist.clear();
}
private:
int m_width;
int m_height;
int m_match_point_x;
int m_match_point_y;
double m_total_weight;
double** m_nodes;
std::map<std::pair<int, int>, double> m_node_pos;
std::map<std::pair<int, int>, double> m_node_dist;
};
int main(int argc, char* argv[]) {
int width = 5;
int height = 5;
if (argc > 1) {
height = atoi(argv[1]);
}
if (argc > 2) {
width = atoi(argv[2]);
}
MatchPointFinder mpf;
//mpf.RandomInit(width, height);
//mpf.Try();
//return 0;
mpf.ValidateEquivalenceOfTwoMethods(1000);
return 0;
/*
( 0) ( 1) ( 2) ( 3) ( 4)
( 0) 1 1 0 0 1
( 1) 1 0 0 1 0
( 2) 1 1 0 1 0
( 3) 0 0 * 1 1 0
( 4) 1 1 1 1 1
*/
mpf.Reset();
mpf.RandomInit(5, 5, false);
mpf.AddNode(0, 0);
mpf.AddNode(0, 1);
mpf.AddNode(0, 4);
mpf.AddNode(1, 0);
mpf.AddNode(1, 3);
mpf.AddNode(2, 0);
mpf.AddNode(2, 1);
mpf.AddNode(2, 3);
mpf.AddNode(3, 2);
mpf.AddNode(3, 3);
mpf.AddNode(4, 0);
mpf.AddNode(4, 1);
mpf.AddNode(4, 2);
mpf.AddNode(4, 3);
mpf.AddNode(4, 4);
mpf.Try();
return 0;
}