One way: define a "wrapper" for special relativistic vectors, and define its behavior under the usual operations with upvalues.
srv[x1_, y1_, z1_, t1_].srv[x2_, y2_, z2_, t2_] ^:=
x1 x2 + y1 y2 + z1 z2 - t1 t2
Check that it has the right properties. Time unit vector:
tu = srv[0, 0, 0, 1];
tu.tu
(* -1 *)
Space unit vector:
xu = srv[1, 0, 0, 0]
(* 1 *)
