Given two $n\times n$ matrices $A$ and $B$, define their Hadamard product $A\circ B$ as the element-wise product, i.e. $$(A\circ B)_{ij} = A_{ij}B_{ij}\,.$$ A well known result is the Schur product theorem, stating that if both $A$ and $B$ are non-negative defined matrices, then $A\circ B$ too is non-negative.
Is it possible to somehow extend this result to integral operators on a Hilbert space?
For instance let us consider two real functions $a,b:[0,1]^2\to \mathbb{R}$. Let's assume that they're both continuous and symmetric (i.e. $a(x,y) = a(y,x)$ and the same for $b$).
Then we can define two compact self-adjoint integral operators $A$ and $B$ on $L^2([0,1])$, by $$A\phi(x) = \int_0^1 a(x,y)\phi(y)dy\,;\qquad B\phi(x) = \int_0^1 b(x,y)\phi(y)dy \,.$$ Let $A\circ B$ be the integral operator given by $$(A\circ B)\phi(x) = \int_0^1 a(x,y)b(x,y)\phi(y)dy\,.$$
Assume that $A$ and $B$ are non-negative, i.e. for all $\phi$ $\langle A\phi,\phi\rangle\geq 0$ and $\langle B\phi,\phi\rangle\geq 0$.
Can we state that $A\circ B$ is non-negative?
At least, is it the case when $A$ and $B$ commute, so that they have a common orthonormal basis?
Yes, we can indeed say that $A \circ B$, as you've defined it, is a non-negative operator. Here's an adapted version of the proof outlined in the introduction to section 7.5 of Horn and Johnson's Matrix Analysis (second edition), which is the section about the Schur product (AKA Hadamard product) theorem for matrices:
Here's an alternative proof:
Define the operator $A \otimes B: L^2([0,1]^2)$ such that for all $f,g \in L^2[0,1]$, the element $f_1 \otimes f_2$ defined by $f_1 \otimes f_2(x_1,x_2) := f_1(x_1)f_2(x_2)$ is mapped as follows: $$ (A\otimes B)(f_1 \otimes f_2)(x_1,x_2) = \int_{0}^1\int_0^1 a(x_1,y_1)b(x_2,y_2)f_1(y_1)f_2(y_2)\,dy_1\,dy_2. $$ In other words, taking $k(x_1,x_2,y_1,y_2) = a(x_1,y_1)b(x_2,y_2)$, $A \otimes B$ is simply the map $$ (A \otimes B) f(x) = \int_{[0,1]^2}k(x,y)f(y)\,dy, \quad x,y \in \Bbb R^2. $$ Because $L^2([0,1]^2)$ is spanned by elements of the form $f_1 \otimes f_2$, we can quickly see that $A \otimes B$ must be non-negative.
With that established: we note that $$ \int_{0}^1 \int_0^1 f(x) a(x,y)b(x,y)f(y)\,dy\,dx = \lim_{n \to \infty} \frac{1}{\mu(D_n)^2} \int_{D_n}\int_{D_n} k(x,y) (f \otimes f)(x)\,dx\,dy, $$ where $D_n \subset [0,1]^2$ is defined by $D_n = \{(x_1,x_2) : |x_1 - x_2| \leq 1/n\}$.