まず「\(S\) から \(T\) への観光ルートは何個あるか」という問題を考えます。
この問題は、巡回セールスマン問題と同様のbitDPにより求めることができます。

元の問題を考えます。観光ルートのスコアの平均を求めるためには、観光ルートのスコアの総和と観光ルートの個数が分かれば十分です。後者は冒頭で述べた問題そのものであり、前者は個数の情報も持ちながら同様のDPを行うことで求めることができます。

具体的には次の通りです。

\(\mathrm{CNT}[v][ss]\) 及び \(\mathrm{SUM}[v][ss]\) を「『スポット \(S\) からスタートし、スポット \(v\) にいて、訪問済みのスポットの集合が \(ss\) であるような単純パス』全ての個数・スコアの総和」と定めます。

このとき、\(v\) の直前に居たスポットがどこであるかを考えることで

\(\mathrm{CNT}[v][ss] = \sum_{u \in ss かつ uv間に辺がある}\mathrm{CNT}[u][ss\setminus\{v\}]\)

\(\mathrm{SUM}[v][ss] =\mathrm{CNT}[v][ss]\times C[v]+ \sum_{u \in ss かつ uv間に辺がある}\mathrm{SUM}[u][ss\setminus\{v\}]\)

という漸化式を得ることができます。

以上より \(O(2^NN^2)\) でこの問題を解くことができました。制約上、登場する数は全て 64 bit符号付き整数の範囲に収まることが示せます。

以下の実装例では、上で説明したもらうDPではなく、配るDPにより実装が行われています。

実装例 (C++)

#include<bits/stdc++.h>
using namespace std;

int main(){
  int n, m, s, t;
  cin >> n >> m >> s >> t;
  s--, t--;
  vector<int>c(n);
  for(int i=0; i<n; i++) cin >> c[i];
  vector<vector<bool>>e(n, vector<bool>(n));
  for(int i=0; i<m; i++){
    int u, v;
    cin >> u >> v;
    u--, v--;
    e[u][v] = e[v][u] = true;
  }

  vector<vector<long long>>sum(n, vector<long long>(1<<n));
  vector<vector<long long>>cnt(n, vector<long long>(1<<n));
  sum[s][1<<s] = c[s];
  cnt[s][1<<s] = 1;

  for(int ss=0; ss<1<<n; ss++){
    for(int k=0; k<n; k++){
      for(int kk=0; kk<n; kk++){
        if( !(ss&(1<<kk)) && e[k][kk]){
          cnt[kk][ss^(1<<kk)] += cnt[k][ss];
          sum[kk][ss^(1<<kk)] += sum[k][ss] + cnt[k][ss] * c[kk];
        }
      }
    }
  }
  long long ans_cnt = 0, ans_sum = 0;
  for(int ss=0; ss<1<<n; ss++){
      ans_cnt += cnt[t][ss];
      ans_sum += sum[t][ss];
  }
  printf("%.10f\n", (double) ans_sum / ans_cnt);
}

実装例 (Python)

N, M, S, T = map(int, input().split())
S -= 1
T -= 1
C = [int(input()) for _ in range(N)]
e = [[False]*N for _ in range(N)]
for _ in range(M):
  u, v = map(int, input().split())
  u -= 1
  v -= 1
  e[u][v] = True
  e[v][u] = True

SUM = [[0]*(1<<N) for _ in range(N)]
CNT = [[0]*(1<<N) for _ in range(N)]
SUM[S][1<<S] = C[S]
CNT[S][1<<S] = 1

for s in range(1<<N):
  for k in range(N):
    for kk in range(N):
      if not s&(1<<kk) and e[k][kk]:
        CNT[kk][s^(1<<kk)] += CNT[k][s]
        SUM[kk][s^(1<<kk)] += SUM[k][s] + CNT[k][s] * C[kk]

ans_cnt = sum(CNT[T])
ans_sum = sum(SUM[T])

print(ans_sum / ans_cnt)