pkglog

C++ 분할정복: 병합 정렬, 퀵소트, 이진 탐색, 가장 가까운 점 쌍, Strassen 행렬 곱

작성 · · 22분 읽기 · 수정 고급 실습

분할정복이 빨라지는 이유

분할정복(Divide and Conquer)은 문제를 작은 하위 문제로 나누고(Divide), 각각을 재귀적으로 풀고(Conquer), 그 결과를 합쳐(Combine) 전체 답을 만드는 설계 기법입니다.

flowchart TD
  A[원본 문제] --> B[Divide: 작은 하위 문제로 분할]
  B --> C[Conquer: 각 하위 문제 해결 - 재귀 또는 기본 케이스]
  C --> D[Combine: 하위 결과를 합쳐 최종 답 생성]
  D --> E[최종 답]

분할정복이 빨라지는 이유는 나눴기 때문이 아니라 합치는 비용이 싸기 때문입니다. 병합정렬에서 크기 n인 문제를 n/2 두 개로 나누면 각 층의 병합 비용은 모두 합쳐 O(n)이고, 층의 수가 log n이므로 전체가 O(n log n)이 됩니다. 식으로 쓰면 T(n) = 2T(n/2) + O(n)이고, 마스터 정리로 바로 풀립니다. 반대로 합치는 단계가 O(n²)이라면 나눠 봐야 전체도 O(n²)에 머뭅니다. 그래서 새로운 문제에 분할정복을 적용할 때 가장 먼저 따져 볼 질문은 “두 부분의 답을 알 때 전체 답을 선형 시간 안에 만들 수 있는가”입니다. 이 글의 최근접 점 쌍 알고리즘이 까다로운 것도 이 Combine 단계를 선형으로 만드는 기하학적 논증 때문입니다.

예제는 구조화된 바인딩과 std::optional을 쓰므로 C++17 이상으로 컴파일합니다.

분할정복이 필요한 전형적인 문제

100만 건의 로그를 타임스탬프로 정렬할 때 버블 정렬 같은 O(n²) 알고리즘은 비교만 약 5천억 번이 필요해 최적화 빌드에서도 분 단위 시간이 걸립니다. 병합정렬은 같은 입력에서 비교 횟수가 약 n log₂ n ≈ 2천만 번에 그칩니다.

// O(n²): 전체를 한 번에 처리
void bubbleSort(std::vector<int>& a) {
    for (size_t i = 0; i < a.size(); ++i)
        for (size_t j = 0; j + 1 < a.size() - i; ++j)
            if (a[j] > a[j + 1]) std::swap(a[j], a[j + 1]);
}

정렬된 100만 개의 ID에서 특정 값을 찾을 때 순차 탐색은 O(n)이지만, 이진 탐색은 중간값과 비교해 범위를 절반씩 좁혀 최대 약 20번의 비교로 끝납니다. 평면 위 10만 개 점에서 가장 가까운 두 점을 찾을 때 모든 쌍을 비교하면 약 50억 쌍을 봐야 하지만, x축으로 분할하고 중앙 밴드만 검사하면 O(n log n)이 됩니다. 1024×1024 행렬 곱셈은 단순 구현으로 곱셈이 약 10억 번 필요한데, Strassen 알고리즘은 2×2 블록 곱셈을 8번 대신 7번의 곱셈으로 처리해 점근적으로 O(n^2.807)로 줄입니다.

문제 유형분할정복 알고리즘시간 복잡도비고
정렬병합정렬O(n log n)안정, 추가 공간 O(n)
정렬퀵소트평균 O(n log n), 최악 O(n²)제자리, 피벗 선택 중요
검색이진 탐색O(log n)정렬된 배열 전제
행렬 곱셈StrassenO(n^2.807)큰 n에서만 유리
최근접 점 쌍Closest PairO(n log n)아래 단순 구현은 O(n log² n)

Divide·Conquer·Combine 템플릿

// 분할정복 일반 패턴
template<typename T>
T divideConquer(const std::vector<T>& data, int lo, int hi) {
    // 1. Base case: 더 이상 분할 불가
    if (lo >= hi) {
        return baseCase(data, lo);
    }
    // 2. Divide: 중간 지점으로 분할
    int mid = lo + (hi - lo) / 2;
    // 3. Conquer: 각 부분 해결
    T leftResult  = divideConquer(data, lo, mid);
    T rightResult = divideConquer(data, mid + 1, hi);
    // 4. Combine: 결과 병합
    return combine(leftResult, rightResult);
}

mid를 (lo + hi) / 2가 아니라 lo + (hi - lo) / 2로 계산하는 것은 lo + hi가 int 범위를 넘는 경우를 피하기 위해서입니다. Base case(lo >= hi)를 가장 먼저 검사하지 않으면 재귀가 끝나지 않습니다.

재귀 vs 반복

// 반복(bottom-up) 병합정렬
void mergesortIterative(std::vector<int>& a) {
    int n = static_cast<int>(a.size());
    std::vector<int> tmp(n);
    for (int sz = 1; sz < n; sz *= 2) {
        for (int lo = 0; lo < n - sz; lo += 2 * sz) {
            int mid = lo + sz - 1;
            int hi = std::min(lo + 2 * sz - 1, n - 1);
            merge(a, lo, mid, hi, tmp);
        }
    }
}

반복 버전은 크기 1짜리 구간끼리 먼저 병합하고, 그다음 크기 2, 4, 8… 순서로 올라갑니다. 재귀 호출이 없으니 스택 깊이 걱정이 없고 함수 호출 오버헤드도 사라집니다. 다만 병합정렬은 재귀 깊이가 log₂n(100만 개여도 약 20)에 불과해서 스택 오버플로우가 실제로 문제가 되는 경우는 드뭅니다. 반복 처리가 정말 필요한 쪽은 최악의 경우 깊이가 n까지 갈 수 있는 퀵소트입니다. 퀵소트에서는 두 부분 중 작은 쪽만 재귀하고 큰 쪽은 반복문으로 처리하는 꼬리 재귀 제거로, 최악의 입력에서도 스택 깊이를 O(log n)으로 제한할 수 있습니다.

병합정렬부터 Strassen까지 분할정복 예제

병합정렬 (Merge Sort)

반으로 나누고, 각각 정렬한 뒤 두 정렬된 배열을 병합합니다. 안정 정렬이며 최악에도 O(n log n)을 보장합니다.

flowchart TD
  A[5,2,8,1,9] --> B["5,2,8 | 1,9"]
  B --> C["5,2 | 8 | 1 | 9"]
  C --> D["2,5,8 | 1,9"]
  D --> E[1,2,5,8,9]
#include <vector>
#include <algorithm>
// 병합: a[lo..mid]와 a[mid+1..hi]는 이미 정렬됨
void merge(std::vector<int>& a, int lo, int mid, int hi,
           std::vector<int>& tmp) {
    int i = lo, j = mid + 1, k = lo;
    while (i <= mid && j <= hi) {
        if (a[i] <= a[j]) tmp[k++] = a[i++];  // <= 로 안정성 유지
        else              tmp[k++] = a[j++];
    }
    while (i <= mid) tmp[k++] = a[i++];
    while (j <= hi)  tmp[k++] = a[j++];
    for (int idx = lo; idx <= hi; ++idx)
        a[idx] = tmp[idx];
}
void mergesort(std::vector<int>& a, int lo, int hi, std::vector<int>& tmp) {
    if (lo >= hi) return;
    int mid = lo + (hi - lo) / 2;
    mergesort(a, lo, mid, tmp);
    mergesort(a, mid + 1, hi, tmp);
    merge(a, lo, mid, hi, tmp);
}
void mergesort(std::vector<int>& a) {
    if (a.empty()) return;
    std::vector<int> tmp(a.size());
    mergesort(a, 0, static_cast<int>(a.size()) - 1, tmp);
}

merge에서 a[i] <= a[j]의 =가 안정성을 결정합니다. 값이 같을 때 왼쪽(원래 앞에 있던) 원소를 먼저 가져가야 같은 키를 가진 원소들의 원래 순서가 유지됩니다. <로 바꾸면 정렬 결과는 여전히 맞지만 안정성이 깨지는데, int만 정렬할 때는 차이가 보이지 않다가 레코드를 타임스탬프로 정렬했을 때 같은 시각 로그의 순서가 뒤섞이는 식으로 나중에 드러납니다. 표준 라이브러리에서 안정성이 필요하면 std::sort가 아니라 std::stable_sort를 써야 하는 이유도 같습니다.

tmp를 바깥에서 한 번만 할당해 모든 재귀 호출이 공유하는 것도 중요한 설계입니다. 병합할 때마다 std::vector를 새로 만들면 할당이 n번 가까이 일어납니다. 이 구현은 tmp의 [lo, hi] 구간만 쓰므로 서로 다른 구간을 처리하는 호출끼리 충돌하지 않고, 뒤에서 볼 병렬 버전도 이 성질 덕분에 tmp를 공유할 수 있습니다. 같은 이유로 tmp를 0번 인덱스부터 채우고(k = 0) 복사할 때는 tmp[idx]로 읽는 식으로 인덱스 기준을 섞으면 lo가 0이 아닌 구간에서 결과가 틀어집니다.

작은 구간에서는 재귀 호출 비용이 실제 작업보다 커지므로, 구간 길이가 일정 값(흔히 16 안팎) 미만이면 삽입 정렬로 처리하는 컷오프를 두는 구현이 많습니다. 적절한 값은 데이터와 하드웨어에 따라 다르므로 측정해서 정합니다.

퀵소트 (Quick Sort)

피벗을 골라 피벗보다 작은 원소는 왼쪽, 큰 원소는 오른쪽에 오도록 분할한 뒤 각 부분을 재귀적으로 정렬합니다. 제자리 정렬이며 평균 O(n log n)입니다.

#include <vector>
#include <algorithm>
// Hoare 분할: 중간 원소 피벗
int partitionHoare(std::vector<int>& a, int lo, int hi) {
    int pivot = a[lo + (hi - lo) / 2];
    int i = lo - 1, j = hi + 1;
    while (true) {
        do ++i; while (a[i] < pivot);
        do --j; while (a[j] > pivot);
        if (i >= j) return j;
        std::swap(a[i], a[j]);
    }
}
void quicksort(std::vector<int>& a, int lo, int hi) {
    if (lo >= hi) return;
    int p = partitionHoare(a, lo, hi);
    quicksort(a, lo, p);
    quicksort(a, p + 1, hi);
}
void quicksort(std::vector<int>& a) {
    if (a.size() <= 1) return;
    quicksort(a, 0, static_cast<int>(a.size()) - 1);
}

Hoare 분할은 Lomuto 분할(마지막 원소를 피벗으로 두고 한 방향으로 훑는 방식)보다 교환 횟수가 적고, 중복 값이 많아도 양쪽이 고르게 나뉜다는 장점이 있습니다. 대신 반환값 j가 피벗의 최종 위치가 아니라 왼쪽 구간의 끝이라서, 재귀를 (lo, p)와 (p + 1, hi)로 나눠야 합니다. Lomuto 방식처럼 (lo, p - 1)과 (p + 1, hi)로 쓰면 원소 하나가 정렬에서 빠지고, 반대로 Lomuto 분할의 반환값으로 (lo, p)를 재귀하면 피벗이 최댓값일 때 구간이 줄지 않아 무한 재귀가 됩니다. 두 방식의 코드 조각을 섞어 쓰면 이런 버그가 생기기 쉽습니다.

중간 원소 피벗은 정렬된 입력에서의 최악 케이스는 피하지만, 최악의 입력이 사라진 것은 아닙니다. 피벗 선택 규칙이 고정되어 있으면 그 규칙을 역이용한 입력(“median-of-3 killer” 같은)을 만들 수 있고, 외부 입력을 정렬하는 서버라면 이것이 서비스 거부 공격 수단이 됩니다. 주요 표준 라이브러리의 std::sort가 순수 퀵소트 대신 인트로소트를 쓰는 이유가 여기에 있습니다. 재귀 깊이가 한도(libstdc++는 2·log₂n)를 넘으면 힙 정렬로 전환해 최악에도 O(n log n)을 보장합니다.

세 원소의 중앙값을 피벗으로 쓰는 방법도 흔합니다. 아래 함수는 a[lo], a[mid], a[hi]를 정렬해 중앙값을 mid에 둡니다.

int medianOfThree(std::vector<int>& a, int lo, int hi) {
    int mid = lo + (hi - lo) / 2;
    if (a[lo] > a[mid]) std::swap(a[lo], a[mid]);
    if (a[lo] > a[hi]) std::swap(a[lo], a[hi]);
    if (a[mid] > a[hi]) std::swap(a[mid], a[hi]);
    return mid;
}

재귀 대신 명시적 스택을 쓰면 호출 스택 오버플로우는 사라지지만, 스택에 쌓이는 구간 수는 순서를 신경 쓰지 않으면 최악의 경우 O(n)이 될 수 있습니다. 두 구간 중 큰 쪽을 먼저 push하고 작은 쪽을 나중에 push해서 작은 쪽을 먼저 처리하면, 스택에 남는 구간 수가 O(log n)으로 제한됩니다.

#include <stack>
void quicksortIterative(std::vector<int>& a) {
    if (a.size() <= 1) return;
    std::stack<std::pair<int,int>> st;
    st.push({0, static_cast<int>(a.size()) - 1});
    while (!st.empty()) {
        auto [lo, hi] = st.top();
        st.pop();
        if (lo >= hi) continue;
        int p = partitionHoare(a, lo, hi);
        // 큰 쪽을 먼저 push → 작은 쪽이 먼저 pop되어 처리됨
        if (p - lo > hi - (p + 1)) {
            st.push({lo, p});
            st.push({p + 1, hi});
        } else {
            st.push({p + 1, hi});
            st.push({lo, p});
        }
    }
}

정렬된 배열에서 중간값과 비교해 범위를 절반씩 좁혀 갑니다.

#include <vector>
// 반복 버전
int binarySearch(const std::vector<int>& a, int target) {
    int lo = 0, hi = static_cast<int>(a.size()) - 1;
    while (lo <= hi) {
        int mid = lo + (hi - lo) / 2;  // 오버플로우 방지
        if (a[mid] == target) return mid;
        if (a[mid] < target) lo = mid + 1;
        else                 hi = mid - 1;
    }
    return -1;
}
// 재귀 버전
int binarySearchRecursive(const std::vector<int>& a, int target, int lo, int hi) {
    if (lo > hi) return -1;
    int mid = lo + (hi - lo) / 2;
    if (a[mid] == target) return mid;
    if (a[mid] < target) return binarySearchRecursive(a, target, mid + 1, hi);
    return binarySearchRecursive(a, target, lo, mid - 1);
}

이 구현은 폐구간 [lo, hi]를 불변식으로 씁니다. while (lo <= hi)와 hi = mid - 1이 짝을 이루는 이유가 그것이고, 반열린 구간 [lo, hi)로 바꾸려면 초기값 hi = size, 조건 lo < hi, 갱신 hi = mid를 함께 바꿔야 합니다. 셋 중 하나만 바꾸면 경계 원소를 놓치거나 무한 루프에 빠집니다. 예를 들어 while (lo < hi)에서 lo = mid로 갱신하면 lo = 3, hi = 4일 때 mid가 계속 3이 되어 범위가 줄지 않습니다. 또 hi를 int로 두는 것에도 이유가 있습니다. size_t로 두면 hi = mid - 1에서 mid가 0일 때 언더플로우가 일어나 아주 큰 값이 되고, 빈 배열에서는 size() - 1 자체가 최댓값이 되어 곧바로 범위 밖을 읽습니다. (lo + hi) / 2의 오버플로우는 무한 루프가 아니라 음수 인덱스 접근으로 나타납니다.

실무에서는 값이 있는지보다 처음으로 target 이상이 되는 위치가 필요한 경우가 많고, 중복 값이 있을 때 위 코드는 여러 위치 중 아무거나 돌려줍니다. 이럴 때는 직접 구현하기보다 std::lower_bound/std::upper_bound/std::equal_range를 쓰는 편이 안전합니다.

#include <optional>
#include <algorithm>
#include <cassert>
std::optional<int> findIndex(const std::vector<int>& a, int target) {
    assert(std::is_sorted(a.begin(), a.end()));  // 디버그 빌드에서만 검사
    auto it = std::lower_bound(a.begin(), a.end(), target);
    if (it == a.end() || *it != target) return std::nullopt;
    return static_cast<int>(std::distance(a.begin(), it));
}

std::is_sorted는 O(n)이라 매 호출마다 실행하면 이진 탐색의 O(log n) 이점이 사라집니다. 그래서 assert로 디버그 빌드에서만 확인하고, 릴리스에서는 정렬된 입력을 함수의 사전 조건으로 문서화하는 것이 일반적인 절충입니다.

최근접 점 쌍 (Closest Pair of Points)

점을 x좌표로 정렬한 뒤 중간에서 나눕니다. 왼쪽의 최소 거리 dL과 오른쪽의 최소 거리 dR을 재귀로 구하고, d = min(dL, dR)일 때 경계선에서 x 거리가 d 미만인 중앙 밴드 안의 점들만 y순으로 검사합니다.

#include <vector>
#include <cmath>
#include <algorithm>
#include <limits>
struct Point {
    double x, y;
};
double dist(const Point& a, const Point& b) {
    return std::hypot(a.x - b.x, a.y - b.y);
}
// 밴드 내 점들만 검사 - y좌표 정렬된 상태
double closestInBand(const std::vector<Point>& band, double d) {
    double best = d;
    for (size_t i = 0; i < band.size(); ++i) {
        for (size_t j = i + 1; j < band.size() && (band[j].y - band[i].y) < best; ++j) {
            best = std::min(best, dist(band[i], band[j]));
        }
    }
    return best;
}
double closestPairRec(const std::vector<Point>& pts, int lo, int hi) {
    if (hi - lo <= 2) {
        double best = std::numeric_limits<double>::max();
        for (int i = lo; i <= hi; ++i)
            for (int j = i + 1; j <= hi; ++j)
                best = std::min(best, dist(pts[i], pts[j]));
        return best;
    }
    int mid = lo + (hi - lo) / 2;
    double midX = pts[mid].x;
    double dL = closestPairRec(pts, lo, mid);
    double dR = closestPairRec(pts, mid + 1, hi);
    double d = std::min(dL, dR);
    // 밴드: 현재 구간 [lo, hi] 안에서 |x - midX| < d 인 점들만
    std::vector<Point> band;
    for (int i = lo; i <= hi; ++i) {
        if (std::abs(pts[i].x - midX) < d) band.push_back(pts[i]);
    }
    std::sort(band.begin(), band.end(), [](const Point& a, const Point& b) {
        return a.y < b.y;
    });
    return std::min(d, closestInBand(band, d));
}
double closestPair(std::vector<Point> pts) {
    if (pts.size() < 2) return std::numeric_limits<double>::max();
    std::sort(pts.begin(), pts.end(), [](const Point& a, const Point& b) {
        return a.x < b.x;
    });
    return closestPairRec(pts, 0, static_cast<int>(pts.size()) - 1);
}

이 알고리즘의 핵심은 closestInBand의 안쪽 루프가 상수 번만 돈다는 사실입니다. 밴드 안의 점을 y순으로 정렬해 두면, 어떤 점에서 y 차이가 d 미만인 이웃은 d×2d 크기의 직사각형 안에 있어야 합니다. 그런데 왼쪽과 오른쪽 각각에서는 이미 모든 점 쌍의 거리가 d 이상이므로, 그 직사각형 안에 들어갈 수 있는 점의 수에는 기하학적인 상한이 있습니다(교과서에서는 7개 이하로 보입니다). 그래서 이중 루프처럼 보여도 밴드 검사는 선형 시간입니다. (band[j].y - band[i].y) < best 조건이 빠지면 이 논증이 무너져 O(n²)으로 돌아갑니다.

위 구현은 이해하기 쉽도록 재귀 단계마다 밴드를 y로 다시 정렬하므로 전체 복잡도가 O(n log² n)입니다. 교과서의 O(n log n) 버전은 병합정렬처럼 재귀가 돌아오면서 y순 정렬 결과를 병합해 올려 보내 이 정렬 비용을 없앱니다. 흔한 실수는 전체 점을 y로 한 번 정렬한 배열을 모든 재귀 단계에서 통째로 훑는 것인데, 그러면 매 호출이 n개를 보게 되어 O(n²)이 됩니다. 결과는 맞게 나오기 때문에 작은 테스트에서는 이 성능 버그를 발견하기 어렵습니다. 또 거리 비교만 필요하다면 hypot 대신 제곱 거리를 비교하면 제곱근 계산을 피할 수 있습니다.

Strassen 행렬 곱셈 (Matrix Multiplication)

2×2 블록 행렬 곱셈은 블록 곱셈 8번이 필요하지만, Strassen은 7번의 곱셈으로 같은 결과를 얻습니다. 큰 행렬을 2×2 블록으로 나눠 이 공식을 재귀적으로 적용합니다. A, B를 네 개의 부분행렬로 나누면 다음과 같습니다.

  • M1 = (A11 + A22)(B11 + B22)
  • M2 = (A21 + A22)B11
  • M3 = A11(B12 - B22)
  • M4 = A22(B21 - B11)
  • M5 = (A11 + A12)B22
  • M6 = (A21 - A11)(B11 + B12)
  • M7 = (A12 - A22)(B21 + B22)

결과: C11 = M1+M4-M5+M7, C12 = M3+M5, C21 = M2+M4, C22 = M1-M2+M3+M6

검산 하나만 해 보면 C12 = M3 + M5 = A11·B12 − A11·B22 + A11·B22 + A12·B22 = A11·B12 + A12·B22로, 일반 행렬 곱의 정의와 일치합니다. 곱셈을 8번에서 7번으로 줄인 대가로 덧셈·뺄셈이 18번으로 늘었는데, 재귀에서는 곱셈 한 번이 곧 크기 n/2 문제 하나이므로 T(n) = 7T(n/2) + O(n²) = O(n^log₂7) ≈ O(n^2.807)이 됩니다. 공식을 옮기다 부호 하나만 틀려도 컴파일은 되고 결과만 틀리므로, 구현 후에는 일반 곱셈 결과와 원소별로 비교하는 테스트를 붙여야 합니다.

#include <vector>
#include <algorithm>
#include <stdexcept>
using Matrix = std::vector<std::vector<double>>;
// 행렬 덧셈
Matrix matrixAdd(const Matrix& a, const Matrix& b) {
    size_t n = a.size();
    Matrix c(n, std::vector<double>(n));
    for (size_t i = 0; i < n; ++i)
        for (size_t j = 0; j < n; ++j)
            c[i][j] = a[i][j] + b[i][j];
    return c;
}
// 행렬 뺄셈
Matrix matrixSub(const Matrix& a, const Matrix& b) {
    size_t n = a.size();
    Matrix c(n, std::vector<double>(n));
    for (size_t i = 0; i < n; ++i)
        for (size_t j = 0; j < n; ++j)
            c[i][j] = a[i][j] - b[i][j];
    return c;
}
// 일반 O(n³) 행렬 곱셈 (기본 케이스)
Matrix matrixMulNaive(const Matrix& a, const Matrix& b) {
    size_t n = a.size();
    Matrix c(n, std::vector<double>(n, 0));
    for (size_t i = 0; i < n; ++i)
        for (size_t j = 0; j < n; ++j)
            for (size_t k = 0; k < n; ++k)
                c[i][j] += a[i][k] * b[k][j];
    return c;
}
// 부분행렬 추출: A에서 [r1:r2)[c1:c2) 구간
Matrix submatrix(const Matrix& a, size_t r1, size_t r2, size_t c1, size_t c2) {
    Matrix sub(r2 - r1, std::vector<double>(c2 - c1));
    for (size_t i = r1; i < r2; ++i)
        for (size_t j = c1; j < c2; ++j)
            sub[i - r1][j - c1] = a[i][j];
    return sub;
}
// 4개 부분행렬을 하나로 합침
Matrix combineQuarters(const Matrix& c11, const Matrix& c12,
                       const Matrix& c21, const Matrix& c22) {
    size_t n = c11.size() * 2;
    Matrix c(n, std::vector<double>(n));
    size_t m = n / 2;
    for (size_t i = 0; i < m; ++i) {
        for (size_t j = 0; j < m; ++j) {
            c[i][j] = c11[i][j];
            c[i][j + m] = c12[i][j];
            c[i + m][j] = c21[i][j];
            c[i + m][j + m] = c22[i][j];
        }
    }
    return c;
}
// Strassen: 7번의 곱셈으로 n×n 행렬 곱셈 (n은 2의 거듭제곱)
Matrix strassen(const Matrix& a, const Matrix& b) {
    size_t n = a.size();
    if (n != b.size() || (n & (n - 1)) != 0)
        throw std::invalid_argument("Strassen: 정방행렬, 2의 거듭제곱 크기 필요");
    if (n <= 64)  // 작은 행렬은 일반 곱셈으로 전환
        return matrixMulNaive(a, b);
    size_t m = n / 2;
    Matrix a11 = submatrix(a, 0, m, 0, m);
    Matrix a12 = submatrix(a, 0, m, m, n);
    Matrix a21 = submatrix(a, m, n, 0, m);
    Matrix a22 = submatrix(a, m, n, m, n);
    Matrix b11 = submatrix(b, 0, m, 0, m);
    Matrix b12 = submatrix(b, 0, m, m, n);
    Matrix b21 = submatrix(b, m, n, 0, m);
    Matrix b22 = submatrix(b, m, n, m, n);
    Matrix m1 = strassen(matrixAdd(a11, a22), matrixAdd(b11, b22));
    Matrix m2 = strassen(matrixAdd(a21, a22), b11);
    Matrix m3 = strassen(a11, matrixSub(b12, b22));
    Matrix m4 = strassen(a22, matrixSub(b21, b11));
    Matrix m5 = strassen(matrixAdd(a11, a12), b22);
    Matrix m6 = strassen(matrixSub(a21, a11), matrixAdd(b11, b12));
    Matrix m7 = strassen(matrixSub(a12, a22), matrixAdd(b21, b22));
    Matrix c11 = matrixAdd(matrixSub(matrixAdd(m1, m4), m5), m7);
    Matrix c12 = matrixAdd(m3, m5);
    Matrix c21 = matrixAdd(m2, m4);
    Matrix c22 = matrixAdd(matrixAdd(matrixSub(m1, m2), m3), m6);
    return combineQuarters(c11, c12, c21, c22);
}
// 2의 거듭제곱이 아닌 행렬용: 0으로 패딩 후 Strassen 적용
Matrix strassenPadded(Matrix a, Matrix b) {
    size_t n = std::max({a.size(), a[0].size(), b.size(), b[0].size()});
    size_t N = 1;
    while (N < n) N *= 2;
    Matrix A(N, std::vector<double>(N, 0));
    Matrix B(N, std::vector<double>(N, 0));
    for (size_t i = 0; i < a.size(); ++i)
        for (size_t j = 0; j < a[0].size(); ++j)
            A[i][j] = a[i][j];
    for (size_t i = 0; i < b.size(); ++i)
        for (size_t j = 0; j < b[0].size(); ++j)
            B[i][j] = b[i][j];
    Matrix C = strassen(A, B);
    Matrix result(a.size(), std::vector<double>(b[0].size()));
    for (size_t i = 0; i < result.size(); ++i)
        for (size_t j = 0; j < result[0].size(); ++j)
            result[i][j] = C[i][j];
    return result;
}

strassen은 2×2 블록 분할을 전제로 하므로 크기가 2의 거듭제곱이 아니면 예외를 던지고, 3×3 같은 행렬은 strassenPadded로 4×4까지 0을 채워 계산합니다. 작은 n에서는 덧셈과 임시 행렬 할당 비용이 곱셈 절감분보다 커서 일반 O(n³) 곱셈이 더 빠르므로, 위 예제처럼 일정 크기(여기서는 64) 이하에서는 일반 곱셈으로 전환합니다. 적절한 전환점은 하드웨어와 구현에 따라 달라지므로 측정해서 정해야 합니다.

이 구현은 알고리즘 구조를 보여 주기 위한 것이고, 실제 성능을 기대할 코드는 아닙니다. std::vector<std::vector<double>>는 행마다 메모리가 흩어져 캐시 효율이 나쁘고, 재귀 한 단계마다 부분행렬 8개와 중간 결과 여러 개를 새로 할당합니다. 캐시 블로킹과 SIMD를 적용한 O(n³) 곱셈(OpenBLAS, Intel MKL, Eigen 등)이 실용적인 크기 대부분에서 이런 Strassen 구현보다 빠르고, 주요 BLAS 라이브러리도 기본 경로에서는 Strassen을 쓰지 않습니다. Strassen은 뺄셈이 섞여 있어 부동소수점 오차가 일반 곱셈보다 커지는 수치 안정성 문제도 있습니다. 행렬 곱이 병목이라면 먼저 BLAS 라이브러리를 쓰고, Strassen은 분할정복 설계 기법을 이해하는 예제로 보는 것이 현실적입니다.

최대 부분 배열 합 (Maximum Subarray)

중간을 기준으로 왼쪽 최대, 오른쪽 최대, 중간을 지나는 최대를 각각 구하고 그중 가장 큰 값을 고릅니다.

#include <vector>
#include <algorithm>
#include <limits>
// 중간을 지나는 최대 부분합 (a[mid]와 a[mid+1]을 반드시 포함)
int maxCrossingSum(const std::vector<int>& a, int lo, int mid, int hi) {
    int leftSum = 0, leftMax = std::numeric_limits<int>::min();
    for (int i = mid; i >= lo; --i) {
        leftSum += a[i];
        leftMax = std::max(leftMax, leftSum);
    }
    int rightSum = 0, rightMax = std::numeric_limits<int>::min();
    for (int i = mid + 1; i <= hi; ++i) {
        rightSum += a[i];
        rightMax = std::max(rightMax, rightSum);
    }
    return leftMax + rightMax;
}
int maxSubarrayDC(const std::vector<int>& a, int lo, int hi) {
    if (lo == hi) return a[lo];
    int mid = lo + (hi - lo) / 2;
    int leftMax  = maxSubarrayDC(a, lo, mid);
    int rightMax = maxSubarrayDC(a, mid + 1, hi);
    int crossMax = maxCrossingSum(a, lo, mid, hi);
    return std::max({leftMax, rightMax, crossMax});
}
int maxSubarray(const std::vector<int>& a) {
    if (a.empty()) return 0;
    return maxSubarrayDC(a, 0, static_cast<int>(a.size()) - 1);
}

leftMax와 rightMax의 초기값이 이 코드의 함정입니다. 흔히 0으로 초기화하는데, 그러면 중간을 지나는 부분 배열이 빈 구간을 허용하게 되어 모든 원소가 음수인 입력(예: {-3, -1})에서 실제 답 -1 대신 0을 돌려줍니다. 중간을 지나는 부분 배열은 a[mid]와 a[mid+1]을 반드시 포함해야 하므로 위처럼 최솟값으로 시작해야 합니다. 루프가 각각 최소 한 번은 돌기 때문에 두 값은 모두 실제 부분합으로 바뀌어 있고, INT_MIN끼리 더해지는 일은 없습니다(부분합 자체가 int 범위를 넘는 입력이라면 long long을 써야 합니다). 이 문제는 카데인(Kadane) 알고리즘으로 O(n)에 풀리므로, 분할정복 버전(O(n log n))은 실전용이라기보다 Combine 단계 설계를 연습하는 예제입니다. 왼쪽과 오른쪽 결과만 비교하고 crossMax를 빠뜨리는 것이 Combine 단계 누락의 전형적인 예로, 경계에 걸친 답을 놓쳐 결과가 틀립니다.

분할정복과 동적 계획법의 경계

분할정복은 하위 문제가 서로 겹치지 않을 때 효율적입니다. 병합정렬과 퀵소트의 하위 구간은 서로 겹치지 않습니다. 반면 같은 하위 문제가 여러 번 등장하면 순수 분할정복은 지수 시간이 됩니다.

// 피보나치: fib(n) = fib(n-1) + fib(n-2)
// 순수 재귀는 같은 fib(k)를 반복 계산해 호출 수가 지수적으로 증가
int fib(int n) {
    if (n <= 1) return n;
    return fib(n - 1) + fib(n - 2);
}
// 메모이제이션이나 상향식 DP로 바꾸면 O(n)

하위 문제가 겹친다면 결과를 저장해 재사용하는 동적 계획법으로 바꿔야 합니다. 피보나치와 LCS가 대표적인 예이고, 자세한 내용은 C++ 동적 계획법: 메모이제이션 vs 타뷸레이션, 2D→1D 공간 최적화, 배낭·LCS·Bitmask DP에서 다룹니다.

std::async로 병렬화하기

#include <future>
#include <algorithm>
void mergesortParallel(std::vector<int>& a, int lo, int hi, std::vector<int>& tmp) {
    if (hi - lo < 10000) {
        mergesort(a, lo, hi, tmp);
        return;
    }
    int mid = lo + (hi - lo) / 2;
    auto leftFuture = std::async(std::launch::async, [&]() {
        mergesortParallel(a, lo, mid, tmp);
    });
    mergesortParallel(a, mid + 1, hi, tmp);
    leftFuture.get();
    merge(a, lo, mid, hi, tmp);
}

이 병렬 버전이 a와 tmp를 두 스레드가 공유해도 데이터 경쟁이 없는 이유는, 왼쪽과 오른쪽 작업이 각각 [lo, mid]와 [mid+1, hi]의 서로 겹치지 않는 구간만 읽고 쓰기 때문입니다. merge는 get()으로 왼쪽 작업이 끝난 것을 확인한 뒤에야 두 구간을 모두 읽습니다. 주의할 점은 스레드 수입니다. std::launch::async는 호출마다 새 스레드를 만들 수 있어서, 임계값 10000으로 100만 개를 정렬하면 재귀 트리가 약 7단계까지 분기해 100개가 넘는 스레드가 생길 수 있습니다. 코어 수보다 훨씬 많은 스레드는 전환 비용만 늘리므로, 재귀 깊이로 병렬 분기를 제한(예: 깊이 log₂(코어 수)까지만 async)하거나 C++17 병렬 알고리즘 std::sort(std::execution::par, ...)를 쓰는 편이 낫습니다.

성능은 직접 측정하기

정렬 알고리즘의 실제 속도는 데이터 분포(무작위, 거의 정렬됨, 중복 많음), 원소 크기, 컴파일러 최적화, 캐시 크기에 따라 크게 달라지므로 자기 데이터로 측정해야 합니다. 일반적인 경향으로는 잘 다듬어진 std::sort(인트로소트 + 작은 구간 삽입 정렬)가 직접 짠 퀵소트·병합정렬보다 빠른 경우가 대부분이고, 안정성이 필요하면 std::stable_sort를 씁니다. 직접 구현이 의미 있는 경우는 병렬화나 외부 정렬처럼 표준 라이브러리가 제공하지 않는 구조가 필요할 때입니다. 측정할 때는 같은 입력을 매 반복마다 새로 복사해 정렬해야 합니다. 이미 정렬된 배열을 다시 정렬하면 퀵소트·삽입 정렬의 결과가 비현실적으로 좋게 나옵니다.

정확성 검증에는 무작위 입력을 만들어 직접 구현한 결과와 std::sort 결과를 비교하는 테스트가 효과적입니다.

#include <vector>
#include <random>
#include <cassert>
std::vector<int> randomVector(size_t n, int minVal = 0, int maxVal = 1000) {
    std::mt19937 gen(std::random_device{}());
    std::uniform_int_distribution<> dist(minVal, maxVal);
    std::vector<int> v(n);
    for (auto& x : v) x = dist(gen);
    return v;
}
void testMergesort() {
    for (size_t n : {0, 1, 2, 17, 1000}) {
        auto v = randomVector(n);
        auto expected = v;
        std::sort(expected.begin(), expected.end());
        mergesort(v);
        assert(v == expected);
    }
}

빈 배열, 원소 하나, 원소 둘처럼 base case 근처의 크기를 반드시 포함해야 경계 버그가 드러납니다. 값 범위를 좁게(예: 0~10) 잡으면 중복이 많은 입력도 함께 시험할 수 있습니다.

참고 자료

같이 보면 좋은 글