It sounds like what you're asking is: how do you construct a representation of SU(2) in terms of 3x3 matrices on a real 3-dimensional vector space? (This representation is also known as the "spin-1" representation, as it's used to describe the spin of spin-1 particles.)
The H you mention, which appears to be some kind of Hamiltonian, is irrelevant to the above question. I assume it is part of a longer homework question which isn't described fully here, so I'll ignore it.
As your teacher mentions, a simple way to construct $S_x$, $S_y$, and $S_z$ is to start with the raising and lowering operators $S_+$ and $S_-$.
If you work in the $S_z$ basis, then you know what the action of $S_z$ is on each of the 3 $S_z$ eigenstates:
$S_z \begin{bmatrix}1 \\ 0 \\ 0 \end{bmatrix} = \hbar \begin{bmatrix}1 \\ 0 \\ 0\end{bmatrix}$
$S_z \begin{bmatrix}0 \\ 1 \\ 0 \end{bmatrix} = 0$
$S_z \begin{bmatrix}0 \\ 0 \\ 1 \end{bmatrix} = -\hbar \begin{bmatrix}0 \\ 0 \\ 1\end{bmatrix}$
So the 3x3 matrix form of $S_z$ in this basis must be:
$S_z = \hbar\begin{bmatrix}1 & 0 & 0 \\ 0 & 0 & 0 \\ 0 & 0 & -1 \end{bmatrix}$
And you also know the actions of the raising and lowering operators on these $S_z$ eigenstates (up to an undetermined constant):
$S_+$ $\begin{bmatrix}0 \\ 0 \\ 1 \end{bmatrix} = c\begin{bmatrix}0 \\ 1 \\ 0 \end{bmatrix}$
$S_+$ $\begin{bmatrix}0 \\ 1 \\ 0 \end{bmatrix} = c\begin{bmatrix}1 \\ 0 \\ 0 \end{bmatrix}$
$S_+$ $\begin{bmatrix}1 \\ 0 \\ 0 \end{bmatrix} = 0$
$S_-$ $\begin{bmatrix}0 \\ 0 \\ 1 \end{bmatrix} = 0$
$S_-$ $\begin{bmatrix}0 \\ 1 \\ 0 \end{bmatrix} = c\begin{bmatrix}0 \\ 0 \\ 1 \end{bmatrix}$
$S_-$ $\begin{bmatrix}1 \\ 0 \\ 0 \end{bmatrix} = c\begin{bmatrix}0 \\ 1 \\ 0 \end{bmatrix}$
If you take those actions and write them in matrix form, you get:
$S_+$ = c$\begin{bmatrix}0 & 1 & 0 \\ 0 & 0 & 1 \\ 0 & 0 & 0\end{bmatrix}$
$S_-$ = c$\begin{bmatrix}0 & 0 & 0 \\ 1 & 0 & 0 \\ 0 & 1 & 0\end{bmatrix}$
Then, you can write down $S_x$ and $S_y$ just by taking the right linear combinations of $S_+$ and $S_-$:
$S_x = \frac{1}{2}(S_+ + S_-)$
$S_y = \frac{1}{2i}(S_+ - S_-)$
The only final step required is to determine the constant c. This can be determined by finding the eigenvalues of the $S_x$ and $S_y$ matrices. You want them to be $-\hbar$, $0$, and $\hbar$. You can accomplish this by setting $c = \hbar\sqrt{2}$.
As for the commutation relations $[J_i,J_j] = i\hbar\epsilon_{ijk} J_k$, it just means that:
$[S_x,S_y] = i\hbar S_z$
$[S_y,S_z] = i\hbar S_x$
$[S_z,S_x] = i\hbar S_y$
You can verify these directly using matrix multiplication, for example by showing that $S_x S_y - S_y S_x = i\hbar S_z$
The commutator of an observable like $S_x$ or $P_y$ (the spin in the x direction and the y component of the momentum respectively) with the hamiltonian will tell you about the time evolution of that observable. It tells you how that observable changes with time.
This is given through the heisenberg equation:
$\frac{dA}{dt} = \frac{i}{\hbar}[H, A]$
where A is an operator (like $S_x$ or $P_y$) and I have assumed that A does not depend explicitly on time, i.e $A = A(x(t),p(t))$ and not $A = A(x(t), p(t), t)$.
Now from the equation above, we can see that if A commutes with the hamiltonian, $[H,A] = 0$, then $\frac{dA}{dt} = 0$ and thus A is constant in time, it is conserved.
So when you showed that $S_z$ commutes with H while $S_x$ and $S_y$ do not, you showed that the z component of a particles spin is constant in time, while the x and y components are not constant, they precess.
Best Answer
Here I basically do what joshphysics has already mentioned, just in a little more detail, and in a bit more intuitive basis (which makes effectively no difference). So definitely not a slick way.
I use the basis $|m_1\rangle\otimes|m_2\rangle$, where $m_i=\pm1/2$. Keeping only the sign we denote the basis as $\{|++\rangle,|+-\rangle,|-+\rangle,|--\rangle\}$.
Now for $S_z=S_{1z}+S_{2z}$ you can check that $$ S_z|++\rangle=\left(\frac{\hbar}{2}+\frac{\hbar}{2}\right)|++\rangle $$ and similarly $$ S_z|--\rangle=\left(-\frac{\hbar}{2}-\frac{\hbar}{2}\right)|--\rangle $$ so that you and up with the matrix $$ S_z\rightarrow\hbar\left[\begin{array}{cccc}1&0&0&0\\0&0&0&0\\0&0&0&0\\0&0&0&-1\end{array}\right] $$
and
$$ S_z^2\rightarrow\hbar^2\left[\begin{array}{cccc}1&0&0&0\\0&0&0&0\\0&0&0&0\\0&0&0&1\end{array}\right] $$
For $S_x$ and $S_y$ you use raising and lowering operators as josh suggested and obtain the matrix elements by acting on each state in the basis and seeing what has a non-zero dot product with the result, e.g.
$$ \begin{align} S_y|+-\rangle=\frac{1}{2i}\left(S_+-S_-\right)|+-\rangle&=\\\frac{1}{2i}\left(S_{1+}+S_{2+}-S_{1-}-S_{2-}\right)|+-\rangle&=\frac{1}{2i}\left(\hbar|++\rangle-\hbar|--\rangle\right) \end{align} $$
thus $\langle++|S_y|+-\rangle=\frac{\hbar}{2i}$ and $\langle--|S_y|+-\rangle=-\frac{\hbar}{2i}$. Once you do this you can get
$$ S_x\rightarrow\hbar\left[\begin{array}{cccc}0&1&1&0\\1&0&0&1\\1&0&0&1\\0&1&1&0\end{array}\right] $$
and
$$ S_y\rightarrow\hbar\left[\begin{array}{cccc}0&1&1&0\\-1&0&0&1\\-1&0&0&1\\0&-1&-1&0\end{array}\right] $$
Finally the Hamiltonian looks like $$ H=\left[\begin{array}{cccc}\hbar^2B&0&0&\hbar^2A\\0&0&0&0\\0&0&0&0\\\hbar^2A&0&0&\hbar^2B\end{array}\right] $$
Therefore $|+-\rangle$ and $|-+\rangle$ are eigenvectors with eigenvalue $0$, and you only need diagonalize $2\times2$ matrix to find that the eigenvalue $\hbar^2\left(B+A\right)$ corresponds to the eigenvector $\left(|++\rangle+|--\rangle\right)/\sqrt{2}$, and $\hbar^2\left(B-A\right)$ corresponds to the eigenvector $\left(|++\rangle-|--\rangle\right)/\sqrt{2}$.