I know that, the sets $\{(x,y)\in I:x\in\mathbb{Q}\}$ and $\{(x,y)\in I:x\not \in\mathbb{Q}\}$ are dense in I. We may required to consider to split the domain to $[0,1/2]\times [0,1]$ and $[1/2,1]\times [0,1]$. But I am unable to prove this. please help.
Consider an arbitrary partition $P$ of $[0,1]\times [0,1]$ with subrectangles $R_{jk}=[x_{j-1},x_j] \times [y_{k-1},y_k]$ for $1 \leqslant j \leqslant n$ and $1 \leqslant k \leqslant m$. Without loss of generality there exists an index $q$ such that $y_{q-1} < 1/2 \leqslant y_q$.
Using the density of rationals and irrationals, we have for all $j = 1,\ldots, n$,
$$\inf_{(x,y) \in R_{jk}}f(x,y) = 2y_{k-1}, \sup_{(x,y) \in R_{jk}}f(x,y) = 1 \,\,\,\text{for} \,\,\,\,1 \leqslant k \leqslant q, \\ \inf_{(x,y) \in R_{jk}}f(x,y) = 2y_{q}, \sup_{(x,y) \in R_{jk}}f(x,y) = \max(1,2 y_{q+1}) \,\,\,\text{for} \,\,\,\,k = q+1,\\ \inf_{(x,y) \in R_{jk}}f(x,y) = 1, \sup_{(x,y) \in R_{jk}}f(x,y) = 2 y_k \,\,\,\text{for} \,\,\,\,q+2 \leqslant k \leqslant m$$
The lower Darboux sum is
$$L(P,f) = \sum_{j=1}^n\sum_{k=1}^{q+1}2y_{k-1}(x_j-x_{j-1})(y_k - y_{k-1})+ \sum_{j=1}^n\sum_{k=q+2}^{m}1\cdot(x_j-x_{j-1})(y_k - y_{k-1})\\ = \sum_{k=1}^{q+1}2y_{k-1}(y_k - y_{k-1}) + 1 - y_{q+1}$$
The sum on the RHS is a lower Darboux sum for the integral of $y \mapsto 2y$ over the interval $[0,y_q]$ and, since $y_q < 1/2 \leqslant y_{q+1}$ , it follows that
$$L(P,f) \leqslant \int_0^{y_q} 2y \, dy + 1 - y_{q+1} \leqslant \int_0^{1/2} 2y \, dy +1 - 1/2 = 3/4$$
The upper Darboux sum is
$$U(P,f) = \sum_{j=1}^n\sum_{k=1}^{q}1 \cdot(x_j-x_{j-1})(y_k - y_{k-1})+ \sum_{j=1}^n\max(1,2y_{q+1})(x_j-x_{j-1})(y_{q+1} - y_{q})\\+\sum_{j=1}^n\sum_{k=q+2}^{m}2y_{k}(x_j-x_{j-1})(y_k - y_{k-1}) \\ = \sum_{k=1}^{q}1 \cdot(y_k - y_{k-1})+ \max(1,2y_{q+1})(y_{q+1} - y_{q})+\sum_{k=q+2}^{m}2y_{k}(y_k - y_{k-1}) \\ = y_q + \max(1,2y_{q+1})(y_{q+1} - y_{q})+\sum_{k=q+2}^{m}2y_{k}(y_k - y_{k-1}) $$
The sum on the RHS is an upper Darboux sum for the integral of $y \mapsto 2y$ over the interval $[y_{q+1},1]$ and, since $y_{q+1} \geqslant 1/2$, it follows that
$$U(P,f) \geqslant y_q + \max(1,2y_{q+1})(y_{q+1} - y_{q})+ \int_{y_{q+1}}^1 2y \, dy \geqslant y_{q+1} + \int_{y_{q+1}}^1 2y \, dy \\ \geqslant 1/2 + \int_{1/2}^1 2y \, dy= 5/4$$
Therefore, $f$ is not Riemann integrable since
$$\sup_P L(P,f) \leqslant 3/4 < 5/4 \leqslant \inf_P U(P,f)$$