Given a lattice $L=\{x = j_1a_1 + j_2a_2 : j_1,j_2 \in \mathbb{Z}\}$ with $a_1,a_2 \in\mathbb{R}^2$ beeing linearly independent.
We can define a torus by $\mathbb{T}=L\,/ (\mathbb{Z}C + \mathbb{Z}T)$, where $C=c_1a_1 + c_2a_2$ and $T=t_1a_1+t_2a_2$, $c_1,c_2,t_1,t_2 \in\mathbb{Z}$, are linearly independent as well.
Thus, we have for $x,y \in \mathbb{T}$: $x=y \Leftrightarrow \exists z_1,z_2 \in\mathbb{Z} $ with $x=y + z_1C+z_2T$.
Let $C^\perp$ and $T^\perp$ such that $$ 0 = \langle C^\perp, C\rangle = \langle T^\perp, T\rangle $$ and $$ 2\pi = \langle C^\perp, T\rangle = \langle T^\perp, C\rangle. $$
I would like to show that the functions $e_k(x) =\frac{1}{\sqrt{|\mathbb{T}|}} e^{i \langle k, x\rangle }$, $e_k:\mathbb{T}\rightarrow \mathbb{C}$, $k=\mathbb{Z}C^\perp + \mathbb{Z}T^\perp$ are orthogonal, i.e.,
$$\langle e_k, e_{\tilde{k}} \rangle = \sum_{x\in\mathbb{T}} e_k(x) \overline{e_{\tilde{k}(x)}} = \sum_{x\in\mathbb{T}} e^{i \langle k - \tilde{k}, x\rangle } = 0 \Leftrightarrow k\neq \tilde{k},$$ but I don't quite know how to do that.
Eventually I would like to show that these functions form an orthonormal basis of $\{f:\mathbb{T}\rightarrow\mathbb{C}\}$ and the orthogonality is the missing part.
Maybe the proof of the cardinality helps in some sense: The number of wavefunctions is equal to the number of lattice points:
Let $b_1,b_2$ the reciprocal lattice vectors, i.e. $\langle a_i,b_j \rangle = 2\pi \delta_{ij}$. Then we have $e_k(x)= e_{k+j_1b_1+j_2b_2}(x)$ for every $x \in\mathbb{T}, j_1,j_2\in\mathbb{Z}$.
Thus, the number $K$ of different wavefunctions $e_k, k=\mathbb{Z}C^\perp + \mathbb{Z}T^\perp$ is equal to $$ K= |\frac{\det\begin{pmatrix}b_1|b_2\end{pmatrix}}{\det{\begin{pmatrix}C^\perp|T^\perp\end{pmatrix}}}|.$$
The number of lattice points in $\mathbb{T}$ is given by $$ |\mathbb{T}| = |\frac{\det{\begin{pmatrix}C|T\end{pmatrix}}}{\det\begin{pmatrix}a_1|a_2\end{pmatrix}}|.$$
It is a straightforward calculation to see that $|\mathbb{T}|=K$ if you write $k$ in terms of the reciprocal lattice vectors $b_1,b_2$.
In matrix/vector notation :
We have a lattice $L = M \mathbb{Z}^2$ for some inversible matrix $M \in \mathbb{R}^{2 \times 2}$
If $L_2= M_2 \mathbb{Z}^2$ is another lattice then $L_2 \supseteq L$ iff $MM_2^{-1} \in \mathbb{Z}^{2 \times 2}$.
The dual lattice is $L^\perp = M^{-1} \mathbb{Z}^2$.
We are looking at the vector space of $L$-periodic functions $L_2 \to \mathbb{C}$, ie. the finite dimensional vector space of functions $L_2/L \to \mathbb{C}$
We have a basis of complex exponentials $$e_k(x) = e^{-2i \pi \langle k, x \rangle}, \qquad k \in L^\perp/L_2^\perp,x \in L_2 $$
They are orthonormal for the inner product
$$\langle e_k,e_m \rangle = \frac{1}{|L_2/L|} \sum_{x \in L_2/L} e_k(x) \overline{e_m(x)} = \frac{1}{|L_2/L|}\sum_{x \in L_2/L} \exp(-2i \pi \langle k-m,x \rangle)\\ =\frac{1}{|L_2/L|}\sum_{x \in L_2/L} \exp(-2i \pi \langle (\underbrace{MM_2^{-1} }_{\in \mathbb{Z}^{2 \times 2}})^{-1}\underbrace{M(k-m)}_{\in \mathbb{Z}^2},\underbrace{M_2^{-1} x}_{\in \mathbb{Z}^2} \rangle)$$
And you can check this is $ = 0$ whenever $M_2(k-m) \not \in \mathbb{Z}^2$ ie. whenever $k \ne m$ in $L_2^\perp/L^\perp$.