I want to find out the intersection of two surfaces that exist in $n$-dimensions. These surfaces are defined by a system of $2$ linear equations.
$$A_1x+B_1y+C_1z+\cdots = D_1$$
$$A_2x+B_2y+C_2z+\cdots = D_2$$
But first, how would these surfaces look like? In $3$-D space, they would be planes, but in dimensions higher than $3$ are these not planes? Can we generalize a way to find the intersection of such linear expressions in some way?
Also can this problem be solved using some tools, preferably Matlab?
Suppose we would like to find the intersection of $2$ hyperplanes in $\mathbb R^n$
$$\begin{array}{rl} a_{11} x_1 + a_{12} x_2 + \cdots + a_{1n} x_n &= b_1\\ a_{21} x_1 + a_{22} x_2 + \cdots + a_{2n} x_n &= b_2\end{array}$$
In matrix form,
$$\begin{bmatrix} — \mathrm a_1^\top —\\ — \mathrm a_2^\top —\end{bmatrix} \mathrm x = \begin{bmatrix} b_1\\ b_2\end{bmatrix}$$
or, more economically, $\rm A x = b$, where $\rm A$ is a $2 \times n$ matrix. Without loss of generality, we assume that both hyperplanes are in Hessian normal form, i.e., $\| \mathrm a_1 \|_2 = \| \mathrm a_2 \|_2 = 1$.
If vectors $\rm a_1$ and $\rm a_2$ are collinear, then the $2$ given hyperplanes are either parallel or overlapping, neither of which are very interesting. Thus, we assume that vectors $\rm a_1$ and $\rm a_2$ are not collinear, i.e., matrix $\rm A$ has full row rank. Let $\theta := \arccos ( \mathrm a_1^\top \mathrm a_2 )$ be the angle formed by vectors $\rm a_1$ and $\rm a_2$.
The solution set of the linear system $\rm A x = b$ is a $(n-2)$-dimensional affine space. To find this affine space, we must find a particular solution and the null space of $\rm A$. One particular solution would be the least-norm solution, i.e., the one closest to the origin
$$\begin{array}{rl} \mathrm x_{\text{LN}} := \mathrm A^\top \left( \mathrm A \mathrm A^\top \right)^{-1} \mathrm b &= \begin{bmatrix} | & | \\ \mathrm a_1 & \mathrm a_2\\ | & |\end{bmatrix} \begin{bmatrix} \| \mathrm a_1 \|_2^2 & \langle \mathrm a_1, \mathrm a_2 \rangle\\ \langle \mathrm a_2, \mathrm a_1 \rangle & \| \mathrm a_2 \|_2^2\end{bmatrix}^{-1} \begin{bmatrix} b_1\\ b_2\end{bmatrix}\\ &= \begin{bmatrix} | & | \\ \mathrm a_1 & \mathrm a_2\\ | & |\end{bmatrix} \begin{bmatrix} 1 & \cos (\theta)\\ \cos (\theta) & 1\end{bmatrix}^{-1} \begin{bmatrix} b_1\\ b_2\end{bmatrix}\\ &= \frac{1}{\sin^2 (\theta)} \begin{bmatrix} | & | \\ \mathrm a_1 & \mathrm a_2\\ | & |\end{bmatrix} \begin{bmatrix} 1 & -\cos (\theta)\\ -\cos (\theta) & 1\end{bmatrix} \begin{bmatrix} b_1\\ b_2\end{bmatrix}\\ &= \left(\frac{b_1 - b_2 \cos (\theta)}{\sin^2 (\theta)}\right) \begin{bmatrix} | \\ \mathrm a_1 \\ | \end{bmatrix} + \left(\frac{b_2 - b_1 \cos (\theta)}{\sin^2 (\theta)}\right) \begin{bmatrix} | \\ \mathrm a_2 \\ | \end{bmatrix}\end{array}$$
where $\sin (\theta) \neq 0$ is ensured by the non-collinearity of vectors $\rm a_1$ and $\rm a_2$.
To find the $(n-2)$-dimensional null space of $\rm A$, we solve the homogeneous linear system $\rm A x = 0_2$. Introducing an $n \times n$ permutation matrix $\rm P$ with certain desired properties, we reorder the columns of $\rm A$, i.e., $\rm A P P^\top x = 0_2$, and obtain a homogeneous linear system in $\rm y:= P^\top x$
$$\rm A P y = 0_2$$
If $\rm P$ is chosen wisely, then there is a $2 \times 2$ matrix $\rm E$, a product of elimination matrices, that puts $\rm A P$ in reduced row echelon form (RREF), i.e.,
$$\rm E A P = \begin{bmatrix} \mathrm I_2 & \mathrm F\end{bmatrix}$$
where $\rm F$ is a $2 \times (n-2)$ matrix. Hence, the null space of $\rm A P$ is given by
$$\left\{ \begin{bmatrix} -\mathrm F \\ \mathrm I_{n-2} \end{bmatrix} \eta : \eta \in \mathbb R^{n-2} \right\}$$
and, thus, the null space of $\rm A$ is given by
$$\left\{ \mathrm P \begin{bmatrix} -\mathrm F \\ \mathrm I_{n-2} \end{bmatrix} \eta : \eta \in \mathbb R^{n-2} \right\}$$
Lastly, the intersection of the $2$ given hyperplanes is given by
$$\left\{ \color{blue}{\mathrm x_{\text{LN}} + \mathrm P \begin{bmatrix} -\mathrm F \\ \mathrm I_{n-2} \end{bmatrix} \eta} : \eta \in \mathbb R^{n-2} \right\}$$
In MATLAB, we can use function rref to find permutation matrix $\rm P$ and matrix $\rm F$. Of course, we can also use function null to find an orthonormal basis for the null space of $\rm A$.