I don't think these objects can be classified in a manner similar to the normed unital division algebras, if you take "algebra" to mean "vector space $V$ equipped with a bilinear map $V \otimes V \to V$". In particular, I suspect you end up with high-dimensional moduli spaces of such structures in all large dimensions.

Here is a naive calculation of degrees of freedom:

Let $a_{i,j}^k$ be the structure constants of our algebra, assembled into matrices $A^k$, and consider a point $x = (x_1, \ldots, x_n)$. We may write $x \ast x = (x^T A^1 x,\ldots, x^T A^n x)$. The length condition becomes $\left(\sum_i x_i^2\right)^2 = \left(\sum_k x^T A^k x \right)^2$. Writing this out in terms of the coordinates of $x$, we obtain an identity of homogeneous polynomials in $x_i$ of total degree 4, with coefficients that are quadratic in the structure constants. In other words, the space of solutions is an intersection of $\binom{n+3}{4}$ quadric hypersurfaces in $n^3$-dimensional space. Subtracting, we get positive dimensions for $2 \leq n \leq 17$ with maximum 377 at $n=13$, but after that, we get more constraints than variables. As it happens, the positivity holds for $n \leq 16$ when we account for the $O(n)$ symmetry of the solution space.

Since the Cayley-Dickson construction exists and provides solutions in arbitrarily large dimension, it is clear that the constraints are highly non-generic. This does not completely eliminate the possibility that in some dimensions there are no solutions, but I think it is at least discouraging as far as classification is concerned.

I don't think these objects can be classified in a manner similar to the normed unital division algebras, if you take "algebra" to mean "vector space $V$ equipped with a bilinear map $V \otimes V \to V$". In particular, I suspect you end up with high-dimensional moduli spaces of such structures in all large dimensions.

Here is a naive calculation of degrees of freedom:

Let $a_{i,j}^k$ be the structure constants of our algebra, assembled into matrices $A^k$, and consider a point $x = (x_1, \ldots, x_n)$. We may write $x \ast x = (x^T A^1 x,\ldots, x^T A^n x)$. The length condition becomes $\left(\sum_i x_i^2\right)^2 = \sum_k \left(x^T A^k x \right)^2$. Writing this out in terms of the coordinates of $x$, we obtain an identity of homogeneous polynomials in $x_i$ of total degree 4, with coefficients that are quadratic in the structure constants. In other words, the space of solutions is an intersection of $\binom{n+3}{4}$ quadric hypersurfaces in $n^3$-dimensional space.

[Revised following YCor's comment:] When we account for the $O(n)$ symmetry of the solution space, we get the formula $$ n^3 - \binom{n+3}{4} - \binom{n}{2}$$ which is positive for $2 \leq n \leq 16$ with maximum 299 at $n=13$. When $n$ is large we get more constraints than variables.

Since the Cayley-Dickson construction exists and provides solutions in arbitrarily large dimension, it is clear that the constraints are highly non-generic. This does not completely eliminate the possibility that in some dimensions there are no solutions, but I think it is at least discouraging as far as classification is concerned.