# Find the solution of a matrix equation

Let $a_k$ and $b_k$ be strictly positive real numbers and $c_{k,q}$ a known symmetric real matrix. All indices range from $1$ to $M$. I have to find $x_{k,q}$ such that

$$c_{k,q} = \sum_{p=1}^{M}\frac{a_k\,a_q\,a^2_p}{b_k\,b_q}\,x_{k,p}\,x_{q,p}+\left(\frac{a_k\,a_{q}^2}{b_k}+\frac{a_q\,a_{k}^2}{b_q}\right)\,x_{k,q}+a_k^2\,\delta_{k,q}$$

My guess is, in the hypotesis that $b_k$ are large and hence the "non-linear" term $\sum_{p=1}^{M}\frac{a_k\,a_q\,a^2_p}{b_k\,b_q}\,x_{k,p}\,x_{q,p}$ is neglectable, to use an iterative alogrithm. I start with a guess where the non-linear term is put to zero: $$x^{\left(0\right)}_{k,q} = \frac{1}{\left(\frac{a_k\,a_{q}^2}{b_k}+\frac{a_q\,a_{k}^2}{b_q}\right)}\,\left(c_{k,q}-a^2_k\,\delta_{k,q}\right)$$

and then I iterate

$$c_{k,q} = \sum_{p=1}^{M}\frac{a_k\,a_q\,a^2_p}{b_k\,b_q}\,x^{(n)}_{k,p}\,x^{(n)}_{q,p}+\left(\frac{a_k\,a_{q}^2}{b_k}+\frac{a_q\,a_{k}^2}{b_q}\right)\,x^{(n+1)}_{k,q}+a_k^2\,\delta_{k,q},$$

solving for $x^{n+1}_{k,q}$, that is

$$x^{(n+1)}_{k,q} = \frac{1}{\left(\frac{a_k\,a_{q}^2}{b_k}+\frac{a_q\,a_{k}^2}{b_q}\right)}\,\left(c_{k,q}-\sum_{p=1}^{M}\frac{a_k\,a_q\,a^2_p}{b_k\,b_q}\,x^{(n)}_{k,p}\,x^{(n)}_{q,p}-a^2_k\,\delta_{k,q}\right).$$

How do I control for convergence? Is there some "standard" approach in facing this kind of problems?

-