  Infinite projected entangled pair states (iPEPS) have emerged as a powerful
tool for studying interacting two-dimensional fermionic systems. In this
review, we discuss the iPEPS construction and some basic properties of this
tensor network (TN) ansatz. Special focus is put on (i) a gentle introduction
of the diagrammatic TN representations forming the basis for deriving the
complex numerical algorithm, and (ii) the technical advance of fully exploiting
non-abelian symmetries for fermionic iPEPS treatments of multi-band lattice
models. The exploitation of non-abelian symmetries substantially increases the
performance of the algorithm, enabling the treatment of fermionic systems up to
a bond dimension $D=24$ on a square lattice. A variety of complex
two-dimensional (2D) models thus become numerically accessible. Here, we
present first promising results for two types of multi-band Hubbard models, one
with $2$ bands of spinful fermions of $\mathrm{SU}(2)_\mathrm{spin} \otimes
\mathrm{SU}(2)_\mathrm{orb}$ symmetry, the other with $3$ flavors of spinless
fermions of $\mathrm{SU}(3)_\mathrm{flavor}$ symmetry.
