I am sorry to disagree with the answer of Peter Walser. Solving the equations gives on my desk:
u = ((bs.y - as.y) * bd.x - (bs.x - as.x) * bd.y) / (bd.x * ad.y - bd.y * ad.x)
v = ((bs.y - as.y) * ad.x - (bs.x - as.x) * ad.y) / (bd.x * ad.y - bd.y * ad.x)
Factoring out the common terms, this comes to:
dx = bs.x - as.x
dy = bs.y - as.y
det = bd.x * ad.y - bd.y * ad.x
u = (dy * bd.x - dx * bd.y) / det
v = (dy * ad.x - dx * ad.y) / det
Five subtractions, six multiplications and two divisions.
If you only need to know if the rays intersect, the signs of u and v are enough, and these two divisons can be replaced by num*denom<0 or (sign(num) != sign(denom)), depending on what is more efficient on your target machine.
Please note that the rare case of det==0 means that the rays do not intersect (one additional comparison).