Three-dimensional rectangular prisms are building blocks for calculating gravity anomalies from irregular 3D mass bodies with spatially variable density contrasts. A 3D vector gravity potential is defined for a 3D rectangular prism with density contrast varying in depth and horizontally. The vertical component of the gravity anomaly equals the flux of the 3D vector gravity potential through the enclosed surface of the prism. Thus, the 3D integral for the gravity anomaly is reduced to a 2D surface integral. In turn, a 2D vector gravity potential is defined. The vertical component of the gravity anomaly equals the net circulation of the 2D vector gravity potential along the enclosed contour bounding the surfaces of the prism. The 3D integral for the gravity anomaly is reduced to 1D line integrals. Further analytical or numerical solutions can then be obtained from the line integrals, depending on the forms of the density contrast functions. If an analytical solution cannot be obtained, the line-integral method is semianalytical, requiring numerical quadratures to be carried out at the final stages. Singularity and discontinuity exist in the algorithm and the method of exclusive infinitesimal sphere or circle is effective to remove them. Then the vector-potential line-integral method can calculate the gravity anomaly resulting from a rectangular prism with density contrast, varying simply in one direction and sophisticatedly in three directions. The advantage of the method is that the constraint to the form of the density contrast is greatly reduced and the numerical calculation for the gravity anomaly is fast.