Solving quadratic equations involving Quaternion variables

I am solving a system of quadratic equations that comes from the 3D rotational kinematics problem expressed by Quaternion variables after simplifying the expressions:

Solve[{-2/dT (-qbr1 qbrp0 + qbr0 qbrp1 - qbr3 qbrp2 + qbr2 qbrp3) + 
     2 constC (-qb1 qbr0 + qb0 qbr1 - qb3 qbr2 + qb2 qbr3) (qb0 qbr0 +
         qb1 qbr1 + qb2 qbr2 + qb3 qbr3) == 
    const1 && -2/
      dT (-qbr2 qbrp0 + qbr3 qbrp1 + qbr0 qbrp2 - qbr1 qbrp3) + 
     2 constC (-qb2 qbr0 + qb3 qbr1 + qb0 qbr2 - qb1 qbr3) (qb0 qbr0 +
         qb1 qbr1 + qb2 qbr2 + qb3 qbr3) == 
    const2 && -2/
      dT (-qbr3 qbrp0 - qbr2 qbrp1 + qbr1 qbrp2 + qbr0 qbrp3) - 
     2 constC (qb3 qbr0 + qb2 qbr1 - qb1 qbr2 - qb0 qbr3) (qb0 qbr0 + 
        qb1 qbr1 + qb2 qbr2 + qb3 qbr3) == const3 && 
   qbr0*qbr0 + qbr1*qbr1 + qbr2*qbr2 + qbr3*qbr3 == 1}, {qbr0, qbr1, 
  qbr2, qbr3}, Cubics -> False, Quartics -> False]

It keeps running but, so far for 5 hours, cannot get the solution yet. This seems possible to solve and not too difficult. Did I make something wrong, or should I solve using other methods or functions? Please advise me.

We can calculate very quickly using a numerical algorithm using NSolve:

dT = 1;
constC = 2;
qbrp0 = 3;
qbrp1 = 4;
qbrp2 = 5;
qbrp3 = 6;
qb0 = 7;
qb1 = 8;
qb2 = 9;
qb3 = 10;
const1 = 11; 
const2 = 12;
const3 = 13;

EQ = {-2/dT *(-qbr1*qbrp0 + qbr0* qbrp1 - qbr3*qbrp2 + qbr2*qbrp3) + 
    2*constC*(-qb1*qbr0 + qb0 *qbr1 - qb3*qbr2 + 
       qb2*qbr3)*(qb0 *qbr0 + qb1* qbr1 + qb2*qbr2 + qb3 *qbr3) == 
   const1, -2/
      dT *(-qbr2 *qbrp0 + qbr3* qbrp1 + qbr0* qbrp2 - qbr1*qbrp3) + 
    2*constC *(-qb2* qbr0 + qb3 *qbr1 + qb0*qbr2 - 
       qb1*qbr3)* (qb0* qbr0 + qb1* qbr1 + qb2* qbr2 + qb3*qbr3) == 
   const2, -2/
      dT* (-qbr3*qbrp0 - qbr2 *qbrp1 + qbr1* qbrp2 + qbr0* qbrp3) - 
    2* constC* (qb3*qbr0 + qb2*qbr1 - qb1*qbr2 - 
       qb0*qbr3)* (qb0* qbr0 + qb1* qbr1 + qb2 *qbr2 + qb3* qbr3) == 
   const3, qbr0*qbr0 + qbr1*qbr1 + qbr2*qbr2 + qbr3*qbr3 == 
   1};
 NSolve[EQ, {qbr0, qbr1, qbr2, qbr3}]
 (*{{qbr0 -> 0.300808 - 0.0888247 I, qbr1 -> 6.31248 - 12.6133 I, 
   qbr2 -> -13.5724 - 0.470986 I, 
   qbr3 -> 6.92771 + 10.5743 I}, {qbr0 -> 0.300808 + 0.0888247 I, 
   qbr1 -> 6.31248 + 12.6133 I, qbr2 -> -13.5724 + 0.470986 I, 
   qbr3 -> 6.92771 - 10.5743 I}, {qbr0 -> -0.393081, qbr1 -> -0.471046,
    qbr2 -> -0.528483, qbr3 -> -0.586778}, {qbr0 -> 0.39098, 
   qbr1 -> 0.469899, qbr2 -> 0.528865, 
   qbr3 -> 0.588754}, {qbr0 -> 0.897994, qbr1 -> -0.21354, 
   qbr2 -> -0.285171, qbr3 -> -0.258236}, {qbr0 -> -0.909058, 
   qbr1 -> 0.243877, qbr2 -> 0.169177, qbr3 -> 0.292432}}*)

 NSolve[EQ, {qbr0, qbr1, qbr2, qbr3}, Reals]
 (* {{qbr0 -> -0.393081, qbr1 -> -0.471046, qbr2 -> -0.528483, 
   qbr3 -> -0.586778}, {qbr0 -> 0.39098, qbr1 -> 0.469899, 
   qbr2 -> 0.528865, qbr3 -> 0.588754}, {qbr0 -> 0.897994, 
   qbr1 -> -0.21354, qbr2 -> -0.285171, 
   qbr3 -> -0.258236}, {qbr0 -> -0.909058, qbr1 -> 0.243877, 
   qbr2 -> 0.169177, qbr3 -> 0.292432}}*)

Regards M.I.

Thanks, M.I. But that’s numerical solution is not what I really look for. I’d like the closed form solution, or the algorithm, so that I can implement in code and solve in real-time.

The numeric solving takes just over 0.1 seconds on my slow laptop. This might be faster than plugging numerical values into a huge complicated set of parametrized Root objects. And given the way the parameters appear in the equations, my guess is the symbolic result will be huge.

Finally, the Solve quit silently with no return or abort. I still believe that there is a way to solve this nice-structured and physically meaningful algebraic system. Maybe, we need some trick to solve this kind of math Olympiad-like problem.