题目：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;
}