#ifndef _GLIBCXX_NO_ASSERT
#include <cassert>
#endif
#include <cctype>
#include <cerrno>
#include <cfloat>
#include <ciso646>
#include <climits>
#include <clocale>
#include <cmath>
#include <csetjmp>
#include <csignal>
#include <cstdarg>
#include <cstddef>
#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <ctime>

#if __cplusplus >= 201103L
#include <ccomplex>
#include <cfenv>
#include <cinttypes>
#include <cstdbool>
#include <cstdint>
#include <ctgmath>
#include <cwchar>
#include <cwctype>
#endif

// C++
#include <algorithm>
#include <bitset>
#include <complex>
#include <deque>
#include <exception>
#include <fstream>
#include <functional>
#include <iomanip>
#include <ios>
#include <iosfwd>
#include <iostream>
#include <istream>
#include <iterator>
#include <limits>
#include <list>
#include <locale>
#include <map>
#include <memory>
#include <new>
#include <numeric>
#include <ostream>
#include <queue>
#include <set>
#include <sstream>
#include <stack>
#include <stdexcept>
#include <streambuf>
#include <string>
#include <typeinfo>
#include <utility>
#include <valarray>
#include <vector>

#if __cplusplus >= 201103L
#include <array>
#include <atomic>
#include <chrono>
#include <condition_variable>
#include <forward_list>
#include <future>
#include <initializer_list>
#include <mutex>
#include <random>
#include <ratio>
#include <regex>
#include <scoped_allocator>
#include <system_error>
#include <thread>
#include <tuple>
#include <typeindex>
#include <type_traits>
#include <unordered_map>
#include <unordered_set>
#endif

using namespace std;
typedef long long ll;
typedef pair<int,int> Pint;
typedef pair<ll, ll> P;
//typedef pair<int, pair<int, int>> P;
//typedef tuple<int,int,int> T;
ll INFL = 1000000000000000010;//10^18 = 2^60
int INF = 2147483600;//10^9
ll MOD  = 1000000007;
//vector<int> dy = {0,0,1,-1};
//vector<int> dx = {1,-1,0,0};
int N;
vector<int> A(110);
double dp[110][110][110];

double solve(int x, int y, int z){
  if(dp[x][y][z] != -1) return dp[x][y][z];
  double exp_num = (N)/(N-z);
  double res = 0;
  if(x != N) res += (solve(x+1, y, z) + exp_num) * (N-x) / (N-z);
  if(y != x) res += (solve(x, y+1, z) + exp_num) * (x-y) / (N-z);
  if(z != y) res += (solve(x, y, z+1) + exp_num) * (y-z) / (N-z);
  return dp[x][y][z] = res;
}

int main(void){
    cin >> N;
    for(int i = 0; i < N; i++) cin >> A[i];

    int a = 0, b = 0, c = 0, d = 0;
    for(int i = 0; i < N; i++){
      if(A[i] >= 1) a++, b++, c++;
      if(A[i] >= 2) b++, c++;
      if(A[i] >= 3) c++;
    }
      for(int i = 0; i < 110; i++){
        for(int j = 0; j < 110; j++){
          for(int k = 0; k < 110; k++){
            dp[i][j][k] = -1;
          }
        }
      }
      dp[N][N][N] = 0;
      cout << solve(a, b, c) << endl;


}