Attempt to find a solution (out of potentially many) to a system of equations with FindInstance runs forever or crashes

Question

I have a system of equations that in principle should have a solution (I would think not a unique one), which I want to solve. I am not interested in all solutions; I just want to locate one, assuming it exists. I was thinking that FindInstance would not take that long to evaluate but after a whole night of running the code (MacBook Pro 2021, 32GB memory), I still do not get an answer. In fact, it seems I run out of memory or something, since the kernel crashed.

I also tried Reduce, Solve and NSolve. I am wondering whether my code is terrible (or there is a mistake) and what I can do to improve it. Before I show the code, let me describe the underlying mathematical problem I am trying to solve.

Mathematical details on the underlying problem

I have 7 unit vectors with complex coefficients, $v_i=(a_i,b_i,c_i,d_i), \,\,\, i=0,\ldots,6$. We define the inner product as $$ (v_i,v_j)\equiv a^*_i a_j+b^*_i b_j+c^*_i c_j+d^*_i d_j, $$ where $^*$ denotes the complex conjugate.

Then, the vectors are constrained to have the following inner products: $$ |(v_0,v_j)| = \sqrt{\frac{5}{6}} \,, \,\, j=1,\ldots,6, \notag \\ |(v_i,v_j)| = \frac{4}{5} \,, \,\, i,j=0, \notag \\ (v_i,v_i) = 1 $$ with the last one being the normalisation condition, i.e. the fact that they should be unit vectors. I can slightly simplify the problem by using the overall freedom to align the first vector with one of the axes, that is, $v_0 = (1,0,0,0)$. Then from the overlap of the other states with state $v_0$, we find that $|a_j|=\sqrt{\frac{5}{6}}, j=1,\ldots,6$. With this assumption, it follows that the condition $(v_i,v_i) = 1$ is equivalent to $$ |b_j|^2+|c_j|^2+|d_j|^2=1-|a_j|^2=\frac{1}{6} $$ or equivalently, in terms of real and complex parts $$ (b_j^{(R)})^2+(b_j^{(I)})^2+(c_j^{(R)})^2+(c_j^{(I)})^2+(d_j^{(R)})^2+(d_j^{(I)})^2=\frac{1}{6}. $$ Now the last equation defines points on a 5-sphere in $\mathbb{R}^6$. It is known that there can exist $n+1$ equidistant points on a sphere in $\mathbb{R}^n$ and thus a solution should exist to my problem.

In the particular implementation, I arbitrarily choose all the $a_j$ to be real and equal to $\sqrt{\frac{5}{6}}$ (this shouldn't be too strong of an assumption). I further assumed vector $v_1$ to be equal to $v_1=(\sqrt{\frac{5}{6}},\frac{1}{\sqrt{6}},0,0)$ and left the rest to Mathematica.

Here is the code:

x[j_] = {Sqrt[5/6], bR[j] + I bI[j], cR[j] + I cI[j], dR[j] + I dI[j]};
x[0] = {1, 0, 0, 0};
bR[1] = 1/Sqrt[6];
bI[1] = 0;
cR[1] = 0;
cI[1] = 0;
dR[1] = 0;
dI[1] = 0;


cs = Table[Abs[x[i]\[Conjugate] . x[j]]^2, {i, 1, 6}, {j, i, 6}] // 
    Flatten // ComplexExpand;
os = Table[
    If[i == j, 1, If[i == 0, (Sqrt[5/6])^2, (4/5)^2]], {i, 1, 6}, {j, 
     i, 6}] // Flatten;

nsol = FindInstance[{cs == os, 
   Element[
    Table[{bR[j], cR[j], dR[j], bI[j], cI[j], dI[j]}, {j, 2, 6}] // 
     Flatten, 
    Reals]}, {Table[{bR[i], cR[i], dR[i], bI[i], cI[i], dI[i]}, {i, 2,
       6}]} // Flatten]

Given that I have 21 equations for 30 unknowns I would expect that it would be easy to locate a solution but my code runs forever. Could you help me improve my code and get a solution to the problem?

Also, I noticed that when I run my code, my CPU usage is something like ~12.5%, even though WolframKernel app runs at a 100%. Is there a way to increase the CPU usage and thus reduce the runtime by using some sort of parallelisation to use more available kernels?