Ilya and Jan, because of floating point arithmetic, it's impossible to find intersections perfectly in all cases. Jan, your example has two curves that touch at a single point; some people might define this as an intersection. Your assertion that (0,0) isn't an intersection is debatable.

It also erroneously finds contact points (no real intersections).
Example: [x0, y0] = intersections([-1,0,-1], [-1,0,1], [1, 0, 1], [-1,0,1], 1);
returns point (0, 0) as intersection point although it isn't.

An excellent function! However, if both intersecting curves already include the intersection point, weird results are possible (see the provided example). This problem was already touched by John Mahoney (the results were different).
My test case: several lines pass through one point and every line must eventually include this point. However, it's not known in advance that we have this situation, it simply may appear.
x1=[-0.49313932739246, -0.02127781500161, 0.450583697389237];
y1=[-10.01, 0, 10.01];
x2=[-1.22073877122679, -0.02127781500161, 1.17818314122357];
y2=[10.01, 0, -10.01];
[x,y,seg1,seg2] = intersections(x1, y1, x2, y2)
The case is published because it may be useful to know about. Otherwise I agree that the function works fine in 99,99999999% of cases.

Comment only