Bdtyktl

How did we derive that the state we get by $n$ qubits is their tensor product? You can use $n=2$ in the explanation for simplicity.

This assertion is an axiom of quantum mechanics. It appears, for example, as Postulate 4 on page 94 of Nielsen and Chuang.

It is considered to be one of the Dirac-von Neuman axioms of quantum mechanics.

However, when the quantum system's Hilbert space is defined as the space of $l^2$ functions on some set, then when the set is a Cartesian product, the corresponding Hilbert space becomes a tensor product. For example, when you put spins on a lattice $Lambda$ as in Kitaev's surface code, then the Systems Hilbert space is:
$$l^2(underline{2}^{Lambda}) = bigotimes_{Lambda} mathbb{C}^2$$
Or, when you quantize a system of two particles moving on the line $mathbb{R}$, the individual particle Hilbert space is $ L^2(mathbb{R})$ and the composite system Hilbert space is $L^2(mathbb{R} times mathbb{R}) $. It is known that
$$L^2(mathbb{R} times mathbb{R}) = L^2(mathbb{R}) bigotimes L^2(mathbb{R}),$$
Even though a specific basis of the tensor product $L^2(mathbb{R}) bigotimes L^2(mathbb{R})$ does not include the entangled states in $L^2(mathbb{R} times mathbb{R})$, these states can be constructed from combinations of this basis.

However, the above explanation doesn't prove the assertions; it only transfers the problem to the classical spaces before quantization. Then why should they be taken as Cartesian products?

Many authors analyzed this problem from the quantum logic point of view, see for example
Aerts and Daubechies. These attempts intend to find a logical equivalent to the above axiom and not to prove it.

Perhaps it helps to take a step back and start with something simpler:

why do we tabulate probability amplitudes for state vectors and unitaries?

For a single quantum system with $d$ distinct states, labelled 0 to $d-1$, we associate a complex number $a_i$ with the probability amplitude for being in state $i$. For a unitary, we associate a probability amplitude $U_{ij}$ with transforming an input $j$ into an output $i$. We choose to tabulate these as
$$
U=left(begin{array}{cccc}
U_{00} & U_{21} & ldots & U_{0,d-1} \
U_{10} & U_{11} & ldots \
vdots & vdots & ddots \
U_{d-1,0} & U_{d-1,1} & ldots & U_{d-1,d-1}
end{array}right)qquad |psirangle=left(begin{array}{c} a_0 \ a_1 \ vdots \ a_{d-1} end{array}right)
$$
for the simple reason that calculating $U|psirangle$ takes are of the two axioms of quantum theory that tell us how to calculate the output probability amplitudes:

• for independent events, the probability amplitude for both events happening is the product of the individual probability amplitudes.


• for mutually exclusive events, the probability amplitude for either of the two events to happen is the sum of their probability amplitudes.

This puts is in a position to answer

Why is the state of multiple qubits given by their tensor product?

If two quantum systems have distinguishable states $i$ and $j$, then we can choose to label the states of the composition system by $i,j$, meaning that the first system is in state $i$ and the second system is in state $j$. Then, if the two systems start in states
$$
|psirangle=left(begin{array}{c} a_0 \ a_1 \ vdots \ a_{d-1} end{array}right) qquad |phirangle=left(begin{array}{c} b_0 \ b_1 \ vdots \ b_{d-1} end{array}right),
$$
then what's the probability amplitude for the overall system to be in $1,0$, for example? By the independence axiom, it must be $a_1b_0$. Hence, if we tabulate the probability amplitudes using the ordering $00,01,02,03,ldots (0,d-1),10,11,ldots,(1,d-1),20,21,ldots,(d-1,d-1)$, the column vector turns out to be exactly $|psirangleotimes|phirangle$, basically by definition of what the tensor product is. If you do similar things with unitaries $U$ and $V$, you find the composite action is $Uotimes V$, and everything is internally consistent, meaning that the outcome of the combined unitary evolution would be
$$
(Uotimes V)(|psirangleotimes|phirangle)=(U|psirangle)otimes(V|phirangle),
$$
as was explicitly checked here.

I'll attempt here to provide a physical justification for why tensor products are the natural way to describe systems comprised of a number of different subsystems (e.g. a number of qubits).
The takeaway message is that tensor products, while not strictly necessary, are a natural choice when dealing with local operations.

First of all, it is important to notice that it is not strictly necessary to use tensor products at all to describe any kind of quantum mechanical system.
It just so happens that not using them makes things much more complicated and inelegant that they could be.

By this, I mean that we can describe, say, the state of an $N$-qubit system, as a $2^N$-dimensional complex vector (a $2^N$-dimensional qudit if you prefer) with no tensor product structure. If you choose to describe things this way, unitary operations, as well as any other operator acting on states, become $2^Ntimes 2^N$ matrices.

However, if we are talking of "$N$ qubits", instead of simply a single system with many degrees of freedom, we are implicitly assuming that we will deal with local operations. Local operations are operations that act only on some subset of the degrees of freedom of the system, leaving the rest untouched.

More generally, if you have a system which is composed of two subsystems, each of which can be in one of a number of different states, it is only natural to describe the whole system in terms of the states of the single subsystems. Therefore, if the first system can be in one of the states $i=1,2,...,N$ and the second system in one of $j=1,2,...,M$, then a natural choice for the basis states of the whole system is the set of $NM$ pairs $(i,j)$.

Consider as a toy example a system with $N=6$ degrees of freedom (a $6$-dimensional qudit), which is obtained by combining a $3$-dimensional system with a $2$-dimensional one. Without tensor products, it is awkward to describe a local operation acting only on the first subsystem: states would be described as single-index vectors $psi_iinmathbb C$, $i=1,...,6$ and operations as matrices $mathcal U_{ij}$ with $i,j=1,...,6$, and a local operation would be some $mathcal U$ such that $mathcal U_{ij}$ is always zero when $i$ and $j$ are indices corresponding to two different states of the second subsystem, assuming some conventional way to arrange the indices.
On the other hand, if you describe the state as an object $psi_{ij}inmathbb C$ with $i=1,2,3$ describing the states of the first subsystem and $j=1,2$ describing the state of the second subsystem, an operation $mathcal U$ acting locally on the first subsystem is naturally written as one satisfying $mathcal U_{ij,kl}=U_{ij}delta_{kl}$ for some unitary $U$.

One can then proceed with building the mathematical infrastructure to reason formally about tensor product spaces etc., but from a physical perspective, this is all there is to it. It is a convenient way to describe systems which are comprised of subsystems in such a way that it is natural to describe the state of the whole system via the tuples of states of the subsystems.