Let's create some example data:

d = {d1, d2, d3};
A = {{B1, B2, B3}};
d = d - PseudoInverse[A].A.d // ComplexExpand // Simplify;
B = {{0, B3, -B2}, {-B3, 0, B1}, {B2, -B1, 0}};
d.s // Simplify

(* {0} *)

Now, we augment the system matrix B and the right hand side d:

BB = ArrayFlatten[{{B, - Transpose[A]}, {A, 0}}];
dd = Join[d, {0}];
x = LinearSolve[BB, dd][[1 ;; 3]]
(* {(-B3 d2 + B2 d3)/(B1^2 + B2^2 + B3^2), (B3 d1 - B1 d3)/(B1^2 + B2^2 + B3^2), (-B2 d1 + B1 d2)/(B1^2 + B2^2 + B3^2)} *)

Let's test the result:

B.x - d // Simplify
A.x // Simplify

(* {0,0,0} *)
(* {0} *)

With the augmentation, we enforced that $B_1 \, \omega_1 + B_2 \, \omega_2 + B_3 \, \omega_3 = 0$ and we added the orthogonal complement of the image of B to the image of BB.

Let's create some example data:

d = {d1, d2, d3};
A = {{B1, B2, B3}};
d = d - PseudoInverse[A].A.d // ComplexExpand // Simplify;
B = {{0, B3, -B2}, {-B3, 0, B1}, {B2, -B1, 0}};
A.d // Simplify

(* {0} *)

Now, we augment the system matrix B and the right hand side d:

BB = ArrayFlatten[{{B, - Transpose[A]}, {A, 0}}];
dd = Join[d, {0}];
x = LinearSolve[BB, dd][[1 ;; 3]]
(* {(-B3 d2 + B2 d3)/(B1^2 + B2^2 + B3^2), (B3 d1 - B1 d3)/(B1^2 + B2^2 + B3^2), (-B2 d1 + B1 d2)/(B1^2 + B2^2 + B3^2)} *)

Let's test the result:

B.x - d // Simplify
A.x // Simplify

(* {0,0,0} *)
(* {0} *)

With the augmentation, we enforced that $B_1 \, \omega_1 + B_2 \, \omega_2 + B_3 \, \omega_3 = 0$ and we added the orthogonal complement of the image of B to the image of BB.

Let's create some example data:

{B1, B2, B3} = RandomReal[{-1, 1}, 3];
A = {{B1, B2, B3}};
\[CapitalDelta] = RandomReal[{-1, 1}, 3];
\[CapitalDelta] = \[CapitalDelta] - PseudoInverse[A].A.\[CapitalDelta];
B = {{0, B3, -B2}, {-B3, 0, B1}, {B2, -B1, 0}};
A.\[CapitalDelta]

(* {0.} *)

Now, we augment the system matrix B and the right hand side \[CapitalDelta]:

BB = ArrayFlatten[{{B, -A\[Transpose]}, {A, 0}}];
\[CapitalDelta]\[CapitalDelta] = Join[\[CapitalDelta], {0}];
\[Omega] = LinearSolve[BB, \[CapitalDelta]\[CapitalDelta]][[1 ;; 3]];
B.\[Omega] - \[CapitalDelta]
A.\[Omega]

(* {1.11022*10^-16, 0., -1.11022*10^-16} *)
(* {-2.77556*10^-17} *)

With the augmentation, we enforced that $B_1 \, \omega_1 + B_2 \, \omega_2 + B_3 \, \omega_3 = 0$.

Let's create some example data:

d = {d1, d2, d3};
A = {{B1, B2, B3}};
d = d - PseudoInverse[A].A.d // ComplexExpand // Simplify;
B = {{0, B3, -B2}, {-B3, 0, B1}, {B2, -B1, 0}};
d.s // Simplify

(* {0} *)

Now, we augment the system matrix B and the right hand side d:

BB = ArrayFlatten[{{B, - Transpose[A]}, {A, 0}}];
dd = Join[d, {0}];
x = LinearSolve[BB, dd][[1 ;; 3]]
(* {(-B3 d2 + B2 d3)/(B1^2 + B2^2 + B3^2), (B3 d1 - B1 d3)/(B1^2 + B2^2 + B3^2), (-B2 d1 + B1 d2)/(B1^2 + B2^2 + B3^2)} *)

Let's test the result:

B.x - d // Simplify
A.x // Simplify

(* {0,0,0} *)
(* {0} *)

With the augmentation, we enforced that $B_1 \, \omega_1 + B_2 \, \omega_2 + B_3 \, \omega_3 = 0$ and we added the orthogonal complement of the image of B to the image of BB.