Here is a complete solution to the main question when $n$ and $q$ are both odd, and a partial solution for the other parities. The partial solution includes a reduction to the case $n=2$.
Let $\text{Tr}_k$ denote the trace map from $\mathbb{F}_{q^n}$ to $\mathbb{F}_{q^k}$, assuming that $k|n$.
The equation is
$$z+z^q + z^{q^2} = x^{q^2+q+1} - xy^q - x^{q^2}y.$$
First, make the change of variables $y \leftarrow -yx^{q^+1},$ so that the equation becomes
$$z+z^q+z^{q^2} = x^{q^2+q+1}(y^q+y+1).$$
Second, when $q$ is odd, we can clarify matters a little with the change of variables $y \leftarrow y - \frac12$ to get rid of the constant. The image of the $q$-linear map $z \mapsto z + z^q + z^{q^2}$, acting on $\mathbb{F}_{q^{3n}}$, consists of those elements $Z'$ such that $\text{Tr}_3(z') = \text{Tr}_1(z')$. Whenever such a $z'$ is reached by the right side, there are $q^2$ solutions for $Z$. On the right side, $y \mapsto y' = y^q+y$ is an $q$-linear isomorphism when $n$ is odd, while when $n$ is even its image is the locus of $\text{Tr}_2(y') = \text{Tr}_1(y')$. Meanwhile $y' \mapsto x^{q^2+q+1}y'$ is a $q$-linear isomorphism unless $x=0$.
Thus when $n$ and $q$ are both odd, there are $(q^{3n}-q^2)(q^{3n}-1)$ solutions with $x,y \ne 0$. There are $q^2(2q^{3n}-1)$ more solutions when one of them is zero.
When $q$ is odd and $n$ is even, then in principle different $x' = x^{q^2 + q + 1}$ could behave differently in the equation $z' = x'y'$. For a fixed $x'$, the set of possible $x'y'$ is a certain $q$-linear hyperplane with codimension $1$, while the set of possible $z'$ is a certain $q$-linear hyperplane with codimension $2$. In the special case that $\{x'y'\}$ contains $\{z'\}$, then there are $q^{3n+1}$ solutions for that value of $x'$. For generic non-zero values of $x'$, there are $q^{3n}$ solutions. Dualize the hyperplanes with respect to $\text{Tr}_1$. The dual of $\{y'\}$ is the line $L$ of trace 0 elements in $\mathbb{F}_{q^2}$, while the dual of $\{z'\}$ is the plane $P$ of trace 0 elements in $\mathbb{F}_{q^3}$. Meanwhile $\{x'\}$ is the subgroup $G$ of $\mathbb{F}_{q^{3n}}$ of index $q^2+q+1$. A special value of $x'$ in this subgroup is one such that $x'L \subset P$. A priori I am not sure that it never happens. What I can say is that if $x'$ is special, then it must lie in $\mathbb{F}_{q^6}$ because both $L$ and $P$ do. So you can reduce the counting problem to the case that $3n = 6$ or $n = 2$.
If $q$ is even, then the equation is
$$z' = x'(y'+1),$$
where as before $z' = z+z^q+z^{q^2}$ and $y' = y + y^q$. In this case $y'$ is any element with zero trace, and the dual line $L$ is just $\mathbb{F}_q$ itself. I have not worked out exactly how it looks, but I suppose that it reduces to the case $n=1$ for similar reasons as above.
Afterthought: I don't feel like changing all of the equations, but I'm wondering now whether there a dual change of variables to put the right side in the form $y''(x^{q^2}+x)$. I think that the map $x \mapsto x^{q^2}+x$ is always non-singular when $q$ is odd.
A remark about where the trace conditions come from. If $a$ is an irreducible element of $\mathbb{F}_{q^n}$, then the map $x \mapsto x^q$ is a cyclic permutation matrix in the basis of conjugates of $a$. A map such as $z \mapsto z+z^q+z^{q^2}$ is then a sum of disjoint permutation matrices and it easy to compute its image and cokernel.
Some remarks about Jared's second, more general question: C.f. the answer to this other mathoverflow question about counting points on varieties. For fixed $q$, the equation of a hypersurface is equivalent to a general Boolean expression, and there may not be much that you can do other than count one by one. There are several strategies that work in the presence of special structure: You can use zeta function information, if you have it, to extrapolate to large values of $q$. You can count the points on a variety if you happen to know that it's linear, or maybe quadratic, or the coset space of a group. And you can use standard combinatorial counting tricks, which in algebraic geometry form amount to looking at fibrations, blowups, inclusion-exclusion for constructible sets, etc.
This particular variety decomposes a lot because it can be made jointly linear in $Y$ and $Z$, and $X$ only enters in a multiplicative form.
