Calculating Hamiltonian matrix elements in a chaotic system

[Pyrus96]

The system in which I tried to calculate the Hamiltonian matrix was a particle in a stadium (Billiard stadium). And I used the principle where we take a rectangle around the stadium in which the parts outside the stadium have a very high potential V₀.
We know the wave function of a rectangular potential well and through this I found the following formula for the Hamiltonian matrix elements:
$$H_{nm} = [(\frac{m_{1}}{a_{1}})^{2} +(\frac{m_{2}}{a_{2}})^{2}]\delta_{nm} + V_{0}v_{nm}$$
where v_nm is given by:
$$v_{nm} = \int_{\Gamma^{'}}\phi_{n}\phi_{m}$$
where $\Gamma^{'}$ denotes the area in the rectangle which isn't part of the stadium. And the wave functions are given by:
$$\phi_{m_{1},m_{2}}(x_{1},x_{2}) = \int sin(\frac{\pi m_{1}x_{1}}{a_{1}})sin(\frac{\pi m_{2}x_{2}}{a_{2}})$$
but here's my problem, the system we work in is two dimensional and thus we get energy levels corresponding to the doublets: $n = (n_{1},n_{2})$ and $m = (m_{1},m_{2})$. How am I supposed to interpret the matrix with this because I can't imagine a matrix from that can correspond to $H_{(2,1)x(1,1)}$

Please note that I didn't bother to write several constants in front of the formula of the hamiltonian and the wave function because they don't matter to the answer.

 

[Misha]

It is difficult to say, what is your final objective with those matrix elements but they will be just a version of your first formula
$$
H_{\mu, \nu} = \epsilon_{\mu} \delta_{\mu, \nu} + V_0 v_{\mu, \nu},
$$
where $\mu$ and $\nu$ are some parametrizations of states in an rectangular box, $\epsilon_\mu$ are respective energies and so on. For instance, $|\mu \rangle = |m_1, m_2 \rangle$ and
$$
\langle \mu | H | \nu \rangle =: H_{\mu, \nu} = H_{m_1, m_2; n_1, n_2} : = \langle m_1, m_2 | H | n_1, n_2 \rangle =
\int dx_1 \ dx_2 \ \phi_{m_1, m_2}^*(x_1, x_2) \widehat{H} \phi_{n_1, n_2}(x_1, x_2).
$$
You deal with such "direct product" matrices in a way very similar to usual matrices. For example,
$$
\langle \mu | \widehat{A} \widehat{B} | \nu \rangle = \sum_\lambda A_{\mu, \lambda} B_{\lambda, \nu} =
\sum_{l_1, l_2} \langle \mu | \widehat{A} | l_1, l_2 \rangle \langle l_1, l_2 | \widehat{B} | \nu \rangle
$$
and so on.