In general, how should we find the pmf given only the moment generating function without comparing its form to that of famous pmf? Background 
It is known that moment generating function generates moments, but does it hold information about the probability of the random variable being realised at a particular value?
Example 
Focusing on only discrete random variable, 
we have 
$$M_Y(t) = \frac{2^t+e^{3t}+4k}{3}$$ for $t \in \mathbb{R}$ and $k$ is a constant. 
Question 
For the moment generating function in the example, how is it possible we know $\text{Pr}(Y=1)$ without comparing to the known form of some famous probability mass function?
 A: As you suggest in your question, the moment generating function holds information on the moments of a distribution. Except for notable examples (e.g. Bernoulli random variable) where the first moment also coincides with the probability of success of the trial, to the best of my knowledge don't hold any direct information on the probability mass.
What you are probably interested in is the probability generating function of a distribution. This on the other hand contains information on the probability mass associated to each value of the random variable's spectrum. 
For a discrete non-negative random variable $X$ with $x \in \mathbb{N}$ and probability mass function $p$, then the probability generating function of $X$ is defined as 
$$G_X(z) = E(z^X) = \sum_{x=0}^\infty p(x)z^x$$
The logic of it is very similar to that of a moment generating function: a clothesline where instead of hanging moments we hang probabilities of points on the spectrum. Taking the appropriate derivative and evaluating the function at 0 we obtain
$$P(X=k) = \frac{d^{(k)}}{dt^k}\frac{G_X(0)}{k!}$$
---Edit---
I'm sorry, I did not read the question thoroughly enough the first time, I'll clarify better off now.
One way you could go by doing this is using the connection between the mgf and the pgf of a discrete random variable. Using your definition of moment generating function and my definition of probability generating function, we can say that 
$$M_Y(t) = E(e^{tY}) = G_Y(e^t)$$
and so a way to obtain probabilities directly from the moment generating function could be the following
$$\frac{1}{k!}\frac{d^{(k)}}{dt^k}M_Y(log(t)) \quad \text{evaluated at} \quad t=0$$
In your example, we would have
$$M_Y(log(t)) = \frac{2^{log(t)} + e^{3log(t) + 4k}}{3}$$
If you wanted to find $P(X=1)$, you would take the first derivative and let $t=0$. 
A: First, for a mgf $\DeclareMathOperator{\E}{\mathbb{E}} M_Y(t)=\E e^ {tY}$ we have that $M_Y(0)=1$. Using that you can find $4k=1$. Using as in the other answer the corresponding probability generating function (assuming that the given function is really the mgf of some random variable with support on the nonnegative integers), we could find the point probabilities by a series expansion. I did so in maple, but that gives a series with many negative coefficients, indicating that the assumption about the support is incorrect.  An example where this approach works.
So we should stay with the mgf. 
If the mgf exists, in the sense that it exist for all arguments in some open interval containing zero, then it uniquely characterizes distributions, so that if two random variables have the same mgf (existing in above sense), then they have the same distribution. So, in principle $\mathrm{Pr}(Y=1)$ can be determined. But in practice that can be difficult. First, if $Y$ is continuous, then $\mathrm{Pr}(Y=1) =0$. But I do not know if/how we can determine from the mgf alone if the distribution is continuous ... In most practical (even theoretical) problems, apart from your quiz-type problem, we know beforehand at least the support of the random variable. Without that knowledge even numerical methods will have problems ... 
A: Note that with a discrete distribution the contribution of an atom of probability at some particular $x$-value $x_0$ has a characteristic effect on the mgf.
For example, the probability $p_0$ at $X=0$ will have a term in the mgf of $p_0$.
A probability of $p_1$ at $X=1$ will have a term in the mgf of $p_1e^t$. 
A probability of $p_2$ at $X=2$ will have a term in the mgf of $p_2e^{2t}$. 
A probabiity of $p_i$ at $X=x_i$ will have a term in the mgf of $p_ie^{x_it}$.
This all follows directly from the definition of the mgf.
The overall mgf will be a sum of such terms. 
Consequently, if you see in the mgf a sum of terms of the form $ae^{bt}$ that means that you have probability $a$ at $X=b$. If you can, simply convert all the terms to be in this form.
Once you have it in the form of a sum of terms like that, you can just write down the probability at any specific value. 
A: Firstly, it is useful to note that the requirement that $M_Y(0) = 1$ means that we must have $k = \tfrac{1}{4}$.  Hence, this constant is fixed, and we can write the characteristic function of $Y$ as:
$$\varphi_Y(t) = \frac{1 + e^{(\ln 2) i t} + e^{3 i t}}{3} = \frac{e^{0 i t} + e^{(\ln 2) i t} + e^{3 i t}}{3}.$$
We can obtain the probability mass function for $Y$ by appropriate inversion of the characteristic function.  As a preliminary matter we note that:
$$\frac{1}{2T} \int \limits_{-T}^{+T} e^{i k t} dt = \text{sinc} (kT) \quad \quad \quad \quad \lim_{T \rightarrow \infty} \frac{1}{2T} \int \limits_{-T}^{+T} e^{i k t} dt = \mathbb{I}(k=0).$$
Hence, using a well-known Fourier inversion formula we obtain:
$$\begin{equation} \begin{aligned}
p_Y(y) 
&= \lim_{T \rightarrow \infty} \frac{1}{2T} \int \limits_{-T}^{+T} e^{-ity} \varphi_Y(t) dt \\[6pt]
&= \lim_{T \rightarrow \infty} \frac{1}{2T} \int \limits_{-T}^{+T} \frac{e^{-iyt} + e^{(\ln 2 - y) i t} + e^{(3-y) i t}}{3}  dt \\[6pt]
&= \frac{\mathbb{I}(y = 0)+ \mathbb{I}(y = \ln 2) + \mathbb{I}(y = 3)}{3}. \\[6pt]
\end{aligned} \end{equation}$$
Hence, we see that $Y$ is uniformly distributed over the three points $Y=0,\ln 2, 3$.  (In fact, as Glen_b points out in his answer, with a bit of familiarity with Fourier transformation, this form can easily be seen from the original characteristic function, without having to do the maths.)
