Go: package main

import (
	"fmt"
	"math"
)

func main(){
	var n, m, k int
	for {
		fmt.Scan(&n, &m, &k)
		if n == 0 && m == 0 && k == 0{
			break
		}

		dp := make([][]int, n+1)
		for i:= 0; i < n+1; i++{
			dp[i] = make([]int, m*n + 1)
		}

		for i := 1; i <=m ; i++{
			dp[1][i] = 1
		}

		for i:= 1; i < n; i++{
			for j := 0; j <= i * m; j++{
				for k := 1; k <= m; k++{
					dp[i+1][j+k] += dp[i][j]
				}
			}
		}

		count := 0.0
		div := math.Pow(float64(m), float64(n))
		for i := 0; i <= n * m; i++{
			bill := 0.0
			if i <= k{
				bill = 1.0
			} else {
				bill = float64(i - k)
			}
			count += bill * float64(dp[n][i]) / div
		}
		fmt.Printf("%.010f\n", count)
	}
}

C++: #include <iostream>
#include <cstring>
#include <algorithm>
#include <cstdio>

using namespace std;

int main(){
    int n,m,k;

    while(cin >> n >> m >> k){
        if((n|m|k) == 0) break;
        double dp[n+1][n*m+1];

        fill((double *)dp, (double *)dp+(n+1)*(n*m+1), 0);

        dp[0][0] = 1;

        for(int i=0;i<n;i++){
            for(int j=0;j<=n*m;j++){
                for(int l=1;l<=m;l++){
                    if(j+l <= n*m){
                        dp[i+1][j+l] += dp[i][j] / m;
                    }
                }
            }
        }
        double ret = 0;

        for(int i=0;i<=n*m;++i){
            ret += dp[n][i] * max(1,(i-k));
        }
        printf("%.8f\n", ret);
    }
}

C++: //#include <bits/stdc++.h>
#include <stack>
#include <cassert>
#include <set>
#include <queue>
#include <iostream>
#include <vector>
#include <complex>
#include <algorithm>
#include <cstring>
#include <cstdio>
#include <map>
#include <cmath>
#include <numeric>

using namespace std;
typedef long long ll;
typedef pair<int,int> PI;
#define EPS (1e-6)
#define rep(i,n) for(int i=0;i<(int)(n);++i)
#define REP(i, n) rep (i, n)
#define F first
#define S second
#define mp(a,b) make_pair(a,b)
#define pb(a) push_back(a)
#define min3(a,b,c) min((a),min((b),(c)))
#define min4(a,b,c,d) min((a),min3((b),(c),(d)))
#define SZ(a) (int)((a).size())
#define ALL(a) (a).begin(),(a).end()
#define RALL(a) a.rbegin(),a.rend()
#define FLL(a,b) memset((a),b,sizeof(a))
#define CLR(a) memset((a),0,sizeof(a))
#define declare(a,it) __typeof(a) it=(a)
#define FOR(it,a) for(declare(a.begin(),it);it!=a.end();++it)
#define FORR(it,a) for(declare(a.rbegin(),it);it!=a.rend();++it)

template<typename T,typename U> ostream& operator<< (ostream& out, const pair<T,U>& val){return out << "(" << val.F << ", " << val.S << ")";}
template<class T> ostream& operator<< (ostream& out, const vector<T>& val){out << "{";rep(i,SZ(val)) out << (i?", ":"") << val[i];return out << "}";}
typedef double FP;
typedef complex<FP> pt;
typedef pt P;
typedef pair<pt,pt> line;
namespace std{
  bool operator<(const P&a,const P&b){
    if(abs(a.real()-b.real())>EPS)
      return a.real()<b.real();
    return a.imag()<b.imag();
  }
}
FP dot(P a,P b){return real(conj(a)*b);}
FP crs(P a,P b){return imag(conj(a)*b);}
P ortho(P a){return P(imag(a),-real(a));}
P ortho(line a){return ortho(a.S-a.F);}
P crspt(P a,P b,P c,P d){b-=a,d-=c;return a+b*crs(d,c-a)/crs(d,b);}
P crspt(line a,line b){return crspt(a.F,a.S,b.F,b.S);}
bool onl(P a1,P a2,P b){return abs(b-a1)+abs(b-a2)<abs(a1-a2)+EPS;}
bool onl(line a,P b){return onl(a.F,a.S,b);}
bool iscrs(line a,line b){P c=crspt(a,b);return onl(a,c)&&onl(b,c);}
void pkuassert(bool t){t=1/t;};
int dx[]={0,1,0,-1,1,1,-1,-1};
int dy[]={1,0,-1,0,-1,1,1,-1};
enum{TOP,BTM,LFT,RGT,FRT,BCK};
int dxdy2ce[]={RGT,FRT,LFT,BCK};
int s2i(string& a){stringstream ss(a);int r;ss>>r;return r;}
template<class T> T shift(T a,int b,int c,int d,int e){
  __typeof(a[0])t=a[b];
  a[b]=a[c];a[c]=a[d];a[d]=a[e];a[e]=t;return a;}
template<class T> T rgt(T a){return shift(a,TOP,LFT,BTM,RGT);}
template<class T> T lft(T a){return shift(a,TOP,RGT,BTM,LFT);}
template<class T> T frt(T a){return shift(a,TOP,BCK,BTM,FRT);}
template<class T> T bck(T a){return shift(a,TOP,FRT,BTM,BCK);}
line mkl(P a,P v){return line(a,a+v);}
FP lpdist(line a,P b){return abs(b-crspt(a,mkl(b,ortho(a))));}
FP spdist(line a,P b){
  P c(crspt(a,mkl(b,ortho(a))));
  return onl(a,c)?abs(b-c):min(abs(a.F-b),abs(a.S-b));
}
FP ssdist(line a,line b){
  return
    iscrs(a,b)?0.:
    min4(spdist(a,b.F),spdist(a,b.S),
         spdist(b,a.F),spdist(b,a.S));
}

int n,m,k;

void solve(){
  double ans = 0;
  if(n==1){
    rep(j,m) ans += max(1,j+1-k)*1./m;
  }else{
    double dp[n*m+1];
    rep(i,n*m+1) dp[i]=0;
    dp[0] = 1;
    rep(jj,n){
      double ndp[n*m+1];
      rep(i,n*m+1) ndp[i]=0;
      
      rep(i,n*m+1) if(dp[i]> 1e-9)
        rep(j,m) ndp[i+j+1]+=dp[i]/m;
      rep(i,n*m+1) dp[i]=ndp[i];
    }
    for(int i =1;i<=n*m;++i)
      ans += max(1,i-k)*dp[i];
  }
  printf("%.13f\n",ans);
}

int main(int argc, char *argv[])
{
  while(cin >> n >> m >> k,n) solve();
  return 0;
}

C++: #include<bits/stdc++.h>
#define range(i,a,b) for(int i = (a); i < (b); i++)
#define rep(i,b) for(int i = 0; i < (b); i++)
#define all(a) (a).begin(), (a).end()
#define debug(x) cout << "debug " << x << endl;
const int INF = 100000000;
using namespace std;


void solve(double dp[100000][2], int n, int m){
    dp[0][0] = 1;
    rep(i,n){
        rep(j, n * m){
            range(k, 1, m + 1){
                dp[j + k][(i + 1) % 2] += dp[j][i % 2];
            }
            dp[j][i % 2] = 0;
        }
    }
}

int main(){
    int n, m;
    double k;
    while(cin >> n >> m >> k, n||m||k){
        double dp[100000][2] = {0};
        solve(dp, n, m);
        //rep(i,n * m + 1){ cout << dp[i][n] << endl; }


        double t = 1;
        rep(i,n){
            t *= m;
        }

        double sum = 0;
        for(int i = n; i <= n * m; i++){
            sum += max(static_cast<double>(i - k), 1.0) * dp[i][n % 2];

        }
        cout << fixed << setprecision(8) << sum / t + 1e-12 << endl;
    }
}

C++: #include <sstream>
#include <string>
#include <vector>
#include <map>
#include <algorithm>
#include <iostream>
#include <utility>
#include <set>
#include <cctype>
#include <queue>
#include <stack>
#include <cstdio>
#include <cstdlib>
#include <cmath>
#include <climits>
using namespace std;
#define INF 100000000
#define MAXN 100
#define MAXM 10000
 
typedef long long ll;
const int dx[] = {1, 0, -1, 0};
const int dy[] = {0, 1, 0, -1};
 
double dp[MAXN][MAXM]; // dp[n][m]はn回サイコロをふったとき合計値がmになる確率
 
int n, m, k;
 
int main(void) {
    while (cin >> n >> m >> k) {
        if (n == 0 && m == 0 && k == 0) break;
        for (int i = 0; i < MAXN; i++) 
            for (int j = 0; j < MAXM; j++) 
                dp[i][j] = 0;
 
        dp[0][0] = 1;
        int ma = n * m;
        for (int i = 1; i <= n; i++) {
            for (int j = 1; j <= ma; j++) {
                for (int k = 1; k <= m; k++) {
                    if (j - k >= 0) {
                        dp[i][j] += dp[i-1][j-k] / m;
                    }
                }
            }
        }
        //for (int i = 0; i <= ma; i++) cout << dp[n][i] << endl;
        double ans = 0;
        for (int i = 1; i <= ma; i++) {
            if (i <= k) ans += 1.0*dp[n][i];
            else ans += (i-k) * dp[n][i];
        }
        printf("%.10f\n", ans);
    }
    return 0;
}

C++: #include <iostream>
#include <cmath>
using namespace std;

int main(){
    int n,m,k;
    cout<<fixed;
    cout.precision(8);
    while(cin>>n>>m>>k,n){
        double dp[2][n*m+1];
        fill(dp[0],dp[2],0.0);
        dp[0][0]=1;
        for(int i=0;i<n;i++){
            for(int j=0;j<=n*m;j++) dp[i&1^1][j]=0;
            for(int j=0;j<=n*m;j++){
                for(int d=1;d<=m&&j+d<=n*m;d++){
                    dp[i&1^1][j+d]+=dp[i&1][j]/m;                
                }
            }
        }
        double ans = 0;
        for(int i=0;i<=n*m;i++){
            ans+=dp[n&1][i]*(max(1,i-k));
        }
        cout<<ans<<endl;
    }
    return 0;
}