this is one of those subjects of Mathematics that physicists learn in Physics courses (because mathematicians do not teach them). Typically, identify the symmetries of the system and choose coordinates "adapted"  to these symmetries  so that eg. the  dependency  of the functions, the components of the vector fields and the limits of the integrals, when transformed to such coordinates, take the simplest form.

No systematic theory is normally used to convey these techniques, but much of it is basically use of invariance under discrete and Lie groups.

1. The task: 

I did not follow your 'physical' thinking, just use IR^2 ~ C and take polar coordinates. Here it goes:

#        IR^2              ->  IR
  f:=   (x,y)              ->  exp(-x^2-y^2);

Compute it using Maple's routines:

  Int(Int(f(x,y),x=-infinity..infinity),y=-infinity...infinity);
  value(%);
                infinity    infinity
               /           /
              |           |                2    2
              |           |          exp(-x  - y ) dx dy
              |           |
             /           /
               -infinity   -infinity


                                  Pi

Now do a change of variables using the formula above:

#        IR * ]-Pi,Pi[ = W ->  V = IR^2
  Phi:= (r,theta)          ->  [r*cos(theta),r*sin(theta)];

  'f( op(Phi(r,theta)) )*detJ(Phi)(r,theta)':    # the change of variables
  'Int(Int(%, r=0..infinity), theta=-Pi .. Pi)'; # according integral
  %;
  simplify(%); value(%);

               f(op(Phi(r, theta))) detJ(Phi)(r, theta)


     Pi    infinity
    /     /
   |     |
   |     |          f(op(Phi(r, theta))) detJ(Phi)(r, theta) dr d
   |     |
  /     /
    -Pi   0

        theta


     Pi    infinity
    /     /
   |     |                2           2    2           2             2                 2
   |     |          exp(-r  cos(theta)  - r  sin(theta) ) (cos(theta)  r + r sin(theta) ) dr dtheta
   |     |
  /     /
    -Pi   0



                  infinity                  Pi
                 /                         /
                |                2        |
                |          exp(-r ) r dr  |    1 dtheta
                |                         |
               /                         /
                 0                         -Pi


                                  Pi


Note that one needs op(phi(...)) to switch from the 'vector notation' to sequences
which is what f expects as input.

--------------------------------------------------------------------

2. I am aware that the practical problem is to find appropriate transforms & coordinates.
But that fits into Alejandro Jakubi's answers above :-)

It's worth thinking about some special cases.  In principle, the general change of variables can be broken down into individual steps, where each step either changes the "inside" variable only or interchanges the order of two integrations.  Thus going from int_R f(x,y,z) dx dy dz to an integral ds dt du might be done in the following steps:

dx dy dz -> du dy dz -> dy du dz -> dt du dz -> dt dz du -> dz dt du -> ds dt du

A few comments on the various posts:

• Jakubi's intuition is absolutely right, to find a 'useful' change of coordinates, symmetries are one of the most useful methods.  Unfortunately, that wasn't at all what I was asking for.  I really meant: given a specific transformation of coordinates, how does one rewrite the bounds?
• What Axel did was right, but he hid exactly the hard part: how to go from -infinity..infinity X -infintiy..infintiy to 0..infinity X -Pi..Pi ?  We all know how to do this, but can this be generalized?
• I don't think the problem is quite as hard as generalized visualization
• Robert definitely is thinking along the same lines I was.

Ok, let's be more specific.  Let's consider something simple: a 2D integral over a domain 0..1 X 0..1 and a linear change of coordinates given by a 2x2 matrix A.  We know exactly how to change the integrand f(x,y) into f(u,v) *|det(A)| via A.  But how does the boundary change?  What conditions on A must be imposed for this to make sense?  This should generalize to n-dimensional integrals and a linear change of coordinates over an arbitrary rectangle R fairly easily.  This should be the first case to be worked out in full detail.

The next case to work out is problably the ``rectifying transformation'', that is a 2D integral like

int(int(f(x,y), x=a(y)..b(y)),y=c..d)

where we try to simultaneously 'straighten' a and b [one definitely needs to impose conditions on a and b for this to be possible!].

Another case to work out in detail is how some particular domains transform under rectangular <-> polar coordinate transforms.

After a day of work I am brainwashed and can look at that differently.

I did not understand Alejandro's suggestion, because in alg Geometry
I would divide out the invariance ... however what he says in his
first contribution should be understood by me in the way his task
was set up and solved, as it is what is shown at the Wiki as well:

For certain invariance one simply knows by experience, that one may
use a different coordinate system (polar, cylinder, sphere) hoping
to achieve a (simplier) problem. Astonishing the Wiki examples do
even turn out to be products (as an experiment one could start with
a more complicated situation from there to get ugly ones in IR^n).

Those transformations are known well enough to determine how the
domains transform.

hmm: would one resort 1-dim Fourier problems or periodics to this
as well?

For the invariance may be one can use Lie stuff from Maple (for DE)
to find that out (never worked with it).

May be also some stuff like Arnold's "Geomterische Methoden" for
DE helps for more ideas (that's what Alejandro may have had in
mind with his response).



The other thing are the domains and the boundarys, seems to be one
of Jacques' central points.

Here I see, that the Wiki suggests to find a parametrization of the
boundary curve (that's what is done at 'normal domains', not sure
whether one could also try to use Stoke's Theorem).

But beyond that I think one would try to use transforms for which
one *knows* which geometry is left invarint (or how it transforms),
like Moebius transforms or whatever.

This also applies to the question 'what happens under linear trans-
forms?' - if non-degenerated (which always is needed) it just maps
quadrilaterals to quadrilaterals, since it respects affine lines
(yes Alejandro already answered that more concretly). And circles
transform to ellipsoids etc.

The second question I think is too general - just think of the
homotopy invariance in complex analysis :-( ...  

x = (-a(y)+b(y))*xi+a(y) for Int(Int(f(x,y), x=a(y)..b(y)),y=0..1)

gives me  Int( (-a(y)+b(y)) * J(y), y = 0 .. 1), 
          J(y) = Int( f((-a(y)+b(y))*xi+a(y)), xi = 0 .. 1)

So it becomes an integral over the unit square, but I think it is
an approach inverse to parametrizing the boundaries.

Hi,

You are all very smart people.I am leanring multi variable integration at the momment. This is right on time and very helpful .

Thank you all.