In your question you impose the condition of pairwise independence, but it is not possible to find the distribution of a linear combination of more than two random variables with this specification. I will therefore assume that you are willing to impose the stronger condition of mutual independence, so that the distribution of the linear combination is fully determined.
On this basis, what you are looking for is the family of stable distributions, which can be obtained by considering the properties of characteristic functions (CFs) (i.e., Fourier transformations of the density functions). If $X_1, ..., X_n$ are mutually independent random variables with characteristic functions denoted by $\varphi_k(t) = \mathbb{E}(\exp(it X_k))$, then any linear combination of these random variables has a CF that is a simple product of the CFs of the underlying random variables:
$$Y = \sum_{k=1}^n \alpha_k X_k \quad \quad \implies \quad \quad \varphi_Y(t) = \prod_{k=1}^n \varphi_k(\alpha_k t).$$
What you are looking for are classes of distributions where the form of the distribution is closed under linear combinations. This property occurs for classes of distributions where the form of the characteristic function is closed under multiplication. This occurs for the normal distribution as a special case, but it also occurs for a broader family of distributions called the Lévy alpha-stable distribution (see below).
Normal distribution: This has characteristic function $\varphi_k(t) = \exp(i \mu_k t - \tfrac{1}{2} \sigma_k^2 t^2)$. We can see that this functional form is closed under multiplication:
$$\begin{equation} \begin{aligned}
\varphi_Y(t) = \prod_{k=1}^n \varphi_k(\alpha_k t)
&= \prod_{k=1}^n \exp(i \alpha_k \mu_k t - \tfrac{1}{2} \alpha_k^2 \sigma_k^2 t^2) \\[8pt]
&= \exp \Big( \sum_{k=1}^n (i \alpha_k \mu_k t - \tfrac{1}{2} \alpha_k^2 \sigma_k^2 t^2) \Big) \\[8pt]
&= \exp \Big( i \mu_Y t - \tfrac{1}{2} \sigma_Y^2 t^2 \Big), \\[8pt]
\end{aligned} \end{equation}$$
where $\mu_Y = \sum_{k=1}^n \alpha_k \mu_k$ and $\sigma_Y^2 = \sum_{k=1}^n \alpha_k \sigma_k^2$. From this characteristic function we see that the linear combination of random variables is also normal, with updated mean and variance that are corresponding linear combinations of the underlying means and variances.
Lévy alpha-stable distribution: This is a generalisation of the normal distribution, allowing for a generalised characteristic function that is still closed under multiplication. It has characteristic function defined by:
$$\varphi_k(t) = \exp(i \mu_k t - \tfrac{1}{\gamma} |\sigma_k t|^{\gamma} (1-ib_k \cdot \text{sgn}(t) \Phi_{\gamma}(t))),$$
where $0 < \gamma \leqslant 2$ is a fixed parameter (common to each of the underlying distributions) and:
$$\Phi_\gamma(t) = \begin{cases}
\tan (\tfrac{\pi \gamma}{2}) & & \text{for } \gamma \neq 1, \\[4pt]
- \tfrac{2}{\pi} \ln |t| & & \text{for } \gamma = 1.
\end{cases}$$
(We also have the parameter restriction $|b_k| \leqslant 1$.)
This characteristic function is also closed under multiplication (when the parameter $\gamma$ is common to all the underlying random variables), and so it defines a broader class of distributions than the normal distribution. The normal distribution occurs in the special case where $\gamma = 2$ and $b_k = 0$. The corresponding density function for this distribution does not exist in closed form, but it is possible to derive many of its properties (which I will leave beyond the scope of the present answer).
