Problem background
I have a 2D array
Map[Height][Width]that stores aboolto represent each cell.There exists two non-overlapping regions
RegionAandRegionB.Each region is a list of unique adjacent integer co-ordinates.
The number of co-ordinates in
RegionAandRegionBaremandnrespectivelyA co-ordinate in
RegionAwill never be horizontally or vertically adjacent to a coordinate inRegionB(i.e. the regions do not touch)The bounding box of
RegionAcould overlap with the bounding box ofRegionBRegionAandRegionBmay have holesRegionAmay or may not surroundRegionB(and vise versa)if a co-ordinate (X, Y) is part of a region then
Map[Y][X] == 1, otherwise its zero.
Problem
I'm looking for an algorithm to determine the two co-ordinates A and B that have the minimum Euclidean distance, where A is a member of RegionA and B is a member of RegionB.
The brute force method has a time complexity of O(mn).
Is there a more time-efficient algorithm (preferably better than O(mn)) for solving this problem?
The name + link, code, or description of one instance of such an algorithm would be greatly appreciated and will be accepted as an answer.
C++ code to create two random regions A and B
The region generation code is immaterial to the problem being asked. and I am only interested in how close the regions are to each other. In the project I am working on I used cellular automata to create the regions. I can include the code for this if necessary. The Region construction procedure exists only to create relevant examples.
#include <cstdlib>
#include <cstring>
#include <vector>
#include <cstdio>
#include <map>
#include <ctime>
#include <cmath>
constexpr int Width = 150;
constexpr int Height = 37;
bool Map[Height][Width];
struct coord {
int X, Y;
coord(int X, int Y): X(X), Y(Y) {}
bool operator==(const coord& Other) {
return Other.X == X && Other.Y == Y;
}
bool IsAdjacentTo(const coord& Other){
int ManhattanDistance = abs(X - Other.X) + abs(Y - Other.Y);
return ManhattanDistance == 1;
}
bool IsOnMapEdge(){
return X == 0 || Y == 0 || X == Width - 1 || Y == Height - 1;
}
};
using region = std::vector<coord>;
void ClearMap(){
std::memset(Map, 0, sizeof Map);
}
bool& GetMapBool(coord Coord) {
return Map[Coord.Y][Coord.X];
}
int RandInt(int LowerBound, int UpperBound){
return std::rand() % (UpperBound - LowerBound + 1) + LowerBound;
}
region CreateRegion(){
std::puts("Creating Region");
ClearMap();
region Region;
int RegionSize = RandInt(50, Width * Height / 2);
Region.reserve(RegionSize);
while (true){
Region.emplace_back(
RandInt(1, Width - 1),
RandInt(1, Height - 1)
);
if (!Region[0].IsOnMapEdge())
break;
else
Region.pop_back();
}
GetMapBool(Region[0]) = 1;
while (Region.size() != RegionSize){
coord Member = Region[RandInt(0, Region.size() - 1)];
coord NeighbourToMakeAMember {
RandInt(-1, 1) + Member.X,
RandInt(-1, 1) + Member.Y
};
if (!Member.IsOnMapEdge() && !GetMapBool(NeighbourToMakeAMember) && Member.IsAdjacentTo(NeighbourToMakeAMember)){
GetMapBool(NeighbourToMakeAMember) = 1;
Region.push_back(NeighbourToMakeAMember);
};
}
std::puts("Created Region");
return Region;
}
std::pair<region, region> CreateTwoRegions(){
bool RegionsOverlap;
std::pair<region, region> Regions;
do {
Regions.first = CreateRegion();
Regions.second = CreateRegion();
ClearMap();
for (coord Member : Regions.first){
GetMapBool(Member) = 1;
}
RegionsOverlap = 0;
for (coord Member : Regions.second){
if (GetMapBool(Member)){
ClearMap();
std::puts("Regions Overlap/Touch");
RegionsOverlap = 1;
break;
} else {
GetMapBool(Member) = 1;
}
}
} while (RegionsOverlap);
return Regions;
}
void DisplayMap(){
for (int Y = 0; Y < Height; ++Y){
for (int X = 0; X < Width; ++X)
std::printf("%c", (Map[Y][X] ? '1' : '-'));
std::puts("");
}
}
int main(){
int Seed = time(NULL);
std::srand(Seed);
ClearMap();
auto[RegionA, RegionB] = CreateTwoRegions();
DisplayMap();
}
Problem illustration
What are the co-ordinates of the points A and B that form the minimum Euclidean distance between the two regions?
Problem background
I have a 2D array
Map[Height][Width]that stores aboolto represent each cell.There exists two non-overlapping regions
RegionAandRegionB.Each region is a list of unique adjacent integer co-ordinates.
The number of co-ordinates in
RegionAandRegionBaremandnrespectivelyA co-ordinate in
RegionAwill never be horizontally or vertically adjacent to a coordinate inRegionB(i.e. the regions do not touch)The bounding box of
RegionAcould overlap with the bounding box ofRegionBRegionAandRegionBmay have holesRegionAmay or may not surroundRegionB(and vise versa)if a co-ordinate (X, Y) is part of a region then
Map[Y][X] == 1, otherwise its zero.
Problem
I'm looking for an algorithm to determine the two co-ordinates A and B that have the minimum Euclidean distance, where A is a member of RegionA and B is a member of RegionB.
The brute force method has a time complexity of O(mn).
Is there a more time-efficient algorithm (preferably better than O(mn)) for solving this problem?
The name + link, code, or description of one instance of such an algorithm would be greatly appreciated and will be accepted as an answer.
C++ code to create two random regions A and B
The region generation code is immaterial to the problem being asked. and I am only interested in how close the regions are to each other. In the project I am working on I used cellular automata to create the regions. I can include the code for this if necessary. The Region construction procedure exists only to create relevant examples.
#include <cstdlib>
#include <cstring>
#include <vector>
#include <cstdio>
#include <map>
#include <ctime>
#include <cmath>
constexpr int Width = 150;
constexpr int Height = 37;
bool Map[Height][Width];
struct coord {
int X, Y;
coord(int X, int Y): X(X), Y(Y) {}
bool operator==(const coord& Other) {
return Other.X == X && Other.Y == Y;
}
bool IsAdjacentTo(const coord& Other){
int ManhattanDistance = abs(X - Other.X) + abs(Y - Other.Y);
return ManhattanDistance == 1;
}
bool IsOnMapEdge(){
return X == 0 || Y == 0 || X == Width - 1 || Y == Height - 1;
}
};
using region = std::vector<coord>;
void ClearMap(){
std::memset(Map, 0, sizeof Map);
}
bool& GetMapBool(coord Coord) {
return Map[Coord.Y][Coord.X];
}
int RandInt(int LowerBound, int UpperBound){
return std::rand() % (UpperBound - LowerBound + 1) + LowerBound;
}
region CreateRegion(){
std::puts("Creating Region");
ClearMap();
region Region;
int RegionSize = RandInt(50, Width * Height / 2);
Region.reserve(RegionSize);
while (true){
Region.emplace_back(
RandInt(1, Width - 1),
RandInt(1, Height - 1)
);
if (!Region[0].IsOnMapEdge())
break;
else
Region.pop_back();
}
GetMapBool(Region[0]) = 1;
while (Region.size() != RegionSize){
coord Member = Region[RandInt(0, Region.size() - 1)];
coord NeighbourToMakeAMember {
RandInt(-1, 1) + Member.X,
RandInt(-1, 1) + Member.Y
};
if (!Member.IsOnMapEdge() && !GetMapBool(NeighbourToMakeAMember) && Member.IsAdjacentTo(NeighbourToMakeAMember)){
GetMapBool(NeighbourToMakeAMember) = 1;
Region.push_back(NeighbourToMakeAMember);
};
}
std::puts("Created Region");
return Region;
}
std::pair<region, region> CreateTwoRegions(){
bool RegionsOverlap;
std::pair<region, region> Regions;
do {
Regions.first = CreateRegion();
Regions.second = CreateRegion();
ClearMap();
for (coord Member : Regions.first){
GetMapBool(Member) = 1;
}
RegionsOverlap = 0;
for (coord Member : Regions.second){
if (GetMapBool(Member)){
ClearMap();
std::puts("Regions Overlap/Touch");
RegionsOverlap = 1;
break;
} else {
GetMapBool(Member) = 1;
}
}
} while (RegionsOverlap);
return Regions;
}
void DisplayMap(){
for (int Y = 0; Y < Height; ++Y){
for (int X = 0; X < Width; ++X)
std::printf("%c", (Map[Y][X] ? '1' : '-'));
std::puts("");
}
}
int main(){
int Seed = time(NULL);
std::srand(Seed);
ClearMap();
auto[RegionA, RegionB] = CreateTwoRegions();
DisplayMap();
}
Problem illustration
What are the co-ordinates of the points A and B that form the minimum Euclidean distance between the two regions?
Share Improve this question edited Mar 30 at 10:20 Mark Rotteveel 110k229 gold badges156 silver badges224 bronze badges asked Mar 18 at 22:34 timmy geetimmy gee 6131 silver badge12 bronze badges 22 | Show 17 more comments4 Answers
Reset to default 11Notes on finding the nearest points between two regions
As others have mentioned, you can improve on finding the closest points between the two regions of n and m points by first filtering each region's points to only those which border unassigned cells and then doing an all-pairs search between the boundaries. For a roughly rectangular region this would take O(m + n + 16 √m √n) time: two linear filters and then the all-pairs search. I've implemented an example of this search below.
But this will not save you any time if your regions are shaped like this:
########
# #
# #### #
# # # #
# # ## #
# # #
# ######
If you are using larger datasets and can tolerate preprocessing you can get further acceleration using a 2D k-d tree.
Construction takes O(n log n) time after which you can search the tree in O(log n) time for m queries. So the total time to find the closest points is then O(m log n + n log n) or O(n log m + m log m).
You can do some benchmarking to determine if it is better to build the k-d tree for the larger or the smaller of the regions. Remember that this answer may change depending on how many times you're able to reuse the same k-d tree for queries against different regions.
For small datasets you'll probably find that the all-pairs comparison is faster because it has much greater cache locality.
Note that implementing a k-d tree is tricky, so you'll probably want to use an existing library or definitely use test cases if you're rolling your own.
Another algorithm
If we happen to know that the regions are within a distance d of each other we can get an even faster algorithm.
For one of the regions we can build a hash table where a reference to each member cell is stored at a hash location (c.y // d) * height + (c.x // d). To find the nearest neighbor to a query point we then do an all-pairs check against the points referenced in the associated "meta-cell" and as well as the meta-cell's neighbors.
Notes on region building
On the off-chance you're trying to grow the regions quickly, but their boundaries can overlap, this should do the trick.
Note that your question can be improved by stating why you are trying to calculate a thing. If you want to know the closest-cell so you have a primitive to build regions more quickly, then reconsidering how region building works might be a better use of your time, and it's what I address here.
If region generation is immaterial and you're only interested in how close they are to each other, then focus on defining, conceptually, what a region is and be very clear that your region construction procedure exists only to create relevant examples.
The fundamental trick we'll use to accelerate region construction is to make a space-time tradeoff. Rather than storing the map as a set of boolean values, we'll store it as a set of 8-bit integers. We reserve a few of these integers for special purposes:
- 0 - An unassigned space that a region may be built on
- 1 - A buffer cell that isn't part of any region, but close to one. This is not available for building on.
- 2, 3, 4, ... - valid region labels.
Thus, the algorithm for building a region becomes simple:
- Find a valid seed cell
- Grow outward from the seed cell via 4-connected cells until the desired number of cells is reached
- A new Region Cell can only be placed in an Unassigned Space
- Iterate over Region Cells marking their 8-connected neighbors as Buffer cells.
- Note that, if we add the new Buffer cells to the list of Region Cells and repeat the last step, we can grow arbitrarily large buffers between regions.
A few other parts of your code needed cleaning. Types are typically defined using Camel Case and variables use snake case. I prefer snake_case for functions as well. See a style guide for more opinions.
Do not use global variables they are bad and, in the fullness of time, they will bring you nothing but sorrow.
You're using C++, avoid using headers such as stdio, stdlib, and cstring. They provide the C way of doing things, which often come with disadvantageous such as diminished type safety.
Any time your if statements have a lot of conditions, ask yourself whether they can be simplified by inverting one or more of the conditions (that is, switching is_good() to !is_good()) and doing an early exit. This can make your code much easier to read and reasonable about.
Nice use of std::vector::reserve(). I have seen companies wasting >$100k/yr doing unnecessary vector reallocations.
Ranged loops such as
for (coord Member : Regions.second){
should be written as
for (const auto& member : Regions.second){
auto is almost always better as a data type because it ensures you do the right thing and it makes it obvious to reviewers of your code that there are no tricks in play. const protects from unintended writes. auto& means "take a reference, don't make a copy". I've seen companies spending >$10M/year because they fot about using & to use references in for-loops.
Source Code
#include <array>
#include <iostream>
#include <optional>
#include <random>
#include <string>
#include <tuple>
#include <vector>
constexpr int REGION_SEED_TRIES = 100;
constexpr int NEW_MEMBER_TRIES = 200;
constexpr std::array<char, 10> REGION_SYMBOLS =
{'0', '1', '2', '3', '4', '5', '6', '7', '8', '9'};
struct coord {
int x;
int y;
coord(int x, int y) : x(x), y(y) {}
bool adjacent(const coord& o) const {
return std::abs(x - o.x) + std::abs(y - o.y) == 1;
}
};
// Methods written in inadvisable one-line format for ease of display on
// StackOverflow
struct Map {
using CellType = uint8_t;
static constexpr int UNASSIGNED = 0;
static constexpr int BUFFER = 1;
int free_label = 2;
std::vector<CellType> data;
int m_width;
int m_height;
Map(int width, int height) {
data.resize(width * height, UNASSIGNED);
m_width = width;
m_height = height;
}
CellType& operator[](const coord& c) {
return data.at(c.y * m_width + c.x);
}
CellType operator[](const coord& c) const {
return data.at(c.y * m_width + c.x);
}
void clear() {
std::fill(data.begin(), data.end(), false);
}
int get_free_label() {
return free_label++;
}
int width() const {
return m_width;
}
int height() const {
return m_height;
}
bool on_edge(const coord& c) const {
return c.x == 0 || c.y == 0 || c.x == m_width - 1 || c.y == m_height - 1;
}
bool in_body(const coord& c) const {
return c.x >= 0 && c.y >= 0 && c.x < m_width && c.y < m_height;
}
bool neighbor_is_buffer(const coord& c) const {
for (int ny = c.y - 1; ny <= c.y + 1; ++ny) {
for (int nx = c.x - 1; nx <= c.x + 1; ++nx) {
const coord n{nx, ny};
if (!in_body(n)) {
continue;
}
if (data[n.y * m_width + n.x] == BUFFER) {
return true;
}
}
}
return false;
}
};
using RegionBoundary = std::vector<coord>;
int rand_int(const int lb, const int ub, std::mt19937& prng) {
return std::uniform_int_distribution<>(lb, ub)(prng);
}
void create_region(
Map& map,
std::mt19937& prng,
const int distance_to_neighboring_regions) {
const auto region_label = map.get_free_label();
const auto region_size = rand_int(50, map.width() * map.height() / 2, prng);
std::vector<coord> region_members;
region_members.reserve(region_size);
// Seed the region
for (int i = 0; i < REGION_SEED_TRIES; i++) {
const int x = rand_int(1, map.width() - 2, prng);
const int y = rand_int(1, map.height() - 2, prng);
if (map[{x, y}] == Map::UNASSIGNED) {
region_members.emplace_back(x, y);
break;
}
}
if (region_members.empty()) {
throw std::runtime_error(
"Unable to find a seed location in a reasonable amount of time");
}
// Grow the region
for (int tries = 0;
tries < NEW_MEMBER_TRIES && region_members.size() < region_size;
tries++) {
// Choose a random member of the region to grow from
const auto& member =
region_members.at(rand_int(0, region_members.size() - 1, prng));
// Choose a random neighbour of the member
const coord neighbour(
member.x + rand_int(-1, 1, prng), member.y + rand_int(-1, 1, prng));
// Early exits so that the logic of the function is easier to follow
if (!neighbour.adjacent(member)) {
continue;
}
if (map.on_edge(neighbour)) {
continue;
}
if (map[neighbour] != Map::UNASSIGNED) {
continue;
}
// Neighbor is a valid candidate for expansion, so add it to the region
map[neighbour] = region_label;
region_members.push_back(neighbour);
tries = 0; // Reset the number of tries
}
// Mark cells adjacent to the region so regions don't grow into each other
// We grow outwards from the region converting UNASSIGNED cells to BUFFER
// cells and we repeat this several times, if necessary, to ensure an adequate
// buffer between regions
for (int buffer = 0; buffer < distance_to_neighboring_regions; buffer++) {
for (const auto& member : region_members) {
// We will include all 8 cells around the focal cell
for (int nx = member.x - 1; nx <= member.x + 1; ++nx) {
for (int ny = member.y - 1; ny <= member.y + 1; ++ny) {
const coord n{nx, ny};
if (!map.in_body(n)) {
continue;
}
if (map[n] == Map::UNASSIGNED) {
map[n] = Map::BUFFER;
// Add buffer cell to the region so that we can continue to grow the
// buffer between regions if we need to
region_members.push_back(n);
}
}
}
}
}
}
RegionBoundary get_region_boundary(const Map& map, const int region_label) {
// Iterate over all cells in the map considering only those that are
// of type `region_label` and are adjacent to a cell of type `Map::BUFFER`
RegionBoundary boundary;
for (int y = 0; y < map.height(); ++y) {
for (int x = 0; x < map.width(); ++x) {
const auto type = map[{x, y}];
if (type == region_label && map.neighbor_is_buffer({x, y})) {
boundary.emplace_back(x, y);
}
}
}
return boundary;
}
void display_map(const Map& map, const std::vector<coord>& highlight) {
for (int y = 0; y < map.height(); ++y) {
for (int x = 0; x < map.width(); ++x) {
const coord c = {x, y};
const auto type = map[c];
const char symbol = [&]() {
if (std::any_of(highlight.begin(), highlight.end(), [&](const auto& x) {
return x.x == c.x && x.y == c.y;
})) {
return '!';
}
if (type == Map::UNASSIGNED || type == Map::BUFFER) {
return ' ';
} else {
return REGION_SYMBOLS.at((type - 2) % REGION_SYMBOLS.size());
}
}();
std::cout << symbol;
}
std::cout << "\n";
}
}
std::optional<std::tuple<coord, coord, double>> closest_points_on_boundary(
const RegionBoundary& boundary1,
const RegionBoundary& boundary2) {
std::optional<std::tuple<coord, coord, double>> closest_points;
for (const auto& p1 : boundary1) {
for (const auto& p2 : boundary2) {
const auto distance = std::pow(p1.x - p2.x, 2) + std::pow(p1.y - p2.y, 2);
if (!closest_points || distance < std::get<2>(*closest_points)) {
closest_points = {p1, p2, distance};
}
}
}
return closest_points;
}
int main(int argc, char** argv) {
std::random_device rd;
std::mt19937 prng(rd());
if (argc != 3) {
std::cerr << "Usage: " << argv[0] << " width height\n";
return 1;
}
const int width = std::stoi(argv[1]);
const int height = std::stoi(argv[2]);
Map map(width, height);
create_region(map, prng, 5);
create_region(map, prng, 5);
const auto boundary2 = get_region_boundary(map, 2);
const auto boundary3 = get_region_boundary(map, 3);
const auto closest = closest_points_on_boundary(boundary2, boundary3);
if (closest) {
display_map(map, {std::get<0>(*closest), std::get<1>(*closest)});
}
return 0;
}
Producing output like this
00
000000000
00 00000000
00000000000
00000000000
0000000000000
000000000000
0000000000
000000 000
000000000
000 0!00
111
111
1111
1111
!11111111111
111111111111
If regions are filled (not just random set of points), minimal distance between regions will be between 2 points on borders of this two regions. So in first step filter points of regions to find only those one on border (O(4*m + 4*n) "if") and then do brute force on borders. In 2d space border points count will be equal to something about sqrt(m). So resulting algorithm will have complexity O(4*sqrt(m)*sqrt(n)).
According to @ufok's answer, one can have a brute force approach only on border of the two regions. So the main purpose is to find these borders.
Here is the main function of a possible solution:
int main()
{
int Seed = time(NULL);
printf ("seed: %d\n", Seed);
std::srand(Seed);
ClearMap();
auto[RegionA, RegionB] = CreateTwoRegions();
std::pair<coord,coord> lookup = std::make_pair (coord(0,0),coord(0,0));
int dmin = std::numeric_limits<int>::max();
// We compute the borders
auto borderA = border(RegionA);
auto borderB = border(RegionB);
// We perform the Cartesian product of the borders of the two regions.
for (auto a : borderA)
{
for (auto b : borderB)
{
int d = distance2(a,b);
if (d<dmin)
{
dmin = d;
lookup = std::make_pair (a,b);
}
}
}
// We tag the found coords in the map.
GetMapBool(lookup.first) = MATCH;
GetMapBool(lookup.second) = MATCH;
std::cout << lookup.first << " " << lookup.second << " " << dmin << "\n";
DisplayMap();
}
So the interesting part is the border function:
auto border (region const& reg)
{
coords_bitset seen;
std::set<coord> result;
border_rec (reg.front(), seen, result);
return result;
}
that calls the recursive function border_rec:
void border_rec (coord const& c, coords_bitset& seen, std::set<coord>& result, int depth=0)
{
// We skip points that 1) are not in a region or 2) have already been processed
// Both checks are done in O(1)
if (GetMapBool(c)==OUT or seen[c.pos()] ) { return; }
// We tag the point as seen.
seen[c.pos()] = true;
size_t count = 0;
size_t total = 0;
auto iter_neighbours = [] <typename FCT> (coord const& c, FCT fct)
{
// We try all possible neighbours of the current point.
for (int dx : {-1,0,1})
{
for (int dy : {-1,0,1})
{
coord o { c.X+dx, c.Y+dy };
if (o.IsInsideMap()) { fct (o); }
}
}
};
// We check whether the current point is an interior point.
iter_neighbours (c, [&] (coord const& o)
{
count += GetMapBool(o)==OUT ? 0 : 1;
total ++;
});
iter_neighbours (c, [&] (coord const& o)
{
// This is not an interior point => we use it as a border point.
if (count<total)
{
result.insert(c);
// uncomment the following for displaying the border
// GetMapBool(c) = BORDER;
}
// We continue recursively.
border_rec (o, seen, result);
});
}
The border points are put in a std::set. We also need to track points already processed; we can use a std::bitset for this, which can work since we define an order between points through the coord::pos function.
Full code here
Here is a sample, where the two required points are tagged with M:
..........|.||....................................................................................................
.......|..||...................................................................................................
.......||...................................................................................................
....|.||M......................|............................................................................
...||.||.......................|.|..........................................................................
...|.....................|..|.........................................................................
..............M|...|.|..........................................................................
..............||.|.|.........................................................................
..||................|........................................................................
...|................|.....................................................................
...................||.....................................................................
......|...............||.....................................................................
...................|...................................................................
.......|...............................................................................
.......|.............|.................................................................
.......|.|.|.............|................................................................
.......|...|..............||..............................................................
............|...............||............................................................
...........................|..............................................................
.............||...............||.............................................................
.............||.|................|...........................................................
...............|.|.|................||.||..........................................................
.............|.|....||...||..................|..|...........................................................
.......................|...|....|..................|...........................................................
...........................|..|...|...................|...........................................................
................................|....|..................................................................................
.....................................|..................................................................................
............................................................|...........................................................
............................................................||.............................................................
.............................................................|.||.............................................................
........................................................................................................................
........................................................................................................................
.............................................................|.|..............................................................
............................................................||....||..||...............................................................
.......................................................................................................................................
.......................................................................|..|.|................................................................
..........................................................................||....|..................................................................
and the same result with border points tagged with + and the interior points with :
..........+ +.+ +....................................................................................................
.......+..+ +++ ++...................................................................................................
.......++++ +...................................................................................................
....+.++ +M......................+............................................................................
...++.+ ++.......................+.+..........................................................................
...++++ +.....................+..++++.........................................................................
..++ ++............M++++...+.++ +..........................................................................
..+ ++............++.+.++++++ ++.........................................................................
..++ +................+++ ++........................................................................
...++ +...............++ ++.+.....................................................................
....+++ +...............+ +++.....................................................................
......+ ++...............+ +.....................................................................
......++++ +.............+.++ ++++..................................................................
.......+.+ +.............+++ ++.................................................................
.......+.+ +++............++ +.................................................................
.......+.++ +.+.............++ ++................................................................
.......+..+++ +.+..............+ +++..............................................................
............+ ++...............++ +++............................................................
............++ +...............++ +..............................................................
.............+ +++ ++...............+ ++.............................................................
.............+ +++.+ +++ +++ +................+ +++...........................................................
.............+++ +..++++.+ ++.+ +++................+ ++.++..........................................................
.............+.++++....++.+++..+ +++ +..................+ +..+...........................................................
.......................+...+...+ +.+++ ++..................+ ++++...........................................................
...........................+..++++...++++...................++ +...........................................................
................................+....+.......................+ ++...........................................................
.....................................+......................++ +............................................................
............................................................+ +++...........................................................
............................................................+ +++++++.............................................................
...........................................................++ +..+.++.............................................................
........................................................++++ ++++................................................................
........................................................++++++++ +................................................................
.............................................................+.++++++++ +++..............................................................
............................................................++....++..++ ++...............................................................
.......................................................................++++ +++ +................................................................
.......................................................................+..++++.++++ ++................................................................
..........................................................................++....+.+ +.................................................................
Update
In function border_rec, I don't know why I put the insertion of a border point in the std::set during the iteration of all neighbours => it can be put outside and it is therefore possible to use a std::vector instead of a std::set, which should be faster:
void border_rec (coord const& c, coords_bitset& seen, std::vector<coord>& result, int depth=0)
{
// We skip points that 1) are not in a region or 2) have already been processed
// Both checks are done in O(1)
if (GetMapBool(c)==OUT or seen[c.pos()] ) { return; }
// We tag the point as seen.
seen[c.pos()] = true;
size_t count = 0;
size_t total = 0;
auto iter_neighbours = [] <typename FCT> (coord const& c, FCT fct)
{
// We try all possible neighbours of the current point.
for (int dx : {-1,0,1})
{
for (int dy : {-1,0,1})
{
coord o { c.X+dx, c.Y+dy };
if (o.IsInsideMap()) { fct (o); }
}
}
};
// We check whether the current point is an interior point.
iter_neighbours (c, [&] (coord const& o)
{
count += GetMapBool(o)==OUT ? 0 : 1;
total ++;
});
// This is not an interior point => we use it as a border point.
if (count<total) { result.push_back (c); }
iter_neighbours (c, [&] (coord const& o)
{
// We continue recursively.
border_rec (o, seen, result);
});
}
Demo
Just for the fun, it is possible to find the solution at compile-time :)
We define for instance our map with a string
constexpr char mapstr[] =
// 0000000000111111111122222222223333333333444444444455555555556666
// 0123456789012345678901234567890123456789012345678901234567890123
/*0*/ " 1 1 1 "
/*1*/ " 1 1 11111111 "
/*2*/ " 11111111 1 111 1111 "
/*3*/ " 111 111111 1 111 11 1 "
/*4*/ " 111111 111111 M111111111 "
/*5*/ " 1111111111111M 1111111 "
/*6*/ " 11111111111111 11111111 "
/*7*/ " 111 11111111 1 1 ";
and at least one point in each region:
constexpr coord pointInRegionA = coord(10,1);
constexpr coord pointInRegionB = coord(44,0);
The hunted points are defined by (and tagged with M in the string just for visualizing them):
constexpr auto truth = std::make_pair (coord(18,5), coord(38,4));
Then you can run the algorithm at compile time with compute_min_distance and check the result with:
static_assert (compute_min_distance() == truth);
Demo
Ok, I know, it's not that useful ;). Note that one may have to help the compiler by setting its recursion level.
Update:
If you like std::views, you can write compute_min_distance this way:
constexpr auto borderA = border (pointInRegionA);
constexpr auto borderB = border (pointInRegionB);
constexpr auto compute_min_distance ()
{
auto v = std::views::cartesian_product (
borderA.first | std::views::take(borderA.second),
borderB.first | std::views::take(borderB.second)
) ;
return *std::ranges::min_element (v, [](auto const& u, auto const& v) { return distance2(u) < distance2(v); });
}
gotoimproves readability and reduces complexity, the time spent explaining and fighting for thegotoexceeds the benefits. – user4581301 Commented Mar 18 at 22:53