Solution (source code)

= Solution

For nonparallel direction vectors, let $\mathbf n=(\mathbf b\times\mathbf d)/|\mathbf b\times\mathbf d|$. Every joining <vector> has the same <normal> component, $\mathbf R\cdot\mathbf n=(\mathbf a-\mathbf c)\cdot\mathbf n$. Its <norm> is therefore at least the absolute value of that component. Minimizing $|\mathbf R|^2$ makes it perpendicular to both directions, and hence parallel to $\mathbf n$, attaining the bound:
$$
\boxed{d_{\min}=\left|(\mathbf a-\mathbf c)\cdot\mathbf n\right|.}
$$
The expression printed without absolute values is an oriented separation, not a nonnegative distance. Reversing either line's parametrization can flip its sign. For parallel nonzero directions, this <normal> is undefined and the distance is instead $| (\mathbf a-\mathbf c)\times\mathbf b|/|\mathbf b|$.