How to draw normal vectors to an ellipse

How do I draw an ellipse with lines of the same length coming out of it?

It's easy to do with a circle, I can just write something like

for (u = 0 ; u < 2*pi ; u += 0.001*pi) {
    drawdot (cos(u), sin(u)) ;
    drawline (cos(u), sin(u), 2*cos(u), 2*sin(u) ;
}

But if I did that for an ellipse, like below, the lines are different lengths.

for (u = 0 ; u < 2*pi ; u += 0.001*pi) {
    drawdot (2*cos(u), sin(u)) ;
    drawline (2*cos(u), sin(u), 4*cos(u), 2*sin(u) ;
}

How do I figure out how to make them the same length?

Solution

There are a few ways of thinking about this.

You can think of an ellipse as a circle that's been stretched in some direction. In this case, you've taken the circle x^2 + y^2 = 1 and applied the transformation to all points on that curve:

x' = 2x
y' = y

You can think of this as multiplying by the matrix:

[ 2  0 ]
[ 0  1 ]

To transform normals, you need to apply the inverse transpose of this matrix (i.e. the inverse of the transpose, or transpose of the inverse; it's the same thing):

[ 1/2  0 ]
[  0   1 ]

(This, by the way, is known as the dual of the previous transformation. This is a very important operation in modern geometry.)

A normal to the circle at the point (x,y) goes in the direction (x,y). So a normal to the ellipse at the point (2x,y) goes in the direction (0.5*x,y). This suggests:

for (u = 0 ; u < 2*pi ; u += 0.001*pi) {
    x = cos(u); y = sin(u);
    drawdot (2*x, y) ;
    drawline (2*x, y, 2*x + 0.5*x, y+y);
}

Or if you need a unit normal:

for (u = 0 ; u < 2*pi ; u += 0.001*pi) {
    x = cos(u); y = sin(u);
    drawdot (2*x, y) ;
    dx = 0.5*x;
    dy = y;
    invm = 1 / sqrt(dx*dx + dy*dy);
    drawline (2*x, y, 2*x + dx * invm, y + dy * invm);
}

Another way to think about it is in terms of an implicit contour. If you define the curve by a function:

f(x,y) = 0

then the normal vector points in the direction:

(df/dx, df/dy)

where the derivatives are partial derivatives. In your case:

f(x,y) = (x/2)^2 + y^2 = 0

df/dx = x/2
df/dy = y

which, you will note, is the same as the dual transformation.