#include <bits/stdc++.h>

typedef long long LL;
const int N = 105;

int n, m, k, mod;
int factor[N], ifactor[N];

void reduce(int &x) { x += x >> 31 & mod; }
void up(int &a, int b, int c) { a = (a + (LL) b * c) % mod; }
int pow(int x, int y, int ans = 1) {
	for (; y; y >>= 1, x = (LL) x * x % mod)
		if (y & 1) ans = (LL) ans * x % mod;
	return ans;
}

int w[N][N * N];
void init(int n) {
	factor[0] = 1;
	for (int i = 1; i <= n; ++i)
		factor[i] = (LL) factor[i - 1] * i % mod;
	ifactor[n] = pow(factor[n], mod - 2);
	for (int i = n; i; --i)
		ifactor[i - 1] = (LL) ifactor[i] * i % mod;
	for (int i = 1; i <= k; ++i) {
		w[i][0] = 1;
		for (int j = 1; j <= n * n; ++j)
			w[i][j] = (LL) w[i][j - 1] * i % mod;
	}
}

int dp[2][N][N], trans[N][N];

int main() {
	std::ios::sync_with_stdio(0), std::cin.tie(0);
	std::cin >> n >> m >> k >> mod, init(std::max(n, m));
	dp[0][0][0] = 1;
	for (int t = 1, now = 0, lst = 1; t <= k; ++t) {
		if (t > 1) {
			now ^= 1, lst ^= 1, std::memset(dp[now], 0, sizeof dp[0]);
			for (int p = 0; p <= n; ++p)
				for (int q = 0; q <= m; ++q)
					trans[p][q] = (LL) w[t - 1][p * (m - q)] * w[k - t + 1][p * q] % mod * (p & 1 ? mod - ifactor[p] : ifactor[p]) % mod;
			for (int i = 0; i <= n; ++i)
				for (int j = 0; j <= m; ++j)
					for (int k = 0; k <= n - i; ++k)
						up(dp[now][i + k][j], dp[lst][i][j], trans[k][j]);
			now ^= 1, lst ^= 1, std::memset(dp[now], 0, sizeof dp[0]);
			for (int p = 0; p <= m; ++p)
				for (int q = 0; q <= n; ++q)
					trans[p][q] = (LL) w[t - 1][p * (n - q)] * w[k - t + 1][p * q] % mod * (p & 1 ? mod - ifactor[p] : ifactor[p]) % mod;
			for (int i = 0; i <= n; ++i)
				for (int j = 0; j <= m; ++j)
					for (int k = 0; k <= m - j; ++k)
						up(dp[now][i][j + k], dp[lst][i][j], trans[k][i]);
		}
		now ^= 1, lst ^= 1, std::memset(dp[now], 0, sizeof dp[0]);
		for (int p = 0; p <= n; ++p)
			for (int q = 0; q <= m; ++q)
				trans[p][q] = (LL) w[t][p * (m - q)] * w[k - t + 1][p * q] % mod * ifactor[p] % mod;
		for (int i = 0; i <= n; ++i)
			for (int j = 0; j <= m; ++j)
				for (int k = 0; k <= n - i; ++k)
					up(dp[now][i + k][j], dp[lst][i][j], trans[k][j]);
		now ^= 1, lst ^= 1, std::memset(dp[now], 0, sizeof dp[0]);
		for (int p = 0; p <= m; ++p)
			for (int q = 0; q <= n; ++q)
				trans[p][q] = (LL) w[t][p * (n - q)] * w[k - t + 1][p * q] % mod * ifactor[p] % mod;
		for (int i = 0; i <= n; ++i)
			for (int j = 0; j <= m; ++j)
				for (int k = 0; k <= m - j; ++k)
					up(dp[now][i][j + k], dp[lst][i][j], trans[k][i]);
	}
	std::cout << (LL) dp[0][n][m] * factor[n] % mod * factor[m] % mod << std::endl;
	return 0;
}