SRM460 div1 medium

40^5でDP。

#include <vector>

using namespace std;

double DPJ[40+10][40+10][40*40+10];
double DPB[40+10][40+10][40*40+10];

class TheFansAndMeetingsDivOne 
{
	public:
		double find(vector <int> minJ, vector <int> maxJ, vector <int> minB, vector <int> maxB, int k) 
		{
			int n=minJ.size();
			for(int city=n;0<=city;city--)
				for(int chance=0;chance<=k;chance++)
					for(int left=0;left<=40*40;left++)
					{
						DPJ[city][chance][left]=0.0;
						if(chance==0 && left==0)
							DPJ[city][chance][left]=1.0;
						else if(city==n || chance==0)
							;
						else
						{
							double h=(double)chance/(double)(n-city);
							double a=maxJ[city]-minJ[city]+1.0;
							for(int i=minJ[city];i<=maxJ[city] && 0<=left-i;i++)
								DPJ[city][chance][left]+=h*DPJ[city+1][chance-1][left-i]/a;
							DPJ[city][chance][left]+=(1-h)*DPJ[city+1][chance][left];
						}

						DPB[city][chance][left]=0.0;
						if(chance==0 && left==0)
							DPB[city][chance][left]=1.0;
						else if(city==n || chance==0)
							;
						else
						{
							double h=(double)chance/(double)(n-city);
							double a=maxB[city]-minB[city]+1.0;
							for(int i=minB[city];i<=maxB[city] && 0<=left-i;i++)
								DPB[city][chance][left]+=h*DPB[city+1][chance-1][left-i]/a;
							DPB[city][chance][left]+=(1-h)*DPB[city+1][chance][left];
						}
					}

			double ans=0.0;
			for(int i=0;i<=40*40;i++)
				ans+=DPJ[0][k][i]*DPB[0][k][i];

			return ans;
		}
};