If one is familiar with the concepts of sufficiency and completeness, then this problem is not too difficult. Note that $f(x; \theta)$ is the density of a $\Gamma(2, \theta)$ random variable. The gamma distribution falls within the class of the exponential family of distributions, which provides rich statements regarding the construction of uniformly minimum variance unbiased estimators via notions of sufficiency and completeness.
The distribution of a random sample of size $n$ from this distribution is
$$
g(x_1,\ldots,x_n; \theta) = \theta^{2n} \exp\Big(-\theta \sum_{i=1}^n x_i + \sum_{i=1}^n \log x_i\Big)
$$
which, again, conforms to the exponential family class.
From this we can conclude that $S_n = \sum_{i=1}^n X_i$ is a complete, sufficient statistic for $\theta$. Operationally, this means that if we can find some function $h(S_n)$ that is unbiased for $\theta$, then we know immediately via the Lehmann-Scheffe theorem that $h(S_n)$ is the unique uniformly minimum variance unbiased (UMVU) estimator.
Now, $S_n$ has distribution $\Gamma(2n, \theta)$ by standard properties of the gamma distribution. (This can be easily checked via the moment-generating function.)
Furthermore, straightforward calculus shows that
$$
\mathbb{E} S_n^{-1} = \int_0^\infty s^{-1} \frac{\theta^{2n} s^{2n - 1}e^{-\theta s}}{\Gamma(2n)} \,\mathrm{d}s = \frac{\theta}{2n - 1} \>.
$$
Hence, $h(S_n) = \frac{2n-1}{S_n}$ is unbiased for $\theta$ and must, therefore, be the UMVU estimator.
Addendum: Using the fact that $\newcommand{\e}{\mathbb{E}}\e S_n^{-2} = \frac{\theta^2}{(2n-1)(2n-2)}$, we conclude that the $\mathbb{V}ar(h(S_n)) = \frac{\theta^2}{2(n-1)}$. On the other hand, the information $I(\theta)$ from a sample of size one is readily computed to be $-\e \frac{\partial^2 \log f}{\partial \theta^2} = 2 \theta^{-2}$ and so the Cramer-Rao lower bound for a sample of size $n$ is
$$
\mathrm{CRLB}(\theta) = \frac{1}{n I(\theta)} = \frac{\theta^2}{2n} \> .
$$
Hence, $h(S_n)$ does not achieve the bound, though it comes close, and indeed, achieves it asymptotically.
However, if we reparametrize the density by taking $\beta = \theta^{-1}$ so that
$$
f(x;\beta) = \beta^{-2} x e^{-x/\beta},\quad x > 0,
$$
then the UMVU estimator for $\beta$ can be shown to be $\tilde{h}(S_n) = \frac{S_n}{2 n}$. (Just check that it's unbiased!) The variance of this estimator is $\mathbb{V}ar(\tilde{h}(S_n)) = \frac{\beta^2}{2n}$ and this coincides with the CRLB for $\beta$.
The point of the addendum is that the ability to achieve (or not) the CRLB depends on the particular parametrization used and even when there is a one-to-one correspondence between two unique parametrizations, an unbiased estimator for one may achieve the Cramer-Rao lower bound while the other one does not.
For the first question, the best unbiased estimator is $\chi\left(\sum_i x_i = n\right)$ as you wrote, because the going probability function for the $n$ observations:
$$
\mathbb{P}\left( X_1=x_1, \ldots, X_n=x_n \right)=p^{x_1}(1-p)^{1-x_1} \cdots p^{x_n} (1-p)^{x_n} = p^{\sum_i x_i} (1-p)^{n - \sum_i x_i}
$$
Thus it factors into $(p^n)^{\chi\left(\sum_i x_i = n\right)} \cdot \left( p^{\sum_i x_i} (1-p)^{n - \sum_i x_i} \right)^{1-\chi\left(\sum_i x_i = n\right)}$.
For the second question $\bar{x}=\frac{1}{n} \sum_{i=1}^n x_i$ is the BUE for $\mu$. The factor of the likelihood that depends on this statistics is $\exp(-\frac{n}{2} \left( \mu - \bar{x} \right)^2 )$.
The variance of $\bar{x}$ is $\mathrm{Var}(\bar{x}) = \frac{1}{n^2} \sum_i \mathrm{Var}(x_i) = \frac{1}{n^2} \cdot n = \frac{1}{n}$, hence the Fisher information is $\mathcal{I}(\mu) = \frac{1}{\mathrm{Var}(\bar{x})} = n$.
For the third question, the joint density for the sample:
$$
f = 2^n \chi_{\theta-\frac{1}{4} \le \min(x_1,\ldots, x_n)} \chi_{\theta+\frac{1}{4} \ge \max(x_1,\ldots,x_n)} = 2^n \chi_{ \max(x_1,\ldots,x_n) -\frac{1}{4} \le \theta \le \min(x_1, \ldots,x_n) + \frac{1}{4} }
$$
Thus $\theta$ is determined by two-component vector statistics consisting of the minimal and maximal element of the sample suitably shifted, and $\theta$ can be anywhere in between. The mean of these two values could be a possible choice for the estimator.
Best Answer
$$E[Y_1]=\frac{n(n+1)}{(2n+1)\theta^n}\int_{\theta}^{2\theta}t(t-\theta)^{n-1}dt=\dots=\theta$$
In fact
$$F_{X_{(n)}}(t)=\Bigg(\frac{t-\theta}{\theta}\Bigg)^n$$
...the rest easy follows
----------‐-------
Your error was that you considered the distribution of an uniform $U(0;\theta)$ and not correctly $U(\theta;2\theta)$