티스토리 뷰
최근 크리스토피데스 알고리즘을 직접 구현해서 사용해볼 기회가 생겼다.
Swift 뿐만 아니라 다른 언어에서도 이 알고리즘을 순수하게 작성한 자료가 없는 것 같아서 따로 글로써 작성해본다.
(*순수하게라는 뜻은 내장 라이브러리를 사용하지 않고 구현한 것을 뜻한다)
다만, 내가 구현한 알고리즘은 TSP문제에서 시작 정점이 정해져있고 다시 시작 정점으로 돌아오지 않는 경우에 적용할 수 있는 알고리즘이다. 만약 기존 TSP 문제의 상황에서 사용하고자 한다면 몇가지 수정을 거쳐야한다. 이 점 참고하면 좋을 듯 하다.
크리스토피데스 알고리즘이란?
TSP(외판원 순회 문제)는 한 세일즈 맨이 자신의 도시에서부터 시작하여 다른 도시들을 방문해가며 출장을 다녀야할 때, 소모해야하는 최소 비용을 구하는 문제이다. 유명한 NP-Hard 문제이고, (아마도) 다항시간에 해결하는 알고리즘이 존재하지 않는다. 그래서 이 문제를 정확하게 해결 하고자 한다면 지수시간 이상의 시간복잡도를 가지고, 이는 도시의 수가 늘어날 때 마다 어마어마하게 연산이 늘어난다는 뜻이다.
그렇다면 우리는 항상 이런 어마어마한 시간을 소요해서 최소 비용을 구해야할까? 문제의 상황마다 다르겠지만, 분명히 정확한 계산이 필요하지 않는 경우도 있다. 이때 적용할 수 있는 알고리즘이 크리스토피데스 알고리즘이다. 이 알고리즘은 최적화에 따라 약 O(N^3)의 다항시간의 시간복잡도를 가진다. 기존 알고리즘에 비하면 실로 어마어마한 개선이다. 하지만 시간 복잡도를 크게 개선할 수 있었던 대신, 정확도는 떨어진다. 최적의 경로로 구한 최소 비용 보다 1.5배까지의 비용 차이가 있을 수 있다. 그래도 충분히 활용가치가 높다고 생각한다.
하지만 크리스토피데스 알고리즘을 사용하기 위한 제약 조건이 있다. triangle inequality를 만족해야한다. 이는 아래와 같은 조건을 만족하는 그래프를 뜻한다.
- 그래프의 지점 집합 V의 내부 정점들이 있을 때, 거리 d에 대해서
- d(u,v) = 0 이라면 u = v
- 임의의 정점 u, v에 대하여, d(u, v) = d(v, u)
- 임의의 정점 u, v, w에 대하여 d(u, w) <= d(u, v) + d(u, w)
이 조건이 만족되지 않는다면, 1.5배의 최대 오차가 더 커질 수 있다.
전체적인 크리스토피데스 알고리즘의 동작
1. 최소 신장 트리를 찾는다
2. 해밀턴 경로를 찾기 위해 잘못된 차수를 갖는 정점을 서로 매칭한다. 이때 매칭 비용이 큰 정점 하나를 제외한다.
- 잘못된 차수는 시작 정점이 짝수 차수인 경우, 시작 정점 외에 정점이 홀수 차수인 경우이다.
- 최소 비용을 갖도록 정점을 서로 매칭하는 것을 minimum-cost perfect matching이라 부른다. 이것을 Hungarian Algorithm이나 blossom algorithm을 사용해서 해결한다. (나는 Hungarian Algorithm을 사용했다)
3. 만약 시작 정점이 매칭에 참가했지만, 매칭된 정점에서 시작 정점이 없다면 시작 정점과 인접한 간선을 하나 제거한다.
4. 오일러 경로를 구한다.
5. 오일러 경로를 따라 방문하며, 이미 방문한 정점이라면 생략하고 다음 정점에 방문한다.
자세한 설명과 원리에 대해서는 크리스토피데스 알고리즘으로 경로 외판원 문제 풀기 (Path TSP by Christofides' Algorithm) 이곳에서 보는게 훨씬 좋을 것 같다.
내가 구현한 코드는 O(N^4)의 시간복잡도를 가진다. N이 30이고 클럭 속도가 3.2GHz일 때, 0.1초 정도의 수행시간이 소요되었다. 하나의 정점을 제외해나가며 Minimum-cost matching을 얻는 것과 정점의 매칭이 중복되지 않도록 재매칭하는 것이 수행시간을 엄청 늘리는 것 같다. 필요하다면 이들을 조율해서 정확도를 포기하고 수행속도를 늘리는 방안도 고려해보면 좋을 것 같다. 그럼 더 이상 크리스토피데스 알고리즘이 아니겠지만...
코드
christofidesTSP
func christofidesTSP(matrix: [[Double]], start: Int) -> ([Int], Double) {
let n = matrix.count
var mstEdges = kruskalMST(matrix: matrix, n: n)
var degree = Array(repeating: 0, count: n)
for (u, v) in mstEdges {
degree[u] += 1
degree[v] += 1
}
var matchVertices = [Int]()
// 시작 정점이 짝수인 경우 minimumCostMatching에 참여시킴
// 오일러 경로를 그릴 때 시작 지점은 홀수 차수여야하기 때문
var isStartVertexJoinMatching = false
if degree[start]%2 == 0 {
matchVertices.append(start)
isStartVertexJoinMatching = true
}
for i in 0..<n {
guard i != start else { continue }
if degree[i] % 2 == 1 {
matchVertices.append(i)
}
}
var minMatchedEdges: [(Int, Int)] = []
var minMatchedCost: Double = Double.infinity
for excludeVertex in matchVertices {
let excludedVertices = matchVertices.filter{ $0 != excludeVertex }
let (matchedEdges, matchedCost) = minimumCostMatching(vertices: excludedVertices, matrix: matrix)
if matchedCost < minMatchedCost {
minMatchedEdges = matchedEdges
minMatchedCost = matchedCost
}
}
// 시작 정점이 짝수라서 minimumCostMatching에 참가 시켰는데, 매칭되지 않은 경우 start와 연결된 임의의 간선을 제거
// 이를 통해 시작 정점의 차수를 홀수로 변경 가능
if isStartVertexJoinMatching, minMatchedEdges.filter({ $0.0 == start || $0.1 == start }).isEmpty {
let arbitraryStartEdgeIndex = mstEdges.firstIndex(where: { $0.0 == start || $0.1 == start })!
mstEdges.remove(at: arbitraryStartEdgeIndex)
}
let combinedEdges = mstEdges + minMatchedEdges
let eulerianPath = eulerianPath(from: combinedEdges, start: start, n: n)
var visited = Array(repeating: false, count: n)
var tspPath: [Int] = []
var cost: Double = 0
var prevVertex: Int?
for v in eulerianPath {
if !visited[v] {
tspPath.append(v)
visited[v] = true
if let prevVertex = prevVertex {
cost += matrix[prevVertex][v]
}
prevVertex = v
}
}
return (tspPath, cost)
}
Kruskal - MST 구하기
func kruskalMST(matrix: [[Double]], n: Int) -> [(Int, Int)] {
var mstEdges: [(Int, Int)] = []
var edges: [(u: Int, v: Int, cost: Double)] = []
for u in 0..<n {
for v in u+1..<n {
if matrix[u][v] != 0 {
edges.append((u, v, matrix[u][v]))
}
}
}
edges.sort { $0.cost < $1.cost }
let unionFind = UnionFind(size: n)
for edge in edges {
if unionFind.find(edge.u) != unionFind.find(edge.v) {
mstEdges.append((edge.u, edge.v))
unionFind.union(edge.u, edge.v)
}
if mstEdges.count == n - 1 { break }
}
return mstEdges
}
class UnionFind {
private var parent: [Int]
private var rank: [Int]
init(size: Int) {
parent = Array(0..<size)
rank = Array(repeating: 0, count: size)
}
func find(_ node: Int) -> Int {
if parent[node] != node {
parent[node] = find(parent[node])
}
return parent[node]
}
func union(_ u: Int, _ v: Int) {
let rootU = find(u)
let rootV = find(v)
if rootU != rootV {
if rank[rootU] > rank[rootV] {
parent[rootV] = rootU
} else if rank[rootU] < rank[rootV] {
parent[rootU] = rootV
} else {
parent[rootV] = rootU
rank[rootU] += 1
}
}
}
}
Minimum-cost matching (Hungarian Algorithm)
func minimumCostMatching(vertices: [Int], matrix: [[Double]]) -> (matchedEdges: [(Int, Int)], cost: Double) {
let n = vertices.count
guard n > 1 else { return ([], 0) }
var originCost: [[Double]] = Array(repeating: Array(repeating: Double.infinity, count: n), count: n)
for i in 0..<vertices.count {
let u = vertices[i]
for j in i+1..<vertices.count {
let v = vertices[j]
originCost[i][j] = matrix[u][v]
originCost[j][i] = matrix[v][u]
}
}
var cost = originCost
// 각 행에서 최소값을 뺌
for i in 0..<n {
let rowMin = cost[i].min() ?? 0
for j in 0..<n {
cost[i][j] -= rowMin
}
}
// 각 열에서 최소값을 뺌
for j in 0..<n {
let colMin = (0..<n).map { cost[$0][j] }.min() ?? 0
for i in 0..<n {
cost[i][j] -= colMin
}
}
let rows = cost.count
let cols = cost[0].count
var rowCovered = Array(repeating: false, count: rows)
var colCovered = Array(repeating: false, count: cols)
var matchedEdges: [(Int, Int)] = []
func isMaximumMatching() -> Bool {
var costCopy = cost
rowCovered = Array(repeating: false, count: rows)
colCovered = Array(repeating: false, count: cols)
matchedEdges = []
var totalLine = 0
var zeroVertex: [(row: Int, col: Int)] = []
for row in 0..<rows {
for col in 0..<cols {
if costCopy[row][col] == 0 {
zeroVertex.append((row, col))
}
}
}
let deletedZero: Double = -1
while !zeroVertex.isEmpty {
// rowScanning
for row in 0..<rows where !rowCovered[row] {
var zeroCount = 0
var zeroColIndex = 0
for col in 0..<cols {
if costCopy[row][col] == 0 {
zeroCount += 1
zeroColIndex = col
}
}
// zero가 하나라면, column으로 Line 긋기
if zeroCount == 1 {
totalLine += 1
for deleteRow in 0..<rows {
if costCopy[deleteRow][zeroColIndex] == 0 {
let zeroVertexIndex = zeroVertex.firstIndex(where: { $0.0 == deleteRow && $0.1 == zeroColIndex })!
zeroVertex.remove(at: zeroVertexIndex)
}
costCopy[deleteRow][zeroColIndex] = deletedZero
}
colCovered[zeroColIndex] = true
matchedEdges.append((row, zeroColIndex))
}
}
// columnScanning
for col in 0..<cols where !colCovered[col] {
var zeroCount = 0
var zeroRowIndex = 0
for row in 0..<rows {
if costCopy[row][col] == 0 {
zeroCount += 1
zeroRowIndex = row
}
}
if zeroCount == 1 {
totalLine += 1
for deleteCols in 0..<cols {
if costCopy[zeroRowIndex][deleteCols] == 0 {
let zeroVertexIndex = zeroVertex.firstIndex(where: { $0.0 == zeroRowIndex && $0.1 == deleteCols })!
zeroVertex.remove(at: zeroVertexIndex)
}
costCopy[zeroRowIndex][deleteCols] = deletedZero
}
rowCovered[zeroRowIndex] = true
matchedEdges.append((zeroRowIndex, col))
}
}
// 모든 0이 cover 안됐다면, digonal selection 실행
guard !zeroVertex.isEmpty else { break }
for row in 0..<rows where !rowCovered[row] {
for col in 0..<cols where !colCovered[col] {
guard costCopy[row][col] == 0 else { continue }
guard let adjacentZeroVertex = zeroVertex.first(where: { $0.0 > row && $0.1 != col }) else { continue }
for (coverRow, coverCol) in [(row, col), adjacentZeroVertex] {
totalLine += 1
colCovered[coverCol] = true
for deleteRow in 0..<rows {
if costCopy[deleteRow][coverCol] == 0 {
let zeroVertexIndex = zeroVertex.firstIndex(where: { $0.0 == deleteRow && $0.1 == coverCol })!
zeroVertex.remove(at: zeroVertexIndex)
}
costCopy[deleteRow][coverCol] = deletedZero
}
matchedEdges.append((coverRow, coverCol))
}
}
}
}
return totalLine == n
}
func findUncoveredMinNumber() -> Double {
var minNumber = Double.infinity
for i in 0..<rows {
for j in 0..<cols {
if !rowCovered[i] && !colCovered[j] {
minNumber = min(minNumber, cost[i][j])
}
}
}
return minNumber
}
while !isMaximumMatching() {
let minUncoveredNumber = findUncoveredMinNumber()
for row in 0..<rows {
for col in 0..<cols {
if rowCovered[row], colCovered[col] {
cost[row][col] += minUncoveredNumber
}
if !rowCovered[row], !colCovered[col] {
cost[row][col] -= minUncoveredNumber
}
}
}
}
var nonDuplicatedMatching: [(Int, Int, Double)] = []
var visitedVertex: Set<Int> = []
let matchedEdgesWithCost = matchedEdges.map{ ($0.0, $0.1, originCost[$0.0][$0.1]) }
let sortedMatchedVertex = matchedEdgesWithCost.sorted{ $0.2 < $1.2 }
for (u, v, cost) in sortedMatchedVertex {
guard !visitedVertex.contains(u), !visitedVertex.contains(v) else { continue }
visitedVertex.insert(u)
visitedVertex.insert(v)
nonDuplicatedMatching.append((u, v, cost))
}
var finalMatchedEdges = nonDuplicatedMatching.map{ (vertices[$0.0], vertices[$0.1]) }
var finalCost = nonDuplicatedMatching.reduce(0){ $0 + $1.2 }
var remainVertices = vertices.filter { vertex in
!finalMatchedEdges.contains(where: { matched in
matched.0 == vertex || matched.1 == vertex
})
}
while !remainVertices.isEmpty {
let (remainMatchedEdges, remainCost) = minimumCostMatching(vertices: remainVertices, matrix: matrix)
remainVertices = remainVertices.filter { vertex in
!remainMatchedEdges.contains(where: { matched in
matched.0 == vertex || matched.1 == vertex
})
}
finalMatchedEdges.append(contentsOf: remainMatchedEdges)
finalCost += remainCost
}
return (finalMatchedEdges, finalCost)
}
오일러 경로
func eulerianPath(from edges: [(Int, Int)], start: Int, n: Int) -> [Int] {
var adjList = Array(repeating: [Int](), count: n)
for (u, v) in edges {
adjList[u].append(v)
adjList[v].append(u)
}
var path = [Int]()
func findEulerianPath(_ path: inout [Int], _ u: Int) {
while !adjList[u].isEmpty {
let v = adjList[u].removeLast()
if let index = adjList[v].firstIndex(of: u) {
adjList[v].remove(at: index)
}
findEulerianPath(&path, v)
}
path.append(u)
}
findEulerianPath(&path, start)
return path.reversed()
}
테스트
실제로 최적해의 1.5배 근사해를 갖는지 확인하고자 임의의 matrix를 생성하여 브루트포스 알고리즘 풀이와 비교해봤다. 3~13개의 정점을 갖는 matrix를 랜덤하게 생성하여 10000개의 케이스에 대해 테스트했다.
모두 다 잘 통과했고, 간략하게 분포를 봤는데 1.0 ~ 1.3의 근사해를 가는 결과가 엄청 많았다. 1.4의 근사해는 잘 안나오는 것 같아서 실제 효용성이 생각했던 것 보다 훨씬 좋은 것 같다.
for _ in 0..<10000 {
let n = Int.random(in: 3...13)
let graph = generateCompleteGraph(n: n)
let (christofidesPath, christofidesCost) = christofidesTSP(matrix: graph, start: 0)
let (shortestPath, shortestPathsCost) = shortestPaths(startPath: 0, distanceMatrix: graph)
print("christofidesPath: ", christofidesPath)
print("christofidesCost: ", christofidesCost)
print("shortestPath: ", shortestPath)
print("shortestCost: ", shortestPathsCost)
let ratio = christofidesCost / shortestPathsCost
if ratio > 1.5 {
assertionFailure()
}
}
func generateCompleteGraph(n: Int) -> [[Double]] {
var graph = Array(repeating: Array(repeating: 0.0, count: n), count: n)
for i in 0..<n {
for j in i+1..<n {
let weight = Double(Int.random(in: 1...1000))
graph[i][j] = weight
graph[j][i] = weight
}
}
for k in 0..<n {
for i in 0..<n {
for j in 0..<n {
if i != j && i != k && j != k {
let sum = graph[i][k] + graph[k][j]
if graph[i][j] > sum {
graph[i][j] = sum - 1
graph[j][i] = graph[i][j]
}
}
}
}
}
return graph
}
마치며
특수하고 희귀한(?) 알고리즘을 구현해보면서 많은 공부와 성취감을 느꼈다. 글에는 안썼지만 구현하기까지 많은 시행착오와 시도가 있었다. 관련한 자료가 부족하다는게 가장 큰 난관이였는데, 퍼즐 조립하듯이 하나하나 공부하고 테스트해내가면서 해결했다. 이미 구현된 라이브러리(예를 들면 Minimum-cost matching를 파이썬 코드에서는 라이브러리를 사용하고 있었다)를 대부분 재사용했기 때문에 나는 시행착오를 겪으며 코드를 처음부터 새롭게 짜야했다. 😢😢😢
그래도 이렇게 끝까지 해서 구현해내고, 잘 동작하니까 뿌듯하다! 어렵게 시간들여서 열심히 구현한 코드가 혹시나 누군가에게 도움이 된다면 좋겠다. 😊
'Swift 알고리즘' 카테고리의 다른 글
| [백준 - 21940] 가운데에서 만나기 - Swift (0) | 2023.10.31 |
|---|