$\newcommand{\De}{\Delta} \newcommand{\ep}{\epsilon} \newcommand{\ga}{\gamma} \newcommand{\Ga}{\Gamma} \newcommand{\la}{\lambda} \newcommand{\Si}{\Sigma} \renewcommand{\th}{\theta} \newcommand{\R}{\mathbb{R}} \newcommand{\F}{\mathcal{F}} \newcommand{\E}{\operatorname{\mathsf E}} \newcommand{\PP}{\operatorname{\mathsf P}} \newcommand{\ii}[1]{\operatorname{\mathsf I}\{#1\}}$ Knowing only the 2nd and 3rd moments of the (presumable centered) summand random variables (r.v.'s) is hardly enough for good approximation. Also, it is somewhat unclear what you mean by "3 moments". In the Berry--Esseen bound, the absolute 3rd moments are used, whereas the Edgeworth expansion (which is in a sense more accurate) is given in terms of initial moments. So, given only the limited information you have, here is a cheap and fast way to approximate and simulate the sum $S_n$ of independent centered r.v.'s $X_1,\dots,X_n$ with given $\E X_i^2$ and $\E X_i^3$: just match the 2nd and 3rd (initial) moments of $S_n$ with the corresponding moments of a r.v. of the form $cZ$, where $c$ is a real number, $Z=Y-\E Y=Y-\la$, and $Y$ has the Poisson distribution with some parameter $\la>0$. Note that $\E S_n^2=B_2:=\sum_1^n \E X_i^2$, $\E S_n^3=B_3:=\sum_1^n \E X_i^3$, and $\E Z^2=\la=\E Z^3$. So, the "matching" system of equations is \begin{equation} c^2\la=B_2,\quad c^3\la=B_3. \end{equation} Solving this system, we get \begin{equation} c=\frac{B_3}{B_2},\quad \la=\frac{B_2^3}{B_3^2}. \end{equation} With these values of $c$ and $\la$, we can approximately simulate $S_n$ as $c(Y-\la)$ with $Y\sim\text{Poisson}(\la)$. [1]: https://en.wikipedia.org/wiki/Edgeworth_series