The following used to contain a mistake, found by Hans. I have found a fix, and have incorporated into the following proof.
Here is a sketch of an argument that the log likelihood function has at most one critical point. It should be read in conjunction with Hans's partial answer.
The idea is to use the "variation diminishing property of the Laplace transform", according to which the Laplace transform of a function with $k$ sign changes cannot have more than $k$ sign-changing zero-crossings. For a function $\phi:\mathbb R^+\to\mathbb R$, let $S(\phi)$ be the maximal $k$ for which there exist $0<x_0<x_1<x_2<\cdots < x_{k}$ for which $\phi(x_i)\phi(x_{i+1})<0$ for all $0\le i < k$. Then the Laplace transform $g(s)=\int_0^\infty e^{-sx} G(x)dx$ of $G$ obeys $S(g)\le S(G)$. This topic is not well explained in Wikipedia articles, but the result used here is in chap V, paragraph 80 in vol 2 of
Pólya and Szegő's Problems and Theorems in Analysis (p. 225 in my copy), is discussed at length in Karlin's book Total Positivity (see Theorem 3.1, page 21, and pages 233 and 237), in papers by I.J. Schoenberg, etc. One can think of it as a continuous analogue of Descartes' Rule of Signs. I used it in answering this MSE problem.
If the logarithm of the likelihood function had two or more local maxima, its derivative would have three or more roots, since between every two local maxima lies a local minimum. So it suffices, by the variation diminishing property of the LT, to show that what the OP, in his draft answer, calls $\tilde f$ has at most two sign changes.
This seems evident numerically, but deserves proof just as much as the original problem does. Here is one way of seeing this, using another application of the variation diminishing property of the Laplace transform.
Here the argument. First, I will change notation, using $s$ instead of $t$ and setting $y=x^2$. The claim is that, for fixed real $y\ge0$,
$$\tilde f(s) = \frac{1-e^{-ys}}s +(1-y)e^{-ys}-\frac 2{1+e^{-s}}$$
has at most two sign changes as a function of $s\in\mathbb R^+$. Let $g(s)=\dfrac{1+e^{-s}}{s^2}\tilde f(s)$; clearly $g$ has as many sign changes as $\tilde f$ does. But $g$ is itself a Laplace transform:
\begin{align}g(s)&=\frac{1+e^{-s}-e^{-ys}-e^{-(y+1)s}}{s^3}+(1-y)\frac{e^{-ys}+e^{-(y+1)s}}{s^2} - \frac2{s^2}\\
&=\int_0^\infty e^{-sx} G(x) dx,\end{align}
from which one reads off
\begin{align}G(x)&=\frac 1 2\left((x)_+^2 - (x-y)_+^2 +(x-1)_+^2 - (x-(y+1))_+^2\right) \\&+ (1-y)\left((x-y)_++(x-y-1)_+\right)-2x.\end{align}
Here $(x)_+=\max(x,0)$. Since $x\mapsto (x)_+$ is continuous, so is $G$.
If $y<1$ the function $G$ is piece-wise quadratic on each of the intervals $(0,y)$, $(y,1)$, $(1,y+1)$, $(y+1,\infty)$; if $y>1$ then $G$ is piecewise quadratic on the intervals $(0,1)$, $(1,y)$, $(y,y+1)$, $(y+1,\infty)$, so verification of the lemma is in principle easy in a case-by case manner. In practice, tedious and error prone.
If $y<1$ the formula for $G(x)$ reduces to
$$G(x)=\begin{cases} x^2/2 -2x&0\le x<y\\
y^2/2-x-y&y\le x<1\\
x^2/2-2x+(y-1)^2/2&1\le x<1+y\\
y^2-2y-1&1+y\le x\end{cases}$$ and if
$y>1$, the formula reduces to
$$G(x)=\begin{cases} x^2/2 -2x&0\le x<1\\
x^2-3x+1/2&1\le x<y\\
x^2/2-2x+(y-1)^2/2&y\le x<1+y\\
y^2-2y-1&1+y\le x.\end{cases}$$
These can be merged into the following, where the cases are referred to below:
$$
G(x)=\begin{cases}
x^2/2 -2x&\text{A: if }0\le x<\min(1,y)\\
y^2/2-x-y&\text{B: if }y\le x< 1\\
x^2-3x+1/2&\text{C: if }1\le x< y\\
x^2/2-2x+(y-1)^2/2&\text{D: if }\max(1,y)\le x< 1+y\\
y^2-2y-1&\text{E: if }1+y< x
\end{cases}
$$
Note that cases B and C are mutually exclusive.
Computations show that for fixed $y$ the function $G(x)$ has at most one sign change; I sketch an argument for this below. (Omitting an analysis of the possibility of sign changes at the case boundaries.)
$G$ has no sign changes in cases A or E (in A, the only possibilities are $x=0$ or $x=4$, the former is not a sign change, and $x=4$ does not obey $0\le x<\min(1,y)$. Constant functions, as in case E, do not have sign changes.) Case B has no sign changes, for the value $x=y^2/2-y$ violates $y<x<1$. In case C, a sign change could only occur at $x=(3\pm\sqrt 7)/2$, and then $1<x<y$ implies $x=(3+\sqrt7)/2$ and $y>(3+\sqrt7)/2$. In case D, a sign change can only occur at $x=2\pm\sqrt{3+2y-y^2}$, and $\max(0,x)<x<y$ is only possible if $x=2+\sqrt{3+2y-y^2}$ and $1+\sqrt 2<y<(3+\sqrt 7)/2$. Putting these together: if $y<1$ then there can be no sign changes in the relevant cases A,B,D,E. If $y>1$ there might be at most one sign change in each of C, D (out of the relevant A,C,D,E), but not actually both, since that would violate $(3+\sqrt7)/2 <y<(3+\sqrt 7)/2$. Hence, $G$ has at most one sign change among A,B,C,D,E.
Finally, since $G(0)=0, G'(0)=-1, G(\infty)=y^2-2y-1$, we see $G$ has exactly one sign change if $y^2-2y-1>0$ and none if $
y^2-2y-1\le0$.
The meta-motivation is to shoehorn the original question into an application of the variation diminishing machinery given in my first paragraphs. The micro-motivation for my choice of $g$ (and hence of $G$) comes from the realization that $\tilde f$ is the Laplace transform of the signed measure $$\mu = \lambda_{[0,y]} + (1-y)\delta_y -2\sum_{k\ge0}(-1)^k \delta_k,$$ where $\lambda_{[0,y]}$ is Lebesgue measure restricted to $[0,y]$ and $\delta_k$ represents the unit point mass at $k$. The signed measure $\mu$ has infinitely many sign changes, but the telescoping series $\mu*(\delta_0+\delta_1)$ does not, where $*$ denotes convolution of measures, so $1+e^{-s}$ times $\tilde f$ is a better candidate for the variation diminishing trick sketched above.
Dividing by a power of $s$ has the effect of smoothing the signed measure, and eliminating some small oscillations that create their own extraneous sign changes. The mistake Hans found in an earlier version of this answer was to divide by $s$, which allowed for 3 sign changes for a certain range of $y$. Dividing by $s^2$ fixed this problem, at the price of making $G$ piecewise quadratic instead of piecewise linear.
