explicit representation of creation/annihilation operators & its fourier transform (matrix form) (tight-binding hamiltonian, graphene) While I'm a mathematician/computer-scientist myself, I've some problems trying to understand some paper about the electronic states of graphene nanoribbons modelled by the tight-binding Hamiltonian (http://iopscience.iop.org/article/10.1088/1468-6996/11/5/054504/meta).

The nearest neighbour tight-binding Hamiltonian is defined as
  $$H = -t \sum_{R_B}\sum_{l=1}^3 a_{R_B+\tau_l}^\dagger b_{R_B} + \hbox {H.c.}\, ,$$ where
   H.c. stands for the Hermitian conjugate, $R_B+\tau_l$ are the neighbouring atoms of $R_B$  and $a_{R_A}^\dagger$ $(b_{R_B}^\dagger)$ create and $a_{R_A}$ $(b_{R_B})$ annihilate an electron at $R_A$ $(R_B)$ of the $A$ $(B)$ sublattice. The operators satisfy the Fermion anticommutation relations
  $$ 
\{a_{R_A},a_{R_A'}^\dagger\} = \delta_{R_A,R_A'}\, ,\quad  \{b_{R_B},b_{R_B'}^\dagger\} = \delta_{R_B,R_B'}\, ,\quad 
\{a_{R_A},a_{R_A'}\} =0\, ,\quad \{b_{R_B},b_{R_B'}\} =0\, .
$$

For me, it was sufficient to interpret the Hamiltonian for some finite graphene sample 
\begin{align}
\mathbb{L}&=\{j_1 a_1 +j_2a_2 + \{0,1\}\tau : \tau=\frac{1}{3}(a_1 + a_2), \\
a_1&=\left(\frac{3}{2}, +\frac{\sqrt{3}}{2}\right)\, , \quad a_2=
\left(\frac{3}{2}, -\frac{\sqrt{3}}{2}\right)\, ,\quad j_i \in \{0,1,\ldots\} \}
\end{align}
as a Adjacency matrix, which has 3 entries per row (all equal to $-t$). With this interpretation of $H$ I am able to verify numerically and understand the extraction of the eigenvalues and eigenfunctions to some extend: F.e.: I can block-diagonalize the matrix H (if permuted correctly) into $2\times2$ blocks using the orthonormal basis $(v_i)_{i=1,\ldots,N}$, where the vector $v_i$ has the entries $\left(\begin{matrix} e^{ikR_A} \\ \pm e^{ik(R_A+\tau)}\end{matrix}\right)$ where the first row corresponds to the values for the atoms of type A at  position $R_A$ and the second row to atoms of type B at the position $(R_A+\tau)$ and $k\in L^*=$reciprocal Lattice of $L$.
What I don't understand is: what would be an explicit representation of the creation and annihilation operators? I interpreted them just as vectors (or functionals) $a_{R_A}=\left(\begin{array}{c} 0 & \cdots & 0 & 1 & 0 & \cdots & 0\end{array}\right)$, where the 1 is at the position which corresponds to the atom $R_A$, since with this interpretation i end up with the correct matrix representation of $H$. BUT, obviously, something like
$ \{a_{R_A},a_{R_A'}\} = a_{R_A}a_{R_A'} =$ row-vector $\cdot$ row-vector
does not make any sense at all. So, is there any matrix representation of these operators? 
And as a follow-up question: How is the Fourier transform defined for these operators? In the paper it says

We apply the following Fourier transformation to the above Hamiltonian $$a_{R_A} = \frac{1}{\sqrt{L_xL_y}} \sum_k e^{ikR_A} \alpha_k\, ,\quad 
b_{R_B} = \frac{1}{\sqrt{L_xL_y}} \sum_k e^{ikR_B} \beta_k.$$ 
  Here $k=(k_x,k_y),$ and  $L_x(L_y)$ denotes the number of unit cells in $x(y)$ direction.
  $H_k = -t \sum_k \sum_{l=1}^3 e^{-ik\tau_l} \alpha_k^\dagger \beta_k + $H.c.

I know that the Fourier transform is defined for vectors or functionals, for example if I interpret the creation/annihilation operators as functionals, i.e.:
$$a_{R_A} : \mathbb{L}=\hbox{graphene Lattice}=\{
\hbox{positions of the atoms of the graphene lattice}\}\rightarrow \mathbb{C}$$ with 
$a_{R_A}(x) = \delta_{R_A,x} $, the Fourier transform can be applied to $a_{R_A}$.
The Fourier transform is defined by 
$$F:\{f:\mathbb{L} \rightarrow \mathbb{C}\}\rightarrow \{f:\mathbb{L^*} \rightarrow \mathbb{C}\}$$ with 
$(Ff)(k) = \sum_x e^{ikx} f(x)$, where $L^*$ is the reciprocal Lattice. Thus, we have
$$(Fa_{R_A})(k) = \sum_x e^{ikx} a_{R_A}(x) = e^{ikR_A},$$ 
which is a result which does not seem to be completely off...
But how can the Fourier transform be applied to an operator, such that everything is consistent?
 A: I think you are confusing the space in which the hopping matrix acts, which has dimension $N$ where $N$ is the number of atomic sites, which the much bigger Fock space on which the $a$ and $a^\dagger$ act.
The Fock   Hilbert space is $\bigotimes_{i=1}^N (V_2)_i $ where $N$ is the number of sites and $(V_2)_i$ is a two dimensional complex vector space, one for each atomic  site labelled a sequence by $i=1,\ldots N$. So the total Hilbert space is $2^N$ dimensional. Acting on this space we can first try 
$$
a_i = {\mathbb I}_1\otimes {\mathbb I}_2\otimes \cdots \otimes (\sigma_-)_i \otimes {\mathbb I}_{i+1}\otimes \cdots\otimes {\mathbb I}_N
$$
where the 
$$
\sigma_-= \left(\matrix{0&0\\ 1&0}\right)
$$
factor acts on the $(V_2)_i$ at site $i$. Similarly
$$
a_i^\dagger  = {\mathbb I}_1\otimes {\mathbb I}_2\otimes \cdots \otimes (\sigma_+)_i \otimes \cdots \otimes {\mathbb I}_N
$$
where 
$$
\sigma_+= \left(\matrix{0&1\\ 0&0}\right)
$$
(I am not going to distinguish between $a$'s and $b$'s here.)
Unifortunately, whle $\{\sigma_+,\sigma_-\}= {\mathbb I}$ these $a_i$ commute at different sites rather than anticommute.  To get the anticommutation we use a Klein transformation where we replace the ${\mathbb I}$ to the left of $\sigma_{\pm}$ in the product by 
$$
\sigma_3= \left(\matrix{1&0\\ 0&-1}\right).
$$
Physicists don't usually explain it in these terms. They just take linear combinations (your fourier series) to reduce the Hamiltonian acting on the big Fock space to the form
$$
\hat H= \sum E_n c^\dagger_n c_n
$$ 
