#include <bits/stdc++.h>
using namespace std;
long long n, m, k, nr;
unordered_map<long long, long long> h;
struct cmp {
  bool operator()(long long a, long long b) { return h[a] < h[b]; }
};
priority_queue<long long, vector<long long>, cmp> q;
struct punct {
  long long l, c;
};
punct p;
long long nroffishes(long long l, long long c) {
  long long lst = max(1LL, l - k + 1);
  long long cst = max(1LL, c - k + 1);
  long long ldr = min(n - k + 1, l);
  long long cdr = min(m - k + 1, c);
  long long sol = 1LL * (cdr - cst + 1) * (ldr - lst + 1);
  return sol;
}
punct pozitie(long long nr) {
  punct sol;
  sol.c = nr % m;
  if (!sol.c) {
    sol.c = m;
    sol.l = nr / m;
  } else
    sol.l = (nr / m) + 1;
  return sol;
}
int main() {
  ios_base::sync_with_stdio(false);
  cin.tie(NULL);
  cout.tie(NULL);
  cin >> n >> m >> k >> nr;
  long long midl = (n + 1) / 2;
  long long midc = (m + 1) / 2;
  long long nnr = (midl - 1) * m + midc;
  h[nnr] = nroffishes(midl, midc);
  long long sol = 0;
  q.push(nnr);
  while (nr) {
    long long val = q.top();
    sol += h[val];
    punct poz = pozitie(val);
    q.pop();
    long long val1 = val - 1;
    poz = pozitie(val1);
    if (!h[val1] && poz.l >= 1 && poz.l <= n && poz.c >= 1 && poz.c <= m) {
      long long x = nroffishes(poz.l, poz.c);
      h[val1] = x;
      q.push(val1);
    }
    val1 = val + 1;
    poz = pozitie(val1);
    if (!h[val1] && poz.l >= 1 && poz.l <= n && poz.c >= 1 && poz.c <= m) {
      long long x = nroffishes(poz.l, poz.c);
      h[val1] = x;
      q.push(val1);
    }
    val1 = val + m;
    poz = pozitie(val1);
    if (!h[val1] && poz.l >= 1 && poz.l <= n && poz.c >= 1 && poz.c <= m) {
      long long x = nroffishes(poz.l, poz.c);
      h[val1] = x;
      q.push(val1);
    }
    val1 = val - m;
    poz = pozitie(val1);
    if (!h[val1] && poz.l >= 1 && poz.l <= n && poz.c >= 1 && poz.c <= m) {
      long long x = nroffishes(poz.l, poz.c);
      h[val1] = x;
      q.push(val1);
    }
    --nr;
  }
  long long hmv = (n - k + 1) * (m - k + 1);
  double avg = (0.0 + sol) / (0.0 + hmv);
  cout << fixed << setprecision(12) << avg << '\n';
  return 0;
}
