백준 1538 공 칠하기

https://www.acmicpc.net/problem/1538

1538번: 공 칠하기

순서를 고려하지 않고 $N$을 자연수의 합으로 분할하는 방법 한 가지는 가방 속 공의 상태와 대응된다. 이 '상태'를 적절히 표현해보자. 예를 들어 $N=5$일 때, 아래와 같은 7가지 상태가 가능하다.

$$(5), (3, 2), (4, 1), (2, 2, 1), (3, 1, 1), (2, 1, 1, 1), (1, 1, 1, 1, 1)$$

어떤 상태 $\sigma$에서 시작했을 때의 기댓값을 $D_\sigma$라고 하면, 이들의 관계식을 세울 수 있다. 예를 들어 $D_{(3, 1, 1)}$은 아래와 같다.

$$D_{(3, 1, 1)} = {6 \over 20}D_{(3, 1, 1)}+{3 \over 20}D_{(2, 2, 1)}+{3 \over 20}D_{(2, 2, 1)}+{3 \over 20}D_{(4, 1)}+{3 \over 20}D_{(4, 1)}+{2 \over 20}D_{(3, 2)}+1$$

이를 정리하면 다음과 같은 관계식을 얻는다.

$${7 \over 10}D_{(3, 1, 1)}-{3 \over 10}D_{(2, 2, 1)}-{3 \over 10}D_{(4, 1)}-{1 \over 10}D_{(3, 2)} = 1$$

모든 상태에 대해 이런 관계식을 세울 수 있으므로 상태의 총 개수를 $p(N)$이라 하면 식 $p(N)$개, 미지수 $p(N)$개짜리 연립방정식을 얻는다. 그러나 한 가지 문제는 $p(24)>1500$이라는 건데, 가우스 소거법은 $O(N^{3})$의 시간복잡도를 가지므로 시간 내에 통과하지 못할 가능성이 있다.

한 가지 관찰을 추가하자. 자명한 방법으로 상태의 길이를 정의할 수 있는데, 이때 $D_\sigma$의 관계식에는 $\sigma$보다 길이가 작거나 같은 상태들만이 나온다. 즉 길이가 1인 상태부터 시작해 상태의 길이를 기준으로 순차적으로 구하는 것이 가능하다. $N=5$일 때의 예시를 들면 다음과 같다. (각 관계식을 유도하는 과정은 생략한다.)

$$D_{(5)} = 0$$

$$\rightarrow {3 \over 10}D_{(3,2)}-{3 \over 10}D_{(4, 1)}=1, {2 \over 5}D_{(4, 1)}-{1 \over 5}D_{(3, 2)}={1 \over 5}D_{(5)}+1$$

$$\rightarrow {3 \over 5}D_{(2, 2, 1)}-{2 \over 5}D_{(3, 1, 1)}={1 \over 5}D_{(3, 2)}+1, {7 \over 10}D_{(3, 1, 1)}-{3 \over 10}D_{(2, 2, 1)}={1 \over 10}D_{(3, 2)}+{3 \over 10}D_{(4, 1)}+1$$

$$\rightarrow {3 \over 5}D_{(2, 1, 1, 1)}={3 \over 10}D_{(2, 2, 1)}+{3 \over 10}D_{(3, 1, 1)}+1$$

$$\rightarrow D_{(1, 1, 1, 1, 1)}=D_{(2, 1, 1, 1)}+1$$

식이 복잡해보일 수 있지만, 등호를 기준으로 왼쪽에는 해당 단계에서 구해야 할 항들이, 오른쪽에는 이미 구해 놓은 항들이 위치한다는 사실만 확인하면 된다. 각 단계에서 식의 개수와 미지수의 개수가 같으므로, 가우스 소거법을 이용하면 모든 값을 구할 수 있다.

 

코드는 다음과 같다.

#include<stdio.h>
#include<string.h>
#include<vector>
#include<algorithm>
#include<map>
using namespace std;
typedef pair<int, int> pi;

int n;
char s[30];
int tmp[30];
vector<int> all[2000];
int t = 0;
map<vector<int>, int> mp;
double dp[2000];
double mat[2000][2000];
pi ran[30];

void gen(int a, int prv, int idx){
    if(a==0){
        for(int i=0;i<idx;i++) all[t].push_back(tmp[i]);
        ++t;
        return;
    }
    for(int i=1;i<=a&&i<=prv;i++){
        tmp[idx] = i;
        gen(a-i, i, idx+1);
    }
}

bool cmp1(vector<int> a, vector<int> b){
    if(a.size()==b.size()) return a<b;
    return a.size()<b.size();
}

bool cmp2(int a, int b){
    return a>b;
}

int main(){
    scanf("%d", &n);
    scanf("%s", &s);
    gen(n, n, 0);
    sort(all, all+t, cmp1);
    ran[1] = pi(0, 0);
    for(int i=1;i<t;i++){
        if(all[i-1].size()!=all[i].size()) ran[all[i].size()].first = i;
        if(all[i].size()!=all[i+1].size()) ran[all[i].size()].second = i;
    }
    for(int i=0;i<t;i++){
        while(all[i].size()<n) all[i].push_back(0);
    }
    for(int i=0;i<t;i++) mp[all[i]] = i;
    vector<int> vec;
    for(int sze=2;sze<=n;sze++){
        int l = ran[sze].second-ran[sze].first+1;
        for(int i=0;i<l;i++){
            for(int j=0;j<=l;j++) mat[i][j] = 0;
        }
        for(int i=ran[sze].first;i<=ran[sze].second;i++){
            int id = i-ran[sze].first;
            mat[id][id] = 1;
            mat[id][l] = 1;
            for(int u=0;u<sze;u++){
                for(int v=0;v<sze;v++){
                    vec = all[i];
                    double x = (double)vec[u]*vec[v]/((double)(n)*(n-1));
                    if(u==v) x = (double)vec[u]*(vec[v]-1)/((double)(n)*(n-1));
                    vec[u]++;
                    vec[v]--;
                    sort(vec.begin(), vec.end(), cmp2);
                    int o = mp[vec];
                    if(o<ran[sze].first) mat[id][l] += dp[o]*x;
                    else mat[id][o-ran[sze].first] -= x;
                }
            }
        }
        for(int i=0;i<l;i++){
            for(int j=0;j<l;j++){
                if(i==j) continue;
                double q = mat[j][i]/mat[i][i];
                for(int k=0;k<=l;k++) mat[j][k] -= mat[i][k]*q;
            }
        }
        for(int i=0;i<l;i++){
            double p = mat[i][i];
            for(int j=0;j<=l;j++) mat[i][j] /= p;
        }
        for(int i=0;i<l;i++) dp[i+ran[sze].first] = mat[i][l];
    }
    vector<int> inp;
    for(int i='A';i<='Z';i++) inp.push_back(0);
    for(int i=0;i<n;i++) inp[s[i]-'A']++;
    sort(inp.begin(), inp.end(), cmp2);
    while(!inp.empty()&&inp.size()>n) inp.pop_back();
    printf("%.9lf", dp[mp[inp]]);
}