Is $dx\,dy$ really a multiplication of $dx$ and $dy$?

In a double integral, you are actually integrating a differential two-form:

$$\int_R \mathrm{f}(x,y) \ \mathrm{d}x \wedge \mathrm{d}y$$

Here, $\mathrm{d}x$ and $\mathrm{d}y$ are the basis differential one-forms and $\mathrm{d}x \wedge \mathrm{d}y$ is their exterior product.

Just as one can think of the derivative in Robinson's framework as a true ratio $\frac{\Delta y}{\Delta x}$ modulo an infinitesimal error (eliminated by applying the shadow), so also one can think of a single-variable integral as an infinite sum of infinitesimal terms of type $dx$ (again up to applying shadow). Double integrals can naturally be viewed as double (infinite) sums, where $dx\,dy$ is most decidedly an ordinary product. And of course this generalizes to multiple integrals as the OP suggested. If one is in Euclidean space, talking about differential forms is an unnecessary obfuscation.

Edit 1: For finite Riemann sums approxiating the double integral, it is obvious that the term $\Delta x \Delta y$ is a product; it seems nobody in his right mind would deny this. The difference is that one cannot deduce the value of the integral from a finite Riemann sum. On the other hand, with an infinite Riemann sum when $\Delta x$ is replaced by $dx$, etc., the value of the integral is deduced from the value of the Riemann sum by taking the shadow (see above). That's the advantage of having the richer syntax of the hyperreal approach.

Edit 2: the OP's question is in fact equivalent to a question about single-variable integrals, namely: does $f(x)dx$ denote multiplication of $f(x)$ by $dx$? Perhaps the right answer is that it denotes a memory of multiplication. Namely, the multiplication is still there at the level of the hyperfinite Riemann sum. To pass from this to the integral one applies the standard part function, after which we have only a "memory" left. Similarly, one can form the differential quotient Δy/Δx which is still a ratio, but one doesn't get the derivative until one applies the standard part function. Here also there is only a memory of a division left. The advantage of the hyperreal framework is that one has a direct procedure for passing from the ratio to the derivative which isn't the case in the traditional real-based framework where one must appeal to an indirect notion of an epsilon, delta limit.

A survey of various approaches to Robinson's framework is due to appear in Real Analysis Exchange.

As others have said $dx\, dy$ does not represent a product of differentials. But it represents a product of measures. We have the "natural" Lebesgue measure $\lambda$ on the $x$-axis, and integration with respect to this measure is signalled by writing ${\rm d}x$ as right parenthesis of the integral. Similarly we have the "natural" Lebesgue measure $\lambda$ on the $y$-axis, and integration with respect to this measure is signalled by writing ${\rm d}y$ as right parenthesis of an integral involving the variable $y$. The individual measures $\lambda$ on the $x$-axis ${\mathbb R}$ and the $y$-axis ${\mathbb R}$ define a product measure $\lambda\otimes\lambda$ on the cartesian product ${\mathbb R}^2$, again called Lebesgue measure on ${\mathbb R}^2$. Integration with respect to this product measure is is signalled by writing ${\rm d}(x,y)$, ${\rm d}x\otimes {\rm d}y$, or simply $dx\,dy$, as right parenthesis of an integral over some subset $A\subset{\mathbb R}^2$. Fubini's theorem then tells us that $$\int\nolimits_A f(x,y)\>{\rm d}(x,y)=\int\nolimits_{A'}\left(\int\nolimits_{A_x} f(x,y)\> {\rm d}y\right)\ {\rm d}x\ ,$$ where $A'$ denotes the projection of $A$ onto the $x$-axis and $A_x:=\{y\mid (x,y)\in A\}$ collects the $y$-values to be weighted in for given $x\in A'$.