Skip navigation

Tag Archives: math

Simple example

f := (x) -> x^2

Is the function f(x)=x^2. You call it by writing


To animate the function, just add another parameter and give it a range in a plot command.

plot( a*f(x), x=-2..2, a=0.1..5 )


There’s this point in Woo’s paper that confused me before, and it just confused me again now.

The paper says to always pick the largest t.

So the situations in red and pink confused me. It would seem at first glance that “further t” (in the pink diagram, the one intersecting min.y) would be the “greater t” but the _wrong one_ to pick.

But it turns out when you cull the back faces _based on ray direction_, the smaller t would be at the max.y plane, even though the ray starts AFTER the max.y plane.

Note Woo’s method as presented _doesn’t work for when the ray starts totally inside the box_.

If the ray starts inside the box, if you only need the boolean result then RETURN TRUE immediately. But if you need the intersection point, you need to turn around some of the rules of the test:

If you’re using numerical indexing to access the xyz components of a vector, (with 0=x, 1=y, 2=z) a good tip is to remember that

axis = 0 ; // for example.
otherAxis1 = ( axis + 1 ) % 3 ; // =1
otherAxis2 = ( axis + 2 ) % 3 ; // =2

This works for all values of axis to get you “the other 2” axes

You could compute them using factorials, but this is far more efficient.

// Computes n CHOOSE k
// ( n )
// ( k )
int binomial( int n, int k )
  // must accumulate nums/dens separately
  // to avoid roundoff error
  int nums = 1;
  int dens = 1;
  for( int i = 1 ; i <= k ; i++ )
    nums *= ( n - k + i ) ;
    dens *= i ;
  return nums/dens ;

int p=7;
for( int i = 0  ; i <= p ; i++ )
  printf( "%dx^%d +", binomial(p, i), p-i ) ;

Arbitrary axis angle rotation C++ code, ROW MAJOR

/// Arbitrary axis angle rotation
Matrix Matrix::Rotation( const Vector & u, real radians )
  real c = cos( radians ) ;
  real l_c = 1 - c ;

  real s = sin( radians ) ;
  return Matrix(
    u.x*u.x + (1 - u.x*u.x)*c,      u.x*u.y*l_c + u.z*s,        u.x*u.z*l_c - u.y*s,  0,
          u.x*u.y*l_c - u.z*s,  u.y*u.y+(1 - u.y*u.y)*c,        u.y*u.z*l_c + u.x*s,  0,
          u.x*u.z*l_c + u.y*s,      u.y*u.z*l_c - u.x*s,  u.z*u.z + (1 - u.z*u.z)*c,  0,
                            0,                        0,                          0,  1
  ) ;

Arbitrary axis angle rotation C++ code, __COLUMN__ MAJOR

/// Arbitrary axis angle rotation
Matrix Matrix::Rotation( const Vector & u, real radians )
  real c = cos( radians ) ;
  real l_c = 1 - c ;

  real s = sin( radians ) ;
  return Matrix(
    u.x*u.x + (1 - u.x*u.x)*c,      u.x*u.y*l_c - u.z*s,        u.x*u.z*l_c + u.y*s,  0,
          u.x*u.y*l_c + u.z*s,  u.y*u.y+(1 - u.y*u.y)*c,        u.y*u.z*l_c - u.x*s,  0,
          u.x*u.z*l_c - u.y*s,      u.y*u.z*l_c + u.x*s,  u.z*u.z + (1 - u.z*u.z)*c,  0,
                            0,                        0,                          0,  1
  ) ;

So, there’s always this question of “what’s the difference between rotate THEN translate VS translate THEN rotate”?



If you do ROTATE THEN TRANSLATE, the effect is __spinning in place__, because ROTATE about axis first (we are still @ origin), then translate OUT to some point in space.


If you do TRANSLATE THEN ROTATE, the effect is __ORBITING__ around an axis, because it translates you OUT FIRST, away from the origin, then you rotate that new further point about the origin:


We all know d/dx sin( x ) = cos( x ).

Also noted that constitutes a +90 degree phase shift in sin( x )

Page 283 of Ashdown’s Radiosity text and page 9 of ppsloan’s cube projection both use this.

The reason r is sqrt( 1 + u^2 + v^2 ) is:

This is easy but usually tricky to get the hang of at first


If you have an incident vector hitting a surface (that has a normal defined for it!), then you can calc the angle of the reflected vector like so:



reflected = incident – 2 × (normal ⋅ incident) normal



All you’ve gots to do is look at the angle the incident makes with the normal. That’s what the formula I gave says anyway.

But also consider this. If the incident is perpendicular to the normal, then the reflected angle has pretty much no change from the incident angle. IE AT the case where the incident vector has an angle of 90 degrees with the normal of the surface, the “reflected” vector (if you can call it that!) is in the same direction as the “incident” vector. That is because the dot product (normal ⋅ incident) goes to 0 when the angle between incident and normal is 90 degrees, so:

reflected = incident – 2 × (normal ⋅ incident) normal
reflected = incident – 2 × (0) normal
reflected = incident

So that’s interesting. Now consider the other case where the angle between the normal and incident is 180 degrees (ie the ray hits the surface head-on):


reflected = incident – 2 × (normal ⋅ incident) normal
reflected = incident – 2 × (-1) normal
reflected = incident + 2 × normal

IE reflected turns completely around.

For everything else in between, there’s cosine.

So there is LESS change in the direction of the reflected angle as the incident-normal angle grows. This is encoded by the dot product/cosine in the normal reflection.

Ok, here’s the math for it.

Using this diagram:


The angle THETA there is going to be the cosine of the NORMAL and the INCIDENT

The SMALLER ANGLE GAMMA is going to be (THETA – 90) because the normal is at an angle of 90 degrees with the surface. That’s why its a normal, duh!

Now then some simple vector math will give you:

The length of the pink side X (perp distance to surface from INCIDENT vector:)

X = sin( GAMMA ) (remember incident is a unit vector)

So, now


This makes sense because WHEN (sin(GAMMA)) is small, (angle between INCIDENT and NORMAL == nearly 90 degrees) we don’t need to add a large “displacement” in the direction of the normal to get the reflected vector (we add 0 x normal in the extreme case of a 90 degree INCIDENT angle with NORMAL). As the angle between NORMAL and INCIDENT INcreases, so does the displacement we need to use to the reflected vector, until we get to the other extreme case where the REFLECTED = INCIDENT + 2 * NORMAL.

Ok then, simplify.

REFLECTED = INCIDENT + 2 (sin (THETA)cos(90) – sin(90)cos(THETA)) (NORMAL)

cos(90)=0, sin(90)=1




And there you have it.

Well, this post was looking for the same thing I was, and nobody seemed to know of any.

The problem with algebra tricks is they’re so context dependent. There’s really only a few things that math books and tests use to make “wicked” problems (list being expanded).

  • 1. Multiplying by the conjugate
  • 2. Factoring using the +C -C trick (complete the square method)

Really the “hard” problems require an algebraic twist, a keen eye, insight – stuff you only get after having practiced with hundreds of problems (and checking the solutions!).

So first,
sosmath: Derivative of sin, cos, tan, csc, sec, cot

Mathworld: Inverse trig functions

trig identities

Now, the rest of this post is just examples of the trick (by name) and a couple of examples that illustrate the trick in action (list being expanded!) (Bear with the non-mathmlness for now!)

1. Multiplying by the conjugate

INTEGRAL( dx / (1 + cos x) )

To solve, you have to multiply the numerator and denominator by ( 1 – cos x ) to get sin2x in the denominator..

2. Factoring using the +C -C trick

INTEGRAL( dx / (x2 + 10x + 30) )

You have to factor (wactor! as my friend used to say) the denominator into ((x+5)2 + 5) by:

x2 + 10x + 25 - 25 + 30
= (x + 5)2 + 5

You know to do +25, -25 BECAUSE half of 10 is 5, and 52 is 25. That method in general is called completing the square and you have to do it all the time.

Sufficient condition

An event A is said to be a sufficient condition for an event B iff

The truth/occurrence of event A guarantees the truth/occurrence of event B

Necessary condition

An event A is said to be necessary for an event B
iff the falsity/nonoccurrence of event A GUARANTEES
the falsity/nonoccurrence of event B.

Rewriting something so that it is written differently, but ultimately “means” the same thing. E.g. in Math, you can refactor x + xy to x( 1 + y ). These two statements (x+xy) and (x(1+y)) mean exactly the same thing mathematically, but they’re just written differently. e.g. 2 in programming, refactor is used more loosely. You can refactor a bunch of code which means you rewrite the code to do the same intended thing, but in a better way (more efficiently, or with more adherence to OOP principles, or rewrite it to remove bugs).

It would seem that

Ax + By + C = 0

doesn’t really have a direction. Its just an infinite line, right?



-x + y = 0

This is the line through the origin, forming the line y = x.

1. Does this line have a direction?

2. Is there another line with a different equation that DESCRIBES THE SAME LINE?

Yes there is, to the second question, and YES IT DOES, to the first question.

The answer to 2.: The equivalent line is x – y = 0. That line is EXACTLY the same as -x + y = 0, only it “points” in the opposite direction.

So which “way” is -x + y pointing, and which way is x – y = 0 pointing?

Well, I’ll tell you.

x – y = 0 points “down and to the left”

-x + y = 0 points “up and to the right”

the reason for this is clear when you use vectors and the equations

A = y1 – y2
B = x2 – x1
C = x1y2 – x2y1

Then, if you use the points

P1 = (1,1)
P2 = (0,0)

You get the equation

A = 1 – 0 = 1
B = 0 – 1 = -1
C = 0


x – y = 0

And so, the opposite equation
-x + y = 0

On the “left side” of a line is +
On the “right side” of a line is –

Taking the one that goes up and right:
plug in ( 2, 0 )
-x + y = 0
-2 + 0 = -2

So, (2,0) is negative in that line equation, meaning (2,0) is “right” of -x + y = 0

Taking (-2,0)
2 + 0 = 2
so, (-2,0) is “left” of -x + y = 0

The dot product


Another way to find the dot product is:


Projections of a vector u onto a vector v are found in two ways I can see so far:


Which is 100% equal to:


WHERE is a unit vector pointing in the same direction as v.

SO we can also say:


Where v has become a unit vector, totally.

Keep in mind the difference between equations 4 and 5 is in 4, v isn’t necessarily normalized when its used in the dot product operation, so you need to divided by |v| to “cancel its effect” on the final result. notice in all of these, we divide by v enough times so that the actual magnitude of v PLAYS NO ROLE in projvu, only the direction of v matters.

Finally you’ll see the following definition usually in physics books:


Where t is the angle between vectors u and v.

This 6th formula follows if you study formula 3 and move things around a bit..

taking formula 3:

Notice that from formula 1, we have:

So, re-arranging formula 3 a bit now..

All I did there was add a |u| to the denominator, and separate the |v|2 into |v||v|


Next, sub in

Then we get:

Now notice that , (v as a unit vector), so

where is a unit vector pointing in the direction of the vector v you are trying to project u onto.

On Rotation Matrices

A common problem you come across is GIVEN a vector describing a “look” direction, HOW CAN WE MAKE A ROTATION MATRIX THAT WILL “LINE UP” WITH IT? You often want to do something like this when you have a character model that you need to rotate to line up with the direction that the character is currently facing. For example, you are writing a 3d space simulator game and you need the ship’s model to be lined up precisely with the direction the ship is facing.

It turns out, this is easy to do. So easy to do in fact, you’ll kick yourself for not having known this. But it isn’t easy to understand how it works.

So, first, let’s understand the principles of what we’re doing here in 2d, then we’ll extend our answer to 3d.

Rotating something to line up with an object in 2d

Assume you work in a world of vectors. In your 2d game, you need to have the model of your character point in the direction that he is looking.

Have: (red arrow is direction the Triangle character is looking).

But the black triangle and black arrow indicate the stock direction the model faces.


The most obvious way to solve this problem is just a 2d rotation matrix.

Matrix that rotates a point t degrees counter-clockwise:

The Triangle starts with its Heading vector pointing (1, 0). What we need to do is rotate the Heading vector 60 degrees counter-clockwise.

To make a 60 degree CCW rotation, let’s use the CCW rotation matrix, plugging in 60 degrees for the rotation we want:

So we end up with this set of numbers that promise to rotate a set of points 60 degrees counterclockwise:

Now here’s another way to get that same matrix. Are you ready?

The red vector (which points in the direction we want to face) is called the forward vector. The blue vector, we’ll simply call the side vector. ASSUME BOTH THE forward AND side vectors are UNIT VECTORS (THIS IS VERY IMPORTANT!!)

Since we know the forward vector is a unit vector and it has a 60 degree angle with the x-axis, we can express it as:
xfwd = cos( 60 )
= 0.5

yfwd = sin( 60 )
= 0.86

So a unit vector in the plane, when broken down into its x, y components, is always going to be x=cos(t), y=sin(t) where t is the angle that the vector forms with the +x axis.

The blue vector makes (90+60=150) degree angle with the +x axis. So it can be expressed as:
xperpleft = cos( 150 )
= -0.86

yperpleft = sin( 150 )
= 0.5

So, speaking in terms of the orientation we want, we have two vectors:

forward (describes direction we want to look in)=( 0.5, 0.86 )
perpleft (90 degrees counter-clockwise perpendicular to the forward vector) = ( -0.86, 0.5 )

We could also derive the perpleft vector by noticing that ( -y, x ) is always going to be a 90 degree counter-clockwise rotated perpendicular to any vector you have.

side note:
You are sure that (-y,x) is a perpendicular to (x,y) because the dot product:
(x,y)*(-y,x) = -yx + xy = 0

and so, (-y,x) is always perpendicular to (x,y).

The other way to derive this is to construct the 90 degree CCW rotation matrix and notice that it always comes out to

[ cos(90)   -sin(90) ] [ x ]
[ sin(90)    cos(90) ] [ y ]

[ 0   -1 ] [ x ] = [ -y ]
[ 1    0 ] [ y ]   [  x ]

So, look at those again:

forward: ( 0.5, 0.86 )
perpleft: ( -0.86, 0.5 )

Does that look familiar?

So we can do this:

Assuming you use column vectors and not row vectors. If you’re using row vectors, you have to transpose it.

Isn’t that neat?

This extends to 3D:

At this point, I’ll direct you to Diana Gruber’s article and leave it at that.

Say you come across

What you do is
break it up

This post really helped. In this post, this guy david says:

Hi Scott,

Here is one possibility:  If the old axis is _a_ and the new axis
is _a'_, use

theta = acos(_a_ dot _a'_)    (where 'dot' means dot product)
_axis_ = _a_ cross _a'_       (cross product)



will provide the desired rotation.  Don't forget to do the 
radians-to-degrees conversion before calling RotateWXYZ.

A much more efficient way to it if the original axis _a_ is
always the x axis is to build the matrix yourself.  First,
you need to find two vectors _b'_ and _c'_ which are
perpendicular to _a'_.  You can use vtkMath::Perpendiculars()
for this if you have a recent version of VTK.

Then, build the matrix by inserting the vectors into the matrix 

   a'[0]  b'[0]  c'[0]  t[0]  
   a'[1]  b'[1]  c'[1]  t[1]
   a'[2]  b'[2]  c'[2]  t[2]
   0      0      0      1

where _t_ is the translation you wish to apply after the rotation.

 - David

I wanna personally thank the dude.

But that wasn’t enough. There was also this post posted by none other than shawn hargreaves. He says:

When you’re working in 3D it’s generally best to avoid thinking in terms of rotations as angles. Trying to extract angle rotations lands you in a world of pain working out complex trig equations, which are not at all my idea of fun! This sort of thing can generally be done much more easily using vector and matrix math.

Ultimately, you want a matrix that describes the rotation of your object, but you can get there from a couple of vectors. One vector alone is not enough to describe a full orientation, though, because that still leaves one axis of freedom that the object could rotate to any angle around that vector.

Let’s start out by choosing a front vector for our orientation: you can get this by normalizing your movement vector.

Now we need vectors pointing right and up. To avoid squashing the object, these need to be perpendicular (at right angles to) to the front vector. We could guess that up might be roughly Vector3.Up, but that might need some adjustment to make it stay perpendicular to our front vector. Here’s how to do that:

– Compute a right vector by taking the cross product of your front vector and Vector3.Up, and normalizing that result

– Compute a properly perpendicular up vector by taking the cross product of the right vector that you just calculated and your front vector.

Now you have three vectors, all perpendicular to each other, that describe the orientation of your object. Turning these into a matrix is simple: start off with Matrix.Identity, then assign your new values to its Forward, Right, and Up properties.

I might have got some of those cross products the wrong way round (I’m terrible at remembering which way they go) so if your model comes out inverted, just swap the order of the arguments to the cross product calls.

The framework actually has a built in method that does almost the same thing as the above, in the Matrix.CreateLookAt method. You could pass Vector3.Zero as cameraPosition, your movement vector as cameraTarget, and Vector3.Up as cameraUpVector, and it will give back almost the same thing as the matrix I described above, only backward. Call Matrix.Transpose on the result, and you can use that for your object orientation.

very helpful indeed. thanks shawn!

Finally for understanding this crap more, the math of viewing matrices, including an explanation of what/how do it.

The XNA docs say:

d: The distance of the Plane along its normal from the origin

However, that is wrong.

The right description is:

d is the distance you must translate the plane along its normal to get to the origin.

For example, the plane

rightBoundaryPlane = new Plane( 1, 0, 0, 1000 );

Is a plane that has a normal pointing in +x. Using the xna vague (and inaccurate) description of d, you might think that d is “the distance of the plane along its normal from the origin..” i.e. this plane should be off 1000 units in the direction of +x.

In fact, that is WRONG. This plane intersects the x axis at ( -1000, 0, 0 ). Why -1000? Because d represents the number of units you must move this plane, along its normal, to get to the origin.

Meaning, this plane must intersect ( -1000, 0, 0 ) because moving the plane 1000 units in the direction of ( 1, 0, 0 ) would take you to the origin.

So how do you get a plane at the right (+X) side?


rightBoundaryPlane = new Plane( 1, 0, 0, -1000 );

This plane intersects ( 1000, 0, 0 ) because you’d have to move it -1000 units ALONG ITS NORMAL ( in the direction of (1, 0, 0) ) to get it to the origin