\mathbb R^2, \\
\\
\int_{\mathbb R^2}e^{\tilde U_i}<\infty, \quad \tilde U_i \mbox{ is radial},\,\, \forall i\in I.
\end{array}
\right.$$

One may expect the limiting entire solution to be different around each blowup point, however the authors proved that $\tilde U$ is independent of blowup points, and only depends on the ratio of $\rho_1^k-\rho_1:\rho_2^k-\rho_2:..\rho_n^k-\rho_n$ (see [@linzhang1]). Naturally it leads to the question: how to construct bubbling solutions with the help of this information?

In this paper and subsequent ones, we are devoted to study the bubbling phenomenon of Liouville systems: how to accurately estimate the bubbling solutions of (\[mainsys\]) and how to construct them. These are quite challenging analytic problems. In general, blowup analysis for a system of equations is much harder than that for the single equation. One reason is that the Pohozaev identity, a balancing condition, is no longer so powerful as in the scalar case. Another reason is that there are too many entire solutions: the parameter $\sigma=(\sigma_1,...,\sigma_n)$ ($\displaystyle{\sigma_i=\frac 1{2\pi}\int_{\mathbb R^2}e^{u_i}}$), which represents the energy of the entire solution, forms a submanifold of $n-1$ dimension. However, for the Liouville equation, the energy is just one number: $\int_{\mathbb R^2}e^u=8\pi$.

In this article we consider the case of one blowup point, 