The equation of the hyperboloid to be identified is of the form
$ (r - r_0)^T Q (r - r_0) = 1 \hspace{15pt}(1) $
where $r$ is the position coordinate vector of a point on the hyperboloid, $r = [x, y, z]^T $, and $r_0 = [x_0, y_0, z_0] $ is the center of the hyperboloid, and symmetric $3 \times 3 $ matrix $Q$ has two positive eigenvalues and one negative eigenvalue.
Setting $r = r_i(t) = P_i + t d_i $ then we have
$ ( P_i - r_0 + t d_i )^T Q (P_i - r_0 + t d_i) = 1 \hspace{5pt} \forall \ t \in \mathbb{R},\ i=1,2,3 \hspace{15pt} (2) $
Since $P_i$ is on the line, and the line is contained in the hyperboloid, then $P_i$ satisfies its equation, i.e.
$ (P_i - r_0)^T Q (P_i - r_0) = 1 \hspace{15pt} (3) $
This reduces the equation $(2)$ to
$ t^2 d_i^T Q d_i + 2 t d_i^T Q (P_i - r_0) = 0 \hspace{15pt} (4)$
Dividing by $t$
$ t d_i^T Q d_i + 2 d_i^T Q (P_i - r_0) = 0 \hspace{15pt} (5)$
Since this is an identity that holds for any $t$ then we must have
$ d_i^T Q (P_i - r_0) = 0 \hspace{15pt} (6)$
and
$ d_i^T Q d_i = 0 \hspace{15pt} (7)$
Now matrix $Q$ is symmetric, and is therefore determined by $6$ unknowns,
$ Q = \begin{bmatrix} a && d && e \\ d && b && f \\ e && f && c \end{bmatrix} \hspace{15pt} (8) $
From equation $(7), (8)$ we can write $3$ linear homogeneous equations in the variables $a,b,c,d,e,f$. The solution will give result in matrix $Q$ being expressed in terms of three basis matrices $Q_1, Q_2, Q_3$ that are determined by solving the $3 \times 6$ homogenous system.
i.e.
$ Q = \alpha Q_1 + \beta Q_2 + \gamma Q_3 $
where $Q_1, Q_2, Q_3 $ are now known, and $\alpha, \beta, \gamma $ are yet to be determined.
Now, from equation (6),
$ d_i^T Q P_i = d_i^T Q r_0 \hspace{15pt} (9)$
Since the three given lines are skew, then matrix $[d_1, d_2, d_3]$ spans $\mathbf{R}^3$, therefore, we can express $ Q r_0$ as
$ Q r_0 = V u $
where $ V = [d_1, d_2, d_3] $
Let $ \Omega = [\alpha, \beta, \gamma]^T$. Now equation (9) becomes
$ [ d_i^T Q_1 P_i , d_i^T Q_2 P_i, d_i^T Q_3 P_i ]^T \Omega = d_i^T V u , i = 1,2,3 \hspace{15pt} (10) $
Iterating equation $(10)$ for $i=1,2,3$ gives the matrix-vector equation
$ F \Omega = G u $
where $F_{ij} = d_i^T Q_j P_i $ and $G = V^T V $
From this, it follows that
$ u = G^{-1} F \Omega $
Hence,
$ Q r_0 = V G^{-1} F \Omega \hspace{15pt}(11)$
Now back to equation (3),
$ (P_i - r_0)^T Q (P_i - r_0) = 1, i = 1, 2, 3 \hspace{15pt} (3) $
Subtracting pairs of equations corresponding to $i=1,2$ and $i=1,3$ gives
$ P_1^T Q P_1 - P_2^T Q P_2 - 2 (P_1 - P_2)^T Q r_0 = 0 \hspace{15pt} (12) $
and
$ P_1^T Q P_1 - P_3^T Q P_3 - 2 (P_1 - P_3)^T Q r_0 = 0 \hspace{15pt} (13) $
Using $(11)$ in $(12), (13)$:
$ P_1^T Q P_1 - P_2^T Q P_2 - 2 (P_1 - P_2)^T V G^{-1} F \Omega = 0 \hspace{15pt} (14) $
and
$ P_1^T Q P_1 - P_3^T Q P_3 - 2 (P_1 - P_3)^T V G^{-1} F \Omega = 0 \hspace{15pt} (15) $
Note $ P_i^T Q P_i = P_i^T [\alpha Q_1 + \beta Q_2 + \gamma Q_3] P_i $
Thus, equations $(14)$ and $(15)$ are two homogeneous linear equations in $ \Omega = [\alpha, \beta, \gamma]^T $
Solving them using Gauss-Jordan elimination results in a solution of the form
$ \Omega = \Omega_0 t $
where $\Omega_0 = [\alpha_0, \beta_0, \gamma_0] $ is a constant vector, and $ t \in \mathbb{R} $ is yet to be found.
Now we have
$ Q = t Q_0$ , where $ Q_0 = \alpha_0 Q_1 + \beta_0 Q_2 + \gamma_0 Q_3$
Now we can find $r_0$ from equation $(11)$,
$ t Q_0 r_0 = V G^{-1} F \Omega_0 t \hspace{15pt}(11')$
Hence
$ r_0 = Q_0^{-1} V G^{-1} F \Omega_0 $
Finally, the value of scalar $t$ can found from equation $(3)$
$ t = \dfrac{1}{ (P_i - r_0)^T Q_0 (P_i - r_0) } $
Now the hyperboloid that contains the three lines is fully specified.