The challenge is creating the intended correlations while, at the same time, achieving the intended non-normal marginals that you mention. Here are two options with references and links:
1. Headrick (2002)
This approach is very fast, but it cannot create all possible distribution shapes. Its usefulness for you will depend on how extreme your power law distribution is (if $\alpha$ is too low, I doubt this method will work). You can find the details in the references below, and I'll summarize it briefly. The basic idea is to simulate standard normal variables, rotate them using intermediate correlations, and then construct your intended variables each as a polynomial of these z-scores. For example, $Y=c_0+c_1z+c_2z^2+c_3z^3+c_4z^4+c_5z^5$. The most difficult part is finding the constants, which is typically done in Mathematica. A tutorial can be found in Headrick et al. (2007).
Headrick, T. C. (2002). Fast fifth-order polynomial transforms for generating univariate and multivariate nonnormal distributions. Computational Statistics & Data Analysis, 40(4), 685-711.
Headrick, T. C., Sheng, Y., & Hodis, F. A. (2007). Numerical Computing and Graphics for the Power Method Transformation Using Mathematica. Journal of Statistical Software, 19(3).
2. Ruscio & Kaczetow (2008)
This is a slower, iterative algorithm. Their article is behind a paywall, but their R code is freely available here. I'll summarize the bivariate case for the sake of simplicity, though it can work with multivariate problems, too. Uncorrelated $X_o$ and $Y_o$ are generated with any shapes (e.g. bimodal). Then, $X_1$ and $Y_1$ are generated as bivariate normal with an intermediate correlation. $X_1$ and $Y_1$ are replaced by $X_0$ and $Y_0$ in a rank-preserving fashion. Adjust the intermediate correlation higher or lower depending on the observed correlation between $X_1$ and $Y_1$. Rinse and repeat until the observed correlation is sufficiently close to the target correlation.
Ruscio, J., & Kaczetow, W. (2008). Simulating multivariate nonnormal data using an iterative algorithm. Multivariate Behavioral Research, 43, 355–381. doi:10.1080/00273170802285693