Fdtxjr

Let $ds:=sigma_x dx+ sigma_y dy + sigma_z dz$. Then squaring $ds$, we get

$$
ds^2=sigma_x^2dx^2+ sigma_y^2dy^2+sigma_z^2dz^2 + (sigma_xsigma_y+ sigma_ysigma_x)dxdy + (sigma_xsigma_z+sigma_zsigma_x)dxdz+(sigma_ysigma_z+sigma_zsigma_y)dydz
$$

We get a lot of cross-terms, but if $sigma_x, sigma_y $ and $sigma_z$ are Pauli matrices with these properties:

$$
sigma_x^2=sigma_y^2=sigma_z^2=1\
sigma_xsigma_y=-sigma_ysigma_x\
sigma_xsigma_z=-sigma_zsigma_x\
sigma_zsigma_y=-sigma_ysigma_z
$$

Then, the metric is now the usual Euclidian metric.

$$
ds^2=dx^2+dy^2+dz^2
$$

Therefore,

$$
sigma_x dx+ sigma_y dy + sigma_z dz=sqrtdx^2+dy^2+dz^2
$$

This is not just a fluke result, because one can repeat the same exercise with

$$
ds:=alpha_x dx+ alpha_y dy + alpha_z dz-alpha_tdt
$$

and, squaring $ds$ one will get

$$
ds^2:=dx^2+ dy^2 + dz^2-dt^2
$$

provided that $alpha_x,alpha_y, alpha_y$ and $alpha_t$ and the Dirac matrices.

What is the source of this connection between eliminating the cross-terms within metrics and matrices relevant to quantum mechanics?

In the case of a generalized metric

$$
ds^2=g_munu dx^mu dx^nu
$$

with a symmetry requirement $g_mu nu=g_nu mu$ (the case of general relativity), the cross terms commute and this seems to be precisely the condition which avoid the need for matrix.

EDIT

Thinking about it some more, it seems general relativity does have a matrix-like construction. The matrices must have the following properties:

Start with

$$
ds:=A_xdx+A_ydy+A_zdz+A_tdt
$$

Squaring $ds$ produces a GR compliant tensor if the matrices have these properties:

$$
A_x^2=g_xx\
A_y^2=g_yy\
A_z^2=g_zz\
A_t^2=-g_tt\
A_xA_y=A_xA_y=frac12g_xy=frac12g_yx\
A_xA_z=A_xA_z=frac12g_xz=frac12g_zx\
A_yA_z=A_zA_y=frac12g_yz=frac12g_zy\
A_xA_t=A_tA_x=frac12g_xt=frac12g_tx\
A_yA_t=A_tA_y=frac12g_yt=frac12g_ty\
A_zA_t=A_tA_z=frac12g_zt=frac12g_tz
$$

I am not sure if this gives quantum mechanically relevant matrices or not.

But we get

$$
A_xdx+A_ydy+A_zdz+A_tdt:=sqrtg_munu dx^mu dx^nu
$$

Your work with infinitesimals in $Bbb R^3$ is a special case of $(vecacdotvecsigma)times(vecbcdotvecsigma)=(vecacdotvecb)I_2+i(vecatimesvecb)cdotvecsigma$; your work in $Bbb R^4$ follows from an analogous identity with gamma matrices, $not a,,not b=2acdot b$. If you want a connection of such identities to your work with $ds^2$, it's simply that these metrics have the same symmetries, be they rotational or Lorentz, as the matrices you're working with. Hence second-order PDEs such as Laplace's equation $nabla^2phi=0$ and the massless Klein-Gordon equation $squarephi=0$ require cross terms' cancellation if they're to be derived from a first-order PDE such as the Dirac equation, and this is why such equations have to use matrices, spinors etc.