0
$\begingroup$

I try to explain what skew lines are, and how to present them, little by little, calculate different intermediate values, this code works, but it is a frankenstein code, even with comments it's hard for myself to follow it after few days. I wouls like to receive advice to avoid so much code .... Sorry this is the monster:

twolns = RandomReal[{-10, 10}, {2, 2, 
   3}];      (* generates two lines point and vector *)
line1 = twolns[[1, 1]] + 
  m  twolns[[1, 2]];          (* produces equation of the lines  *)
line2 = twolns[[2, 1]] + n  twolns[[2, 2]];
PQ = line2 - 
  line1;                 (*   vector that goes from line1 to line2    \
*) 
eq1 = Dot[twolns[[1, 2]], 
  PQ];            (* escalar products of lines to vector from line1 \
to line1 *)
eq2 = Dot[twolns[[2, 2]], PQ];
mn = Solve[{eq1 == 0, eq2 == 0}, {m, 
   n}];      (* founds the escalars to apply to lines to make dot \
products equal to 0 so they are perpendicular *)
dist = Norm[
  PQ /. mn];            (*   calculates de distance between line1 and \
line2   *)
P = First[ 
  line1 /. mn];            (*   real vaslues of points P and Q using \
solution of scalars m and n  *)
Q = First[line2 /. mn];
dihedralangle = 
 DMSList[ArcCos[Dot[twolns[[1, 2]], twolns[[2, 2]]]/(
    Norm[twolns[[1, 2]]] Norm[
      twolns[[2, 
       2]]])]*180/\[Pi]];        (*  calculates angle between planes \
that contain lines and PQ  *)
{P - Q, a = (P - Q)[[1]], b = (P - Q)[[2]], 
 c = (P - Q)[[3]]};    (*  calculates the components of vector PQ  *)
\

{plano11 = {0, 0, Dot[(P - Q), P]/(P - Q)[[3]]}, 
 plano12 = {0, Dot[(P - Q), P]/(P - Q)[[2]], 0}, 
 plano13 =  {Dot[(P - Q), P]/(P - Q)[[1]], 0, 0}, 
 plano21 = {0, 0, Dot[(P - Q), Q]/(P - Q)[[3]]}, 
 plano22 = {0, Dot[(P - Q), Q]/(P - Q)[[2]], 0}, 
 plano23 =  {Dot[(P - Q), Q]/(P - Q)[[1]], 0, 
   0}};    (*  defining three points for each plane normal to PQ that \
contains lines *)

{xp = P[[1]], xq = Q[[1]], yp = P[[2]], yq = Q[[2]], zp = P[[3]], 
 zq = Q[[3]]};  (* calculates coordinates od P and Q *)
ms = Max[{sx = Max[xp, xq] - Min[xp, xq], 
    sy = Max[yp, yq] - Min[yp, yq], 
    sz = Max[zp, zq] - 
      Min[zp, zq]}]*1.25;    (* help calculate range to plot *)
pv = {{cx = (xp + xq)/2, cy = (yp + yq)/2, cz = (zp + zq)/2} - 
   ms, {cx = (xp + xq)/2, cy = (yp + yq)/2, cz = (zp + zq)/2} + 
   ms};        (* help calculate range to plot *)
Line1 = InfiniteLine[twolns[[1, 1]], 
  twolns[[1, 
    2]]];           (* set the ploting information of lines  *)
Line2 = InfiniteLine[twolns[[2, 1]], twolns[[2, 2]]];
Line3 = InfiniteLine[{P, Q}];
plane1 = InfinitePlane[{twolns[[1, 1]], P, 
   Q}];              (*  set informatión to plot planes   *)
plane2 =  InfinitePlane[{twolns[[2, 1]], P, Q}];
planeP = InfinitePlane[{plano11, plano12, plano13}];
planeQ =  InfinitePlane[{plano21, plano22, plano23}];
str = StringTemplate["Line1 = `` + m ``\nLine2 = ``+ n `` "][
   twolns[[1, 1]], twolns[[1, 2]], twolns[[2, 1]], twolns[[2, 2]]] <>


  "\ndist " <> TextString[dist] <> "\ndihedral angle " <> 
  TextString[dihedralangle] <> "\n P" <> TextString[ P] <> "\n Q" <> 
  TextString[Q];   (* Label setting  *)
Graphics3D[{Thick, Red, Line1, Blue, Line2, Green, Line3, 
  PointSize[0.03], Black, Point[P], Point[Q], Opacity[0.1, Red], 
  EdgeForm[], plane1, Opacity[0.08, Blue], plane2, 
  Opacity[0.08, Yellow], planeQ, planeP}, 
 PlotRange -> {{Floor[pv[[1, 1]]], 
    Ceiling[pv[[2, 1]]]}, {Floor[pv[[1, 2]]], 
    Ceiling[pv[[2, 2]]]}, {Floor[pv[[1, 3]]], Ceiling[pv[[2, 3]]]}}, 
 Axes -> True, AxesLabel -> {x, y, z}, 
 PlotLabel -> Style[Framed[str], Bold], 
 FaceGrids -> {{-1, 0, 0}, {0, 1, 0}, {0, 0, -1}}, 
 FaceGridsStyle -> Directive[LightGray]]
(* try to look from different angles so that the normal (green) line \
disapears or the blu or red then you would look that green line is \
normal to both red or blue... also can apreciate the planes which \
contain the blue or red lines have green line as normal   or how the \
yellow planes are parallel ... enjoy it ... always can return to \
default plot rightclicking over the plot *)
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.