conformal mappings between the flat torus and the embedded torus (I just might solve this one fairly shortly and post an answer here if no one else does.  And maybe even if someone else does.)
$$
\begin{align}
& R > r > 0 \\[6pt]
x & =  (R + r \cos v) \cos u \\[4pt]
y & =  (R + r \cos v) \sin u \\[4pt]
z & =  r \sin v
\end{align}
$$
The point $(u,v)$ above is in $\mathbb{R}^2/(2\pi\mathbb{Z})^2$, which is a flat torus.  The point $(x,y,z)$ in on a torus embedded in $\mathbb{R}^3$.  The mapping $(u,v)\mapsto (x,y,z)$ is not conformal.  Putting something else in place of $\cos v$ and $\sin v$ should make it conformal (for now I prefer to leave $\cos u$ and $\sin u$ undisturbed).  What functions should those be?
 A: If you replace $v$ by $f(v)$ where $f$ satisfies
$$
f' = \frac{R}{r} + \cos(f)
$$
then your mapping becomes locally conformal. This equation has a closed form solution, which I happily leave to you (or WA) to figure out.  ;-) This closed form also suggest that it is a globally conformal map on your "flat" torus only if
$$
r = \frac{R}{\sqrt{1 + k^2}}
$$
for some positive integer $k$. (This because $2 \pi$ must be a period of $\cos(f)$ .)
A: This will not be a complete answer, but it's a different point of view from other stuff posted here.  Recall the conventional paramatrization of the torus:
\begin{align}
& R > r > 0 \\[6pt]
x & =  (R + r \cos v) \cos u \\[4pt]
y & =  (R + r \cos v) \sin u \\[4pt]
z & =  r \sin v
\end{align}
As a mapping from the flat torus $\mathbb{R}^2/(2\pi\mathbb{Z})^2)$ to the embedded torus, this is not conformal.  I was suggesting replacing it by
\begin{align}
& R > r > 0 \\[6pt]
x & =  (R + r k(v)) \cos u \\[4pt]
y & =  (R + r \ell(v)) \sin u \\[4pt]
z & =  r \sin v
\end{align}
where $v\mapsto(k(v),\ell(v))$ is some other parametrization of the unit circle.
Let's call $u$ the longitude and $v$ the latitude.  Then the meridians of longitude are closer together on the inside of the embedded torus than on the outside, but the parallels of latitude are still just as far apart for equal values of $\Delta v$.  That tells us the mapping is not conformal.  So the idea would be to make $(x,y,z)$ move only as fast with respect to $v$ as it moves with respect to $u$.  That rate would be proportional to the distance from the center of the embedded torus.
So look at the circle
\begin{align}
x & = R + r \cos v \\
z & = r \sin v
\end{align}
i.e. $(x-R)^2 + z^2 = r^2$.  Draw any line through $(x,z)=(0,0)$ that intersects that circle twice: once on the "inside", i.e. closer to $(x,z)=(0,0)$ and thus closer to the center of the torus, and once on the "outside".  As that line moves, so that both of those two points of intersection move, what is the ratio of their two rates of motion along the arc?  Here's an easy result I proved: That ratio of rates of motion is the same as the ratio of their $x$-coordinates!  I rashly conjectured this based only on the fact that it's obviously true in the case where the aforementioned line through the origin happens to be the $x$-axis in the $xz$-plane.  There didn't seem to be any obvious geometric reason why it should happen, and yet one can verify that it does.
What we would like, then, is that points $v$ and corresponding points $\pi-v$ would be mapped by
$$
v\mapsto
\begin{cases}
x = R + r k(v) \\
z = r \ell(v)
\end{cases}
$$
to two points on the circle such that if you draw the line through them it goes straight through the origin in the $xz$-plane.
A corollary is that the point $v=\pi/2$ would get mapped to the point on the circle at which the tangent line goes through the origin.
Although this tells us exactly the mapping $(k(v),\ell(v))\mapsto(k(\pi-v),\ell(\pi-v))$, it stops short of identifying $v\mapsto(k(v),\ell(v))$, except at four special points ($(\pi/2)\mathbb{Z}\bmod 2\pi\mathbb{Z}$).  Hence this is not a complete answer.
A: Very belated, but in the hope it's useful to posterity:
WimC's answer (especially as expanded in the comments) can be located geometrically rather than by ODEs. For any two positive real numbers $a$ and $b$ with $a^{2} + b^{2} = 1$, the flat (Clifford) torus
$$
(u, v) \mapsto (a\cos(u/a), a\sin(u/a), b\cos(v/b), b\sin(v/b))
$$
in the three-sphere is clearly conformal. Since stereographic projection from the round three-sphere (with the point $(0, 0, 0, 1)$ removed) to Euclidean three-space is also conformal, the composition
$$
(u, v) \mapsto \frac{(a\cos(u/a), a\sin(u/a), b\cos(v/b))}{1 - b\sin(v/b)}
$$
conformally parametrizes a circular torus.

Part of the fun of this viewpoint is, it's easy to rotate the Clifford torus before projecting. Rotating the real plane spanned by the second and fourth coordinates, for example, yields the torus
$$
(u, v) \mapsto (a\cos(u/a), a\sin(u/a)\cos t - b\sin(v/b)\sin t, b\cos(v/b), b\sin(v/b)\cos t + a\sin(u/a)\cos t),
$$
and conformal parametrizations of circular and parabolic cyclides given by
$$
(u, v) \mapsto \frac{(a\cos(u/a), a\sin(u/a)\cos t - b\sin(v/b)\sin t, b\cos(v/b))}{1 - b\sin(v/b)\cos t + a\sin(u/a)\cos t}.
$$
A: Another direction my thinking about this took, is one that WimC says is in some sense equivalent to his own answer.  He approaches the problem by asking what differential equations should be satisfied by the mapping, and that seems sensible since the property of conformality is a differential equation---it's about how fast things are changing at each point.
The image point $(x,y,z)$ should move just as fast as $v$ changes, as it does as $u$ changes.  But the rate of change of $(x,y,z)$ with respect to $u$ is proportional to the distance from $(x,y,z)$ to the central axis $x=y=0$.  So look at just the $xz$-plane and look at the point moving around the circle of radius $r$ as $v$ changes.  It should change fast at $x=R+r$ and slowly at $x=R-r$.  Instead of $v$, think of the angle $\theta$, the "latitude".  Latitudes $0$ and $\pi$ should map to latitudes $0$ and $\pi$ while points in between move closer to the point where $x=R-r$.  So we have a mapping that fixes two points on the circle while leaving the circle invariant.  But we've all seen functions like that before, namely
$$
f(z) = \frac{az+b}{bz+a}.
$$
This leaves $\pm1$ pointwise fixed and leaves the circle $|z|=1$ invariant.  I looked at the case $a=3$ and $b=1$.  It maps $i$ to $(-3/5)+(4/5)i$.  The tangent line to the circle at that point meets the real axis at $-5/3$.
And lo and behold: the rates of change of $|f(z)|$ are proportional to the differences between the real parts of $f(z)$ and $-5/3$.
A: I was looking at this question because of the tokamak fields and the use of Poisson's equations to solve for the electric and magnetic fields inside a region. By mapping the boundaries of the torus to a rectangular grid, the classical second derivative of the potential of the charge can be found analytically and the boundary conditions become simpler. Hypothetically, one shape can be mapped to another. This is seen in the mapping of the stereographic projection, the semicircle to the rectangular plane, and it is proposed now that the torus can be mapped to a rectangular plane. If this is so, imagine a practice problem (E and M, Jackson 2.14) and subsequently revolving the two dimensional cylindrical shape to create a surface of revolution, a torus with four charged plates. Each plate has a voltage and then this creates a potential landscape. The curvature of the plate quadrupoles shields stray fields, but I imagine the stable or unstable equilibria within the Poisson equations solution (using the Lorentz force EOM) might indicate the preference of concave curved plates over convex plates. This problem seems to be confusing. In essence, what is the best way to solve for the Poisson equations fields inside a torus with hemispherical bosses located on the inner surface. These bosses create ELMs in real world tokamaks. I guess the questions I have are about super positioning Poisson's equations solutions to boss problems, cylinder-toroid quadrupole problems and similar ones to generate the actual field shape inside the tokamak. 
To summarize: I create a torus, I map a torus to a grid (x,y,z), I place boundary conditions for Poisson's equation on my grid, I solve for the field in the grid, I map to the toroid. Then, I perturb the boundaries with hemispheres and superimpose the new fields onto the grid and subsequently onto the torus. This generates my field. Then I do the perturbation method to generate neutron magnetic moment contributions, and the poloidal field.
I worked on one of the Jackson problems in grad school and was going to do the toroid in the next few weeks. I have some code in matlab for what I am going to turn into the Poisson equation but I have yet to write the PDE ...
 a=5; %radius of circle
 c=15; % radius of torus
 [u,v]=meshgrid(0:10:360);% this generates the entire torus, but I need to generate a quadrupole like section
 x=(c+a*cosd(v)).*cosd(u); 
 y=(c+a*cosd(v)).*sind(u);
 z=a*sind(v);
 % This where I plug in the second derivative

My brain is full.
Steve 
