I have a problem which I don't know how to attack. Actually I am not even sure there is a way to do it.
Is it possible to solve an equation of this form $$A²-\frac{tr(A)²}{4}Id_{4\times 4}=B$$ where A, B are (4*4) matrices, A is the unknown. Is it possible to get A as function of B.
I shall assume that $B$ is diagonalizable, because this gets messy if it's not (plus, the approach would be too numerically unstable).
Notice that $X$ and $X - \xi {\rm I}$ are similar, with $\lambda_i(X) = \lambda_i(X - \xi {\rm I}) + \xi$. Denoting a Jordan decomposition of $B$ as $B = S J S^{-1}$ and $\xi := (\operatorname{tr} A)^2 / 4$, we have
$$J = S^{-1} B S = S^{-1} (A^2 - \xi {\rm I}) S = S^{-1} A^2 S - \xi {\rm I},$$
i.e.,
$$S^{-1} A^2 S = J + \xi {\rm I},$$
which means that $A^2 = S (J + \xi {\rm I}) S^{-1}$, i.e., $A$ has the same similarity matrix in its Jordan normal form (diagonalization) as $B$. (This part wouldn't be true if $A$ could have been nilpotent)
Since we have assumed that $B$ is diagonalizable, $J$ is a diagonal matrix, and so is $J + \xi {\rm I}$. Also, if $A = S K S^{-1}$, then $A^2 = S K^2 S^{-1}$, $K$ is diagonal, and $\operatorname{tr} A = \operatorname{tr} K$. This means that we're solving a diagonal system
$$4J = 4K^2 - (\operatorname{tr} K)^2,$$
or, if we write it elementwise, with $J = \operatorname{diag}(j_1,\dots,j_n)$, $K = \operatorname{diag}(k_1,\dots,k_n)$,
$$4j_i = 4k_i^2 - \left( \sum_{i=1}^n k_i \right)^2.$$
That's a system of $n$ quadratic equations with $n$ unknowns (in your case, $n = 4$), solvable via numerical methods.
Once you solve that for $K$, you have $A = S K S^{-1}$.
Note: Jordan decomposition is not a numerically stable one, so be careful here. Given that your matrices are small, this might work, but do take some care here. If $B$ is real symmetric or complex Hermitian, you can use eigendecomposition, which is stable and should produce fine results.