I am trying to find a function $f(x)$ that minimizes the following cost function
$$E = \left(\int_{-p}^p{{x^2\mathcal{N}(x)}f(x)}dx-\epsilon\right)^2$$
With $\epsilon\geq 0$ and $\mathcal{N}(x)$ a normal distribution with mean $0$ and standard deviation $\sigma$
$$\mathcal{N}(x) = \frac{e^{-\frac{x^2}{2\sigma^2}}}{\sqrt{2\pi\sigma^2}}$$
and subject to the following constraints
$$\int_{-p}^p{f(x)\mathcal{N}(x)}dx = \beta$$
$$0\leq f(x) \leq 1$$
So far, I tried to solve this problem in a similar approach of a calculus of variations problem, by slghtly varying $f(x)\rightarrow f(x)+\delta\eta(x)$, with $\delta<<1$ and $\eta(x)$ an arbitrary continuous function. However, I am not sure that this is the right approach, nor am I sure that the final equation I am getting if correct.
I proposed the following functional $S = E[f] + \lambda(\int_{-p}^p{f(x)\mathcal{N}(x)}dx -\beta)$, with $\lambda$ a Lagrange multiplier, and tried to extreme it (e.g. first order variations in $f$ should verify $S[f+\delta\eta] - S[f] = \Delta S[f] + o(\delta^2)\rightarrow\Delta S[f]=0$). Then
$$S[f+\delta\eta] - S[f] = E[ f+\delta\eta] + \lambda\left(\int_{-p}^p{[f(x)+\delta\eta(x)]\mathcal{N}(x)}dx - \beta\right) - S[f]\\ =2\delta\int_{-p}^p{{x^2\mathcal{N}(x)}f(x)}\int_{-p}^p{{x^2\mathcal{N}(x)}\eta(x)}dx -2\delta\epsilon\int_{-p}^p{x^2\mathcal{N}(x)\eta(x)}dx+\delta\lambda\int_{-p}^p{\mathcal{N}(x)\eta(x)}dx + o(\delta^2)$$
Which leads to
$$\Delta S=2\int_{-p}^p{{x^2\mathcal{N}(x)}\eta(x)}dx\left[\int_{-p}^p{{x^2\mathcal{N}(x)}f(x)dx} -\epsilon\right]+\lambda\int_{-p}^p{\mathcal{N}(x)\eta(x)}dx = 0$$
But here is where I get stuck. This final equation seems far more complex than Euler-Lagrange equations. I really have no idea how to continue. Does this make any sense or is there another way to proceed from scratch?