Update to the above problem at the request of gAr to see my solution.

I was able to eliminate the numerical methods and come up with a quasi analytical method.

Using ac = bd

Using Lagrangian multipliers:

Adding the constraint ( 2 ) to the above equations:

Solving the simultaneous set of non linear equations:

Using b as a free parameter we can generate my original solution and all the rotations. Of course a is generated using ac = bd.