[proofplan] View the inclusion $G \hookrightarrow \operatorname{GL}_n(\mathbb{C})$ as a complex representation of $G$ on $\mathbb{C}^n$ and apply Weyl's unitary trick to produce a $G$-invariant Hermitian inner product $\langle \cdot, \cdot \rangle$ on $\mathbb{C}^n$. Choose a basis orthonormal for $\langle \cdot, \cdot \rangle$ and let $P$ be the change-of-basis matrix from the standard basis to this basis. The dictionary identity $\langle \mathbf{x}, \mathbf{y} \rangle = (P^{-1}\mathbf{x},\, P^{-1}\mathbf{y})_0$, where $(\cdot, \cdot)_0$ is the standard Hermitian inner product, converts the $G$-invariance of $\langle \cdot, \cdot \rangle$ into the statement that $P^{-1}gP$ preserves $(\cdot, \cdot)_0$ for every $g \in G$, i.e., $P^{-1}GP \leq \operatorname{U}(n)$. [/proofplan]

[step:Cast the inclusion $G \hookrightarrow \operatorname{GL}_n(\mathbb{C})$ as a representation on $\mathbb{C}^n$] Let $G \leq \operatorname{GL}_n(\mathbb{C})$ be a finite subgroup. Define \begin{align*} \rho: G &\to \operatorname{GL}(\mathbb{C}^n) \\ g &\mapsto g, \end{align*} where the matrix $g \in \operatorname{GL}_n(\mathbb{C})$ is identified with the linear map $\mathbf{v} \mapsto g\mathbf{v}$ acting by matrix multiplication in the standard basis. Matrix multiplication respects products, so $\rho$ is a group homomorphism, hence a complex representation of $G$ on $\mathbb{C}^n$. [/step]

[step:Apply Weyl's unitary trick to produce a $G$-invariant Hermitian inner product] Since $G$ is finite, [Weyl's Unitary Trick](/theorems/2411) supplies a Hermitian inner product \begin{align*} \langle \cdot, \cdot \rangle: \mathbb{C}^n \times \mathbb{C}^n \to \mathbb{C} \end{align*} that is $G$-invariant: $\langle g\mathbf{v}, g\mathbf{w} \rangle = \langle \mathbf{v}, \mathbf{w} \rangle$ for all $g \in G$ and $\mathbf{v}, \mathbf{w} \in \mathbb{C}^n$. [/step]

[step:Choose an orthonormal basis $(\mathbf{f}_i)$ for $\langle \cdot, \cdot \rangle$ and form the change-of-basis matrix $P$] By Gram–Schmidt orthogonalisation applied to the standard basis with respect to $\langle \cdot, \cdot \rangle$, there exists a basis $\mathbf{f}_1, \dots, \mathbf{f}_n$ of $\mathbb{C}^n$ orthonormal for $\langle \cdot, \cdot \rangle$: \begin{align*} \langle \mathbf{f}_i, \mathbf{f}_j \rangle = \delta_{ij}, \qquad 1 \leq i, j \leq n. \end{align*} Let $\mathbf{e}_1, \dots, \mathbf{e}_n$ be the standard basis. Define $P \in \operatorname{GL}_n(\mathbb{C})$ to be the unique matrix satisfying \begin{align*} P\mathbf{e}_i = \mathbf{f}_i, \qquad 1 \leq i \leq n. \end{align*} Equivalently, $P^{-1}\mathbf{f}_i = \mathbf{e}_i$, and the $i$-th column of $P$ is the coordinate vector of $\mathbf{f}_i$ in the standard basis. [/step]

[step:Establish the dictionary identity $\langle \mathbf{x}, \mathbf{y} \rangle = (P^{-1}\mathbf{x},\, P^{-1}\mathbf{y})_0$] Let $(\cdot, \cdot)_0: \mathbb{C}^n \times \mathbb{C}^n \to \mathbb{C}$ denote the standard Hermitian inner product, defined by sesquilinear extension of $(\mathbf{e}_i, \mathbf{e}_j)_0 = \delta_{ij}$. [claim:For all $\mathbf{x}, \mathbf{y} \in \mathbb{C}^n$, $\langle \mathbf{x}, \mathbf{y} \rangle = (P^{-1}\mathbf{x},\, P^{-1}\mathbf{y})_0$] [/claim] [proof] Both sides are sesquilinear forms on $\mathbb{C}^n \times \mathbb{C}^n$, so it suffices to check the identity on a basis. We use the basis $(\mathbf{f}_k)_{k=1}^n$. For any $1 \leq i, j \leq n$, \begin{align*} \langle \mathbf{f}_i, \mathbf{f}_j \rangle = \delta_{ij} \end{align*} by orthonormality of $(\mathbf{f}_k)$ for $\langle \cdot, \cdot \rangle$, and \begin{align*} (P^{-1}\mathbf{f}_i,\, P^{-1}\mathbf{f}_j)_0 = (\mathbf{e}_i, \mathbf{e}_j)_0 = \delta_{ij} \end{align*} using $P^{-1}\mathbf{f}_i = \mathbf{e}_i$ and orthonormality of $(\mathbf{e}_k)$ for $(\cdot, \cdot)_0$. The two sesquilinear forms agree on the basis $(\mathbf{f}_k)$, so they agree as forms on $\mathbb{C}^n \times \mathbb{C}^n$. [/proof] Substituting $\mathbf{x} \to P\mathbf{x}$, $\mathbf{y} \to P\mathbf{y}$ (valid since $P$ is bijective) yields the equivalent dictionary identity that we will use in the calculation: \begin{align*} \langle P\mathbf{x},\, P\mathbf{y} \rangle = (\mathbf{x}, \mathbf{y})_0 \qquad \text{for all } \mathbf{x}, \mathbf{y} \in \mathbb{C}^n. \quad (\diamond) \end{align*} [/step]

[step:Show that $P^{-1}gP \in \operatorname{U}(n)$ for every $g \in G$]Let $g \in G$ and $\mathbf{u}, \mathbf{u}' \in \mathbb{C}^n$. Compute: \begin{align*} (P^{-1}gP\mathbf{u},\; P^{-1}gP\mathbf{u}')_0 &= \langle P \cdot P^{-1}gP\mathbf{u},\; P \cdot P^{-1}gP\mathbf{u}' \rangle && (\diamond) \text{ with } \mathbf{x} := P^{-1}gP\mathbf{u},\, \mathbf{y} := P^{-1}gP\mathbf{u}' \\ &= \langle gP\mathbf{u},\; gP\mathbf{u}' \rangle && (P \cdot P^{-1} = I_n) \\ &= \langle P\mathbf{u},\, P\mathbf{u}' \rangle && (G\text{-invariance of } \langle \cdot, \cdot \rangle) \\ &= (\mathbf{u}, \mathbf{u}')_0 && (\diamond) \text{ with } \mathbf{x} := \mathbf{u},\, \mathbf{y} := \mathbf{u}'. \end{align*} So $P^{-1}gP$ preserves $(\cdot, \cdot)_0$ for every $g \in G$. By the matrix characterisation of unitarity — a matrix $A \in \operatorname{GL}_n(\mathbb{C})$ preserves $(\cdot, \cdot)_0$ iff $A^* A = I_n$ iff $A \in \operatorname{U}(n)$ — we conclude $P^{-1}gP \in \operatorname{U}(n)$. Therefore \begin{align*} P^{-1} G P \leq \operatorname{U}(n). \end{align*} Setting $Q := P^{-1}$, we have $QGQ^{-1} = P^{-1} G P \leq \operatorname{U}(n)$, exhibiting $G$ as conjugate to a subgroup of $\operatorname{U}(n)$.[/step]

[guided]The whole proof is a single idea: **a $G$-invariant inner product becomes the standard inner product after a change of basis**. The unitary group $\operatorname{U}(n)$ is by definition the matrices that preserve the standard Hermitian inner product $(\cdot, \cdot)_0$. So if our $G$-invariant inner product $\langle \cdot, \cdot \rangle$ were the standard one, every $g \in G$ would already be unitary by virtue of the $G$-invariance condition $\langle g\mathbf{v}, g\mathbf{w}\rangle = \langle \mathbf{v}, \mathbf{w}\rangle$. The change of basis $P$ is the dictionary translating between the two inner products. **Why is the dictionary identity $\langle \mathbf{x}, \mathbf{y} \rangle = (P^{-1}\mathbf{x}, P^{-1}\mathbf{y})_0$ true?** Expand any $\mathbf{x} \in \mathbb{C}^n$ in the basis $(\mathbf{f}_i)$: $\mathbf{x} = \sum_i x_i \mathbf{f}_i$. The coordinate vector of $\mathbf{x}$ in the basis $(\mathbf{f}_i)$ is $P^{-1}\mathbf{x}$ (because $P$ takes the standard basis to $(\mathbf{f}_i)$, so $P^{-1}$ takes vectors to their $(\mathbf{f}_i)$-coordinates). Hence $x_i = (P^{-1}\mathbf{x})_i$. By orthonormality of $(\mathbf{f}_i)$ for $\langle \cdot, \cdot \rangle$, \begin{align*} \langle \mathbf{x}, \mathbf{y} \rangle = \sum_i x_i \overline{y_i} = \sum_i (P^{-1}\mathbf{x})_i \overline{(P^{-1}\mathbf{y})_i} = (P^{-1}\mathbf{x},\, P^{-1}\mathbf{y})_0. \end{align*} **Why does $P^{-1}gP$ work, not $PgP^{-1}$?** Because the dictionary $(\diamond)$ has the form $\langle P\mathbf{x}, P\mathbf{y}\rangle = (\mathbf{x},\mathbf{y})_0$ — that is, $P$ pulls $\langle \cdot, \cdot\rangle$ back to $(\cdot, \cdot)_0$. To make a matrix $A$ preserve $(\cdot, \cdot)_0$, we must arrange that $P A$ preserves $\langle \cdot, \cdot\rangle$ (so that the dictionary translates back to $(\cdot, \cdot)_0$ on both sides). The simplest $A$ arising from $g$ that has this property is $A = P^{-1} g P$: then $P A = g P$, and $g$ preserves $\langle \cdot, \cdot\rangle$ by hypothesis, so $\langle g P \mathbf{x}, g P \mathbf{y}\rangle = \langle P\mathbf{x}, P\mathbf{y}\rangle = (\mathbf{x}, \mathbf{y})_0$. Tracing back, $(A\mathbf{x}, A\mathbf{y})_0 = \langle P A \mathbf{x}, P A \mathbf{y}\rangle = (\mathbf{x}, \mathbf{y})_0$. The hypothesis "$G$ finite" is consumed entirely inside Weyl's unitary trick — its averaging argument requires a finite group. The remainder of the proof is general linear algebra, valid for any finite-dimensional complex inner-product space.[/guided]