Sorry, there cannot be a simple solution. For example,
taking $n=3$ and $B_{ij} = (i+j-1)/6$, we compute numerically
(by iterating the contraction mapping as you suggest)
$$
(x_1,x_2,x_3) = (1.26922421\ldots, 1.54095434\ldots, 1.77148256\ldots)
$$
and then (using algdep in gp) that these satisfy irreducible
sextic equations, respectively
$$
3888 x^6 - 2592 x^5 - 3888 x^4 + 936 x^3 + 210 x^2 + 90 x + 7 = 0,
$$
$$
3888 x^6 - 7776 x^5 + 1080 x^4 + 2592 x^3 - 129 x^2 + 162 x - 10 = 0,
$$
$$
3888 x^6 - 12960 x^5 + 12960 x^4 - 4176 x^3 + 678 x^2 - 342 x + 1 = 0
$$
with maximal Galois group $S_6$, and thus not solvable in radicals because
the group $S_6$ is not solvable.
Added later: some further notes . . .
1) The calculation reported above is not quite a complete proof:
to finish it, we can choose one of the variables, say $x_1$,
which we guessed is a zero of the polynomial $P$,
then write the others as polynomials in that variable
(whose coefficients we can surmise in gp using lindep),
and verify that the equations $x_j^2 = \sum_j B_{ij} x_i$ all hold
modulo $P$. Here we find
$$
x_2 = \frac1{9013}
(81648 X^5 - 19764 X^4 - 101628 X^3 - 4182 X^2 + 7570 X + 1170),
$$
$$
x_3 = \frac1{9013}
(-54432 X^5 + 13176 X^4 + 67752 X^3 + 20814 X^2 - 8051 X - 780).
$$
2) The formulas for our example become somewhat simpler when written as
polynomials satisfied by $6x_i$ (to counteract the denominators of $6$
in $B_{ij} = (i+j-1)/6$); for instance $6x_1$ is a root of
$x^6 - 4 x^5 - 36 x^4 + 52 x^3 + 70 x^2 + 180 x + 84 = 0$.
3) Usually one expects that an $n$-variable system of quadratic
equations should yield solutions satisfying equations of degree
$2^n-1$ (subtracting $1$ for the trivial solution at the origin),
so $7$ for $n=3$; we got $6$ because our $B_{ij}$ form a
degeenrate matrix, making the origin a solution of multiplicity $2$.
4) The contraction mapping argument is not actually correct:
extracting square roots yields an expansion mapping for small $x_i$.
The existence of a positive solution can nevertheless be proved
using the fixed-point theorem on the $(n-1)$-simplex
$\Delta = \{ x_i \geq 0, \sum_i x_i = 1 \}$:
the composite of the linear transformation $B$,
the coordinatewise square root, and the scaling
$(x_i) \mapsto (x_i) / \sum_i x_i$ gives a continuous map $\Delta \to \Delta$;
any fixed point has all coordinates $x_i$ positive, and can be scaled to
a solution $c(x_i)$ of the original system. Uniqueness might still hold
(it does for $n=1$ (trivially) and $n=2$, and I've not been able to find
a counterexample for $n=3$ either), but would have to be proved some other way.
