Ion trap is one of the most promising candidates for quantum computing. Current schemes mainly focus on a linear chain of up to about one hundred ions in a Paul trap. To further scale up the qubit number, one possible direction is to use 2D or 3D ion crystals (Wigner crystals). In these systems, ions are generally subjected to large micromotion due to the strong fast-oscillating electric field, which can significantly influence the performance of entangling gates. In this work, we develop an efficient numerical method to design high-fidelity entangling gates in a general 3D ion crystal. We present numerical algorithms to solve the equilibrium configuration of the ions and their collective normal modes. We then give a mathematical description of the micromotion and use it to generalize the gate scheme for linear ion chains into a general 3D crystal. The involved time integral of highly oscillatory functions is expanded into a fast-converging series for accurate and efficient evaluation and optimization. As a numerical example, we show a high-fidelity entangling gate design between two ions in a 100-ion crystal, with a theoretical fidelity of 99.9%.