For convenient comparison and clear physical meaning, the gravity on the surface of a homogeneous cube and on spheres inside, outside, and intersecting about it is calculated by polyhedral or harmonic expansion methods. In addition, the gravity coefficients of a rectangle’s spherical harmonics are both derived analytically and evaluated numerically, where only five terms are nontrivial up to the order of 4, which is somewhat unexpected when we first obtained them. There are some similarities of these coefficients to an ellipsoid for the terms C20,C22,C42, but they are much different for the terms C40,C44. Thence, a few special gravity characteristics are here revealed or visualized. For example, it is shown as expected that the maximum gravity appears at the sphere intersecting the cube, but maximum surface gravity at the center of the mid-plane of a rectangle’s surface is different from the gravity on an ellipsoid at the end of its short axis. Based on these results, an orbit around a cube is integrated by a polyhedral method, and its secular motion analysis by averaging theory is investigated where the numerical and analytic results fit very well. Finally, a few special trajectories on a surface plane of a cube are simulated; the physical meaning is quite clear, and some insights are shown, such as why a natural celestial body in the shape of a rectangle with sharp corners is rarely found due to its surface gravity distribution. All gravity calculations are visualized on 3D figures both for cubes or rectangles. Additionally, examples of an asteroid and an ellipsoid are shown so that the techniques discussed here can be adopted to directly analyzing the gravity of other shapes.