We have
[tex]\mathrm{div}(\mathbf f)=\nabla\cdot\mathbf f=\dfrac{\partial(2x)}{\partial x}+\dfrac{\partial(xy)}{\partial y}+\dfrac{\partial(2xz)}{\partial z}=2+3x[/tex]
By the divergence theorem, the flux of [tex]\mathbf f[/tex] over the region [tex]\mathcal S[/tex] that bounds the region [tex]\mathcal E[/tex] is
[tex]\displaystyle\iint_{\mathcal S}\mathbf f\cdot\mathrm d\mathbf S=\iiint_{\mathcal E}\nabla\cdot\mathbf f\,\mathrm dV[/tex]
[tex]=\displaystyle\int_{x=0}^{x=3}\int_{y=0}^{y=3}\int_{z=0}^{z=3}(3+2x)\,\mathrm dz\,\mathrm dy\,\mathrm dx[/tex]
[tex]=\displaystyle3^2\int_{x=0}^{x=3}(3+2x)\,\mathrm dx=162[/tex]