Submission #102510

---

Source Code Expand

#include<cstdio>
#include<vector>
#include<iostream>
#include<algorithm>
#include<cassert>
#include<cstdlib>
#include<cmath>
#include<cstdio>
#include<cstring>
#include<complex>

using namespace std;
typedef long long LL;
//typedef double LDouble;
//typedef long double LDouble;
typedef __float128 LDouble;
//typedef long long double LDouble;

//typedef complex<LDouble> Complex;


const double INF = 1e12;

typedef vector<LDouble> Poly;

// a*b where x^n = c
Poly mulPoly(const Poly&a, const Poly&b, const Poly&c) {
    int n=a.size();
    Poly r(n*2);
    for (int i=0; i<n; i++) {
        for (int j=0; j<n; j++) {
            r[i+j] = (r[i+j]+a[i]*b[j]);
        }
    }
    for (int i=n*2; i-->n; ) {
        for (int j=1; j<=n; j++) {
            r[i-j] = (r[i-j]+r[i]*c[n-j]);
        }
    }
    r.resize(n);
    //for (int i=0; i<n; i++) assert(abs((double)r[i]) < INF);
    return r;
}

Poly powPoly(Poly a, LL y, const Poly&c) {
    int n = a.size();
    Poly r(n);
    r[0]=1;
    for (;y; y>>=1) {
        if (y&1) r = mulPoly(r, a, c);
        a = mulPoly(a, a, c);
    }
    return r;
}


int main() {
    //printf("double %d\n__float128 %d\nlong double %d\n", sizeof(double), sizeof(__float128), sizeof(long double));
    double _p, _q, _r;
    LDouble p, q, r;
    LL N, K, M;
    cin >> _p >> _q >> _r;
    cin >> N >> K >> M;
    p=_p; q=_q; r=_r;

    vector<LDouble> v(N*K+1);
    LDouble base = p+q+r;
    v[0] += -r / base;
    v[1] += r / base;
    v[N*(K-1)] += -p /base;
    v[N*(K-1) + 1] += p / base;
    v[N*K-1] += -q / base;
    v[N*K] += (p+2*q+r) / base;


    vector<LDouble>alpha(N*K+1);
    for (int i=1; i<=N*K; i++) {
        alpha[i] = (i>=N ? p*alpha[i-N] : 0) 
            + (i>=1 ? q*alpha[i-1] : 0)
            + (i>=N*K ? r*alpha[i-N*K] : 0)
            + 1;
        alpha[i] /= base;
    }

    //for (int i=0; i<N*K+1; i++) { printf("%d %f\n", i, double(alpha[i])); } puts("");
    Poly beta(N*K+1);
    beta[1]=1;
    beta = powPoly(beta, M*N, v);
    //for (int i=0; i<N*K+1; i++) { printf("%d %f\n", i, double(beta[i])); } puts("");

    LDouble ans=0;
    for (int i=0; i<N*K+1; i++) ans += beta[i] * alpha[i];

    printf("%.9f\n", double(ans));
    //cout << ans << endl;

    return 0;
}

Submission Info

Judge Result