Consider the finite dimensional unitary representations $\alpha,\beta,\gamma$ of the given compact group $G$ on corresponding vector spaces $V_1,V_2,V_3$. Let $|i>_j,i=1,\dots,n_j$ be an orthnormal basis of $V_j$ where $dim V_j=n_j$. Then
$\{|i>_1\otimes|j>_2\otimes|k>_3\}$ forms an orthonormal basis of $V=V_1\otimes V_2 \otimes V_3$. An element $g\in G$ acts on the tensor product space $V$ as

$$|i>_1\otimes|j>_2\otimes|k>_3 \to \alpha(g) |i>_1\otimes\beta(g)|j>_2\otimes \gamma(g)|k>_3 \tag 1$$

We can also find an orthogonal basis $e_{\mu},\mu=1,\dots,N$ (where $N=n_1n_2n_3$) of $V$ relative to which all elements $g \in G$ act as block diagonal matrices. More precisely, suppose with respect to basis $\{e_{\mu}\}$, the action of an element $g\in G$ on $V$ be denoted as
$$e_{\mu}\to \rho (g)e_\mu \tag 2$$
then $\rho(g)^{\nu}_{\mu}=<\nu|\rho(g)|\mu>$ is a block diagonal matrix where the dimensions of different blocks$^1$ are independent of $g$. Let

$$|i>_1\otimes|j>_2\otimes|k>_3=\sum_{\mu}\epsilon ^{\mu}_{ijk}e_{\mu}\tag 3$$

Acting with $g\in G$ on both sides of this equation gives

$$\alpha (g)|i>_1\otimes \beta (g)|j>_2\otimes \gamma(g)|k>_3=\sum_{\mu}\epsilon ^{\mu}_{ijk}\rho (g)e_{\mu} \tag 4$$

Taking the scalar product with $|i'>_1\otimes |j'>_2 \otimes |k'>_3$ and using (3) we get

$$\alpha(g)^{i'}_{i}\beta (g)^{j'}_{j}\gamma(g)^{k'}_{k}=\sum _{\mu,\nu} \rho^{\nu}_{\mu}(g)\epsilon^{*\nu}_{i'j'k'}{\epsilon}^{\mu}_{ijk}\tag 5$$

Now, according to the Peter-Weyl theorem (part 2), matrix elements of the irreducible representations of $G$ form an orthogonal basis of the space of square integrable functions on $G$ wrt the inner product

$$(A,B)=\int_G dg\; A(g)^{*}B(g) \tag 6$$

where $dg$ is the Haar measure. So, if we integrate both sides of (5), the nonzero contribution on RHS will only come from the part of $\rho$ which is direct sum of identity representations. In other words, let $W\subseteq V$ be the subspace of $V$ on which $G$ acts trivially, and let $\{e_1,\dots e_m\}\subseteq $ $\{e_1,\dots e_m,\dots,e_N\}$ be the basis of $W$, then the integration of (5) gives

$$\int_G dg\;\alpha(g)^{i'}_{i}\beta (g)^{j'}_{j}\gamma(g)^{k'}_{k}=\sum _{\mu,\nu =1}^{m} \delta^{\nu}_{\mu}\epsilon^{*\nu}_{i'j'k'}{\epsilon}^{\mu}_{ijk}=\sum _{\mu =1}^{m}\epsilon^{*\mu}_{i'j'k'}{\epsilon}^{\mu}_{ijk} \tag 7$$

where we have assumed that $Vol(G)=\displaystyle\int_G\; dg =1$

For the second part of your question, I would recommend these lecture notes. The basic idea for computing Wilson loop averages is following -

For a surface with boundary, the partition function of two dimensional Yang-Mills theory depends on the holonomy along the boundary. Let the partition function of a surface of genus $h$ and one boundary be denoted as $Z_h(U,ag^2)$ where $U$ is the fixed holonomy along the given boundary, $a$ is the area of the surface and $g$ is the Yang-Mills coupling constant. Now, consider the simplest situation in which a Wilson loop $W$ in representation $R_W$ is inserted along a contractible loop $C$ on a closed surface of genus $h$, and area $a$. To compute the Wilson loop average, we first cut the surface along $C$, which gives a disc $D$ of area (say) $b$ and another surface $S$ of area $c=a-b$, genus $h$ and one boundary. Now the Wilson loop average is given by integrating over $G$ the product of i) the partition functions of $D$ ii) partition function of $S$ and iii) the trace of the Wilson loop in representation $R_W$ -

$$<W>=\frac {1}{Z_h(ag^2)}\int \: dU \:Z_h(U,(a-b)g^2)\chi_{R_W}(U)Z_0(U^{-1},bg^2) \tag 8$$

The case of a self-intersecting Wilson loop too is not very different.

$^1$ The smallest blocks form irreducible representations of $G$; Exactly which irreducible representations show up will depend on $\alpha,\beta,\gamma$; The same irreducible representations may also appear more than once

This post imported from StackExchange Physics at 2014-06-17 22:45 (UCT), posted by SE-user user10001