Regulated function-higher dimensions 
Let $\mathbb{X}$ be a Banach space and let the function $f:[0,1]^2 \rightarrow \mathbb{X}$ be continuous.
(a) Show that $\forall\ \epsilon>0\ \exists\ \delta>0$ s.t. $\|f(x)-f(y)\|<\epsilon$ for every $x,y \in [0,1]^2 $ with $\|x-y\|_{\infty}<\delta$.

Suppose $f$ is continuous but not uniformly continuous. Then $\exists\ \epsilon>0$ s.t. $\forall\ \delta>0$, $\exists\ x,y \in [0,1]^2$ s.t. $\|x-y\|_{\infty}<\delta$ but $\|f(x)-f(y)\| \geq \epsilon$. Take $(f_n)=\frac{1}{n} \to 0$ as $n \to \infty$. Then take two sequences and use Bolzano theorem? Is this right? I will reach some sort of contradiction.

(b) Show that $\forall\ \epsilon>0\ \exists\ g:[0,1]^2 \rightarrow \mathbb{X}$
s.t. $\|f-g\|_{\infty} \leq \epsilon$ and $g$ is a 2-dimensional step function. (For some reals, the function $g$ is constant on each open rectangle)

I have no idea about this bit.

(c) Use part (b) to show that the order of integration does not matter.

 A: For a), I would suggest a constructive approach. You method may be better though.
Since $[0,1]^2$ is compact, then from each cover of it by open subsets, one can find a finite subcover. 
Let $\epsilon >0$, since $f$ is continuous, $\forall x,\exists \delta_x, \parallel x- y \parallel \le \delta_x \Rightarrow \parallel f(x)-f(y) \parallel \le \epsilon/2$.
Consider : $$[0,1]^2 \subset \cup_{x \in [0,1]^2} B(x,\delta_x/2) $$ where $B(x,\delta_x)$ is the ball centered on $x$ of radius $\delta_x$, which is the $\delta_x$ that occurs in the definition of continuity, so that $\parallel f(x)-f(y) \parallel \le \epsilon /2$. From the union we can extract a finite subcover : $$[0,1]^2 \subset \cup_{k=1}^{k=n} B(x_k,\delta_{x_k}/2)$$
So now, if we take $\delta=\min \delta_{x_k}/2$, then I let you see that the definition of uniform continuity is verified.
The pro of that method is that it can be easily generalized for any compact of $\mathbb{R}^n$.
For b, we can construct the step function $g$ with the help of indicator functions. For this one I think there is way easier than what I will write now, let me hear your critics.
Let $$A_{m,i,j}=\{(x,y) \in \left[ \frac i m, \frac {i+1} m \right[\times\left[ \frac j m, \frac {j+1} m \right[\}$$
And : $$g(x,y)=\sum_{i=0}^{m-1} \sum_{j=0}^{m-1} f\left(\frac i m, \frac j m\right) 1_{A_{m,i,j}}$$
Hemce : $$\|f-g\|_{\infty} = \sup_{(x,y) \in [0,1]^2} \parallel f(x,y)-g(x,y)\parallel=\sup_{i,j\in\{0,...,m\}} \sup_{(x,y) \in A_{m,i,j}} \parallel f(x,y)-f\left(\frac i m, \frac j m\right)\parallel $$
Then, for $m$ big enough, for every $(i,j) \in \{0,...,n\}$ and for every $(x,y) \in A_{m,i,j}$ $$\left[ \parallel (x,y) - \left(\frac i m, \frac j m\right) \parallel \le \delta \right] \Rightarrow \left[ \sup_{(x,y) \in A_{m,i,j}} \parallel f(x,y)-f\left(\frac i m, \frac j m\right)\parallel \le \epsilon \right]$$
Then for that $m$ : $$\|f-g\|_{\infty} \le \epsilon$$
For c),
Since $g$ is made by a linear combination indicators functions on cartesian products of intervals (all distincts), it's easy to show that the integration order doesn't matter.
For one product of intervals :
$$\iint_{I\times J}1dxdy=L(I)L(J)$$
Where $L(.)$ gives the length of the interval. Then $$\int_0^1\int_0^1g(x,y)dxdy=\sum_{i=0}^{m-1}\sum_{j=0}^{m-1}\frac {f(\frac i m,\frac j m)} {m^2}=\int_0^1\int_0^1g(x,y)dydx$$
Furthermore, for m big enough, we had : $\|f-g\|_{\infty} \le \epsilon$, let :
$$\|\int_0^1\int^1_0 f(x,y) dxdy-\int^1_0\int_0^1 f(x,y) dydx \| = I$$
Hence :
$$\begin{align*}
 I&=\|\int_0^1\int^1_0 f dxdy-\int_0^1\int_0^1gdxdy+\int_0^1\int_0^1gdydx-\int_0^1\int^1_0 f dydx\| \\ 
 &\le \|\int_0^1\int^1_0 f dxdy-\int_0^1\int_0^1gdxdy\|+\|\int_0^1\int_0^1gdydx-\int_0^1\int^1_0 f dydx\| \\ 
 &\le \int_0^1\int^1_0 \|f-g\|dxdy+\int_0^1\int^1_0 \|f-g\|dydx\\ 
 &\le \int_0^1\int^1_0 \|f-g\|_\infty dxdy+\int_0^1\int^1_0 \|f-g\|_\infty dydx\\
 &\le 2 \epsilon
\end{align*}
$$ 
So we can make the difference between the two order of integration as small as we want, so it's won.
