SCPC 2019 2차예선 4번 구현 소스 코드

SCPC 2019 2차예선 4번. 폭격의 소스 코드입니다.

휴리스틱 + Simulated Annealing을 써서 해결했습니다.

 

구현 아이디어 설명은 https://paido.tistory.com/16을 참고해주세요.

 

#define _CRT_SECURE_NO_WARNINGS
 
#include <stdio.h>
#include <stdlib.h>
#include <algorithm>
#include <vector>
#include <math.h>
#include <memory.h>
 
///////////////////////////////////////////
 
int M, N;
char Targets[50][501];
 
void ReadInput() {
    scanf("%d %d", &M, &N);
    for (int i = 0; i < M; i++) scanf("%s", Targets[i]);
}
 
///////////////////////////////////////////
 
struct Attack {
    int R;
    int C;
    int CoverageCount;
};
 
bool NecessaryMap[50][500];
 
int CoverageMap[50][500]; // counts excess coverage =>
                          // valid solution should not have any negative entry
bool AttackMap[50][500];
int PCount;
 
double RandDouble() {
    int val = rand();
    return val / (double)RAND_MAX;
}
 
struct Event {
    int r;
    int c;
    bool add;
};
Event Journal[25000];
int JournalCount;
 
void MutateSolution(bool cleanup) {
    JournalCount = 0;
 
    bool removed = false;
    if (cleanup || RandDouble() < 0.5) {
        /* remove */
        int r = rand() % (M - 2);
        int c = rand() % (N - 2);
        for (int i = 0; i < M - 2 && (cleanup || !removed); i++) {
            for (int j = 0; j < N - 2 && (cleanup || !removed); j++) {
                int rr = (i + r) % (M - 2);
                int cc = (j + c) % (N - 2);
 
                if (!AttackMap[rr][cc]) continue;
 
                bool removable = true;
                for (int m = rr; m <= rr + 2; m++) {
                    for (int n = cc; n <= cc + 2; n++) {
                        if (CoverageMap[m][n] <= 0) {
                            removable = false;
                            break;
                        }
                    }
                }
 
                if (removable) {
                    AttackMap[rr][cc] = false;
                    for (int m = rr; m <= rr + 2; m++) {
                        for (int n = cc; n <= cc + 2; n++) {
                            CoverageMap[m][n]--;
                        }
                    }
                    PCount--;
                    auto &J = Journal[JournalCount++];
                    J.r = rr;
                    J.c = cc;
                    J.add = false;
                    removed = true;
                }
            }
        }
    }
    if (cleanup || removed) return;
 
    /* add */
    auto try_add = [](int r, int c) -> bool {
        if (NecessaryMap[r][c] && !AttackMap[r][c]) {
            AttackMap[r][c] = true;
            for (int i = r; i <= r + 2; i++) {
                for (int j = c; j <= c + 2; j++) {
                    CoverageMap[i][j]++;
                }
            }
            PCount++;
            auto &J = Journal[JournalCount++];
            J.r = r;
            J.c = c;
            J.add = true;
            return true;
        }
        return false;
    };
    for (int t = 0; t < 50; t++) { // if we do not limit # of tries, we get TLE
        int r = rand() % (M - 2);
        int c = rand() % (N - 2);
        if (try_add(r, c)) break;
    }
}
 
void Undo() {
    for (int i = 0; i < JournalCount; i++) {
        auto &J = Journal[i];
        int delta = J.add ? -1 : +1;
        for (int m = J.r; m <= J.r + 2; m++) {
            for (int n = J.c; n <= J.c + 2; n++) {
                CoverageMap[m][n] += delta;
            }
        }
        AttackMap[J.r][J.c] = !J.add;
        PCount += delta;
    }
}
 
double EvaluateEnergy() {
    return PCount;
}
 
 
void RunSA() {
    /* obtain approximate greedy solution */
    std::vector<Attack> attacks;
    attacks.reserve((M - 2)*(N - 2));
    memset(NecessaryMap, 0, sizeof(NecessaryMap));
    for (int i = 0; i < M - 2; i++) {
        for (int j = 0; j < N - 2; j++) {
            Attack a;
            a.R = i;
            a.C = j;
            a.CoverageCount = 0;
            for (int m = i; m <= i + 2; m++) {
                for (int n = j; n <= j + 2; n++) {
                    a.CoverageCount += (Targets[m][n] == '1') ? 1 : 0;
                }
            }
            if (a.CoverageCount > 0) {
                NecessaryMap[i][j] = true;
                attacks.push_back(a);
            }
        }
    }
    std::sort(attacks.begin(), attacks.end(),
      [](const Attack &a, const Attack &b) -> bool { return a.CoverageCount > b.CoverageCount; });
 
    memset(CoverageMap, 0, sizeof(CoverageMap));
    memset(AttackMap, 0, sizeof(AttackMap));
    for (int i = 0; i < M; i++) {
        for (int j = 0; j < N; j++) {
            CoverageMap[i][j] = (Targets[i][j] == '1') ? -1 : 0;
        }
    }
 
    PCount = 0;
    for (auto &a : attacks) {
        /* if this attack is needed, apply it */
        bool needed = false;
        for (int i = a.R; i <= a.R + 2; i++) {
            for (int j = a.C; j <= a.C + 2; j++) {
                if (CoverageMap[i][j] < 0) needed = true;
            }
        }
 
        if (needed) {
            AttackMap[a.R][a.C] = true;
            for (int i = a.R; i <= a.R + 2; i++) {
                for (int j = a.C; j <= a.C + 2; j++) {
                    CoverageMap[i][j]++;
                }
            }
            PCount++;
        }
    }
 
    /* simulated annealing */
    static const int cool_iter = 50;
    double E = EvaluateEnergy();
    for (int i = 0; i < cool_iter; i++) {
        MutateSolution(false);
        double NewE = EvaluateEnergy();
 
        bool replace;
        if (NewE < E) replace = true;
        else {
            double T = (1 / (double)(i + 1) - 1 / (double)(cool_iter + 1)) * cool_iter;
            double boltzmann_factor = std::exp(-(NewE - E) / T);
            replace = RandDouble() < boltzmann_factor;
        }
 
        if (replace) E = NewE; else Undo();
    }
    MutateSolution(true);
}
 
///////////////////////////////////////////
 
bool AttackMapOpt[50][500];
int PCountOpt;
 
void Solve() {
    PCountOpt = 1000000;
    for (int i = 0; i < 10; i++) {
        RunSA();
        if (PCount < PCountOpt) {
            memcpy(AttackMapOpt, AttackMap, sizeof(AttackMap));
            PCountOpt = PCount;
        }
    }
 
    /* print answer */
    printf("%d\n", PCountOpt);
    for (int i = 0; i < M - 2; i++) {
        for (int j = 0; j < N - 2; j++) {
            if (AttackMapOpt[i][j]) {
                printf("%d %d\n", i + 1, j + 1);
            }
        }
    }
}
 
///////////////////////////////////////////
 
int main() {
    srand(70938491);
 
    int TT;
    scanf("%d", &TT);
    for (int test_case = 1; test_case <= TT; test_case++) {
        ReadInput();
        printf("Case #%d\n", test_case);
        Solve();
    }
    return 0;
}