Let $A$ be an $n \times n$ complex matrix and $0\leq q \leq 1$. The boundary of the $q$-numerical range of $A$ is the orthogonal projection of a hypersurface defined by the dual surface of the homogeneous polynomial\[F(t, x, y, z)= {\rm det}(t\, I_n +x(A +A^*)/2+ y(A -A^*)/(2i) +z \, A^* A).\]We construct different types of cubic surfaces $S_F$ corresponding to the homogeneous polynomial $F(t, x, y, z)$ induced by some $3\times 3$ matrices. The degree of the boundary of the Davis-Wielandt shell of a $3 \times 3$ upper triangular matrix is determined by the cubic surface~$S_F$.