There is increasing interest in improving trapped ions lifetime and the controllability of their equilibrium positions for large scale quantum computation and simulation. Thus a well-designed ion trap is especially important. Here, we present a numerical toolset to find optimal trap parameters. Trap potential field is solved via boundary element method with improved precision. Trapped ion system is chaotic because of the nonlinear Coulomb interaction, thus any computational error will accumulate exponentially with the evolution time. We separate the energy of the system from such errors by simulating molecular dynamics via position extended Forest-Ruth integration method, which preserve time-reversal symmetry. System characteristics such as equilibrium positions and Coulomb crystal phase transitions are calculated after considering background collisions and laser interactions. The relationship between RF heating rate and geometrical parameters of electrodes, the surface roughness and the ambient electromagnetic noise is estimated through simulation. Finally, suitable design parameters are found through the optimization algorithm, and a series of high performance ion traps are designed and fabricated.