MOSFET drain current is:
$$ I_D = k(V_{GS} - V_{TH})^2 (1+\lambda V_{DS}) $$
for \$V_{GS} > V_{TH}\$, where \$\lambda\$ is the fractional change in current due to channel length modulation, and:
$$ k = \frac{1}{2}\mu C_{ox}\frac{W}{L} $$
For a CMOS inverter, the condition \$V_{OUT} = V_{IN}\$ can be produced by manually adjusting gate potential, or by simply joining OUT to IN. Both will produce the same conditions:
simulate this circuit – Schematic created using CircuitLab
In both cases and for both transistors, \$V_{DS} = V_{GS}\$, and the expression for \$I_D\$ becomes:
$$ I_D = k(V_{DS} - V_{TH})^2 (1+\lambda V_{DS}) $$
Each transistor has its own distinct properties \$k\$, \$\lambda\$ and \$V_{TH}\$, and its own \$V_{DS}\$. For the N-channel one I'll give those names a "1" suffix, referring to M1:
$$ I_{D1} = k_1(V_{DS1} - V_{TH1})^2 (1+\lambda_1 V_{DS1}) $$
For the upper P-channel MOSFET the variable suffix will be "2":
$$ I_{D2} = k_2(V_{DS2} - V_{TH2})^2 (1+\lambda_2 V_{DS2}) $$
Note: for the P-channel transistor, technically \$I_D\$, \$V_{TH2}\$ and \$V_{DS2}\$ should be negative, but I can't be bothered to be technically correct. You should still check that the application of KCL and KVL is correct, though.
KCL is easy:
$$
\begin{aligned}
I_{D1} &= I_{D2} \\ \\
k_1(V_{DS1} - V_{TH1})^2 (1+\lambda_1 V_{DS1}) &= k_2(V_{DS2} - V_{TH2})^2 (1+\lambda_2 V_{DS2}) \\ \\
\end{aligned}
$$
A second constraint is set by KVL, for which the two transistors' \$V_{DS}\$ should add up to equal the supply voltage \$V_P\$:
$$
\begin{aligned}
V_{DS1} + V_{DS2} &= V_P \\ \\
V_{DS2} &= V_P - V_{DS1} \\ \\
\end{aligned}
$$
Substituting that expression for \$V_{DS2}\$ into the KCL equation yields an equation in terms of only one unknown, \$V_{DS1}\$, and a few known constants:
$$
\begin{aligned}
k_1(V_{DS1} - V_{TH1})^2 (1+\lambda_1 V_{DS1}) &= k_2((V_P - V_{DS1}) - V_{TH2})^2 (1+\lambda_2 (V_P - V_{DS1})) \\ \\
\end{aligned}
$$
While I won't do it here, I think that you can demonstrate that \$V_{DS1}\$ has a single real solution lying between 0V and \$V_P\$.
As for the question of saturation, the criteria that both transistors must meet (to be saturated) are:
$$
\begin{aligned}
V_{GS} &> V_{TH} \\ \\
V_{DS} &> V_{GS} - V_{TH}
\end{aligned}
$$
Take another look at the self-biased push-pull pair, with some annotations to visualise what's coming:
simulate this circuit
Given that for both transistors \$V_{DS} = V_{GS}\$ (because \$V_{OUT}=V_{IN}\$), then if the first criterion is met, then the second is too. We need only focus on the first, that \$V_{GS}\$ exceeds threshold \$V_{TH}\$:
$$
\begin{aligned}
V_{GS1} &> V_{TH1} \\ \\
V_{GS2} &> V_{TH2} \\ \\
\end{aligned}
$$
That's only possible if supply \$V_P\$ is large enough to permit it:
$$ V_P > V_{TH1} + V_{TH2} $$
One things is certain, if that condition is not met, then at least one of the transistors will not be saturated.
Even if this condition for \$V_P\$ is met, at first glance that may not be a guarantee that both transistors will always be saturated. My first intuition is that it might be possible for one to be saturated while the other isn't, and I would imagine it depends on how similar the two MOSFETs are in terms of their characteristic parameters. If they are similar, then we would find \$V_{GS1} \approx V_{GS2}\$, with \$V_{IN} = V_{OUT} \approx \frac{V_P}{2}\$. Then, if \$V_P\$ is large enough, both conditions \$V_{GS1} > V_{TH1}\$ and \$V_{GS2} > V_{TH2}\$ will be true.
In general though, for the case where \$V_P > V_{TH1} + V_{TH2}\$, permitting the situation where both transistors saturate, and \$V_{OUT}=V_{IN}\$, the question remains is this a guarantee of saturation of both devices?
Saturation of both requires that \$V_{GS1}\$ and \$V_{GS2}\$ must exceed their corresponding thresholds \$V_{TH}\$, but since we know their sum \$V_{GS1} + V_{GS2} = V_P\$, we can write a general constraint in terms of \$V_{GS1}\$ alone:
$$ V_{TH1} < V_{GS1} < (V_P - V_{TH2}) $$
Since \$V_{DS1} = V_{GS1}\$ this is also the bounds for \$V_{DS1}\$:
$$ V_{TH1} < V_{DS1} < (V_P - V_{TH2}) $$
To determine if a given pair of MOSFETs with different characteristics will settle at \$V_{OUT}=V_{IN}\$ with both saturated, it will be necessary to combine these conditions with the long expression for \$V_{DS1}\$ from before:
$$
\begin{aligned}
k_1(V_{DS1} - V_{TH1})^2 (1+\lambda_1 V_{DS1}) &= k_2((V_P - V_{DS1}) - V_{TH2})^2 (1+\lambda_2 (V_P - V_{DS1})) \\ \\
V_{DS1} &> V_{TH1} \\ \\
V_{DS1} &< (V_P - V_{TH2}) \\ \\
\end{aligned}
$$
If you can show that there's a (unique) solution for \$V_{DS1}\$, and that it always falls in the interval \$(V_{TH1}, V_P - V_{TH2})\$, then you've proved that the two transistors will always be saturated when \$V_{OUT}=V_{IN}\$.
In the meantime, though, it's safe to assume that if
The MOSFETs are sufficiently close complements of each other, with very similar characteristics and behaviour
\$V_P > V_{TH1} + V_{TH2}\$
then both will saturate, and \$V_{OUT}=V_{IN}=\frac{V_P}{2}\$
