Prove that this expression is greater than 1/2
Let $$f(x,y):=4x^{2}+4y^{2}-4xy-4y+1 + \frac{4}{\pi^2}\Bigl( \sin^{2}(\pi x)+ \sin^{2}(\pi y) + \sin^{2}(\pi y-\pi x) \Bigr).$$ I will show that $$\min_{0\leq x\leq y\leq 1}f(x,y)=\min_{0\leq x\leq 1/3}f(x,2x).\tag{$\ast$}$$ This suffices, because the minimum of the one-variable function $f(x,2x)$ is easy to analyze numerically: it occurs around $x_0\approx 0.146625$, and has value $f(x_0,2x_0)\approx 0.502187$.
In order to prove $(\ast)$, observe that $f(x,y)$ is invariant under the following two permutations of the closed triangle $T:=\{(x,y):0\leq x\leq y\leq 1\}$: $$(x,y)\mapsto (y-x,1-x)\qquad\text{and}\qquad (x,y)\mapsto (x,1-y+x).$$ These two bijections generate a group $G\cong S_3$ of permutations of $T$. The elements of $G$ are: \begin{align*} &(x,y)\mapsto(x,y),&&(x,y)\mapsto(y-x,1-x),&&(x,y)\mapsto(1-y,1-y+x),\\ &(x,y)\mapsto(x,1-y+x),&&(x,y)\mapsto(1-y,1-x),&&(x,y)\mapsto(y-x,y). \end{align*} The closed triangle $U$ with vertices $$A:=(0,0),\qquad B:=(1/2,1/2), \qquad C:=(1/3,2/3)$$ is a fundamental domain for the action of $G$ on $T$, hence $$\min_T f(x,y)=\min_U f(x,y).$$ For any interior point $(x,y)\in\mathrm{int\,} U$, we have $0<y<2x<1$ and $0<y<2/3$, therefore $$\frac{\pi}{4}\cdot\frac{\partial f}{\partial x}=(2\pi x-\pi y)+2\cos(\pi y)\sin(2\pi x-\pi y)>(2\pi x-\pi y)-\sin(2\pi x-\pi y)>0.$$ So there is no local extremum on the interior of $U$, which implies that $$\min_T f(x,y)=\min_{\partial U} f(x,y).$$ It is straightforward to verify that the minimum over the side $AB$ (resp. $BC$) exceeds $0.6$ (resp. $0.57$). Hence the minimum over $\partial U$ equals the minimum over the side $CA$, and $(\ast)$ follows.
Let $F(x,y)$ denote the left-hand side of your inequality. It is easy to see that $|\nabla F(x,y)|\sqrt2/n<0.002$ if $0<x<y<1$, where $n:=6600$. A direct calculation shows that $F(i/n,j/n)>0.5021\dots$ for all integers $i,j$ such that $0\le i\le j\le n$. It follows that $F(x,y)>0.502-0.002=1/2$ if $0<x<y<1$, as desired.
Details of the calculations can be seen in the image of a Mathematica notebook below. Note that Mathematica can compute any expression in elementary functions with any degree of accuracy.
One can see that the calculation of the minimum of $F(i/n,j/n)$ over all integers $i,j$ such that $0\le i\le j\le n=6600$ took about 422 sec (on 12 cores working in parallel). Without parallelization, the execution time was about 2290 sec.
It is also seen that the infimum of $F(x,y)$ over all $x,y$ such that $0<x<y<1$ is $<0.5022$, which is pretty close to $1/2$.
$\newcommand{\R}{\mathbb{R}}$ An advantage of my previous answer was that, while the computer calculations were pretty heavy there, the logic was extremely simple; virtually no thinking or ingenuity was needed.
On the other hand, one can use a bit of thinking in order to greatly reduce the amount of calculations. More specifically, one can use second-order partial derivatives instead of the first-order ones, based on the following simple
Lemma 1 Suppose that a function $F\colon[x_1,x_2]\times[y_1,y_2]\to\R$ is such that for some real $c$ and $M$ we have $F(x_i,y_j)\ge c$ for $i,j$ in $\{1,2\}$, and the second-order partial derivatives $F''_{xx}$ and $F''_{yy}$ are bounded from above by $M$. Suppose also that $x_2-x_1=y_2-y_1=h>0$. Then $$ F\ge c-Mh^2/4 $$ (on $[x_1,x_2]\times[y_1,y_2]$).
It is easy to see that for $F(x,y)$ denoting the left-hand side of your inequality we have $F''_{xx}\le M$ and $F''_{yy}\le M$, where $M:=32$. Hence, for $n=64$ and $h:=1/n$ we have $Mh^2/4<0.002$. One the other hand, a direct calculation shows that $F(i/n,j/n)>c:=0.5023\dots$ for all integers $i,j$ such that $0\le i\le j\le n$. (This calculation now takes about 0.1 sec -- which may be compared with 2290 sec for the previous calculation, based on first-order partial derivatives, with $n=6600$.)
Thus, by Lemma 1, $F(x,y)>c-Mh^2/4=0.5023-0.002>1/2$ if $0<x<y<1$, as desired.
It remains to prove Lemma 1. Let us first establish its one-dimensional analogue:
Lemma 2 Suppose that a function $f\colon[x_1,x_2]\to\R$ is such that for some real $c$ and $M$ we have $f(x_i)\ge c$ for $i\in\{1,2\}$, and $f''\le M$. Suppose also that $x_2-x_1=h>0$. Then $$ f\ge c-Mh^2/8 $$ (on $[x_1,x_2]$).
Proof of Lemma 2 Let $g(x):=f(x)+M(x-x_1)(x_2-x)/2$. Then $g(x_i)=f(x_i)\ge c$ for $i\in\{1,2\}$ and $g''=f''-M\le0$, so that $g$ is concave and hence $g\ge c$ (on $[x_1,x_2]$). Thus, $f(x)=g(x)-M(x-x_1)(x_2-x)/2\ge c-Mh^2/8$, as claimed. $\Box$
Now we are ready for
Proof of Lemma 1 Take any $(x_*,y_*)\in[x_1,x_2]\times[y_1,y_2]$. For each $j\in\{1,2\}$, applying Lemma 1 to the function $x\mapsto f_j(x):=F(x,y_j)$ in place of $f$, we get $F(x_*,y_j)\ge c_1:=c-Mh^2/8$. Applying now Lemma 1 to the function $y\mapsto g(y):=F(x_*,y)$ in place of $f$, we get $F(x_*,y_*)\ge c_1-Mh^2/8=c-Mh^2/4$, as claimed. $\Box$
Remark In the answer by GH from MO, the minimization of the function $F$ of two arguments was reduced to the minimization of a function of one argument. Using that reduction together with Lemma 2 above, it should be possible to further reduce the execution time from 0.1 sec to something like 0.1/20=0.05 sec.