There is nothing wrong with your mathematica code. The problem is with the volume element. The volume element in spherical coordinates is $$r^2\sin(\theta ) \;dr\,d\theta \,d\varphi $$ not $$r^2\sin(\varphi ) \;dr\,d\theta \,d\varphi $$ Simply replace $$Sin[phi]$$ with $$Sin[Theta]$$ in your code, and you get the answer.

Integrate[Integrate[Integrate[r *Cos[theta]*r^2*Sin[theta], {r, 0, 4}], {theta, 0, Pi/2}], {phi, 0, 2 Pi}]

There is nothing wrong with your mathematica code. The problem is with the volume element. The volume element in spherical coordinates is $$r^2\sin(\theta ) \;dr\,d\theta \,d\varphi $$ not $$r^2\sin(\varphi ) \;dr\,d\theta \,d\varphi $$ Simply replace $$Sin[phi]$$ with $$Sin[theta]$$ in your code, and you get the answer.

Integrate[Integrate[Integrate[r *Cos[theta]*r^2*Sin[theta], {r, 0, 4}], {theta, 0, Pi/2}], {phi, 0, 2 Pi}]