A Monte Carlo approach is used to estimate hole mobilities in molecular β-copper phthalocyanine (CuPc) crystal for different applied electric field directions. Due to the crystal symmetry, the twelve neighboring molecules in the three-dimensional crystal are selected in the hopping rate calculation. Density functional theory is employed to derive the molecular interaction between the central and neighboring molecules for various applied electric fields. The derived molecular hopping rate is applied to 80 × 80 × 80 lattice sites under periodic boundary conditions. In order to achieve accurate statistics, each calculation includes 6561 particles with more than 10000 hopping steps under an applied electric field of 0.5-3.5 MV/cm. The results indicate that the molecular hopping strongly depends on the molecular orientation and neighboring sites related to the applied electric field direction. The estimated carrier mobility can be described by the percentage occupation in each neighboring site and the obtained hole mobility value is in the same range of the measured values of single crystal CuPc. The calculated mobility for applied electric field along the c crystal axis exhibits the highest values while the mobility along the b axis has the smallest value.