Honestly, the argument you're making here is a mess - the question is based on bad premises. So let me show you how to do it properly, and hopefully that will resolve your confusion.
You're right that in order for a wavefunction to be normalized, it must satisfy
$$\int_\text{all space} P(x)\mathrm{d}x
= \int_\text{all space} \psi^*(x)\psi(x)\mathrm{d}x = 1\tag{1}$$
But this statement:
On top of that, any $\psi$ can also be expressed as $\psi\psi^∗$
is not correct. Given a function $\psi(x)$, you can write $\psi(x)\psi^*(x)$, but that's a different function.
Anyway, given that your wavefunction can be written
$$\psi(x) = a\phi_1(x) + b\phi_2(x)$$
then you just plug that into the normalization condition (1) and get
$$\int_0^L \bigl(a^* \phi_1^*(x) + b^* \phi_2^*(x)\bigr)\bigl(a \phi_1(x) + b \phi_2(x)\bigr)\mathrm{d}x = 1$$
which expands to
$$\begin{multline}
a^*a \int_0^L \phi_1^*(x)\phi_1(x)\mathrm{d}x
+ a^*b \int_0^L \phi_1^*(x)\phi_2(x)\mathrm{d}x \\
+ b^*a \int_0^L \phi_2^*(x)\phi_1(x)\mathrm{d}x
+ b^*b \int_0^L \phi_2^*(x)\phi_2(x)\mathrm{d}x
= 1
\end{multline}\tag{2}$$
Now you can use the identity
$$\int_0^L\phi_1^*(x)\phi_2(x)\mathrm{d}x = \int_0^L\phi_2^*(x)\phi_1(x)\mathrm{d}x = 0$$
which follows from the fact that $\phi_1$ and $\phi_2$ are orthogonal functions (it's not enough that they are eigenfunctions of an operator, they have to be orthogonal), and the identity
$$\int_0^L\phi_1^*(x)\phi_1(x)\mathrm{d}x = \int_0^L\phi_2^*(x)\phi_2(x)\mathrm{d}x = 1$$
which simply reflects the fact that $\phi_1$ and $\phi_2$ are normalized. (Check for yourself that this is the same as the normalization condition, equation (1).) With these two identities, equation (2) reduces to
$$\lvert a\rvert^2 + \lvert b\rvert^2 = 1$$
The conceptual question I had was that if we have the probability squared here, is it that or the square root of that probability that is your normalization constant?
That all depends, how do you define your normalization constant? It depends on what you're normalizing and how exactly you express it as a function. However you do it, the end requirement for normalization is just that $\lvert a\rvert^2 + \lvert b\rvert^2 = 1$.
As far as using the specific sinusoidal form for the $\phi_i$, you can do that in this case, because you're given enough information to figure out that the eigenfunctions are in fact sinusoidal. But you don't really need to know that they are sinusoidal for the preceding argument to work; all you need to know is that the $\phi_i$s are orthonormal.
Best Answer
Yes, the Slater determinant should have a factor of $1/\sqrt{N!}$ for the state to be properly normalized.