How to derive the reduced density matrix in a mathematically precise and correct way

Consinder a system $$\mathcal{AB}$$, consisting of two subsystems $$\mathcal{A}$$ and $$\mathcal{B}$$ (with hilbert-spaces $$\mathcal{H_A}$$ and $$\mathcal{H_B}$$). I won't be introducing all formalisms of quantum mechanics here - I assume it is clear what is used in the following.

I am now trying to proof the following (=derive the reduced density matrix of the system $$\mathcal{A}$$):

If an observer only has only access to $$\mathcal{A}$$, there is a simpler mathematical object, that contains all information the observer can learn, than the density matrix $$\rho$$, namely the reduced density matrix.

Schlosshauer (ISBN: 978-3540357735) derives it on page 46 in the following way:

Following the above request, one must be able to calculate the expectation value of an Observable $$O_{\mathcal{A}}\otimes I$$ from this 'reduced density matrix'. We calculate with $$\{\psi_k\}$$ and $$\{\phi_l\}$$ orthonormal bases of the hilbert spaces introduced above:

$$\langle O\rangle=\text{Tr}(\rho O)=\\=\sum_{kl}\langle\phi_l|\langle\psi_k|\rho(O_{\mathcal{A}}\otimes I)|\psi_k\rangle|\phi_k\rangle=\\= \sum_k\langle\psi_k|\left(\sum_k\langle\phi_l|\rho|\phi_k\rangle\right)O_{\mathcal{A}}|\psi_k\rangle=\\= \sum_k\langle\psi_k|(\text{Tr}_{\mathcal{B}}\rho)O_{\mathcal{A}}|\psi_k\rangle=\\=\text{Tr}_{\mathcal{A}}(\rho_{\mathcal{A}}O_{\mathcal{A}})$$

with $$\rho_{\mathcal{A}}$$ the (as in the above calculation) defined reduced density matrix of the system $$\mathcal{A}$$.

What concerns do I have?

1. Wrong definition of the reduced density matrix?: How can one calculate $$\rho|\psi_k\rangle$$? The dimensions of those elements (if seen as matrix and vector) don't match - in linear algebra this operation is not defined / would be forbidden!
2. If the trace is defined by $$\text{Tr}(\gamma)=\langle\alpha_i|\gamma|\alpha_i\rangle$$ and $$\langle\phi_l|\langle\psi_k|$$ above means $$\langle\phi_l|\otimes\langle\psi_k|$$ then the definition of the trace is used in the wrong way, because $$(|\phi_l\rangle\otimes|\psi_k\rangle)^{\dagger}=\langle\phi_l|\otimes\langle\psi_k|$$ - so the second line above should rather read $$\sum_{kl}\langle\psi_k|\langle\phi_l|\rho(O_{\mathcal{A}}\otimes I)|\psi_k\rangle|\phi_k\rangle.$$

What I'd like to know

So ultimately I fail to derive the reduced density matrix. I would gladly appreciate any hint on books that derive it properly.

I would like to ask anyone to explain to me, if my concerns are correct or if I am making a mistake there. If somebody is able to clear up Schlosshauers derivation for me, I'd also appreciate if you gave some references to the used mathematics. Thank you!

In this post all manipulations are at a rather formal level and no deep mathematical discussion takes places. For that see the references at the end.

To start, consider a bipartite Hilbert space $$H=H_A\otimes H_B$$ and define the operator $$\mathbb I_A \otimes |\psi\rangle : H_A\longrightarrow H$$ for some $$|\psi\rangle \in H_B$$ and the identity operator $$\mathbb I_A$$ on $$H_A$$, such that for $$|\varphi\rangle \in H_A$$ it holds that

$$\left(\mathbb I_A \otimes |\psi\rangle \right) |\varphi\rangle := |\varphi\rangle \otimes |\psi\rangle \quad . \tag{1}$$

Its adjoint is the following operator $$\left(\mathbb I_A \otimes \langle \psi| \right) : H \longrightarrow H_A$$, where for $$|\varphi\rangle \in H_A, |\phi\rangle \in H_B$$ we have

$$\left(\mathbb I_A \otimes \langle \psi| \right) (|\varphi\rangle \otimes|\phi\rangle) = |\varphi\rangle \langle \psi|\phi \rangle \quad . \tag{2}$$ For a density operator $$\rho$$ on $$H$$ we define the partial trace of $$\rho$$ as $$\rho_A:=\mathrm{Tr_B}\,\rho := \sum\limits_{k \in K} \left(\mathbb I_A \otimes \langle \psi_k| \right) \rho\left(\mathbb I_A \otimes |\psi_k\rangle \right) \quad , \tag{3}$$ where $$\{|\psi_k\rangle\}_{k\in K}$$ for an index set $$K$$ is an orthonormal basis in $$H_B$$. Note that $$\rho_A$$ is a linear operator on $$H_A$$. We can similarly define an operator $$\rho_B:=\mathrm{Tr}_A\,\rho$$.

Instead of $$(3)$$ people often write, as a shorthand, something like

$$\rho_A= \sum\limits_{k \in K} \langle \psi_k| \rho |\psi_k\rangle \quad , \tag{4}$$ which can indeed be confusing. Regarding your second concern: I agree (cf. the calculations below), but I guess this could also be a matter of definition and convention, although I've not seen this before. In the end, something like $$\langle \varphi|\otimes \langle \psi|$$ should denote an element of $$H^*$$ with $$\langle \varphi|\otimes \langle \psi| \left(|\alpha\rangle \otimes |\beta\rangle\right) = \langle \varphi|\alpha\rangle \langle \psi|\beta\rangle$$, for $$|\varphi\rangle,|\alpha\rangle \in H_A$$ and $$|\psi\rangle, |\beta\rangle \in H_B$$. So if we put elements of $$H_A$$ in the first slot of $$\otimes$$, i.e. on the left of side the tensor product symbol, then it seems natural to me to write $$\langle \varphi|\otimes \langle \psi| \in H^*$$ instead of $$\langle \psi|\otimes \langle \varphi|$$.

To proceed, let us now do the explicit calculations: For an operator $$O_A$$ on $$H_A$$ we compute \begin{align} \mathrm{Tr^{(A)}} \rho_A\,O_A &= \sum\limits_{j\in J} \langle \varphi_j | \rho_A \, O_A|\varphi_j\rangle\\ &= \sum\limits_{j\in J} \langle \varphi_j | \sum\limits_{k \in K} \left(\mathbb I_A \otimes \langle \psi_k| \right) \rho\left(\mathbb I_A \otimes |\psi_k\rangle \right) \, O_A |\varphi_j\rangle \\ &=\sum\limits_{j\in J}\langle \varphi_j | \sum\limits_{k \in K} \left(\mathbb I_A \otimes \langle \psi_k| \right) \rho\, O_A |\varphi_j\rangle \otimes |\psi_k\rangle\tag{5}\\ &=\sum\limits_{j\in J}\sum\limits_{k\in K}\langle \varphi_j|\otimes \langle \psi_k| \,\rho\,O \,|\varphi_j\rangle \otimes |\psi_k\rangle \\ &= \mathrm{Tr}\rho\, O \quad . \end{align} Here, $$\{|\varphi_j\rangle\}_{j\in J}$$ denotes an orthonormal basis in $$H_A$$, $$\mathrm{Tr}^{(A)}$$ the trace operation on $$H_A$$ and $$O:= O_A\otimes \mathbb I_B$$.

For a more rigorous treatment, see for example S. Attal. Tensor products and partial traces. Lecture Notes, especially section $$2.3$$. or Michael M. Wolf. Mathematical Introduction to Quantum Information Processing. Lecture notes, especially theorem 1.35. You can find a pdf for the first reference here and for the second here.

2: This is not correct. Both vectors act on different subspaces of the Hilbert space, you can commute them as you like. The order does not matter. However the hermitian conjugate of a vector will always act on its own dual space. Hence,

$$\langle \psi_k | \langle \phi_l | \psi_k \rangle | \phi_l \rangle = \langle \psi_k | \psi_k \rangle \langle \phi_l | \phi_l \rangle$$