Let $H \in \mathbb{F}^{m \times n}$ for $\mathbb{F} = \mathbb{R}$ or $\mathbb{F} = \mathbb{C}$. If I understood you right, you want $P \in \mathbb{F}^{n \times n}$ of a maximum rank such that
Let $H = U \Sigma V^*$ be an SVD of $H$, where $U$ and $V$ are unitary or orthogonal, depending on $\mathbb{F}$, of order $m$ and $n$, respectively, and
$$\Sigma = \mathop{\rm diag}(\sigma_1, \dots, \sigma_k) \in \mathbb{R}^{m \times n}, \quad \sigma_1 \ge \cdots \ge \sigma_p > \sigma_{p+1} = \cdots \sigma_k = 0, \quad k = \min\{m,n\}.$$
Let $\Sigma' := \mathop{\rm diag}(\sigma_1', \dots, \sigma_n') \in \mathbb{R}^{n \times n}$ such that
$$\sigma_i' = \begin{cases}
0, & i \le p, \\
1, & i > p.
\end{cases}.$$
Note that $\Sigma \Sigma' = 0$ and $(\Sigma')^2 = \Sigma'$.
We define $P = V \Sigma' V^*$ and check the desired properties:
- $HP = U \Sigma V^* V \Sigma' V^* = U \Sigma \Sigma' V^* = U 0 V^* = 0$,
- $P^2 = V \Sigma' V^* V \Sigma' V^* = V \Sigma' \Sigma' V^* = V (\Sigma')^2 V^* = V \Sigma' V^* = P$.
Also, if $P$ had a higher rank, than $\Sigma \Sigma' \ne 0$, so $P$ is the desired projection.
Of course, if you're doing this in Matlab (or anything else on a computer), be careful what you treat as a zero.
if the columns of matrix $A$ are linearly independent,
the projection of a vector, $b$, onto the column space of A can be computed
as
$$P=A(A^TA)^{-1}A^T$$
From here.
Wiki seems to say the same. It also says here that
The column space of $A$ is equal to the row space of $A^T$.
I'm guessing that
if the rows of matrix $A$ are linearly independent,
the projection of a vector, $b$, onto the row space of A can be computed
as
$$P=A^T(AA^T)^{-1}A$$
Best Answer
We can mimic Householder transformation. Let $y = x_1 + Ax_2$. Define: $P = {\rm I} - yy^T / y^Ty$ (Householder would have factor $2$ in the $y$ part of the expression).
Check:
Your condition: $$Px_1 + PAx_2 = Py = ({\rm I} - yy^T / y^Ty) y = y - yy^Ty/y^Ty = y - y = 0,$$
$P$ is a projection: \begin{align} P^2 &= ({\rm I} - yy^T / y^Ty) ({\rm I} - yy^T / y^Ty) = {\rm I} - yy^T / y^Ty - yy^T / y^Ty + yy^Tyy^T / y^Tyy^Ty \\ &= {\rm I} - 2yy^T / y^Ty + yy^T / y^Ty = {\rm I} - yy^T / y^Ty = P. \end{align}
(if needed) $P$ is an orthogonal projection (condition explained on the previous link): $$P^T = ({\rm I} - yy^T / y^Ty)^T = {\rm I} - yy^T / y^Ty = P.$$
You sure that these are the only conditions?