#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;
int A;
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 = 1.0 * N / (N-z);
  double res = 0.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;
    int a = 0, b = 0, c = 0;
    for(int i = 0; i < N; i++){
      int A; cin >> A;
      if(A >= 3) a++, b++, c++;
      else if(A >= 2) a++, b++;
      else if(A >= 1) a++;
    }
      for(int i = 0; i <= N; i++){
        for(int j = 0; j <= N; j++){
          for(int k = 0; k <= N; k++){
            dp[i][j][k] = -1;
          }
        }
      }
      dp[N+1][N+1][N+1] = 0;
      cout << fixed << setprecision(10) << solve(a, b, c) << endl;

}